Ignore:
Timestamp:
Jun 24, 2020, 11:27:16 AM (4 years ago)
Author:
dubos
Message:

simple_physics : turn zc, zd, capcal, fluxgrd into local temporaries

File:
1 edited

Legend:

Unmodified
Added
Removed

  • dynamico_lmdz/simple_physics/phyparam/physics/surface.F90

    r4242 r4245  
    1616  REAL :: lambda
    1717  REAL, ALLOCATABLE :: dz1(:),dz2(:)
    18   !$OMP THREADPRIVATE(dz1,dz2,zc)
     18  !$OMP THREADPRIVATE(dz1,dz2)
    1919  REAL, ALLOCATABLE :: rnatur(:), albedo(:),emissiv(:), z0(:), inertie(:)
    2020  !$OMP THREADPRIVATE( rnatur, albedo, emissiv, z0, inertie)
    2121
    2222  ! internal state, written to / read from disk at checkpoint / restart
    23   REAL, ALLOCATABLE :: zc(:,:),zd(:,:)
    24   !$OMP THREADPRIVATE(zc,zd)
    25   REAL, ALLOCATABLE :: tsurf(:), tsoil(:,:), capcal(:), fluxgrd(:)
    26   !$OMP THREADPRIVATE(tsurf, tsoil, capcal, fluxgrd)
    27 
    28   PUBLIC :: soil, zc, zd, &
     23  REAL, ALLOCATABLE :: tsurf(:), tsoil(:,:)
     24  !$OMP THREADPRIVATE(tsurf, tsoil)
     25  ! variables below should be temporary arrays, not persistent
     26  REAL, ALLOCATABLE :: zc(:,:),zd(:,:), capcal(:), fluxgrd(:)
     27  !$OMP THREADPRIVATE(zc,zd, capcal, fluxgrd)
     28
     29  PUBLIC :: init_soil, &
     30       soil, soil_new, soil_forward, soil_backward, &
     31       zc, zd, &
    2932       rnatur, albedo, emissiv, z0, inertie, &
    3033       tsurf, tsoil, capcal, fluxgrd
    … …  
    3235CONTAINS
    3336
    34   SUBROUTINE init_soil(ngrid,nsoil)
    35     INTEGER, INTENT(IN) :: ngrid, nsoil
     37  SUBROUTINE init_soil(nsoil)
     38    INTEGER, INTENT(IN) :: nsoil
    3639    REAL :: min_period,dalph_soil, rk,fz1,rk1,rk2
    3740    INTEGER :: jk
    … …  
    4245    !   --------------------------------------------------------
    4346
    44     WRITELOG(*,*) 'nsoil,ngrid,firstcall=',nsoil,ngrid, .TRUE.
     47    WRITELOG(*,*) 'nsoil,firstcall=',nsoil, .TRUE.
    4548
    4649    ALLOCATE(dz1(nsoil),dz2(nsoil))
    … …  
    8285
    8386  END SUBROUTINE init_soil
     87
     88  PURE SUBROUTINE soil_backward(ngrid,nsoil, zc,zd, ptsrf,ptsoil)
     89    INTEGER, INTENT(IN) :: ngrid, nsoil         ! number of columns, of soil layers
     90    REAL, INTENT(IN)    :: zc(ngrid, nsoil), zd(ngrid, nsoil) ! LU factorization
     91    REAL, INTENT(IN)    :: ptsrf(ngrid)         ! new surface temperature
     92    REAL, INTENT(INOUT) :: ptsoil(ngrid,nsoil)  ! soil temperature
     93    INTEGER :: ig, jk
     94
     95    !-----------------------------------------------------------------------
     96    !  Computation of the soil temperatures using the Cgrd and Dgrd
     97    !  coefficient computed during the forward sweep
     98    !  -----------------------------------------------
     99   
     100    !  surface temperature => temperature in first soil layer
     101    DO ig=1,ngrid
     102       ptsoil(ig,1)=(lambda*zc(ig,1)+ptsrf(ig))/                   &
     103            &      (lambda*(1.-zd(ig,1))+1.)
     104    ENDDO
     105   
     106    !   other temperatures
     107    DO jk=1,nsoil-1
     108       DO ig=1,ngrid
     109          ptsoil(ig,jk+1)=zc(ig,jk)+zd(ig,jk)*ptsoil(ig,jk)
     110       ENDDO
     111    ENDDO
     112  END SUBROUTINE Soil_backward
     113
     114  PURE SUBROUTINE soil_forward(ngrid, nsoil, ptimestep, ptherm_i, ptsrf, ptsoil, &
     115       &                       zc, zd, pcapcal, pfluxgrd)
     116    INTEGER, INTENT(IN) :: ngrid, nsoil         ! number of columns, of soil layers
     117    REAL, INTENT(IN)    :: ptimestep,         & ! time step
     118         &                 ptherm_i(ngrid),   & ! thermal inertia ??
     119         &                 ptsrf(ngrid),      & ! surface temperature before heat conduction
     120         &                 ptsoil(ngrid, nsoil) ! soil temperature before heat conduction
     121    REAL, INTENT(OUT)   :: zc(ngrid,nsoil),   &
     122         &                 zd(ngrid, nsoil),  & ! LU factorization for backward sweep
     123         &                 pcapcal(ngrid),    & ! effective calorific capacity
     124         &                 pfluxgrd(ngrid)      ! conductive heat flux at the ground
     125    REAL :: z1, zdz2(ngrid)
     126    INTEGER :: jk, ig
     127
     128    !-----------------------------------------------------------------------
     129    !   Computation of the Cgrd and Dgrd coefficients the backward sweep :
     130    !   ---------------------------------------------------------------
     131   
     132    DO jk=1,nsoil
     133       zdz2(jk)=dz2(jk)/ptimestep
     134    ENDDO
     135   
     136    DO ig=1,ngrid
     137       z1=zdz2(nsoil)+dz1(nsoil-1)
     138       zc(ig,nsoil-1)=zdz2(nsoil)*ptsoil(ig,nsoil)/z1
     139       zd(ig,nsoil-1)=dz1(nsoil-1)/z1
     140    ENDDO
     141   
     142    DO jk=nsoil-1,2,-1
     143       DO ig=1,ngrid
     144          z1=1./(zdz2(jk)+dz1(jk-1)+dz1(jk)*(1.-zd(ig,jk)))
     145          zc(ig,jk-1)=                                                &
     146               &      (ptsoil(ig,jk)*zdz2(jk)+dz1(jk)*zc(ig,jk))*z1
     147          zd(ig,jk-1)=dz1(jk-1)*z1
     148       ENDDO
     149    ENDDO
     150   
     151    !-----------------------------------------------------------------------
     152    !   computation of the surface diffusive flux from ground and
     153    !   calorific capacity of the ground:
     154    !   ---------------------------------
     155   
     156    DO ig=1,ngrid
     157       pfluxgrd(ig)=ptherm_i(ig)*dz1(1)*                              &
     158            &   (zc(ig,1)+(zd(ig,1)-1.)*ptsoil(ig,1))
     159       z1=lambda*(1.-zd(ig,1))+1.
     160       pcapcal(ig)=ptherm_i(ig)*                                      &
     161            &   ptimestep*(zdz2(1)+(1.-zd(ig,1))*dz1(1))/z1
     162       pfluxgrd(ig)=pfluxgrd(ig)                                      &
     163            &   +pcapcal(ig)*(ptsoil(ig,1)*z1-lambda*zc(ig,1)-ptsrf(ig))   &
     164            &   /ptimestep
     165    ENDDO
     166  END SUBROUTINE soil_forward
     167
     168  SUBROUTINE soil_new(ngrid,nsoil,ptimestep,ptherm_i, ptsrf,ptsoil, pcapcal,pfluxgrd)
     169    INTEGER, INTENT(IN) :: ngrid, nsoil         ! number of columns, of soil layers
     170    REAL, INTENT(IN)    :: ptimestep,         & ! time step
     171         &                 ptherm_i(ngrid)      ! thermal inertia ??
     172    REAL, INTENT(INOUT) :: ptsrf(ngrid),      & ! surface temperature
     173         &                 ptsoil(ngrid,nsoil)  ! soil temperature
     174    REAL, INTENT(OUT)   :: pcapcal(ngrid),    & ! effective calorific capacity
     175         &                 pfluxgrd(ngrid)      ! conductive heat flux at the ground
     176    CALL soil_backward(ngrid,nsoil, zc,zd, ptsrf,ptsoil)
     177    CALL soil_forward(ngrid, nsoil, ptimestep, ptherm_i, ptsrf, ptsoil, &
     178            &            zc, zd, pcapcal, pfluxgrd)
     179  END SUBROUTINE soil_new
    84180
    85181  SUBROUTINE soil(ngrid,nsoil,firstcall,ptherm_i,          &
    … …  
    149245
    150246    IF (firstcall) THEN
    151        CALL init_soil(ngrid, nsoil)
     247       ! init_soil is now called by iniphyparam
     248       ! CALL init_soil(ngrid, nsoil)
    152249    ELSE
     250       IF(.FALSE.) THEN
     251          !-----------------------------------------------------------------------
     252          !   Computation of the soil temperatures using the Cgrd and Dgrd
     253          !  coefficient computed at the previous time-step:
     254          !  -----------------------------------------------
     255         
     256          !    surface temperature
     257          DO ig=1,ngrid
     258             ptsoil(ig,1)=(lambda*zc(ig,1)+ptsrf(ig))/                   &
     259                  &      (lambda*(1.-zd(ig,1))+1.)
     260          ENDDO
     261         
     262          !   other temperatures
     263          DO jk=1,nsoil-1
     264             DO ig=1,ngrid
     265                ptsoil(ig,jk+1)=zc(ig,jk)+zd(ig,jk)*ptsoil(ig,jk)
     266             ENDDO
     267          ENDDO
     268       ELSE
     269          CALL soil_backward(ngrid,nsoil, zc,zd, ptsrf,ptsoil)
     270       END IF
     271       
     272    ENDIF
     273
     274    IF(.FALSE.) THEN
    153275       !-----------------------------------------------------------------------
    154        !   Computation of the soil temperatures using the Cgrd and Dgrd
    155        !  coefficient computed at the previous time-step:
    156        !  -----------------------------------------------
    157 
    158        !    surface temperature
     276       !   Computation of the Cgrd and Dgrd coefficient for the next step:
     277       !   ---------------------------------------------------------------
     278       
     279       DO jk=1,nsoil
     280          zdz2(jk)=dz2(jk)/ptimestep
     281       ENDDO
     282       
    159283       DO ig=1,ngrid
    160           ptsoil(ig,1)=(lambda*zc(ig,1)+ptsrf(ig))/                   &
    161                &      (lambda*(1.-zd(ig,1))+1.)
    162        ENDDO
    163 
    164        !   other temperatures
    165        DO jk=1,nsoil-1
     284          z1(ig)=zdz2(nsoil)+dz1(nsoil-1)
     285          zc(ig,nsoil-1)=zdz2(nsoil)*ptsoil(ig,nsoil)/z1(ig)
     286          zd(ig,nsoil-1)=dz1(nsoil-1)/z1(ig)
     287       ENDDO
     288       
     289       DO jk=nsoil-1,2,-1
    166290          DO ig=1,ngrid
    167              ptsoil(ig,jk+1)=zc(ig,jk)+zd(ig,jk)*ptsoil(ig,jk)
     291             z1(ig)=1./(zdz2(jk)+dz1(jk-1)+dz1(jk)*(1.-zd(ig,jk)))
     292             zc(ig,jk-1)=                                                &
     293                  &      (ptsoil(ig,jk)*zdz2(jk)+dz1(jk)*zc(ig,jk))*z1(ig)
     294             zd(ig,jk-1)=dz1(jk-1)*z1(ig)
    168295          ENDDO
    169296       ENDDO
    170 
    171     ENDIF
    172 
    173     !-----------------------------------------------------------------------
    174     !   Computation of the Cgrd and Dgrd coefficient for the next step:
    175     !   ---------------------------------------------------------------
    176 
    177     DO jk=1,nsoil
    178        zdz2(jk)=dz2(jk)/ptimestep
    179     ENDDO
    180 
    181     DO ig=1,ngrid
    182        z1(ig)=zdz2(nsoil)+dz1(nsoil-1)
    183        zc(ig,nsoil-1)=zdz2(nsoil)*ptsoil(ig,nsoil)/z1(ig)
    184        zd(ig,nsoil-1)=dz1(nsoil-1)/z1(ig)
    185     ENDDO
    186 
    187     DO jk=nsoil-1,2,-1
     297       
     298       !-----------------------------------------------------------------------
     299       !   computation of the surface diffusive flux from ground and
     300       !   calorific capacity of the ground:
     301       !   ---------------------------------
     302       
    188303       DO ig=1,ngrid
    189           z1(ig)=1./(zdz2(jk)+dz1(jk-1)+dz1(jk)*(1.-zd(ig,jk)))
    190           zc(ig,jk-1)=                                                &
    191                &      (ptsoil(ig,jk)*zdz2(jk)+dz1(jk)*zc(ig,jk))*z1(ig)
    192           zd(ig,jk-1)=dz1(jk-1)*z1(ig)
    193        ENDDO
    194     ENDDO
    195 
    196     !-----------------------------------------------------------------------
    197     !   computation of the surface diffusive flux from ground and
    198     !   calorific capacity of the ground:
    199     !   ---------------------------------
    200 
    201     DO ig=1,ngrid
    202        pfluxgrd(ig)=ptherm_i(ig)*dz1(1)*                              &
    203             &   (zc(ig,1)+(zd(ig,1)-1.)*ptsoil(ig,1))
    204        z1(ig)=lambda*(1.-zd(ig,1))+1.
    205        pcapcal(ig)=ptherm_i(ig)*                                      &
    206             &   ptimestep*(zdz2(1)+(1.-zd(ig,1))*dz1(1))/z1(ig)
    207        pfluxgrd(ig)=pfluxgrd(ig)                                      &
    208             &   +pcapcal(ig)*(ptsoil(ig,1)*z1(ig)-lambda*zc(ig,1)-ptsrf(ig))   &
    209             &   /ptimestep
    210     ENDDO
    211 
     304          pfluxgrd(ig)=ptherm_i(ig)*dz1(1)*                              &
     305               &   (zc(ig,1)+(zd(ig,1)-1.)*ptsoil(ig,1))
     306          z1(ig)=lambda*(1.-zd(ig,1))+1.
     307          pcapcal(ig)=ptherm_i(ig)*                                      &
     308               &   ptimestep*(zdz2(1)+(1.-zd(ig,1))*dz1(1))/z1(ig)
     309          pfluxgrd(ig)=pfluxgrd(ig)                                      &
     310               &   +pcapcal(ig)*(ptsoil(ig,1)*z1(ig)-lambda*zc(ig,1)-ptsrf(ig))   &
     311               &   /ptimestep
     312       ENDDO
     313    ELSE
     314       CALL soil_forward(ngrid, nsoil, ptimestep, ptherm_i, ptsrf, ptsoil, &
     315            &            zc, zd, pcapcal, pfluxgrd)
     316    END IF
    212317  END SUBROUTINE soil
    213 
    214 END MODULE surface
     318     
     319   END MODULE surface
     320   
Note: See TracChangeset for help on using the changeset viewer.