SURFEX v7.3
General documentation of Surfex
 All Classes Files Functions Variables Typedefs
/home/dasprezs/EXPORT_v7_3/src/SURFEX/mode_hydro_dif.F90
Go to the documentation of this file.
00001 !     ######spl
00002       MODULE MODE_HYDRO_DIF 
00003 !     ################
00004 !
00005 !!****  *MODE_HYDRO_DIF * - pedo-transfert functions
00006 !!
00007 !!    PURPOSE
00008 !!    -------
00009 !
00010 !!**  METHOD
00011 !!    ------
00012 !!    
00013 !!
00014 !!    EXTERNAL
00015 !!    --------
00016 !!       
00017 !!    IMPLICIT ARGUMENTS
00018 !!    ------------------ 
00019 !!
00020 !!    REFERENCE
00021 !!    ---------
00022 !!      
00023 !!
00024 !!    AUTHOR
00025 !!    ------
00026 !!      B. Decharme       * Meteo France *
00027 !!
00028 !!    MODIFICATIONS
00029 !!    -------------
00030 !!      Original        11/2010
00031 !-----------------------------------------------------------------------------
00032 !
00033 !*       0.    DECLARATIONS
00034 !
00035 !
00036 INTERFACE VAPCONDCF
00037   MODULE PROCEDURE VAPCONDCF       
00038 END INTERFACE
00039 !
00040 INTERFACE INFMAX_FUNC
00041   MODULE PROCEDURE INFMAX_FUNC
00042 END INTERFACE
00043 !
00044 INTERFACE TRIDIAG_DIF
00045   MODULE PROCEDURE TRIDIAG_DIF
00046 END INTERFACE
00047 !
00048 !-------------------------------------------------------------------------------
00049 CONTAINS
00050 !-------------------------------------------------------------------------------
00051 !
00052 !-------------------------------------------------------------------------------
00053 !vapor conductivity (m s-1)
00054 !-------------------------------------------------------------------------------
00055 !
00056 FUNCTION VAPCONDCF(PTG,PPS,PWG,PWGI,PPSIA,PWSAT,PWFC,PQSAT,PQSATI,KWG_LAYER,KNL) RESULT(PVAPCOND)
00057 !
00058 ! Uses method of Braud et al. (1993) for
00059 !
00060 USE MODD_CSTS,       ONLY : XMV, XMD, XTT, XP00, XG, XRV, XRHOLW
00061 USE MODD_ISBA_PAR,   ONLY : XWGMIN
00062 USE MODD_SURF_PAR,   ONLY : XUNDEF
00063 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
00064 USE PARKIND1  ,ONLY : JPRB
00065 !
00066 IMPLICIT NONE
00067 !
00068 !
00069 !*      0.1    declarations of arguments
00070 !
00071 REAL, DIMENSION(:  ), INTENT(IN)         :: PPS
00072 REAL, DIMENSION(:,:), INTENT(IN)         :: PWG,PWGI,PPSIA,PWSAT, 
00073                                             PWFC,PTG,PQSAT,PQSATI
00074 INTEGER, DIMENSION(:), INTENT(IN)        :: KWG_LAYER       !Moisture layer
00075 INTEGER,               INTENT(IN)        :: KNL             ! number of vertical levels
00076 !
00077 REAL, DIMENSION(SIZE(PWG,1),SIZE(PWG,2)) :: PVAPCOND
00078 !
00079 !
00080 !*      0.2    declarations of local variables
00081 !
00082 REAL    :: ZDVA, ZFVA, ZCHI, ZHUM, ZWORK,  
00083            ZPV, ZESAT, ZESATI, ZWG, ZVC
00084 !
00085 INTEGER :: INI, JJ, JL, IDEPTH
00086 !
00087 ! Parameters:
00088 !
00089 REAL, PARAMETER                     :: ZTORTY = 0.66         ! (-)
00090 REAL, PARAMETER                     :: ZNV    = 1.88         ! (-)
00091 REAL, PARAMETER                     :: ZCV    = 2.17e-5      ! (m2/s)
00092 REAL, PARAMETER                     :: ZWK    = 0.05         ! (m3 m-3)
00093 REAL, PARAMETER                     :: ZLIM   = TINY(1.0)
00094 !
00095 REAL(KIND=JPRB) :: ZHOOK_HANDLE
00096 !-------------------------------------------------------------------------------
00097 !
00098 IF (LHOOK) CALL DR_HOOK('MODE_HYDRO_DIF:VAPCONDCF',0,ZHOOK_HANDLE)
00099 !
00100 INI = SIZE(PWG,1)
00101 !
00102 PVAPCOND(:,:) = 0.0
00103 !
00104 ! Only perform this computation if the soil is sufficiently
00105 ! dry (as otherwise the hydraulic conductivity dominates
00106 ! the diffusion coefficient). Arbitrarily base threshold on field
00107 ! capacity water content:
00108 !
00109 DO JL=1,KNL
00110    DO JJ=1,INI
00111 !
00112       IDEPTH = KWG_LAYER(JJ)
00113       ZWG    = PWG (JJ,JL) + PWGI(JJ,JL)
00114 !      
00115       IF(JL<=IDEPTH .AND. ZWG < PWFC(JJ,JL) .AND. ZWG > XWGMIN)THEN
00116 !
00117 !        Vapor pressure over liquid and solid water surfaces (Pa), respectively:
00118 !
00119          ZESAT  = PQSAT(JJ,JL)* PPS(JJ)/((XMV/XMD)+PQSAT(JJ,JL) *(1.-(XMV/XMD)))
00120 !
00121          ZESATI = PQSATI(JJ,JL)*PPS(JJ)/((XMV/XMD)+PQSATI(JJ,JL)*(1.-(XMV/XMD)))
00122 !
00123 !        molecular diffusivity of water vapor (m2 s-1):
00124 !
00125          ZWORK  = ZNV*LOG(PTG(JJ,JL)/XTT)
00126          ZDVA   = ZCV*(XP00/PPS(JJ))*EXP(ZWORK)
00127 !
00128 !        function of pore space: 
00129 !
00130          ZFVA   = (PWSAT(JJ,JL) - ZWG)*(1.+(ZWG/(PWSAT(JJ,JL)-ZWK)))
00131          ZFVA   = MIN(ZFVA,PWSAT(JJ,JL))
00132 !
00133 !        relative humidity of air in soil pores:
00134 !
00135          ZHUM   = MAX(ZLIM,EXP(PPSIA(JJ,JL)*XG/(XRV*PTG(JJ,JL))))
00136 !
00137 !        fraction of frozen water:
00138 !
00139          ZCHI   = PWGI(JJ,JL)/ZWG
00140 !
00141 !        vapor pressure within pore space (Pa):
00142 !
00143          ZPV    = ZHUM*(ZCHI*ZESAT + (1.-ZCHI)*ZESATI)
00144 !
00145 !        vapor conductivity (kg m-2 s-1)
00146 !
00147          ZVC    = ZTORTY*PPS(JJ)*ZDVA*ZFVA*XG*ZPV/                  &
00148                   ((PPS(JJ)-ZPV)*(XRV*XRV*PTG(JJ,JL)*PTG(JJ,JL)))
00149 !
00150 !        vapor conductivity (m s-1)
00151 !
00152          PVAPCOND(JJ,JL) = ZVC/XRHOLW
00153 !
00154       ENDIF
00155 !
00156    ENDDO
00157 ENDDO
00158 !
00159 IF (LHOOK) CALL DR_HOOK('MODE_HYDRO_DIF:VAPCONDCF',1,ZHOOK_HANDLE)
00160 !
00161 END FUNCTION VAPCONDCF
00162 !
00163 !-------------------------------------------------------------------------------
00164 !Green-Ampt approximation (derived form) for maximum infiltration
00165 !-------------------------------------------------------------------------------
00166 !
00167 FUNCTION INFMAX_FUNC(PWG,PWSAT,PFRZ,PCONDSAT,PMPOTSAT,PBCOEF,PDZG,PDG,KLAYER_HORT)
00168 USE YOMHOOK      ,ONLY : LHOOK,   DR_HOOK
00169 USE MODD_SGH_PAR, ONLY : XHORT_DEPTH
00170 USE PARKIND1     ,ONLY : JPRB
00171 IMPLICIT NONE
00172 REAL, DIMENSION(:,:), INTENT(IN) :: PWG
00173 REAL, DIMENSION(:,:), INTENT(IN) :: PWSAT           
00174 REAL, DIMENSION(:,:), INTENT(IN) :: PFRZ
00175 REAL, DIMENSION(:,:), INTENT(IN) :: PCONDSAT            
00176 REAL, DIMENSION(:,:), INTENT(IN) :: PMPOTSAT    
00177 REAL, DIMENSION(:,:), INTENT(IN) :: PBCOEF 
00178 REAL, DIMENSION(:,:), INTENT(IN) :: PDZG
00179 REAL, DIMENSION(:,:), INTENT(IN) :: PDG
00180 INTEGER,              INTENT(IN) :: KLAYER_HORT   
00181 !
00182 REAL, DIMENSION(SIZE(PWG,1)) :: ZGREEN_AMPT, ZDEPTH
00183 REAL                         :: ZS, ZCOEF
00184 INTEGER                      :: JJ,JL,INI
00185 !
00186 REAL, DIMENSION(SIZE(PWG,1)) :: INFMAX_FUNC
00187 !
00188 REAL(KIND=JPRB) :: ZHOOK_HANDLE
00189 !
00190 IF (LHOOK) CALL DR_HOOK('MODE_HYDRO_DIF:INFMAX_FUNC',0,ZHOOK_HANDLE)
00191 !
00192 INI   =SIZE(PWG,1)
00193 !
00194 ZGREEN_AMPT(:) = 0.0
00195 ZDEPTH     (:) = 0.0
00196 !
00197 DO JL=1,KLAYER_HORT
00198    DO JJ=1,INI  
00199       IF(ZDEPTH(JJ)<XHORT_DEPTH)THEN
00200          ZS              = MIN(1.0,PWG(JJ,JL)/PWSAT(JJ,JL))
00201          ZCOEF           = PBCOEF(JJ,JL)*PMPOTSAT(JJ,JL)*(ZS-1.0)/PDZG(JJ,JL)        
00202          ZGREEN_AMPT(JJ) = ZGREEN_AMPT(JJ)+PDZG(JJ,JL)*PFRZ(JJ,JL)*PCONDSAT(JJ,JL)*(ZCOEF+1.0)
00203          ZDEPTH     (JJ) = PDG(JJ,JL)
00204       ENDIF      
00205    ENDDO
00206 ENDDO
00207 !
00208 INFMAX_FUNC(:) = ZGREEN_AMPT(:)/ZDEPTH(:)
00209 !
00210 IF (LHOOK) CALL DR_HOOK('MODE_HYDRO_DIF:INFMAX_FUNC',1,ZHOOK_HANDLE)
00211 END FUNCTION INFMAX_FUNC
00212 !
00213 !-------------------------------------------------------------------------------
00214 !Solve tridiagonal matrix (for method see tridiag_ground.F90)
00215 !-------------------------------------------------------------------------------
00216 !
00217 SUBROUTINE TRIDIAG_DIF(PAMTRX,PBMTRX,PCMTRX,PFRC,KWG_LAYER,KNL,PSOL)
00218 USE MODD_SURF_PAR, ONLY : XUNDEF
00219 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
00220 USE PARKIND1  ,ONLY : JPRB
00221 IMPLICIT NONE
00222 REAL,    DIMENSION(:,:), INTENT(IN)  :: PAMTRX    ! lower diag. elements of A matrix
00223 REAL,    DIMENSION(:,:), INTENT(IN)  :: PBMTRX    ! main  diag. elements of A matrix
00224 REAL,    DIMENSION(:,:), INTENT(IN)  :: PCMTRX    ! upper diag. elements of A matrix
00225 REAL,    DIMENSION(:,:), INTENT(IN)  :: PFRC      ! Forcing term
00226 INTEGER, DIMENSION(:),   INTENT(IN)  :: KWG_LAYER !Moisture layer
00227 INTEGER,                 INTENT(IN)  :: KNL       ! number of vertical levels
00228 REAL,    DIMENSION(:,:), INTENT(OUT) :: PSOL      ! solution of A.SOL = FRC
00229 !
00230 REAL, DIMENSION(SIZE(PFRC,1),SIZE(PFRC,2)) :: ZWORK! work array
00231 REAL, DIMENSION(SIZE(PFRC,1))              :: ZDET ! work array
00232 INTEGER                                    :: JL   ! vertical loop control
00233 INTEGER :: JJ, INI, IDEPTH
00234 !
00235 REAL(KIND=JPRB) :: ZHOOK_HANDLE
00236 !
00237 IF (LHOOK) CALL DR_HOOK('MODE_HYDRO_DIF:TRIDIAG_DIF',0,ZHOOK_HANDLE)
00238 !
00239 INI=SIZE(PFRC,1)
00240 !
00241 PSOL (:,:)=XUNDEF
00242 ZWORK(:,:)=XUNDEF
00243 !
00244 !first level
00245 ZDET(:)   = PBMTRX(:,1)
00246 PSOL(:,1) = PFRC(:,1) / ZDET(:)
00247 !
00248 !other levels
00249 DO JL=2,KNL
00250    DO JJ=1,INI
00251       IDEPTH=KWG_LAYER(JJ)
00252       IF(JL<=IDEPTH)THEN
00253         ZWORK(JJ,JL) = PCMTRX(JJ,JL-1)/ZDET(JJ)
00254         ZDET (JJ)    = PBMTRX(JJ,JL) - PAMTRX(JJ,JL)*ZWORK(JJ,JL)
00255         PSOL (JJ,JL) = (PFRC (JJ,JL) - PAMTRX(JJ,JL)*PSOL(JJ,JL-1))/ZDET(JJ)  
00256       ENDIF
00257    ENDDO 
00258 ENDDO        
00259 !
00260 !levels going down
00261 DO JL=KNL-1,1,-1
00262    DO JJ=1,INI
00263       IDEPTH=KWG_LAYER(JJ)
00264       IF(JL<IDEPTH)THEN
00265          PSOL(JJ,JL) = PSOL(JJ,JL)-ZWORK(JJ,JL+1)*PSOL(JJ,JL+1)
00266       ENDIF
00267    ENDDO
00268 ENDDO
00269 !
00270 IF (LHOOK) CALL DR_HOOK('MODE_HYDRO_DIF:TRIDIAG_DIF',1,ZHOOK_HANDLE)
00271 !
00272 END SUBROUTINE TRIDIAG_DIF
00273 !
00274 !-------------------------------------------------------------------------------
00275 !
00276 END MODULE MODE_HYDRO_DIF