Changeset 7206


Ignore:
Timestamp:
05/31/21 16:50:43 (8 weeks ago)
Author:
agnes.ducharne
Message:

Complement of r6950: mcs and so_capa_ice defined as global arrays, and mcs provided from hydrol. Checked step by 5d simulations on jean-zay: no change.

Location:
branches/ORCHIDEE_2_2/ORCHIDEE/src_sechiba
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • branches/ORCHIDEE_2_2/ORCHIDEE/src_sechiba/sechiba.f90

    r6954 r7206  
    497497 
    498498    !! 1.10 Initialization of soil thermodynamics 
    499     CALL thermosoil_initialize (kjit, kjpindex, rest_id, & 
     499    CALL thermosoil_initialize (kjit, kjpindex, rest_id, mcs, & 
    500500         temp_sol_new, snow,       shumdiag_perma,        & 
    501501         soilcap,      soilflx,    stempdiag,             & 
     
    771771    !! 7. Compute soil thermodynamics 
    772772    CALL thermosoil_main (kjit, kjpindex, & 
    773          index, indexgrnd, & 
     773         index, indexgrnd, mcs, & 
    774774         temp_sol_new, snow, soilcap, soilflx, & 
    775775         shumdiag_perma, stempdiag, ptnlev1, rest_id, hist_id, hist2_id, & 
  • branches/ORCHIDEE_2_2/ORCHIDEE/src_sechiba/thermosoil.f90

    r7199 r7206  
    152152  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)   :: dz5                      !! Used for numerical calculation [-] 
    153153!$OMP THREADPRIVATE(dz5) 
    154   REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)   :: mcs                      !! Saturation humidity [m3/m3] 
    155 !$OMP THREADPRIVATE(mcs) 
    156154  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)   :: QZ                       !! quartz content [-] 
    157155!$OMP THREADPRIVATE(QZ) 
     
    231229  !! \n 
    232230  !_ ============================================================================================================================== 
    233   SUBROUTINE thermosoil_initialize (kjit,          kjpindex,   rest_id,                   & 
     231  SUBROUTINE thermosoil_initialize (kjit,          kjpindex,   rest_id,          mcs,     & 
    234232                                    temp_sol_new,  snow,       shumdiag_perma,            & 
    235233                                    soilcap,       soilflx,    stempdiag,                 & 
     
    244242    INTEGER(i_std), INTENT(in)                            :: kjpindex         !! Domain size (unitless) 
    245243    INTEGER(i_std),INTENT (in)                            :: rest_id          !! Restart file identifier (unitless) 
     244    REAL(r_std),DIMENSION (kjpindex), INTENT (in)         :: mcs              !! Saturated moisture content (m3/m3) 
    246245    REAL(r_std),DIMENSION (kjpindex), INTENT (in)         :: temp_sol_new     !! Surface temperature at the present time-step, 
    247246    REAL(r_std),DIMENSION (kjpindex), INTENT (in)         :: snow             !! Snow mass (kg) 
     
    370369    IF (ier /= 0) CALL ipslerr_p(3,'thermosoil_initialize', 'Error in allocation of tmc_layt','','') 
    371370 
    372     ALLOCATE (mcs(nscm),stat=ier) 
    373     IF (ier /= 0) CALL ipslerr_p(3,'thermosoil_initialize', 'Error in allocation of mcs','','') 
    374  
    375371    ALLOCATE (QZ(nscm),stat=ier) 
    376372    IF (ier /= 0) CALL ipslerr_p(3,'thermosoil_initialize', 'Error in allocation of QZ','','') 
     
    379375    IF (ier /= 0) CALL ipslerr_p(3,'thermosoil_initialize', 'Error in allocation of so_capa_dry_ns','','') 
    380376     
    381     ALLOCATE (so_capa_ice(nscm),stat=ier) 
     377    ALLOCATE (so_capa_ice(kjpindex),stat=ier) 
    382378    IF (ier /= 0) CALL ipslerr_p(3,'thermosoil_initialize', 'Error in allocation of so_capa_ice','','') 
    383379 
     
    388384       QZ(:) = QZ_fao(:) 
    389385       so_capa_dry_ns(:) = so_capa_dry_ns_fao(:) 
    390        mcs(:) = mcs_fao(:) 
    391386    CASE (13) !Salma changed from 12 to 13 for the new class Oxisols 
    392387       QZ(:) = QZ_usda(:) 
    393388       so_capa_dry_ns(:) = so_capa_dry_ns_usda(:) 
    394        mcs(:) = mcs_usda(:) 
    395389    CASE DEFAULT 
    396390       WRITE (numout,*) 'Unsupported soil type classification. Choose between zobler, fao and usda according to the map' 
     
    504498 
    505499    ! Calculate so_capa_ice 
    506     DO jsc = 1, nscm   
    507        so_capa_ice(jsc) = so_capa_dry + mcs(jsc)*capa_ice*rho_ice  
    508     END DO 
     500    so_capa_ice(:) = so_capa_dry + mcs(:)*capa_ice*rho_ice  
    509501    IF (printlev>=2) WRITE(numout,*) 'Calculation of so_capa_ice(:)=', so_capa_ice(:),' and capa_ice=',capa_ice 
    510502     
     
    530522       CALL thermosoil_coef (& 
    531523            kjpindex,      temp_sol_new,    snow,           njsc, & 
    532             frac_snow_veg, frac_snow_nobio, totfrac_nobio,        & 
     524            mcs, frac_snow_veg, frac_snow_nobio, totfrac_nobio,        & 
    533525            snowdz,        snowrho,         snowtemp,       pb,   & 
    534526            ptn,                                                  & 
     
    582574 
    583575  SUBROUTINE thermosoil_main (kjit, kjpindex, & 
    584        index, indexgrnd, & 
     576       index, indexgrnd, mcs, & 
    585577       temp_sol_new, snow, soilcap, soilflx, & 
    586578       shumdiag_perma, stempdiag, ptnlev1, rest_id, hist_id, hist2_id, & 
     
    601593    INTEGER(i_std),DIMENSION (kjpindex*ngrnd), INTENT (in):: indexgrnd        !! Indeces of the points on the 3D map (vertical  
    602594                                                                              !! dimension towards the ground) (unitless) 
     595    REAL(i_std),DIMENSION (kjpindex), INTENT (in)         :: mcs              !! Saturated moisture content (m3/m3) 
    603596    REAL(r_std),DIMENSION (kjpindex), INTENT (inout)      :: temp_sol_new     !! Surface temperature at the present time-step, 
    604597                                                                              !! Ts @tex ($K$) @endtex 
     
    756749    CALL thermosoil_coef (& 
    757750         kjpindex,      temp_sol_new,    snow,         njsc, & 
    758          frac_snow_veg, frac_snow_nobio, totfrac_nobio,      & 
     751         mcs, frac_snow_veg, frac_snow_nobio, totfrac_nobio,      & 
    759752         snowdz,        snowrho,         snowtemp,     pb,   & 
    760753         ptn,                                                & 
     
    918911 
    919912  SUBROUTINE thermosoil_coef (kjpindex,      temp_sol_new,    snow,           njsc, & 
    920                               frac_snow_veg, frac_snow_nobio, totfrac_nobio,        & 
     913                              mcs,           frac_snow_veg,   frac_snow_nobio,totfrac_nobio, & 
    921914                              snowdz,        snowrho,         snowtemp,       pb,   & 
    922915                              ptn,                                                  & 
     
    932925    REAL(r_std), DIMENSION (kjpindex), INTENT (in)         :: snow         !! snow mass @tex ($Kg$) @endtex 
    933926    INTEGER(i_std),DIMENSION (kjpindex), INTENT (in)       :: njsc         !! Index of the dominant soil textural class in the grid cell (1-nscm, unitless) 
     927    REAL(r_std), DIMENSION (kjpindex), INTENT (in)         :: mcs          !! Saturated moisture content (m3/m3) 
    934928    REAL(r_std),DIMENSION (kjpindex), INTENT(in)           :: frac_snow_veg   !! Snow cover fraction on vegeted area 
    935929    REAL(r_std),DIMENSION (kjpindex,nnobio), INTENT(in)    :: frac_snow_nobio !! Snow cover fraction on non-vegeted area 
     
    994988    ! Computation of the soil thermal properties; snow properties are also accounted for 
    995989    IF (ok_freeze_thermix) THEN 
    996        CALL thermosoil_getdiff( kjpindex, snow, ptn, njsc, snowrho, snowtemp, pb ) 
     990       CALL thermosoil_getdiff( kjpindex, snow, ptn, mcs, njsc, snowrho, snowtemp, pb ) 
    997991    ELSE 
    998992       ! Special case without soil freezing 
    999        CALL thermosoil_getdiff_old_thermix_without_snow( kjpindex, njsc, snowrho, snowtemp, pb ) 
     993       CALL thermosoil_getdiff_old_thermix_without_snow( kjpindex, mcs, njsc, snowrho, snowtemp, pb ) 
    1000994    ENDIF 
    1001995 
     
    12581252!================================================================================================================================ 
    12591253 
    1260   SUBROUTINE thermosoil_cond (kjpindex, njsc, smc, qz, sh2o, cnd) 
     1254  SUBROUTINE thermosoil_cond (kjpindex, njsc, mcs, smc, qz, sh2o, cnd) 
    12611255 
    12621256    !! 0. Variables and parameter declaration 
     
    12651259    INTEGER(i_std), INTENT(in)                                 :: kjpindex      !! Domain size (unitless) 
    12661260    INTEGER(i_std), DIMENSION (kjpindex), INTENT (in)          :: njsc          !! Index of the dominant soil textural class in the grid cell (1-nscm, unitless) 
     1261    REAL(r_std), DIMENSION (kjpindex), INTENT(IN)              :: mcs           !! Saturated moisture content (m3/m3) 
    12671262    REAL(r_std), DIMENSION (kjpindex,ngrnd), INTENT(IN)        :: smc           !! Volumetric Soil Moisture Content (m3/m3) 
    12681263    REAL(r_std), DIMENSION (nscm), INTENT(IN)                  :: qz            !! Quartz Content (Soil Type Dependent) (0-1) 
     
    12971292 
    12981293      !! 1.1. Dry density (Kg/m3) and Dry thermal conductivity (W.M-1.K-1) 
    1299       gammd = (1. - mcs(jst))*2700. 
     1294      gammd = (1. - mcs(ji))*2700. 
    13001295      thkdry = (0.135* gammd+ 64.7)/ (2700. - 0.947* gammd) 
    13011296 
     
    13121307      DO jg = 1,ngrnd       
    13131308        !! 1.4. saturation ratio 
    1314         satratio(ji,jg) = smc(ji,jg) / mcs(jst) 
     1309        satratio(ji,jg) = smc(ji,jg) / mcs(ji) 
    13151310     
    13161311        !! 1.5. Saturated Thermal Conductivity (thksat) 
    13171312        IF ( smc(ji,jg) > min_sechiba ) THEN 
    13181313           xunfroz(ji,jg) = sh2o(ji,jg) / smc(ji,jg)   ! Unfrozen Fraction (From i.e., 100%Liquid, to 0. (100% Frozen)) 
    1319            xu(ji,jg) = xunfroz(ji,jg) * mcs(jst)  ! Unfrozen volume for saturation (porosity*xunfroz) 
    1320            thksat(ji,jg) = thks ** (1. - mcs(jst))* THKICE ** (mcs(jst) - xu(ji,jg))* THKW ** (xu(ji,jg)) 
     1314           xu(ji,jg) = xunfroz(ji,jg) * mcs(ji)  ! Unfrozen volume for saturation (porosity*xunfroz) 
     1315           thksat(ji,jg) = thks ** (1. - mcs(ji))* THKICE ** (mcs(ji) - xu(ji,jg))* THKW ** (xu(ji,jg)) 
    13211316        ELSE 
    13221317           ! this value will not be used since ake=0 for this case 
     
    16061601!_ ================================================================================================================================ 
    16071602 
    1608   SUBROUTINE thermosoil_getdiff( kjpindex, snow, ptn, njsc, snowrho, snowtemp, pb ) 
     1603  SUBROUTINE thermosoil_getdiff( kjpindex, snow, ptn, mcs, njsc, snowrho, snowtemp, pb ) 
    16091604 
    16101605   !! 0. Variables and parameter declaration 
     
    16141609    REAL(r_std),DIMENSION(kjpindex),INTENT (in)         :: snow       !! Snow mass 
    16151610    INTEGER(i_std),DIMENSION (kjpindex), INTENT (in)    :: njsc       !! Index of the dominant soil textural class in the grid cell (1-nscm, unitless) 
     1611    REAL(r_std),DIMENSION (kjpindex), INTENT (in)       :: mcs        !! Saturated moisture content (m3/m3) 
    16161612    REAL(r_std), DIMENSION (kjpindex,nsnow), INTENT(in) :: snowrho    !! Snow density 
    16171613    REAL(r_std), DIMENSION (kjpindex,nsnow), INTENT(in) :: snowtemp   !! Snow temperature (K) 
    1618     REAL(r_std),DIMENSION (kjpindex), INTENT (in)       :: pb         !! Surface presure (hPa) 
     1614    REAL(r_std),DIMENSION (kjpindex), INTENT (in)       :: pb         !! Surface pressure (hPa) 
    16191615    REAL(r_std),DIMENSION(kjpindex,ngrnd),INTENT(in)    :: ptn        !! Soil temperature profile 
    16201616 
     
    16401636      DO jg = 1, ngrnd 
    16411637         jst = njsc(ji) 
    1642          pcapa_tmp(ji, jg) = so_capa_dry_ns(jst) * (1-mcs(jst)) + water_capa * tmc_layt(ji,jg)/mille/dlt(jg) 
     1638         pcapa_tmp(ji, jg) = so_capa_dry_ns(jst) * (1-mcs(ji)) + water_capa * tmc_layt(ji,jg)/mille/dlt(jg) 
    16431639         ! 
    16441640         ! 2. Calculate volumetric heat capacity with allowance for permafrost 
     
    16511647            profil_froz(ji,jg) = 1. 
    16521648            pcappa_supp(ji,jg)= 0. 
    1653             pcapa(ji, jg) = so_capa_dry_ns(jst) * (1-mcs(jst)) + so_capa_ice(jst) * tmc_layt(ji,jg) / mille / dlt(jg) 
    1654             rho_tot = rho_soil * (1-mcs(jst)) + rho_ice * tmc_layt(ji,jg) / mille / dlt(jg)  
     1649            pcapa(ji, jg) = so_capa_dry_ns(jst) * (1-mcs(ji)) + so_capa_ice(ji) * tmc_layt(ji,jg) / mille / dlt(jg) 
     1650            rho_tot = rho_soil * (1-mcs(ji)) + rho_ice * tmc_layt(ji,jg) / mille / dlt(jg)  
    16551651            pcapa_spec(ji, jg) = pcapa(ji, jg) / rho_tot 
    16561652 
     
    16601656            profil_froz(ji,jg) = 0. 
    16611657            pcappa_supp(ji,jg)= 0. 
    1662             rho_tot = rho_soil * (1-mcs(jst)) + rho_water * tmc_layt(ji,jg)/mille/dlt(jg) 
     1658            rho_tot = rho_soil * (1-mcs(ji)) + rho_water * tmc_layt(ji,jg)/mille/dlt(jg) 
    16631659            pcapa_spec(ji, jg) = pcapa(ji, jg) / rho_tot 
    16641660         ELSE 
     
    16681664 
    16691665           IF (ok_freeze_thaw_latent_heat) THEN 
    1670               pcapa(ji, jg) = so_capa_dry_ns(jst) * (1-mcs(jst)) + & 
     1666              pcapa(ji, jg) = so_capa_dry_ns(jst) * (1-mcs(ji)) + & 
    16711667                water_capa * tmc_layt(ji,jg)/mille / dlt(jg) * xx + & 
    1672                 so_capa_ice(jst) * tmc_layt(ji,jg) / mille/dlt(jg) * (1.-xx) + & 
    1673                 shum_ngrnd_perma(ji,jg)*mcs(jst)*lhf*rho_water/fr_dT 
     1668                so_capa_ice(ji) * tmc_layt(ji,jg) / mille/dlt(jg) * (1.-xx) + & 
     1669                shum_ngrnd_perma(ji,jg)*mcs(ji)*lhf*rho_water/fr_dT 
    16741670           ELSE 
    1675               pcapa(ji, jg) = so_capa_dry_ns(jst) * (1-mcs(jst)) + & 
     1671              pcapa(ji, jg) = so_capa_dry_ns(jst) * (1-mcs(ji)) + & 
    16761672                water_capa * tmc_layt(ji,jg)/mille / dlt(jg) * xx + & 
    1677                 so_capa_ice(jst) * tmc_layt(ji,jg) / mille/dlt(jg) * (1.-xx) 
     1673                so_capa_ice(ji) * tmc_layt(ji,jg) / mille/dlt(jg) * (1.-xx) 
    16781674           ENDIF 
    16791675 
    1680            rho_tot =  rho_soil* (1-mcs(jst)) + & 
     1676           rho_tot =  rho_soil* (1-mcs(ji)) + & 
    16811677                rho_water * tmc_layt(ji,jg)/mille / dlt(jg) * xx + & 
    16821678                rho_ice * tmc_layt(ji,jg) / mille/dlt(jg) * (1.-xx) 
    16831679           pcapa_spec(ji, jg) = pcapa(ji, jg) / rho_tot 
    16841680 
    1685            pcappa_supp(ji,jg)= shum_ngrnd_perma(ji,jg)*mcs(jst)*lhf*rho_water/fr_dT*zx2(ji,jg)*dlt(jg) 
     1681           pcappa_supp(ji,jg)= shum_ngrnd_perma(ji,jg)*mcs(ji)*lhf*rho_water/fr_dT*zx2(ji,jg)*dlt(jg) 
    16861682 
    16871683         ENDIF 
     
    17101706    ! 
    17111707    IF (ok_freeze_thaw_latent_heat) THEN 
    1712         CALL thermosoil_cond (kjpindex, njsc, mc_layt, QZ, mcl_layt*(1-profil_froz), pkappa) 
     1708        CALL thermosoil_cond (kjpindex, njsc, mcs, mc_layt, QZ, mcl_layt*(1-profil_froz), pkappa) 
    17131709    ELSE 
    1714         CALL thermosoil_cond (kjpindex, njsc, mc_layt, QZ, mcl_layt, pkappa) 
     1710        CALL thermosoil_cond (kjpindex, njsc, mcs, mc_layt, QZ, mcl_layt, pkappa) 
    17151711    ENDIF 
    17161712 
     
    17441740!_ ================================================================================================================================ 
    17451741 
    1746     SUBROUTINE thermosoil_getdiff_old_thermix_without_snow( kjpindex, njsc, snowrho, snowtemp, pb ) 
     1742    SUBROUTINE thermosoil_getdiff_old_thermix_without_snow( kjpindex, mcs, njsc, snowrho, snowtemp, pb ) 
    17471743 
    17481744   !! 0. Variables and parameter declaration 
     
    17511747      INTEGER(i_std), INTENT(in) :: kjpindex 
    17521748      INTEGER(i_std),DIMENSION (kjpindex), INTENT (in)    :: njsc     !! Index of the dominant soil textural class in the grid cell (1-nscm, unitless) 
     1749      REAL(r_std),DIMENSION (kjpindex), INTENT (in)       :: mcs      !! Saturated moisture content (m3/m3) 
    17531750      REAL(r_std), DIMENSION (kjpindex,nsnow), INTENT(in) :: snowrho  !! Snow density 
    17541751      REAL(r_std), DIMENSION (kjpindex,nsnow), INTENT(in) :: snowtemp !! Snow temperature (K) 
    1755       REAL(r_std),DIMENSION (kjpindex), INTENT (in)       :: pb       !! Surface presure (hPa) 
     1752      REAL(r_std),DIMENSION (kjpindex), INTENT (in)       :: pb       !! Surface pressure (hPa) 
    17561753 
    17571754 
     
    17641761         DO ji = 1,kjpindex 
    17651762            jst = njsc(ji) 
    1766             pcapa_tmp(ji, jg) = so_capa_dry_ns(jst) * (1-mcs(jst)) + water_capa * tmc_layt(ji,jg)/mille/dlt(jg) 
     1763            pcapa_tmp(ji, jg) = so_capa_dry_ns(jst) * (1-mcs(ji)) + water_capa * tmc_layt(ji,jg)/mille/dlt(jg) 
    17671764            pcapa(ji,jg) = pcapa_tmp(ji, jg) 
    17681765            pcapa_en(ji,jg) = pcapa_tmp(ji, jg) 
     
    17701767      ENDDO 
    17711768 
    1772       CALL thermosoil_cond (kjpindex, njsc, mc_layt, QZ, mcl_layt, pkappa) 
     1769      CALL thermosoil_cond (kjpindex, njsc, mcs, mc_layt, QZ, mcl_layt, pkappa) 
    17731770 
    17741771      IF (brk_flag == 1) THEN 
Note: See TracChangeset for help on using the changeset viewer.