1 ! TO DO: 2 ! Check the formulation based on MCp 3 ! 4 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC 5 ! C 6 ! Subroutine: CALC_e_e C 7 ! Purpose: Calculate coefficients linking velocity correction to C 8 ! solids volume fraction correction -- East C 9 ! C 10 ! Author: M. Syamlal Date: 1-JUL-96 C 11 ! Reviewer: Date: C 12 ! C 13 ! Literature/Document References: C 14 ! Variables referenced: C 15 ! Variables modified: C 16 ! Local variables: C 17 ! C 18 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C 19 20 SUBROUTINE CALC_E_E(A_M, MCP, E_E) 21 22 !----------------------------------------------- 23 ! Modules 24 !----------------------------------------------- 25 USE param 26 USE param1 27 USE parallel 28 USE fldvar 29 USE geometry 30 USE indices 31 USE physprop 32 USE run 33 USE constant 34 USE compar 35 USE sendrecv 36 USE fun_avg 37 USE functions 38 IMPLICIT NONE 39 !----------------------------------------------- 40 ! Dummy arguments 41 !----------------------------------------------- 42 ! Septadiagonal matrix A_m 43 DOUBLE PRECISION, INTENT(IN) :: A_m(DIMENSION_3, -3:3, 0:DIMENSION_M) 44 ! Index of close packed solids phase. note mcp is the lowest index of 45 ! those solids phases that have close_packed=t and of the solids 46 ! phase that is used for the solids correction equation (when mmax=1). 47 INTEGER, INTENT(IN) :: Mcp 48 ! Coefficient for solids correction equation. These coefficients are 49 ! initialized to zero in the subroutine time_march before the time loop. 50 ! Note that the E_e coefficients are only used when mcp is assigned 51 ! (when any solids phase has close_packed=T). 52 DOUBLE PRECISION, INTENT(INOUT) :: e_e(DIMENSION_3) 53 !----------------------------------------------- 54 ! Local variables 55 !----------------------------------------------- 56 ! Indices 57 INTEGER :: IJK 58 !----------------------------------------------- 59 60 ! returning if mcp is not defined. 61 IF (MCP == UNDEFINED_I) RETURN 62 63 ! returning if the x-momentum equation of this phase is not to be 64 ! solved. 65 IF (.NOT.MOMENTUM_X_EQ(MCP)) RETURN 66 67 !!$omp parallel do private(ijk) 68 DO IJK = ijkstart3, ijkend3 69 IF (SIP_AT_E(IJK) .OR. MFLOW_AT_E(IJK)) THEN 70 E_E(IJK) = ZERO 71 ELSE 72 IF (A_M(IJK,0,MCP) /= ZERO) THEN 73 ! calculating the correction coefficient 74 E_E(IJK) = AYZ(IJK)/(-A_M(IJK,0,MCP)) 75 ELSE 76 E_E(IJK) = LARGE_NUMBER 77 ENDIF 78 ENDIF 79 ENDDO 80 81 RETURN 82 END SUBROUTINE CALC_E_E 83 84 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC 85 ! C 86 ! Subroutine: CALC_e_n C 87 ! Purpose: Calculate coefficients linking velocity correction to C 88 ! solids volume fraction correction -- North C 89 ! C 90 ! Author: M. Syamlal Date: 1-JUL-96 C 91 ! Reviewer: Date: C 92 ! C 93 ! Literature/Document References: C 94 ! Variables referenced: C 95 ! Variables modified: C 96 ! Local variables: C 97 ! C 98 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C 99 100 SUBROUTINE CALC_E_N(A_M, MCP, E_N) 101 102 !----------------------------------------------- 103 ! Modules 104 !----------------------------------------------- 105 USE param 106 USE param1 107 USE parallel 108 USE fldvar 109 USE geometry 110 USE indices 111 USE physprop 112 USE run 113 USE constant 114 USE compar 115 USE sendrecv 116 USE fun_avg 117 USE functions 118 IMPLICIT NONE 119 !----------------------------------------------- 120 ! Dummy arguments 121 !----------------------------------------------- 122 ! Septadiagonal matrix A_m 123 DOUBLE PRECISION, INTENT(IN) :: A_m(DIMENSION_3, -3:3, 0:DIMENSION_M) 124 ! Index of close packed solids phase. note mcp is the lowest index of 125 ! those solids phases that have close_packed=t and of the solids 126 ! phase that is used for the solids correction equation (when mmax=1). 127 INTEGER, INTENT(IN) :: Mcp 128 ! Coefficient for solids correction equation. These coefficients are 129 ! initialized to zero in the subroutine time_march before the time loop. 130 ! Note that the E_n coefficients are only used when mcp is assigned 131 ! (when the solids phase has close_packed=T) 132 DOUBLE PRECISION, INTENT(INOUT) :: e_n(DIMENSION_3) 133 !----------------------------------------------- 134 ! Local variables 135 !----------------------------------------------- 136 ! Indices 137 INTEGER :: IJK 138 !----------------------------------------------- 139 140 ! returning if mcp is not defined. 141 IF (MCP == UNDEFINED_I) RETURN 142 ! returning if the y-momentum equation of this phase is not to be 143 ! solved. 144 IF (.NOT.MOMENTUM_Y_EQ(MCP)) RETURN 145 146 !!$omp parallel do private(IJK) 147 DO IJK = ijkstart3, ijkend3 148 IF (SIP_AT_N(IJK) .OR. MFLOW_AT_N(IJK)) THEN 149 E_N(IJK) = ZERO 150 ELSE 151 IF ((-A_M(IJK,0,MCP)) /= ZERO) THEN 152 ! calculating the correction coefficient 153 E_N(IJK) = AXZ(IJK)/(-A_M(IJK,0,MCP)) 154 ELSE 155 E_N(IJK) = LARGE_NUMBER 156 ENDIF 157 ENDIF 158 ENDDO 159 160 RETURN 161 END SUBROUTINE CALC_E_N 162 163 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC 164 ! C 165 ! Subroutine: CALC_e_t C 166 ! Purpose: Calculate coefficients linking velocity correction to C 167 ! solids volume fraction correction -- Top C 168 ! C 169 ! Author: M. Syamlal Date: 1-JUL-96 C 170 ! Reviewer: Date: C 171 ! C 172 ! Literature/Document References: C 173 ! Variables referenced: C 174 ! Variables modified: C 175 ! Local variables: C 176 ! C 177 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C 178 179 SUBROUTINE CALC_E_T(A_M, MCP, E_T) 180 181 !----------------------------------------------- 182 ! Modules 183 !----------------------------------------------- 184 USE param 185 USE param1 186 USE parallel 187 USE fldvar 188 USE geometry 189 USE indices 190 USE physprop 191 USE run 192 USE constant 193 USE compar 194 USE sendrecv 195 USE fun_avg 196 USE functions 197 IMPLICIT NONE 198 !----------------------------------------------- 199 ! Dummy arguments 200 !----------------------------------------------- 201 ! Septadiagonal matrix A_m 202 DOUBLE PRECISION, INTENT(IN) :: A_m(DIMENSION_3, -3:3, 0:DIMENSION_M) 203 ! Index of close packed solids phase. note mcp is the lowest index of 204 ! those solids phases that have close_packed=t and of the solids 205 ! phase that is used for the solids correction equation (when mmax=1). 206 INTEGER, INTENT(IN) :: Mcp 207 ! Coefficient for solids correction equation. These coefficients are 208 ! initialized to zero in the subroutine time_march before the time loop. 209 ! Note that the E_t coefficients are only used when mcp is assigned 210 ! (when the solids phase has close_packed=T) 211 DOUBLE PRECISION, INTENT(INOUT) :: e_t(DIMENSION_3) 212 !----------------------------------------------- 213 ! Local variables 214 !----------------------------------------------- 215 ! Indices 216 INTEGER :: IJK 217 !----------------------------------------------- 218 219 ! returning if mcp is not defined. 220 IF (MCP == UNDEFINED_I) RETURN 221 ! returning if the z-momentum equation of this phase is not to be 222 ! solved 223 IF (.NOT.MOMENTUM_Z_EQ(MCP)) RETURN 224 225 !!$omp parallel do private(IJK) 226 DO IJK = ijkstart3, ijkend3 227 IF (SIP_AT_T(IJK) .OR. MFLOW_AT_T(IJK)) THEN 228 E_T(IJK) = ZERO 229 ELSE 230 IF ((-A_M(IJK,0,MCP)) /= ZERO) THEN 231 ! calculating the correction coefficient 232 E_T(IJK) = AXY(IJK)/(-A_M(IJK,0,MCP)) 233 ELSE 234 E_T(IJK) = LARGE_NUMBER 235 ENDIF 236 ENDIF 237 ENDDO 238 239 RETURN 240 END SUBROUTINE CALC_E_T 241 242