Changeset 2245

Timestamp:

Feb 24, 2020, 6:55:50 PM (5 years ago)

Author:

jvatant

Message:

+ Add a 'versH2H2cia' int key in callphys that allows two values (2011 or 2018) to

deal with updated HITRAN file (for interpolateH2H2.F90) from 2018 that includes the
H2H2 dimer from Fletcher et al. 2018, useful for giant planets.
Retrocompatibility is ok, default value to 2011.

--JVO

Location:

trunk

Files:

• LMDZ.GENERIC/README (modified) (1 diff)
• LMDZ.GENERIC/libf/phystd/callkeys_mod.F90 (modified) (1 diff)
• LMDZ.GENERIC/libf/phystd/inifis_mod.F90 (modified) (1 diff)
• LMDZ.GENERIC/libf/phystd/interpolateH2H2.F90 (modified) (5 diffs)
• LMDZ.TITAN/README (modified) (1 diff)
• LMDZ.TITAN/libf/phytitan/callkeys_mod.F90 (modified) (1 diff)
• LMDZ.TITAN/libf/phytitan/inifis_mod.F90 (modified) (1 diff)
• LMDZ.TITAN/libf/phytitan/interpolateH2H2.F90 (modified) (5 diffs)

trunk/LMDZ.GENERIC/README

@@ -1511 +1511 @@
 - Get rid of the old 'sparadrap' in iniorbit.F assuming Mkm instead of AU if periastre or apostre gt 25, (otherwise Neptune is a pb)
   Everybody should be using AU, that's all.
+
+== 24/02/2020 (r2244) == JVO
++ Add a 'versH2H2cia' int key in callphys that allows two values (2011 or 2018) to
+  deal with updated HITRAN file (for interpolateH2H2.F90) from 2018 that includes the
+  H2H2 dimer from Fletcher et al. 2018, useful for giant planets.
+  Retrocompatibility is ok, default value to 2011.

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

@@ -62 +62 @@
 !$OMP THREADPRIVATE(photochem)
 
+      integer,save :: versH2H2cia
       integer,save :: iddist
       integer,save :: iaervar
       integer,save :: iradia
       integer,save :: startype
-!$OMP THREADPRIVATE(iddist,iaervar,iradia,startype)
+!$OMP THREADPRIVATE(versH2H2cia,iddist,iaervar,iradia,startype)
 
       real,save :: topdustref

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

@@ -250 +250 @@
      write(*,*) " H2Ocont_simple = ",H2Ocont_simple
  
+     write(*,*) "version for H2H2 CIA file ?"
+     versH2H2cia=2011 ! default value (should be 2018 but retrocompatibility first)
+     call getin_p("versH2H2cia",versH2H2cia)
+     write(*,*) " versH2H2cia = ",versH2H2cia
+     ! Sanity check
+     if (versH2H2cia.ne.2011 .and. versH2H2cia.ne.2018) then
+        print*,'Error: Choose a correct value (2011 or 2018) for versH2H2cia !'
+        call abort
+     endif 
+
      write(*,*) "call turbulent vertical diffusion ?"
      calldifv=.true. ! default value

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

@@ -6 +6 @@
 !     -------
 !     Calculates the H2-H2 CIA opacity, using a lookup table from
-!     HITRAN (2011)
+!     HITRAN (2011 or later)
 !
 !     Authors
 !     -------
 !     R. Wordsworth (2011)
-!     
+!
+!     + J.Vatant d'Ollone (2019) 
+!        - Enable updated HITRAN file (Karman2019,Fletcher2018)
+!           (2018 one should be default for giant planets)
 !==================================================================
 
+      use callkeys_mod, only: versH2H2cia
       use datafile_mod, only: datadir
 
@@ -27 +31 @@
 
       integer nS,nT
-      parameter(nS=2428)
       parameter(nT=10)
 
@@ -36 +39 @@
 
       double precision amagat
-      double precision wn_arr(nS)
       double precision temp_arr(nT)
-      double precision abs_arr(nS,nT)
+      
+      double precision, dimension(:),   allocatable :: wn_arr
+      double precision, dimension(:,:), allocatable :: abs_arr
 
       integer k,iT
       logical firstcall
 
-      save wn_arr, temp_arr, abs_arr !read by master
+      save nS, wn_arr, temp_arr, abs_arr !read by master
 
       character*100 dt_file
-      integer strlen,ios
+      integer ios
 
