Ignore:
Timestamp:
Jan 24, 2013, 4:12:48 PM (12 years ago)
Author:
aslmd
Message:

24/01/2013 == AS + JL

A more robust way to refer to gas type.

  • Gas names with an arbitrary number of characters (<20) can be used This is good for C2H2, C2H6, H2SO4, C17H21NO4, etc... !!! Remember this must be compliant with Q.dat in corrk_data !!!
  • igas_... labels are assigned once for all in su_gases Then using igas_... everywhere instead of gnom (except for kcm stuff)
  • Users can still use e.g. H2_ but H2 also works
  • Simplified condense_cloud so that igas_CO2 is used directly
Location:
trunk/LMDZ.GENERIC/libf/phystd
Files:
9 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.GENERIC/libf/phystd/calc_cpp_mugaz.F90

    r868 r869  
    3939         else
    4040            ! all values at 300 K from Engineering Toolbox
    41             if(gnom(igas).eq.'CO2')then
     41            if(igas.eq.igas_CO2)then
    4242               !cpp_c   = cpp_c   + 0.744*gfrac(igas) ! @ ~210 K (better for Mars conditions)
    4343               cpp_c   = cpp_c   + 0.846*gfrac(igas)
    4444               mugaz_c = mugaz_c + 44.01*gfrac(igas)
    45             elseif(gnom(igas).eq.'N2_')then
     45            elseif(igas.eq.igas_N2)then
    4646               cpp_c   = cpp_c   + 1.040*gfrac(igas)
    4747               mugaz_c = mugaz_c + 28.01*gfrac(igas)
    48             elseif(gnom(igas).eq.'H2_')then
     48            elseif(igas.eq.igas_H2)then
    4949               cpp_c   = cpp_c   + 14.31*gfrac(igas)
    5050               mugaz_c = mugaz_c + 2.01*gfrac(igas)
    51             elseif(gnom(igas).eq.'He_')then
     51            elseif(igas.eq.igas_He)then
    5252               cpp_c   = cpp_c   + 5.19*gfrac(igas)
    5353               mugaz_c = mugaz_c + 4.003*gfrac(igas)
    54             elseif(gnom(igas).eq.'H2O')then
     54            elseif(igas.eq.igas_H2O)then
    5555               cpp_c   = cpp_c   + 1.864*gfrac(igas)
    5656               mugaz_c = mugaz_c + 18.02*gfrac(igas)
    57             elseif(gnom(igas).eq.'SO2')then
     57            elseif(igas.eq.igas_SO2)then
    5858               cpp_c   = cpp_c   + 0.64*gfrac(igas)
    5959               mugaz_c = mugaz_c + 64.066*gfrac(igas)
    60             elseif(gnom(igas).eq.'H2S')then
     60            elseif(igas.eq.igas_H2S)then
    6161               cpp_c   = cpp_c   + 1.003*gfrac(igas) ! from wikipedia...
    6262               mugaz_c = mugaz_c + 34.08*gfrac(igas)
    63             elseif(gnom(igas).eq.'CH4')then
     63            elseif(igas.eq.igas_CH4)then
    6464               cpp_c   = cpp_c   + 2.226*gfrac(igas)
    6565               mugaz_c = mugaz_c + 16.04*gfrac(igas)
    66             elseif(gnom(igas).eq.'NH3')then
     66            elseif(igas.eq.igas_NH3)then
    6767               cpp_c   = cpp_c   + 2.175*gfrac(igas)
    6868               mugaz_c = mugaz_c + 17.03*gfrac(igas)
    6969               print*,'WARNING, cpp for NH3 may be for liquid'
    70             elseif(gnom(igas).eq.'2H6')then
     70            elseif(igas.eq.igas_C2H6)then
    7171               ! C2H6 http://encyclopedia.airliquide.com/Encyclopedia.asp?GasID=28
    7272               cpp_c   = cpp_c   + 1.763*gfrac(igas)
    7373               mugaz_c = mugaz_c + 30.07*gfrac(igas)
    74             elseif(gnom(igas).eq.'2H2')then
     74            elseif(igas.eq.igas_C2H2)then
    7575               ! C2H2 http://encyclopedia.airliquide.com/Encyclopedia.asp?GasID=1
    7676               cpp_c   = cpp_c   + 1.575*gfrac(igas)

  • trunk/LMDZ.GENERIC/libf/phystd/calc_rayleigh.F90

    r861 r869  
    5454      do igas=1,ngasmx
    5555         if(igas.eq.vgas)then
    56             print*,'Ignoring ',gnom(igas),' in Rayleigh scattering '// &
     56            print*,'Ignoring ',trim(gnom(igas)),' in Rayleigh scattering '// &
    5757            'as it is variable.'
    5858         elseif(gfrac(igas).lt.5.e-2)then
    59             print*,'Ignoring ',gnom(igas),' in Rayleigh scattering '// &
     59            print*,'Ignoring ',trim(gnom(igas)),' in Rayleigh scattering '// &
    6060            'as its mixing ratio is less than 0.05.'
    6161            ! ignore variable gas in Rayleigh calculation
    … …  
    6363            tauconsti(igas) = 0.0
    6464         else
    65             if(gnom(igas).eq.'CO2')then
     65            if(igas.eq.igas_CO2) then
    6666               tauconsti(igas) = (8.7/g)*1.527*scalep/P0
    67             elseif(gnom(igas).eq.'N2_')then
     67            elseif(igas.eq.igas_N2)then
    6868               tauconsti(igas) = (9.81/g)*8.569E-3*scalep/(P0/93.0)
    69             elseif(gnom(igas).eq.'H2O')then
     69            elseif(igas.eq.igas_H2O)then
    7070               tauconsti(igas) = (10.0/g)*9.22E-3*scalep/101325.0
    71             elseif(gnom(igas).eq.'H2_')then
     71            elseif(igas.eq.igas_H2)then
    7272               tauconsti(igas) = (10.0/g)*0.0241*scalep/101325.0
    7373               !tauconsti(igas) = (10.0/g)*0.0487*scalep/(101325.0)
    7474               ! uses n=1.000132 from Optics, Fourth Edition.
    7575               ! was out by a factor 2!
    76             elseif(gnom(igas).eq.'He_')then
     76            elseif(igas.eq.igas_He)then
    7777               print*,'Helium not ready yet!'
    7878               tauconsti(igas) = 0.0
    … …  
    8484
    8585            if(gfrac(igas).eq.1.0)then
    86                print*,'Rayleigh scattering is for a pure ',gnom(igas),' atmosphere.'
     86               print*,'Rayleigh scattering is for a pure ',trim(gnom(igas)),' atmosphere.'
    8787               typeknown=.true.
    8888            endif
    … …  
    113113                  tauvari(igas)   = 0.0
    114114               else
    115                   if(gnom(igas).eq.'CO2')then
     115                  if(igas.eq.igas_CO2)then
    116116                     tauvari(igas) = (1.0+0.013/wl**2)/wl**4
    117                   elseif(gnom(igas).eq.'N2_')then
     117                  elseif(igas.eq.igas_N2)then
    118118                     tauvari(igas) = (1.0+0.0113/wl**2+0.00013/wl**4)/wl**4
    119                   elseif(gnom(igas).eq.'H2O')then
     119                  elseif(igas.eq.igas_H2O)then
    120120                     tauvari(igas) = 1.0/wl**4 ! to be improved...
    121                   elseif(gnom(igas).eq.'H2_')then
     121                  elseif(igas.eq.igas_H2)then
    122122                     tauvari(igas) = 1.0/wl**4
    123                   elseif(gnom(igas).eq.'He_')then
     123                  elseif(igas.eq.igas_He)then
    124124                     print*,'Helium not ready yet!'
    125125                     tauvari(igas) = 0.0

  • trunk/LMDZ.GENERIC/libf/phystd/condense_cloud.F90

    r858 r869  
    143143
    144144      integer igas
    145       integer,save :: igasco2=0
    146       character(len=3) :: gasname
    147145
    148146      real ppco2
    … …  
    169167         enddo
    170168         
    171          ! find CO2 gas
    172          do igas=1,ngasmx
    173             gasname=gnom(igas)
    174 !            gasname=noms(igas) ! was a bug
    175             if (gasname.eq."CO2") then
    176                igasco2=igas
    177             endif
    178          enddo
    179 
    180169        write(*,*) "condense_cloud: i_co2ice=",i_co2ice       
    181170
    … …  
    265254      DO l=1,nlayer
    266255         DO ig=1,ngrid
    267             ppco2=gfrac(igasco2)*pplay(ig,l)
     256            ppco2=gfrac(igas_CO2)*pplay(ig,l)
    268257            call get_tcond_co2(ppco2,ztcond(ig,l))
    269258            call get_tnuc_co2(ppco2,ztnuc(ig,l))
    … …  
    397386!     forecast of ground temperature ztsrf and frost temperature ztcondsol
    398387      DO ig=1,ngrid
    399          ppco2=gfrac(igasco2)*pplay(ig,1)
     388         ppco2=gfrac(igas_CO2)*pplay(ig,1)
    400389         call get_tcond_co2(ppco2,ztcondsol(ig))
    401390         

  • trunk/LMDZ.GENERIC/libf/phystd/cp_neutral.F90

    r716 r869  
    99
    1010
    11   ! this function has been disabled in gradients_kcm.F90 because it dont
     11  ! this function has been disabled in gradients_kcm.F90 because it doesnt
    1212  ! work if you have gaseous mixtures. need to decide whether to generalise
    1313  ! it or simply remove entirely...
    … …  
    1616  !cp_neutral = (22.26+5.981d-2*T-3.501d-5*T**2+7.469d-9*T**3)/m_n
    1717
    18   if(gnom(1).eq.'N2_')then
     18  if(trim(gnom(1)).eq.'N2_')then
    1919     cp_neutral = 1040.0
    20   elseif(gnom(1).eq.'H2_')then
     20  elseif(trim(gnom(1)).eq.'H2_')then
    2121     cp_neutral = 14310.0
    2222  else

  • trunk/LMDZ.GENERIC/libf/phystd/gases_h.F90

    r868 r869  
    1414      integer :: ngasmx
    1515      integer :: vgas
    16       character*3,allocatable,DIMENSION(:) :: gnom ! name of the gas
     16      character*20,allocatable,DIMENSION(:) :: gnom ! name of the gas
    1717      real,allocatable,DIMENSION(:) :: gfrac
    1818
    … …  
    2929      integer :: igas_CH4
    3030      integer :: igas_NH3
    31       integer :: igas_2H2
    32       integer :: igas_2H6
     31      integer :: igas_C2H2
     32      integer :: igas_C2H6
    3333
    3434      end module gases_h

  • trunk/LMDZ.GENERIC/libf/phystd/gradients_kcm.F90

    r716 r869  
    4646     endif
    4747
    48      if(gnom(ngasmx).eq.'H2O')then
     48     if(trim(gnom(ngasmx)).eq.'H2O')then
    4949
    5050        call psat_H2O(T-2d-1,psat_minus)
    … …  
    7272        endif
    7373
    74      elseif(gnom(ngasmx).eq.'NH3')then
     74     elseif(trim(gnom(ngasmx)).eq.'NH3')then
    7575
    7676        call psat_NH3(T-2d-1,psat_minus)
    … …  
    109109
    110110     cp_v=0.0
    111      if(gnom(ngasmx).eq.'H2O')then
     111     if(trim(gnom(ngasmx)).eq.'H2O')then
    112112        cp_v = (32.24+1.923d-3*T+1.055d-5*T**2-3.511d-9*T**3)/m_v
    113      elseif(gnom(ngasmx).eq.'NH3')then
     113     elseif(trim(gnom(ngasmx)).eq.'NH3')then
    114114        cp_v = 2.058d3
    115      elseif(gnom(ngasmx).eq.'CH4')then
     115     elseif(trim(gnom(ngasmx)).eq.'CH4')then
    116116        cp_v = 2.226d3
    117117     endif

  • trunk/LMDZ.GENERIC/libf/phystd/kcmprof_fn.F90

    r716 r869  
    9999     m_v   = 1.0
    100100     tcrit = 1000.0
    101   elseif(gnom(vgas).eq.'H2O')then
     101  elseif(trim(gnom(vgas)).eq.'H2O')then
    102102     m_v   = dble(mH2O/1000.)
    103103     tcrit = 6.47d2
    104   elseif(gnom(vgas).eq.'NH3')then
     104  elseif(trim(gnom(vgas)).eq.'NH3')then
    105105     m_v   = 17.031/1000.
    106106     tcrit = 4.06d2
    107   elseif(gnom(vgas).eq.'CH4')then
     107  elseif(trim(gnom(vgas)).eq.'CH4')then
    108108     m_v   = 16.04/1000.
    109109     tcrit = 1.91d2
    … …  
    120120  psat_v  = psat_max
    121121  if(vgas.gt.0)then
    122      if(gnom(vgas).eq.'H2O')then
     122     if(trim(gnom(vgas)).eq.'H2O')then
    123123        call Psat_H2O(tsurf,psat_v)
    124      elseif(gnom(vgas).eq.'NH3')then
     124     elseif(trim(gnom(vgas)).eq.'NH3')then
    125125        call Psat_NH3(tsurf,psat_v)
    126126     endif
    … …  
    229229
    230230     if(vgas.gt.0)then
    231      if(gnom(vgas).eq.'H2O')then
     231     if(trim(gnom(vgas)).eq.'H2O')then
    232232        call Psat_H2O(T(ilay+1),psat_v)
    233      elseif(gnom(vgas).eq.'NH3')then
     233     elseif(trim(gnom(vgas)).eq.'NH3')then
    234234        call Psat_NH3(T(ilay+1),psat_v)
    235235     endif
    … …  
    257257
    258258        if(vgas.gt.0)then
    259         if(gnom(vgas).eq.'H2O')then
     259        if(trim(gnom(vgas)).eq.'H2O')then
    260260           call Psat_H2O(T(ilay+1),psat_v)
    261         elseif(gnom(vgas).eq.'NH3')then
     261        elseif(trim(gnom(vgas)).eq.'NH3')then
    262262           call Psat_NH3(T(ilay+1),psat_v)
    263263        endif

  • trunk/LMDZ.GENERIC/libf/phystd/su_gases.F90

    r868 r869  
    6969     count=0
    7070     do igas=1,ngasmx
    71         if (gnom(igas).eq."H2_") then
     71        if (trim(gnom(igas)).eq."H2_" .or. trim(gnom(igas)).eq."H2") then
    7272           igas_H2=igas
    7373           count=count+1
    74         elseif (gnom(igas).eq."He_") then
     74        elseif (trim(gnom(igas)).eq."He_" .or. trim(gnom(igas)).eq."He") then
    7575           igas_He=igas
    7676           count=count+1
    77         elseif (gnom(igas).eq."H2O") then
     77        elseif (trim(gnom(igas)).eq."H2O") then
    7878           igas_H2O=igas
    7979           count=count+1
    80         elseif (gnom(igas).eq."CO2") then
     80        elseif (trim(gnom(igas)).eq."CO2") then
    8181           igas_CO2=igas
    8282           count=count+1
    83         elseif (gnom(igas).eq."CO_") then
     83        elseif (trim(gnom(igas)).eq."CO_" .or. trim(gnom(igas)).eq."CO") then
    8484           igas_CO=igas
    8585           count=count+1
    86         elseif (gnom(igas).eq."N2_") then
     86        elseif (trim(gnom(igas)).eq."N2_" .or. trim(gnom(igas)).eq."N2") then
    8787           igas_N2=igas
    8888           count=count+1
    89         elseif (gnom(igas).eq."O2_") then
     89        elseif (trim(gnom(igas)).eq."O2_" .or. trim(gnom(igas)).eq."O2") then
    9090           igas_O2=igas
    9191           count=count+1
    92         elseif (gnom(igas).eq."SO2") then
     92        elseif (trim(gnom(igas)).eq."SO2") then
    9393           igas_SO2=igas
    9494           count=count+1
    95         elseif (gnom(igas).eq."H2S") then
     95        elseif (trim(gnom(igas)).eq."H2S") then
    9696           igas_H2S=igas
    9797           count=count+1
    98         elseif (gnom(igas).eq."CH4") then
     98        elseif (trim(gnom(igas)).eq."CH4") then
    9999           igas_CH4=igas
    100100           count=count+1
    101         elseif (gnom(igas).eq."NH3") then
     101        elseif (trim(gnom(igas)).eq."NH3") then
    102102           igas_NH3=igas
    103103           count=count+1
    104         elseif (gnom(igas).eq."2H6") then
    105            ! this is ethane
    106            igas_2H6=igas
     104        elseif (trim(gnom(igas)).eq."C2H6") then
     105           igas_C2H6=igas
    107106           count=count+1
    108         elseif (gnom(igas).eq."2H2") then
    109            ! this is acetylene
    110            igas_2H2=igas
     107        elseif (trim(gnom(igas)).eq."C2H2") then
     108           igas_C2H2=igas
    111109           count=count+1
    112110        endif

  • trunk/LMDZ.GENERIC/libf/phystd/sugas_corrk.F90

    r868 r869  
    4646      ! ALLOCATABLE ARRAYS -- AS 12/2011
    4747      REAL*8, DIMENSION(:,:,:,:,:), ALLOCATABLE :: gasi8, gasv8
    48       character*3,allocatable,DIMENSION(:) :: gastype ! for check with gnom
     48      character*20,allocatable,DIMENSION(:) :: gastype ! for check with gnom
    4949
    5050      real*8 x, xi(4), yi(4), ans, p
    … …  
    119119      ! Check that gastype and gnom match
    120120      do igas=1,ngas
    121          print*,'Gas ',igas,' is ',gnom(igas)
    122          if(gnom(igas).ne.gastype(igas))then
    123             print*,'Name of a gas in radiative transfer data (',gastype(igas),') does not ', &
    124                  'match that in gases.def (',gnom(igas),'), exiting. You should compare ', &
     121         print*,'Gas ',igas,' is ',trim(gnom(igas))
     122         if (trim(gnom(igas)).ne.trim(gastype(igas))) then
     123            print*,'Name of a gas in radiative transfer data (',trim(gastype(igas)),') does not ', &
     124                 'match that in gases.def (',trim(gnom(igas)),'), exiting. You should compare ', &
    125125                 'gases.def with Q.dat in your radiative transfer directory.'
    126             !call abort
     126            call abort
    127127         endif
    128128      enddo
    … …  
    589589      do igas=1,ngasmx
    590590
    591          if(gnom(igas).eq.'N2_')then
     591         if (igas .eq. igas_N2) then
    592592
    593593            call interpolateN2N2(100.D+0,250.D+0,17500.D+0,testcont,.true.)
    594594
    595          elseif(gnom(igas).eq.'H2_')then
     595         elseif (igas .eq. igas_H2) then
    596596
    597597            ! first do self-induced absorption
    … …  
    599599            ! then cross-interactions with other gases
    600600            do jgas=1,ngasmx
    601                if(gnom(jgas).eq.'N2_')then
     601               if (jgas .eq. igas_N2) then
    602602                  call interpolateN2H2(592.D+0,278.15D+0,200000.D+0,10000.D+0,testcont,.true.)
    603                elseif(gnom(jgas).eq.'He_')then
     603               elseif (jgas .eq. igas_He) then
    604604                  call interpolateH2He(500.D+0,250.D+0,200000.D+0,10000.D+0,testcont,.true.)
    605605               endif
    606606            enddo
    607607
    608          elseif(gnom(igas).eq.'H2O')then
     608         elseif (igas .eq. igas_H2O) then
    609609
    610610            ! H2O is special
Note: See TracChangeset for help on using the changeset viewer.