|
SURFEX v7.3
General documentation of Surfex
|
|
SURFEX v7.3
General documentation of Surfex
|
00001 ! ######### 00002 SUBROUTINE SOIL_HEATDIF(PTSTEP,PDZG,PDZDIF,PSOILCONDZ, & 00003 PSOILHCAPZ,PCT,PTERM1,PTERM2,PTDEEP_A,PTDEEP_B,PTG,PDEEP_FLUX) 00004 ! ############################################################################ 00005 ! 00006 !!**** *SOIL_HEATDIF* 00007 ! 00008 !! PURPOSE 00009 !! ------- 00010 ! This subroutine solves the 1-d diffusion of 'PTG' using a 00011 ! layer averaged set of equations which are time differenced 00012 ! using the backward-difference scheme (implicit). 00013 ! The eqs are solved rapidly by taking advantage of the 00014 ! fact that the matrix is tridiagonal. This is a very 00015 ! general routine and can be used for the 1-d diffusion of any 00016 ! quantity as long as the diffusity is not a function of the 00017 ! quantity being diffused. Aaron Boone 8-98. Soln to the eq: 00018 ! 00019 ! dQ d dQ 00020 ! c -- = -- k -- 00021 ! dt dx dx 00022 ! 00023 ! where k = k(x) (thermal conductivity), c = c(x) (heat capacity) 00024 ! Diffusivity is k/c 00025 ! 00026 ! 00027 !!** METHOD 00028 !! ------ 00029 ! 00030 ! Direct calculation 00031 ! 00032 !! EXTERNAL 00033 !! -------- 00034 ! 00035 ! None 00036 !! 00037 !! IMPLICIT ARGUMENTS 00038 !! ------------------ 00039 !! 00040 !! USE MODD_PARAMETERS 00041 !! USE MODI_TRIDIAG_GROUND 00042 !! 00043 !! REFERENCE 00044 !! --------- 00045 !! 00046 !! Boone (2000) 00047 !! Boone et al. (2000) 00048 !! 00049 !! AUTHOR 00050 !! ------ 00051 !! A. Boone * Meteo-France * 00052 !! 00053 !! MODIFICATIONS 00054 !! ------------- 00055 !! Original 16/02/00 Boone 00056 !! Modif 08/2011 B. Decharme : Optimization 00057 ! 00058 !------------------------------------------------------------------------------- 00059 ! 00060 !* 0. DECLARATIONS 00061 ! ------------ 00062 ! 00063 USE MODD_SURF_PAR, ONLY : XUNDEF 00064 ! 00065 USE MODI_TRIDIAG_GROUND 00066 ! 00067 ! 00068 USE YOMHOOK ,ONLY : LHOOK, DR_HOOK 00069 USE PARKIND1 ,ONLY : JPRB 00070 ! 00071 IMPLICIT NONE 00072 ! 00073 ! 00074 REAL, INTENT(IN) :: PTSTEP ! Model time step (s) 00075 ! 00076 REAL, DIMENSION(:), INTENT(IN) :: PCT, PTERM1, PTERM2, PTDEEP_A, PTDEEP_B 00077 ! PCT = thermal inertia [(m2 K)/J] 00078 ! PTERM1 = coefficient of linearization 00079 ! of surface energy budget 00080 ! PTERM2 = coefficient of linearization 00081 ! of surface energy budget 00082 ! PTDEEP_A = Deep soil temperature 00083 ! coefficient depending on flux 00084 ! PTDEEP_B = Deep soil temperature (prescribed) 00085 ! which models heating/cooling from 00086 ! below the diurnal wave penetration 00087 ! (surface temperature) depth. If it 00088 ! is FLAGGED as undefined, then the zero 00089 ! flux lower BC is applied. 00090 ! Tdeep = PTDEEP_B + PTDEEP_A * PDEEP_FLUX 00091 ! (with PDEEP_FLUX in W/m2) 00092 ! 00093 REAL, DIMENSION(:,:), INTENT(IN) :: PSOILCONDZ, PSOILHCAPZ 00094 ! PSOILCONDZ = soil thermal conductivity [W/(K m)] 00095 ! PSOILHCAPZ = soil heat capacity [J/(m3 K)] 00096 REAL, DIMENSION(:,:), INTENT(IN) :: PDZDIF, PDZG 00097 ! PDZDIF = distance between consecuative layer mid-points 00098 ! PDZG = soil layers thicknesses 00099 ! 00100 REAL, DIMENSION(:,:), INTENT(INOUT) :: PTG 00101 ! PTG = soil temperature (K) 00102 ! 00103 REAL, DIMENSION(:), INTENT(OUT) :: PDEEP_FLUX ! Heat flux at bottom of ISBA (W/m2) 00104 ! 00105 ! 00106 !* 0.2 declarations of local variables 00107 ! 00108 INTEGER :: JJ, JL 00109 ! 00110 INTEGER :: INI, INLVLD ! Number of point and grid layers 00111 ! 00112 REAL, DIMENSION(SIZE(PTG,1),SIZE(PTG,2)) :: ZTGM, ZDTERM, ZCTERM, 00113 ZFRCV, ZAMTRX, ZBMTRX, 00114 ZCMTRX 00115 ! 00116 REAL :: ZWORK1, ZWORK2 00117 ! 00118 REAL(KIND=JPRB) :: ZHOOK_HANDLE 00119 ! 00120 !------------------------------------------------------------------------------- 00121 ! 00122 ! Initialize local variables: 00123 ! 00124 IF (LHOOK) CALL DR_HOOK('SOIL_HEATDIF',0,ZHOOK_HANDLE) 00125 ZDTERM(:,:) = 0.0 00126 ZCTERM(:,:) = 0.0 00127 ZFRCV(:,:) = 0.0 00128 ZAMTRX(:,:) = 0.0 00129 ZBMTRX(:,:) = 0.0 00130 ZCMTRX(:,:) = 0.0 00131 ZTGM(:,:) = PTG(:,:) 00132 ! 00133 INI = SIZE(PTG(:,:),1) 00134 INLVLD = SIZE(PTG(:,:),2) 00135 ! 00136 !------------------------------------------------------------------------------- 00137 ! 00138 ! Calculate tri-diagonal matrix coefficients: 00139 ! 00140 DO JL=1,INLVLD 00141 DO JJ=1,INI 00142 ZWORK1 = PDZG(JJ,JL)*PSOILCONDZ(JJ,JL) 00143 IF(JL<INLVLD)THEN 00144 ZWORK2 = PDZG(JJ,JL+1)*PSOILCONDZ(JJ,JL+1) 00145 ELSE 00146 ZWORK2 = 0.0 00147 ENDIF 00148 ZDTERM(JJ,JL)=(ZWORK1+ZWORK2)/(2.0*PDZDIF(JJ,JL)*PDZDIF(JJ,JL)) 00149 ZCTERM(JJ,JL)= PSOILHCAPZ(JJ,JL)*PDZG(JJ,JL)/PTSTEP 00150 ENDDO 00151 ENDDO 00152 ! 00153 ! - - ------------------------------------------------- 00154 ! 00155 ! Upper BC 00156 ! 00157 ZAMTRX(:,1) = 0.0 00158 ZBMTRX(:,1) = 1./(PCT(:)*PTSTEP) 00159 ZCMTRX(:,1) = -PTERM2(:)*ZBMTRX(:,1) 00160 ZFRCV (:,1) = PTERM1(:)*ZBMTRX(:,1) 00161 ! 00162 ! Interior Grid 00163 ! 00164 DO JL=2,INLVLD-1 00165 DO JJ=1,INI 00166 ZAMTRX(JJ,JL) = -ZDTERM(JJ,JL-1) 00167 ZBMTRX(JJ,JL) = ZCTERM(JJ,JL) + ZDTERM(JJ,JL-1) + ZDTERM(JJ,JL) 00168 ZCMTRX(JJ,JL) = -ZDTERM(JJ,JL) 00169 ZFRCV (JJ,JL) = ZCTERM(JJ,JL)*ZTGM(JJ,JL) 00170 ENDDO 00171 ENDDO 00172 ! 00173 ! Lower BC: 2 currently accounted for, Either zero flux 00174 ! or a fixed temperature 'TDEEP' 00175 ! 00176 ZAMTRX(:,INLVLD) = -ZDTERM(:,INLVLD-1) 00177 ZCMTRX(:,INLVLD) = 0.0 00178 ! 00179 00180 !ZDEEP_FLUX=ZDTERM(:,INLVLD)*(ZTGM(:,INLVLD)-ZTDEEP(:)) 00181 !ZTDEEP=PTDEEP_A*ZDEEP_FLUX+PTDEEP_B 00182 ! 00183 !ZDEEP_FLUX=ZDTERM(:,INLVLD)*(ZTGM(:,INLVLD)-PTDEEP_A*ZDEEP_FLUX-PTDEEP_B) 00184 !ZDEEP_FLUX=ZDTERM(:,INLVLD)*(ZTGM(:,INLVLD)-PTDEEP_B)/(1.+ZDTERM(:,INLVLD)*PTDEEP_A) 00185 ! 00186 ! 00187 WHERE(PTDEEP_B(:) /= XUNDEF) 00188 ! ZBMTRX(:,INLVLD) = ZCTERM(:,INLVLD) + ZDTERM(:,INLVLD-1) + ZDTERM(:,INLVLD) 00189 ! ZFRCV(:,INLVLD) = ZCTERM(:,INLVLD)*ZTGM(:,INLVLD) + ZDTERM(:,INLVLD)*PTDEEP(:) 00190 ZBMTRX(:,INLVLD) = ZCTERM(:,INLVLD) + ZDTERM(:,INLVLD-1) & 00191 + ZDTERM(:,INLVLD)/(1.+ZDTERM(:,INLVLD)*PTDEEP_A) 00192 ZFRCV(:,INLVLD) = ZCTERM(:,INLVLD)*ZTGM(:,INLVLD) & 00193 + ZDTERM(:,INLVLD)*PTDEEP_B(:)/(1.+ZDTERM(:,INLVLD)*PTDEEP_A) 00194 ELSEWHERE 00195 ZBMTRX(:,INLVLD) = ZCTERM(:,INLVLD) + ZDTERM(:,INLVLD-1) 00196 ZFRCV(:,INLVLD) = ZCTERM(:,INLVLD)*ZTGM(:,INLVLD) 00197 END WHERE 00198 ! 00199 ! - - ------------------------------------------------- 00200 ! 00201 ! Compute ZTGM (solution vector) 00202 ! used for systems of equations involving tridiagonal 00203 ! matricies. 00204 ! 00205 CALL TRIDIAG_GROUND(ZAMTRX,ZBMTRX,ZCMTRX,ZFRCV,ZTGM) 00206 ! 00207 ! 00208 ! Update values in time: 00209 ! 00210 PTG(:,:) = ZTGM(:,:) 00211 ! 00212 !* Deep soil Flux 00213 ! 00214 PDEEP_FLUX(:) = 0. 00215 WHERE(PTDEEP_B(:) /= XUNDEF) 00216 PDEEP_FLUX=ZDTERM(:,INLVLD)*(ZTGM(:,INLVLD)-PTDEEP_B)/(1.+ZDTERM(:,INLVLD)*PTDEEP_A) 00217 END WHERE 00218 ! 00219 00220 IF (LHOOK) CALL DR_HOOK('SOIL_HEATDIF',1,ZHOOK_HANDLE) 00221 ! 00222 ! 00223 ! 00224 END SUBROUTINE SOIL_HEATDIF
1.8.0
