1 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv! 2 ! Subroutine name: PHYSICAL_PROP ! 3 ! ! 4 ! Purpose: Calculate physical properties that vary with time. ! 5 ! ! 6 ! Author: J.Musser Date: 09-May-11 ! 7 ! ! 8 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^! 9 SUBROUTINE DES_PHYSICAL_PROP(NP, FOCUS) 10 11 Use des_rxns 12 Use des_thermo 13 Use discretelement 14 Use funits 15 Use param 16 Use param1 17 Use physprop 18 Use run 19 20 ! Calculate the specific heat from polynomical data obtained from the 21 ! thermodynamic databases. 22 use read_thermochemical, only: DES_calc_CpoR 23 24 ! Universal gas constant in cal/mol.K 25 use constant, only: RGAS => GAS_CONST_cal 26 27 IMPLICIT NONE 28 29 ! Passed Variables 30 !----------------------------------------------------------------------- 31 ! Index of particle 32 INTEGER, INTENT(IN) :: NP 33 ! Logical indicating to write additional data about the particle for 34 ! debugging purposes. 35 LOGICAL, INTENT(IN) :: FOCUS 36 37 ! Local Variables 38 !----------------------------------------------------------------------- 39 ! Dummy indices 40 INTEGER M, N 41 ! error indicator 42 INTEGER IER 43 44 DOUBLE PRECISION :: lCPoR 45 46 ! Get the solids phase of the particle 47 M = PIJK(NP,5) + SMAX 48 49 ! Specific heat 50 !----------------------------------------------------------------------- 51 ! This only needs calculated when solving the energy equations. 52 IF(ENERGY_EQ) THEN 53 ! If a constant value specific heat has not been defined in the mfix.dat 54 ! file, calculate the temperature dependent value based on data from the 55 ! thermodynamic databases. 56 IF (C_PS0(M) == UNDEFINED) THEN 57 ! Read the thermodynamic database if it has not already been done. 58 DES_C_PS(NP) = ZERO 59 ! Calculate the specific heat based on the species composition of the 60 ! particle and the data from the thermodynamic databases. 61 DO N = 1, NMAX_s(M) 62 lCPoR = DES_calc_CpoR(DES_T_s_NEW(NP), M, N, IER) 63 DES_C_PS(NP) = DES_C_PS(NP) + & 64 lCpoR * DES_X_s(NP,N) * RGAS / MW_s(M,N) 65 ENDDO 66 ! Convert to SI units if needed. 67 IF (UNITS == 'SI') DES_C_PS(NP) = 4183.925D0*DES_C_PS(NP) 68 ELSE 69 ! If a constant value specific heat has been assigned to the particle 70 ! in the mfix.dat file, use this value. 71 DES_C_PS(NP) = C_PS0(M) 72 ENDIF 73 ENDIF 74 75 ! WRITE(*,"(3X,A,I2,A,F10.6)")'DES_C_PS(',NP,'): ',DES_C_PS(NP) 76 77 RETURN 78 79 END SUBROUTINE DES_PHYSICAL_PROP 80