-      character(LEN=*), parameter :: fmat1 = "(A20,F10.3,F10.3,I7,F7.1,E10.3,F5.3)"
+      character(LEN=*), parameter :: fmat11 = "(A20,F10.3,F10.3,I7,F7.1,E10.3,F5.3)"
+      character(LEN=*), parameter :: fmat18 = "(A12,A3,A5,F10.6,F10.4,I7,F7.3,E10.3,F5.3)"
 
       character*20 bleh
@@ -69 +74 @@
          print*,'Initialising H2-H2 continuum from HITRAN database...'
 
-!     1.1 Open the ASCII files
-         dt_file=TRIM(datadir)//'/continuum_data/H2-H2_norm_2011.cia'
+!     1.1 Open the ASCII files and set nS according to version
+         ! Only two possible versions for now : 2011 or 2018 (sanity check in inifis_mod)
+         if (versH2H2cia.eq.2011) then
+           dt_file=TRIM(datadir)//'/continuum_data/H2-H2_norm_2011.cia'
+           nS = 2428
+         else if (versH2H2cia.eq.2018) then
+           dt_file=TRIM(datadir)//'/continuum_data/H2-H2_norm_2018.cia'
+           nS = 9600
+         endif
+
+         if(.not.allocated(wn_arr))  allocate(wn_arr(nS)) 
+         if(.not.allocated(abs_arr)) allocate(abs_arr(nS,nT)) 
 
 !$OMP MASTER
@@ -80 +95 @@
            write(*,*) 'is correct. You can change it in callphys.def with:'
            write(*,*) 'datadir = /absolute/path/to/datagcm'
-           write(*,*) 'Also check that the continuum data continuum_data/H2-H2_norm_2011.cia is there.'
+           write(*,*) 'Also check that the continuum data continuum_data/H2-H2_norm_2011.cia or H2-H2_norm_2018.cia is there.'
            call abort
          else
+         
+         if(versH2H2cia.eq.2011) then
+           write(*,*) '... You are using H2-H2 CIA from 2011 but you should use more recent data available on HITRAN !'
+           write(*,*) '... (Especially if you are running a giant planet atmosphere)'
+           write(*,*) '... Just find out the H2-H2_norm_2018.cia, put it in your datadir and have a look at interpolateH2H2.F90 ! .'
+         endif
 
             do iT=1,nT
+               
+               ! Only two possibles values for now : 2011 or 2018 (sanity check in inifis_mod)
+               if(versH2H2cia.eq.2011) then
+                 read(33,fmat11) bleh,blah,blah,nres,Ttemp
+               else if (versH2H2cia.eq.2018) then
+                 read(33,fmat18) bleh,bleh,bleh,blah,blah,nres,Ttemp
+               endif
 
-               read(33,fmat1) bleh,blah,blah,nres,Ttemp
                if(nS.ne.nres)then
                   print*,'Resolution given in file: ',trim(dt_file)

trunk/LMDZ.TITAN/README

@@ -1576 +1576 @@
 - Given theses pbs with the radii, also update the optics so that it computes the opacities in 
 a simpler way, directly for M3, through look-up tables, M3 being a good proxy.
+
+== 24/02/2020 (r2244) == JVO
++ Add a 'versH2H2cia' int key in callphys that allows two values (2011 or 2018) to
+  deal with updated HITRAN file (for interpolateH2H2.F90) from 2018 that includes the
+  H2H2 dimer from Fletcher et al. 2018, useful for giant planets.
+  Retrocompatibility is ok, default value to 2011.
+

trunk/LMDZ.TITAN/libf/phytitan/callkeys_mod.F90

@@ -48 +48 @@
       integer,save :: ichim
 !$OMP THREADPRIVATE(ichim)
+      integer,save :: versH2H2cia
       integer,save :: iddist
       integer,save :: iradia
       integer,save :: startype
-!$OMP THREADPRIVATE(iddist,iradia,startype)
+!$OMP THREADPRIVATE(versH2H2cia,iddist,iradia,startype)
       
       real,save :: p_prod, tx_prod, rc_prod

trunk/LMDZ.TITAN/libf/phytitan/inifis_mod.F90

@@ -288 +288 @@
      write(*,*) " continuum = ",continuum
  
+     write(*,*) "version for H2H2 CIA file ?"
+     versH2H2cia=2011 ! default value (should be 2018 but retrocompatibility first)
+     call getin_p("versH2H2cia",versH2H2cia)
+     write(*,*) " versH2H2cia = ",versH2H2cia
+     ! Sanity check
+     if (versH2H2cia.ne.2011 .and. versH2H2cia.ne.2018) then
+        print*,'Error: Choose a correct value (2011 or 2018) for versH2H2cia !'
+        call abort
+     endif 
+
      write(*,*) "call turbulent vertical diffusion ?"
      calldifv=.true. ! default value

