Changeset 2297 for trunk/LMDZ.GENERIC/libf/phystd/inifis_mod.F90

Timestamp:

Apr 26, 2020, 7:45:46 PM (5 years ago)

Author:

jvatant

Message:

Add a generic n-layer aerosol scheme to replace the former buggy 2-layer scheme as well as the hard-coded NH3 cloud.

It can be called using 'aeronlay=.true.' in callphys.def, and set the number of layers (up to 4) with 'nlayaero'.
Then, the following parameters are read as arrays of size nlayaero in callphys.def (separated by blank space)

*aeronlay_tauref (Optical depth of aerosol layer at ref wavelenght)
*aeronlay_lamref (Ref wavelenght (m))
*aeronlay_choice (Choice for vertical profile - 1:tau follows atm scale height btwn top and bottom - 2:tau follows it own scale height)
*aeronlay_pbot (Bottom pressure (Pa))
*aeronlay_ptop (Top pressure (Pa) - useful only if choice=1)
*aeronlay_sclhght (Ratio of aerosol layer scale height / atmospheric scale height - useful only if choice=2 )
*aeronlay_size (Particle size (m))
*optprop_aeronlay_vis File for VIS opt properties.
*optprop_aeronlay_ir File for IR opt properties.

+Extra info :

+ In addition of solving the bug from 2-layer it enables different optical properties.
+ The former scheme are left for retrocompatibility (for now) but you should use the new one.
+ See aeropacity.F90 for the calculations

+ Each layer can have different optical properties, size of particle ...
+ You have different choices for vertical profile of the aerosol layers :

• aeronlay_choice = 1 : Layer tau is spread between ptop and pbot following atm scale height.
• aeronlay_choice = 2 : Layer tau follows its own scale height above cloud deck (pbot).

In this case ptop is dummy and sclhght gives the ratio H_cl/H_atm.

+ The reference wavelenght for input optical depth is now read as input (aeronlay_lamref)
+ Following the last point some comment is added in suaer_corrk about the 'not-really-dummy'ness of IR lamref..

File:

• trunk/LMDZ.GENERIC/libf/phystd/inifis_mod.F90 (modified) (5 diffs)

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

@@ -554 +554 @@
      call getin_p("aeroh2so4",aeroh2so4)
      write(*,*)" aeroh2so4 = ",aeroh2so4
-         
+ 
+     write(*,*)"Radiatively active auroral aerosols?"
+     aeroaurora=.false.     ! default value
+     call getin_p("aeroaurora",aeroaurora)
+     write(*,*)" aeroaurora = ",aeroaurora
+
 !=================================
+! TWOLAY scheme and NH3 cloudare left for retrocompatibility only,
+! You should now use N-LAYER scheme (see below).
 
      write(*,*)"Radiatively active two-layer aerosols?"
@@ -562 +569 @@
      write(*,*)" aeroback2lay = ",aeroback2lay
 
+     if (aeroback2lay) then
+       print*,'Warning : The TWOLAY AEROSOL scheme is deprecated and buggy...'
+       print*,'You should use the generic n-layer scheme (see aeronlay).'
+     endif
+
      write(*,*)"Radiatively active ammonia cloud?"
      aeronh3=.false.     ! default value
@@ -567 +579 @@
      write(*,*)" aeronh3 = ",aeronh3
 
-     write(*,*)"Radiatively active auroral aerosols?"
-     aeroaurora=.false.     ! default value
-     call getin_p("aeroaurora",aeroaurora)
-     write(*,*)" aeroaurora = ",aeroaurora
+     if (aeronh3) then
+       print*,'Warning : You are using specific NH3 cloud scheme ...'
+       print*,'You should use the generic n-layer scheme (see aeronlay).'
+     endif
 
      write(*,*)"TWOLAY AEROSOL: total optical depth ", &
@@ -613 +625 @@
        print*,'Error : TWOLAY AEROSOL, Please ensure that in callphys.def'
        print*,'you have pres_top_tropo > pres_bottom_strato !'
-       stop
+       call abort
      endif
 
@@ -647 +659 @@
      call getin_p("size_nh3_cloud",size_nh3_cloud)
      write(*,*)" size_nh3_cloud = ",size_nh3_cloud
