Changeset 2898 for trunk/LMDZ.GENERIC/libf

Timestamp:

Feb 21, 2023, 9:08:49 AM (20 months ago)

Author:

emillour

Message:

Generic PCM:
Some cleanup in aerosol_mod.F90 : integrate "iniaerosol" in the module.
Also adapted kcm1d to follow up on recent code updates.
EM

Location:

trunk/LMDZ.GENERIC/libf/phystd

Files:

• aerosol_mod.F90 (modified) (1 diff)
• dyn1d/kcm1d.F90 (modified) (10 diffs)
• iniaerosol.F (deleted)
• physiq_mod.F90 (modified) (2 diffs)

trunk/LMDZ.GENERIC/libf/phystd/aerosol_mod.F90

@@ -2 +2 @@
 module aerosol_mod
 implicit none
-save
+
 !==================================================================
 
 !  aerosol indexes: these are initialized to be 0 if the
 !                 corresponding aerosol was not activated in callphys.def
-!                 -- otherwise a value is given in iniaerosol
-      integer :: iaero_co2 = 0 
-      integer :: iaero_h2o = 0
-      integer :: iaero_dust = 0
-      integer :: iaero_h2so4 = 0
-      logical :: noaero = .false.
+!                 -- otherwise a value is set via iniaerosol
+      integer, save, protected :: iaero_co2 = 0 
+      integer, save, protected :: iaero_h2o = 0
+      integer, save, protected :: iaero_dust = 0
+      integer, save, protected :: iaero_h2so4 = 0
+      logical, save, protected :: noaero = .false.
 !$OMP THREADPRIVATE(iaero_co2,iaero_h2o,iaero_dust,iaero_h2so4,noaero)
 
 ! two-layer simple aerosol model
-      integer :: iaero_back2lay = 0
+      integer, save, protected :: iaero_back2lay = 0
  ! NH3 cloud
-      integer :: iaero_nh3 = 0
+      integer, save, protected :: iaero_nh3 = 0
 ! N-layer aerosol model (replaces the 2-layer and hard-coded clouds)
-      integer,dimension(:),allocatable,save :: iaero_nlay
+      integer,dimension(:), allocatable, save, protected :: iaero_nlay
 ! Auroral aerosols
-      integer :: iaero_aurora = 0
+      integer, save, protected :: iaero_aurora = 0
 !$OMP THREADPRIVATE(iaero_back2lay,iaero_nh3,iaero_nlay,iaero_aurora)
 
 ! Generic aerosols
-      integer, dimension(:),allocatable,save :: iaero_generic
-      integer,dimension(:),allocatable,save :: i_rgcs_ice
+      integer, dimension(:), allocatable, save, protected :: iaero_generic
+      integer, dimension(:), allocatable, save, protected :: i_rgcs_ice
 !$OMP THREADPRIVATE(iaero_generic,i_rgcs_ice)
 
 ! Venus clouds
-      integer :: iaero_venus1 = 0
-      integer :: iaero_venus2 = 0
-      integer :: iaero_venus2p = 0
-      integer :: iaero_venus3 = 0
-      integer :: iaero_venusUV = 0
+      integer, save, protected :: iaero_venus1 = 0
+      integer, save, protected :: iaero_venus2 = 0
+      integer, save, protected :: iaero_venus2p = 0
+      integer, save, protected :: iaero_venus3 = 0
+      integer, save, protected :: iaero_venusUV = 0
 !$OMP THREADPRIVATE(iaero_venus1,iaero_venus2,iaero_venus2p)
 !$OMP THREADPRIVATE(iaero_venus3,iaero_venusUV)
 
