1 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv! 2 ! ! 3 ! Module name: DES_CONVECTION ! 4 ! ! 5 ! Purpose: This routine is called from the discrete phase. It is used ! 6 ! to determine the rate of heat transfer to the particle from the gas.! 7 ! ! 8 ! Author: J.Musser Date: 16-Jun-10 ! 9 ! ! 10 ! Comments: ! 11 ! ! 12 ! REF: Zhou, Yu, and Zulli, "Particle scale study of heat transfer in ! 13 ! packed and bubbling fluidized beds," AIChE Journal, Vol. 55, ! 14 ! no 4, pp 868-884, 2009. ! 15 ! ! 16 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^! 17 SUBROUTINE DES_CONVECTION(NP, M, IJK, & 18 INTERP_IJK, INTERP_WEIGHTS, FOCUS) 19 20 use constant, only: Pi 21 Use des_thermo 22 Use discretelement 23 Use fldvar 24 Use interpolation 25 Use param1 26 Use physprop 27 Use usr 28 29 IMPLICIT NONE 30 31 ! Passed variables 32 !---------------------------------------------------------------------// 33 ! Index of particle being looped over 34 INTEGER, INTENT(IN) :: NP 35 ! Solid phase index of particle NP 36 INTEGER, INTENT(IN) :: M 37 ! IJK value of cell containing particle NP 38 INTEGER, INTENT(IN) :: IJK 39 ! IJK indicies of fluid cells involved in interpolation 40 INTEGER, INTENT(IN) :: INTERP_IJK(2**3) 41 ! Weights associated with interpolation 42 DOUBLE PRECISION, INTENT(IN) :: INTERP_WEIGHTS(2**3) 43 ! Indicates that debugging information for the particle 44 LOGICAL, INTENT(IN) :: FOCUS 45 46 ! Local variables 47 !---------------------------------------------------------------------// 48 ! Convective heat transfer coefficient 49 DOUBLE PRECISION GAMMA_CP 50 ! Temperature of the gas 51 DOUBLE PRECISION Tg 52 ! Surface area of particle 53 DOUBLE PRECISION Sa 54 ! Convection source 55 DOUBLE PRECISION Qcv 56 57 ! Obtain the temperature of the gas. --> Not interpolated. 58 Tg = T_g(IJK) 59 ! Obtain the convective heat transfer coefficient of the particle 60 CALL DES_CALC_GAMMA(NP, IJK, GAMMA_CP) 61 ! Calculate the surface area of the particle 62 Sa = 4.0d0 * Pi * DES_RADIUS(NP) * DES_RADIUS(NP) 63 ! Calculate the rate of heat transfer to the particle 64 Qcv = GAMMA_CP * Sa * (Tg - DES_T_s_NEW(NP)) 65 ! Store convection source in global energy source array. 66 Q_Source(NP) = Q_Source(NP) + Qcv 67 68 ! Calculate the source term components --> Not interpolated. 69 DES_ENERGY_SOURCE(IJK) = DES_ENERGY_SOURCE(IJK) - Qcv * DTSOLID 70 71 ! Write out the debugging information. 72 ! IF(FOCUS) THEN 73 ! WRITE(*,"(/5X,A)")'From: DES_CONVECTION -' 74 ! WRITE(*,"(8X,A,D12.6)")'Tg: ',Tg 75 ! WRITE(*,"(8X,A,D12.6)")'Tp: ',DES_T_s_NEW(NP) 76 ! WRITE(*,"(8X,A,D12.6)")'Sa: ',Sa 77 ! WRITE(*,"(8X,A,D12.6)")'GAMMA_CP: ',GAMMA_CP 78 ! WRITE(*,"(8X,A,D12.6)")'Qcv: ',Qcv 79 ! WRITE(*,"(5X,50('-')/)") 80 ! ENDIF 81 82 RETURN 83 84 END SUBROUTINE DES_CONVECTION 85 86 87 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv! 88 ! Subroutine: DES_CALC_GAMMA ! 89 ! ! 90 ! Purpose: Calculate the heat transfer coefficient (GAMMA_CP) for ! 91 ! particle-fluid heat transfer. ! 92 ! ! 93 ! Author: J.Musser Date: 16-Jun-10 ! 94 ! ! 95 ! Comments: ! 96 ! ! 97 ! REF: Ranz, W.E. and Marshall, W.R., "Friction and transfer ! 98 ! coefficients for single particles and packed beds," Chemical ! 99 ! Engineering Science, Vol. 48, No. 5, pp 247-253, 1925. ! 100 ! ! 101 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^! 102 SUBROUTINE DES_CALC_GAMMA(NP, IJK, GAMMA_CP) 103 104 USE compar 105 Use constant 106 Use des_thermo 107 Use discretelement 108 Use fldvar 109 USE geometry 110 USE indices 111 USE param1 112 USE physprop 113 USE fun_avg 114 USE functions 115 116 IMPLICIT NONE 117 118 ! Passed variables 119 !---------------------------------------------------------------------// 120 ! Index value of particle 121 INTEGER, INTENT(IN) :: NP 122 ! Index value of fluid cell 123 INTEGER, INTENT(IN) :: IJK 124 ! Convective heat transfer coefficient 125 DOUBLE PRECISION, INTENT(OUT) :: GAMMA_CP 126 127 ! Local variables 128 !---------------------------------------------------------------------// 129 ! Fluid cell indices 130 INTEGER IMJK, IJMK, IJKM 131 ! Double precision value for 1/3 132 DOUBLE PRECISION, PARAMETER :: THIRD = (1.0d0/3.0d0) 133 134 DOUBLE PRECISION N_Pr ! Prandtl Number 135 DOUBLE PRECISION N_Re ! Reynolds Number 136 DOUBLE PRECISION N_Nu ! Nusselt Number 137 138 ! Magnitude of slip velocity 139 DOUBLE PRECISION SLIP 140 ! Fluid velocity 141 DOUBLE PRECISION cUg, cVg, cWg 142 DOUBLE PRECISION Us, Vs, Ws 143 !---------------------------------------------------------------------// 144 145 ! Initialization 146 IMJK = IM_OF(IJK) 147 IJMK = JM_OF(IJK) 148 IJKM = KM_OF(IJK) 149 150 SELECT CASE(TRIM(DES_CONV_CORR)) 151 152 CASE ('RANZ_1952') ! (Ranz and Mrshall, 1952) 153 ! Initialize variables 154 SLIP = ZERO 155 N_Re = ZERO 156 N_Nu = ZERO 157 ! Gas velocity in fluid cell IJK 158 cUg = AVG_X_E(U_g(IMJK), U_g(IJK), 1) 159 cVg = AVG_Y_N(V_g(IJMK), V_g(IJK)) 160 ! Particle Velocity 161 Us = DES_VEL_NEW(1,NP) 162 Vs = DES_VEL_NEW(2,NP) 163 164 ! Calculate the magnitude of the slip velocity 165 IF(NO_K) THEN 166 SLIP = SQRT((cUg-Us)**2 + (cVg-Vs)**2) 167 ELSE 168 cWg = AVG_Z_T(W_g(IJKM), W_g(IJK)) 169 Ws = DES_VEL_NEW(3,NP) 170 SLIP = SQRT((cUg-Us)**2 + (cVg-Vs)**2 + (cWg-Ws)**2) 171 ENDIF 172 173 ! Calculate the Prandtl Number 174 IF(K_G(IJK) > ZERO) THEN 175 N_Pr = (C_PG(IJK)*MU_G(IJK))/K_G(IJK) 176 ELSE 177 N_Pr = LARGE_NUMBER 178 ENDIF 179 180 ! Calculate the particle Reynolds Number 181 IF(MU_G(IJK) > ZERO) THEN 182 N_Re = (2.0d0*DES_RADIUS(NP)*SLIP*RO_g(IJK)) / MU_g(IJK) 183 ELSE 184 N_Re = LARGE_NUMBER 185 ENDIF 186 187 ! Calculate the Nusselt Number 188 N_Nu = 2.0d0 + 0.6d0 *((N_Re)**HALF * (N_Pr)**THIRD) 189 190 ! Calculate the convective heat transfer coefficient 191 GAMMA_CP = (N_Nu * K_G(IJK))/(2.0d0 * DES_RADIUS(NP)) 192 193 CASE DEFAULT 194 WRITE(*,*)'INVALID DES CONVECTION MODEL' 195 STOP 196 CALL MFIX_EXIT(myPE) 197 END SELECT 198 199 RETURN 200 END SUBROUTINE DES_CALC_GAMMA 201 202 203 204 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv! 205 ! Subroutine: DES_Hgm ! 206 ! ! 207 ! Purpose: This routine is called from the continuum phase and ! 208 ! calculates the source term from the particles to the fluid. ! 209 ! ! 210 ! Author: J.Musser Date: 15-Jan-11 ! 211 ! ! 212 ! Comments: ! 213 ! ! 214 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^! 215 SUBROUTINE DES_Hgm(S_C) 216 217 USE compar 218 Use constant 219 Use des_thermo 220 Use discretelement 221 Use fldvar 222 USE geometry 223 USE indices 224 Use interpolation 225 Use param1 226 Use physprop 227 228 use run, only: ODT 229 use functions 230 231 IMPLICIT NONE 232 233 ! Passed Variables 234 !---------------------------------------------------------------------// 235 ! Source term on RHS 236 DOUBLE PRECISION, INTENT(INOUT) :: S_C(DIMENSION_3) 237 238 ! Local variables 239 !---------------------------------------------------------------------// 240 ! IJK value of cell containing particle NP 241 INTEGER :: IJK 242 !---------------------------------------------------------------------// 243 244 ! Loop over fluid cells. 245 IJK_LP: DO IJK = IJKSTART3, IJKEND3 246 IF(.NOT.FLUID_AT(IJK)) CYCLE IJK_LP 247 ! Redistribute the energy over the fluid time step. Note that by the 248 ! time this routine is called, S_C and S_P have already been multiplied 249 ! by the fluid cell volume. Thus, the mapping should result in units 250 ! of energy per time. 251 S_C(IJK) = S_C(IJK) + DES_ENERGY_SOURCE(IJK)*ODT 252 253 ENDDO IJK_LP ! End loop over fluid cells 254 255 RETURN 256 END SUBROUTINE DES_Hgm 257 258 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv! 259 ! Subroutine: DES_Hgm ! 260 ! ! 261 ! Purpose: ZERO out the array that passes energy source terms back to ! 262 ! the continuum model. Additional entries may be needed to include ! 263 ! heat transfer to the hybrid mode. ! 264 ! ! 265 ! Author: J.Musser Date: 15-Jan-11 ! 266 ! ! 267 ! Comments: ! 268 ! ! 269 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^! 270 SUBROUTINE ZERO_ENERGY_SOURCE 271 272 Use des_thermo 273 Use param1 274 275 IMPLICIT NONE 276 277 DES_ENERGY_SOURCE(:) = ZERO 278 279 280 RETURN 281 END SUBROUTINE ZERO_ENERGY_SOURCE 282 283