+
+!=================================
+! Generic N-LAYER aerosol scheme
+
+     write(*,*)"Radiatively active generic n-layer aerosols?"
+     aeronlay=.false.     ! default value
+     call getin_p("aeronlay",aeronlay)
+     write(*,*)" aeronlay = ",aeronlay
+
+     write(*,*)"Number of generic aerosols layers?"
+     nlayaero=1     ! default value
+     call getin_p("nlayaero",nlayaero)
+     ! Avoid to allocate arrays of size 0
+     if (aeronlay .and. nlayaero.lt.1) then
+       print*, " You are trying to set no generic aerosols..."
+       print*, " Set aeronlay=.false. instead ! I abort."
+       call abort
+     endif
+     if (.not. aeronlay) nlayaero=1
+     write(*,*)" nlayaero = ",nlayaero
+
+     ! This is necessary, we just set the number of aerosol layers
+     IF(.NOT.ALLOCATED(aeronlay_tauref))      ALLOCATE(aeronlay_tauref(nlayaero))     
+     IF(.NOT.ALLOCATED(aeronlay_lamref))      ALLOCATE(aeronlay_lamref(nlayaero))     
+     IF(.NOT.ALLOCATED(aeronlay_choice))      ALLOCATE(aeronlay_choice(nlayaero))     
+     IF(.NOT.ALLOCATED(aeronlay_pbot))        ALLOCATE(aeronlay_pbot(nlayaero))     
+     IF(.NOT.ALLOCATED(aeronlay_ptop))        ALLOCATE(aeronlay_ptop(nlayaero))     
+     IF(.NOT.ALLOCATED(aeronlay_sclhght))     ALLOCATE(aeronlay_sclhght(nlayaero))     
+     IF(.NOT.ALLOCATED(aeronlay_size))        ALLOCATE(aeronlay_size(nlayaero))     
+     IF(.NOT.ALLOCATED(optprop_aeronlay_ir))  ALLOCATE(optprop_aeronlay_ir(nlayaero))     
+     IF(.NOT.ALLOCATED(optprop_aeronlay_vis)) ALLOCATE(optprop_aeronlay_vis(nlayaero))     
+
+     write(*,*)"Generic n-layer aerosols: Optical depth at reference wavelenght"
+     aeronlay_tauref=1.0E-1
+     call getin_p("aeronlay_tauref",aeronlay_tauref)
+     write(*,*)" aeronlay_tauref = ",aeronlay_tauref
+
+     write(*,*)"Generic n-layer aerosols: Reference wavelenght for optical depths (m)"
+     aeronlay_lamref=0.6E-6
+     call getin_p("aeronlay_lamref",aeronlay_lamref)
+     write(*,*)" aeronlay_lamref = ",aeronlay_lamref
+
+     write(*,*)"Generic n-layer aerosols: Vertical profile choice : &
+                     (1) Tau btwn ptop and pbot follows atm. scale height & 
+                 or  (2) Tau above pbot follows its own scale height"
+     aeronlay_choice=1
+     call getin_p("aeronlay_choice",aeronlay_choice)
+     write(*,*)" aeronlay_choice = ",aeronlay_choice
+
+     write(*,*)"Generic n-layer aerosols: bottom pressures (Pa)"
+     aeronlay_pbot=2000.0
+     call getin_p("aeronlay_pbot",aeronlay_pbot)
+     write(*,*)" aeronlay_pbot = ",aeronlay_pbot
+
+     write(*,*)"Generic n-layer aerosols: (if choice=1) Top pressures (Pa) "
+     aeronlay_ptop=300000.0
+     call getin_p("aeronlay_ptop",aeronlay_ptop)
+     write(*,*)" aeronlay_ptop = ",aeronlay_ptop
+
+     write(*,*)"Generic n-layer aerosols: (if choice=2) Scale height / atm. scale height"
+     aeronlay_ptop=0.2
+     call getin_p("aeronlay_sclhght",aeronlay_sclhght)
+     write(*,*)" aeronlay_sclhght = ",aeronlay_sclhght
+
+     write(*,*)"Generic n-layer aerosols: particles sizes (m)"
+     aeronlay_size=1.e-6
+     call getin_p("aeronlay_size",aeronlay_size)
+     write(*,*)" aeronlay_size = ",aeronlay_size
+
+     write(*,*)"Generic n-layer aerosols: VIS optical properties file"
+     optprop_aeronlay_vis = 'optprop_saturn_vis_n20.dat'
+     call getin_p("optprop_aeronlay_vis",optprop_aeronlay_vis)
+     write(*,*)" optprop_aeronlay_vis = ",optprop_aeronlay_vis
+
+     write(*,*)"Generic n-layer aerosols: IR optical properties file"
+     optprop_aeronlay_ir = 'optprop_saturn_ir_n20.dat'
+     call getin_p("optprop_aeronlay_ir",optprop_aeronlay_ir)
+     write(*,*)" optprop_aeronlay_ir = ",optprop_aeronlay_ir
+     
 
 !=================================