+!==================================================================
+
+contains
+
+  SUBROUTINE iniaerosol
+
+  use mod_phys_lmdz_para, only : is_master
+  use radinc_h, only: naerkind
+  use tracer_h, only: n_rgcs, nqtot, is_rgcs
+  use callkeys_mod, only: aeroco2, aeroh2o, dusttau, aeroh2so4, &
+                          aeroback2lay, aeronh3, nlayaero, aeronlay, &
+                          aeroaurora, aerogeneric, &
+                          aerovenus1, aerovenus2, aerovenus2p, &
+                          aerovenus3, aerovenusUV
+
+  IMPLICIT NONE
+!=======================================================================
+!   subject:
+!   --------
+!   Initialization related to aerosols 
+!   (CO2 aerosols, dust, water, chemical species, ice...)   
+!
+!   author: Laura Kerber, S. Guerlet
+!   ------
+!        
+!=======================================================================
+
+  integer :: i, ia, iq
+
+  ! Special case, dyn. allocation for n-layer depending on callphys.def
+  IF(.NOT.ALLOCATED(iaero_nlay)) ALLOCATE(iaero_nlay(nlayaero))
+  iaero_nlay(:) = 0
+  ! Do the same for iaero_generic and i_rgcs_ice
+  IF (.not. allocated(iaero_generic)) allocate(iaero_generic(aerogeneric))
+  if (.not. allocated(i_rgcs_ice)) allocate(i_rgcs_ice(aerogeneric)) 
+
+  ! Init of i_rgcs_ice
+  i_rgcs_ice(:) =0
+  ia = 1
+  do iq=1,nqtot
+    if (is_rgcs(iq) .eq. 1) then 
+        i_rgcs_ice(ia)=iq
+        ia = ia+1
+     endif
+  enddo
+
+  iaero_generic(:)=0
+  ia=0
+  if (aeroco2) then
+     ia=ia+1
+     iaero_co2=ia
+  endif
+  if (is_master) write(*,*) '--- CO2 aerosol = ', iaero_co2
+
+  if (aeroh2o) then
+     ia=ia+1
+     iaero_h2o=ia
+      endif
+  if (is_master) write(*,*) '--- H2O aerosol = ', iaero_h2o
+
+  if (dusttau.gt.0) then
+     ia=ia+1
+     iaero_dust=ia
+  endif
+  if (is_master) write(*,*) '--- Dust aerosol = ', iaero_dust
+
+  if (aeroh2so4) then
+     ia=ia+1
+     iaero_h2so4=ia
+  endif
+  if (is_master) write(*,*) '--- H2SO4 aerosol = ', iaero_h2so4
       