trunk/LMDZ.TITAN/libf/phytitan/interpolateH2H2.F90

@@ -6 +6 @@
 !     -------
 !     Calculates the H2-H2 CIA opacity, using a lookup table from
-!     HITRAN (2011)
+!     HITRAN (2011 or later)
 !
 !     Authors
 !     -------
 !     R. Wordsworth (2011)
-!     
+!
+!     + J.Vatant d'Ollone (2019) 
+!        - Enable updated HITRAN file (Karman2019,Fletcher2018)
+!           (2018 one should be default for giant planets)
 !==================================================================
 
+      use callkeys_mod, only: versH2H2cia
       use datafile_mod, only: datadir
 
@@ -27 +31 @@
 
       integer nS,nT
-      parameter(nS=2428)
       parameter(nT=10)
 
@@ -36 +39 @@
 
       double precision amagat
-      double precision wn_arr(nS)
       double precision temp_arr(nT)
-      double precision abs_arr(nS,nT)
+      
+      double precision, dimension(:),   allocatable :: wn_arr
+      double precision, dimension(:,:), allocatable :: abs_arr
 
       integer k,iT
       logical firstcall
 
-      save wn_arr, temp_arr, abs_arr !read by master
+      save nS, wn_arr, temp_arr, abs_arr !read by master
 
       character*100 dt_file
-      integer strlen,ios
+      integer ios
 
-      character(LEN=*), parameter :: fmat1 = "(A20,F10.3,F10.3,I7,F7.1,E10.3,F5.3)"
+      character(LEN=*), parameter :: fmat11 = "(A20,F10.3,F10.3,I7,F7.1,E10.3,F5.3)"
+      character(LEN=*), parameter :: fmat18 = "(A12,A3,A5,F10.6,F10.4,I7,F7.3,E10.3,F5.3)"
 
       character*20 bleh
@@ -69 +74 @@
          print*,'Initialising H2-H2 continuum from HITRAN database...'
 
-!     1.1 Open the ASCII files
-         dt_file=TRIM(datadir)//'/continuum_data/H2-H2_norm_2011.cia'
+!     1.1 Open the ASCII files and set nS according to version
+         ! Only two possible versions for now : 2011 or 2018 (sanity check in inifis_mod)
+         if (versH2H2cia.eq.2011) then
+           dt_file=TRIM(datadir)//'/continuum_data/H2-H2_norm_2011.cia'
+           nS = 2428
+         else if (versH2H2cia.eq.2018) then
+           dt_file=TRIM(datadir)//'/continuum_data/H2-H2_norm_2018.cia'
+           nS = 9600
+         endif
+
+         if(.not.allocated(wn_arr))  allocate(wn_arr(nS)) 
+         if(.not.allocated(abs_arr)) allocate(abs_arr(nS,nT)) 
 
 !$OMP MASTER
@@ -80 +95 @@
            write(*,*) 'is correct. You can change it in callphys.def with:'
            write(*,*) 'datadir = /absolute/path/to/datagcm'
-           write(*,*) 'Also check that the continuum data continuum_data/H2-H2_norm_2011.cia is there.'
+           write(*,*) 'Also check that the continuum data continuum_data/H2-H2_norm_2011.cia or H2-H2_norm_2018.cia is there.'
            call abort
          else
+         
+         if(versH2H2cia.eq.2011) then
+           write(*,*) '... You are using H2-H2 CIA from 2011 but you should use more recent data available on HITRAN !'
+           write(*,*) '... (Especially if you are running a giant planet atmosphere)'
+           write(*,*) '... Just find out the H2-H2_norm_2018.cia, put it in your datadir and have a look at interpolateH2H2.F90 ! .'
+         endif
 
             do iT=1,nT
+               
+               ! Only two possibles values for now : 2011 or 2018 (sanity check in inifis_mod)
+               if(versH2H2cia.eq.2011) then
+                 read(33,fmat11) bleh,blah,blah,nres,Ttemp
+               else if (versH2H2cia.eq.2018) then
+                 read(33,fmat18) bleh,bleh,bleh,blah,blah,nres,Ttemp
+               endif
 
-               read(33,fmat1) bleh,blah,blah,nres,Ttemp
                if(nS.ne.nres)then
                   print*,'Resolution given in file: ',trim(dt_file)