-!==================================================================
+  if (aeroback2lay) then
+     ia=ia+1
+     iaero_back2lay=ia
+  endif
+  if (is_master) write(*,*) '--- Two-layer aerosol = ', iaero_back2lay
+
+  if (aeronh3) then
+     ia=ia+1
+     iaero_nh3=ia
+  endif
+  if (is_master) write(*,*) '--- NH3 Cloud = ', iaero_nh3
+
+  if (aeronlay) then
+     do i=1,nlayaero
+       ia=ia+1
+       iaero_nlay(i)=ia
+     enddo
+  endif
+  if (is_master) write(*,*) '--- N-layer aerosol = ', iaero_nlay
+
+  if (aeroaurora) then
+     ia=ia+1
+     iaero_aurora=ia
+  endif
+  if (is_master) write(*,*) '--- Auroral aerosols = ', iaero_aurora
+
+  if (aerovenus1) then
+     ia=ia+1
+     iaero_venus1=ia
+  endif
+  if (is_master) write(*,*) '--- Venus cloud, mode 1 aerosol = ', iaero_venus1
+
+  if (aerovenus2) then
+     ia=ia+1
+     iaero_venus2=ia
+  endif
+  if (is_master) write(*,*) '--- Venus cloud, mode 2 aerosol = ', iaero_venus2
+
+  if (aerovenus2p) then
+     ia=ia+1
+     iaero_venus2p=ia
+  endif
+  if (is_master) write(*,*) '--- Venus cloud, mode 2p aerosol = ', iaero_venus2p
+
+  if (aerovenus3) then
+     ia=ia+1
+     iaero_venus3=ia
+  endif
+  if (is_master) write(*,*) '--- Venus cloud, mode 3 aerosol = ', iaero_venus3
+
+  if (aerovenusUV) then
+     ia=ia+1
+     iaero_venusUV=ia
+  endif
+  if (is_master) write(*,*) '--- Venus cloud, UV absorber = ', iaero_venusUV
+
+  if (aerogeneric .ne. 0) then 
+     do i=1,aerogeneric
+        ia = ia+1
+        iaero_generic(i) = ia
+     enddo
+  endif 
+      
+  if (is_master) then
+    write(*,*)'--- Radiative Generic Condensable Species = ',iaero_generic
+
+    write(*,*) '=== Number of aerosols= ', ia
+  endif ! of is_master
+
+! For the zero aerosol case, we currently make a dummy co2 aerosol which is zero everywhere.
+! (See aeropacity.F90 for how this works). A better solution would be to turn off the 
+! aerosol machinery in the no aerosol case, but this would be complicated. LK
+
+  if (ia.eq.0) then  !For the zero aerosol case. 
+     ia = 1
+     noaero = .true.
+     iaero_co2=ia
+  endif
+
+  if (ia.ne.naerkind) then
+    if (is_master) then
+      print*, 'Aerosols counted not equal to naerkind'
+      print*, 'Compile with tag -s',ia,'to run'
+      print*, 'or change options in callphys.def'
+      print*, 'Abort in iniaerosol'
+    endif
+    call abort_physic("iniaerosl",'wrong number of aerosols',1)
+  endif ! of if (ia.ne.naerkind)
+
+  END SUBROUTINE iniaerosol
+
 end module aerosol_mod
 !==================================================================

trunk/LMDZ.GENERIC/libf/phystd/dyn1d/kcm1d.F90

@@ -10 +10 @@
   use callkeys_mod, only: pceil, tstrat, tracer, global1d
   use inifis_mod, only: inifis
+  use aerosol_mod, only: iniaerosol
+  use callcorrk_mod, only: callcorrk
   use comcstfi_mod
   use mod_grid_phy_lmdz, only : regular_lonlat
@@ -56 +58 @@
   real,allocatable :: vmr(:,:) ! tracer mixing ratio [mol/mol]
   real,allocatable :: qsurf(:)        ! tracer surface budget [kg/kg] ????
-  real psurf,psurf_n,tsurf      
-
-  real emis, albedo, albedo_equivalent
+  real psurf,psurf_n,tsurf(1)      
+
+  real emis(1), albedo(1), albedo_equivalent(1)
   real albedo_wv(1,L_NSPECTV)
 
@@ -65 +67 @@
   real dtsw(llm) ! heating rate (K/s) due to SW
   real dtlw(llm) ! heating rate (K/s) due to LW
-  real fluxsurf_lw   ! incident LW flux to surf (W/m2)
-  real fluxtop_lw    ! outgoing LW flux to space (W/m2)
-  real fluxsurf_sw   ! incident SW flux to surf (W/m2)
-  real fluxsurfabs_sw ! absorbed SW flux by the surf (W/m2)
-  real fluxabs_sw    ! SW flux absorbed by planet (W/m2)
-  real fluxtop_dn    ! incident top of atmosphere SW flux (W/m2)
+  real fluxsurf_lw(1)   ! incident LW flux to surf (W/m2)
+  real fluxtop_lw(1)    ! outgoing LW flux to space (W/m2)
+  real fluxsurf_sw(1)   ! incident SW flux to surf (W/m2)
+  real fluxsurfabs_sw(1) ! absorbed SW flux by the surf (W/m2)
+  real fluxabs_sw(1)    ! SW flux absorbed by planet (W/m2)
+  real fluxtop_dn(1)    ! incident top of atmosphere SW flux (W/m2)
 
   ! not used
-  real cloudfrac(llm)
-  real totcloudfrac
-  real tau_col
+  real cloudfrac(1,llm)
+  real totcloudfrac(1)
+  real tau_col(1)
 
   real dTstrat
@@ -81 +83 @@
   real OLR_nu(1,L_NSPECTI)
   real OSR_nu(1,L_NSPECTV)
+  real GSR_nu(1,L_NSPECTV)
+  real int_dtaui(1,llm,L_NSPECTI)
+  real int_dtauv(1,llm,L_NSPECTV)
   real Eatmtot
 
@@ -100 +105 @@
   pi=2.E+0*asin(1.E+0)
 
-  cloudfrac(:)  = 0.0
-  totcloudfrac  = 0.0
+  cloudfrac(1,:)  = 0.0
+  totcloudfrac(1)  = 0.0
 
 
@@ -254 +259 @@
   write(*,*) " albedo = ",albedo
   do iw=1,L_NSPECTV
-     albedo_wv(1,iw)=albedo
+     albedo_wv(1,iw)=albedo(1)
   enddo
 
@@ -360 +365 @@
          
                
-     call iniaerosol()
+     call iniaerosol
      
 
@@ -372 +377 @@
           dtlw,dtsw,fluxsurf_lw,fluxsurf_sw,              & 
           fluxsurfabs_sw,fluxtop_lw,                      &
-          fluxabs_sw,fluxtop_dn,OLR_nu,OSR_nu,tau_col,    &
-          cloudfrac,totcloudfrac,.false.,firstcall,lastcall)
+          fluxabs_sw,fluxtop_dn,OLR_nu,OSR_nu,GSR_nu,     &
+          int_dtaui,int_dtauv,                            &
+          tau_col,cloudfrac,totcloudfrac,                 &
+          .false.,firstcall,lastcall)
 
      !write(*,*) 'BASE 3'
@@ -397 +404 @@
      ! skin temperature (gray approx.) using analytic pure H2 expression
      !Tstrat  = (fluxabs_sw/(2*sigma))**0.25 ! skin temperature (gray approx.)
-     Tstrat  = (fluxtop_lw/(2*sigma))**0.25 ! skin temperature (gray approx.)
+     Tstrat  = (fluxtop_lw(1)/(2*sigma))**0.25 ! skin temperature (gray approx.)
      dTstrat = dTstrat-Tstrat
 
@@ -422 +429 @@
        tsurf,fract,dist_star,aerosol,muvar,                    &
        dtlw,dtsw,fluxsurf_lw,fluxsurf_sw,fluxsurfabs_sw,       &
-       fluxtop_lw, fluxabs_sw,fluxtop_dn,OLR_nu,OSR_nu,        &
-       tau_col,                                                &
-       cloudfrac,totcloudfrac,.false.,firstcall,lastcall)
+       fluxtop_lw, fluxabs_sw,fluxtop_dn,OLR_nu,OSR_nu,GSR_nu, &
+       int_dtaui,int_dtauv,                                    &
+       tau_col,cloudfrac,totcloudfrac,                         &
+       .false.,firstcall,lastcall)
 
 

trunk/LMDZ.GENERIC/libf/phystd/physiq_mod.F90

@@ -21 +21 @@
       use radcommon_h, only: sigma, glat, grav, BWNV, WNOI, DWNI, DWNV, WNOV
       use radii_mod, only: h2o_reffrad, co2_reffrad
-      use aerosol_mod, only: iaero_co2, iaero_h2o
+      use aerosol_mod, only: iniaerosol, iaero_co2, iaero_h2o
       use surfdat_h, only: phisfi, zmea, zstd, zsig, zgam, zthe, &
                            dryness
@@ -512 +512 @@
 !        Initialize aerosol indexes.
 !        ~~~~~~~~~~~~~~~~~~~~~~~~~~~
-         call iniaerosol()
+         call iniaerosol
 
 #ifdef CPP_XIOS