Changeset 4055 for trunk/LMDZ.GENERIC/libf/phygeneric/optci.F90

Timestamp:

Feb 6, 2026, 10:20:21 PM (3 weeks ago)

Author:

mturbet

Message:

Cleaning of continuum routines following 04-02-26 hackathon. RIP our beloved bilinear and bilinearbig routines. You will be deeply missed...

File:

• trunk/LMDZ.GENERIC/libf/phygeneric/optci.F90 (modified) (5 diffs)

trunk/LMDZ.GENERIC/libf/phygeneric/optci.F90

@@ -11 +11 @@
   use radinc_h, only: L_LEVELS, L_NLAYRAD, L_NSPECTI, L_NGAUSS, &
                       L_NLEVRAD, L_REFVAR, naerkind
-  use radcommon_h, only: gasi,tlimit,wrefVAR,Cmk,tgasref,pfgasref,wnoi,scalep,indi,glat_ig
+  use radcommon_h, only: gasi,tlimit,wrefVAR,Cmk,tgasref,pfgasref,wnoi,scalep,glat_ig
   use gases_h, only: gfrac, ngasmx, igas_N2, igas_He, igas_H2O, igas_H2, &
                      igas_CH4, igas_CO2, igas_O2
   use comcstfi_mod, only: g, r, mugaz
-  use callkeys_mod, only: kastprof,continuum,graybody,varspec, &
-                          generic_continuum_database
+  use callkeys_mod, only: kastprof,continuum,graybody,varspec
   use recombin_corrk_mod, only: corrk_recombin, gasi_recomb
   use tpindex_mod, only: tpindex
@@ -98 +97 @@
 
   integer igas, jgas
-
-  integer interm
   
   logical :: firstcall=.true.
@@ -126 +123 @@
   endif ! of if (firstcall)
 
-  !! AS: to save time in computing continuum (see bilinearbig)
-  IF (.not.ALLOCATED(indi)) THEN
-      ALLOCATE(indi(L_NSPECTI,ngasmx,ngasmx))
-      indi = -9999  ! this initial value means "to be calculated"
-  ENDIF
-
   !=======================================================================
   !     Determine the total gas opacity throughout the column, for each
@@ -191 +182 @@
            ! include continua if necessary
            
-           if(generic_continuum_database)then
              T_cont  = dble(TMID(k))
              do igas=1,ngasmx
@@ -241 +231 @@
                
              enddo ! igas=1,ngasmx
-             
-           else ! generic_continuum_database
-           
-            wn_cont = dble(wnoi(nw))
-            T_cont  = dble(TMID(k))
-            do igas=1,ngasmx
-
-              if(gfrac(igas).eq.-1)then ! variable
-                 p_cont  = dble(PMID(k)*scalep*QVAR(k)) ! qvar = mol/mol
-              elseif(varspec) then
-                 p_cont  = dble(PMID(k)*scalep*FRACVAR(igas,k)*(1.-QVAR(k)))
-              else
-                 p_cont  = dble(PMID(k)*scalep*gfrac(igas)*(1.-QVAR(k)))
-              endif
-
-              dtemp=0.0d0
-              if(igas.eq.igas_N2)then
-
-                 interm = indi(nw,igas,igas)
-                 call interpolateN2N2(wn_cont,T_cont,p_cont,dtemp,.false.,interm)
-                 indi(nw,igas,igas) = interm
-
-              elseif(igas.eq.igas_H2)then
-
-                 ! first do self-induced absorption
-                 interm = indi(nw,igas,igas)
-                 call interpolateH2H2(wn_cont,T_cont,p_cont,dtemp,.false.,interm)
-                 indi(nw,igas,igas) = interm
-
-                 ! then cross-interactions with other gases
-                 do jgas=1,ngasmx
-                    if(varspec) then
-                      p_cross = dble(PMID(k)*scalep*FRACVAR(jgas,k)*(1.-QVAR(k)))
-                    else
-                      p_cross = dble(PMID(k)*scalep*gfrac(jgas)*(1.-QVAR(k)))
-                    endif
-                    dtempc  = 0.0d0
-                    if(jgas.eq.igas_N2)then 
-                       interm = indi(nw,igas,jgas)
-                       call interpolateN2H2(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)
-                       indi(nw,igas,jgas) = interm
-                    elseif(jgas.eq.igas_CO2)then 
-                       interm = indi(nw,igas,jgas)
-                       call interpolateCO2H2(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)
-                       indi(nw,igas,jgas) = interm
-                    elseif(jgas.eq.igas_He)then 
-                       interm = indi(nw,igas,jgas)
-                       call interpolateH2He(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)
-                       indi(nw,igas,jgas) = interm
-                    endif
-                    dtemp = dtemp + dtempc
-                 enddo
-
-              elseif(igas.eq.igas_CH4)then
-
-                 ! first do self-induced absorption
-                 interm = indi(nw,igas,igas)
-                 call interpolateCH4CH4(wn_cont,T_cont,p_cont,dtemp,.false.,interm)
-                 indi(nw,igas,igas) = interm
-
-                 ! then cross-interactions with other gases
-                 do jgas=1,ngasmx
-                    if(varspec) then
-                      p_cross = dble(PMID(k)*scalep*FRACVAR(jgas,k)*(1.-QVAR(k)))
-                    else
-                      p_cross = dble(PMID(k)*scalep*gfrac(jgas)*(1.-QVAR(k)))
-                    endif
-                    dtempc  = 0.0d0
-                    if(jgas.eq.igas_H2)then 
-                       interm = indi(nw,igas,jgas)
-                       call interpolateH2CH4(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)
-                       indi(nw,igas,jgas) = interm
-                    elseif(jgas.eq.igas_CO2)then 
-                       interm = indi(nw,igas,jgas)
-                       call interpolateCO2CH4(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)
-                       indi(nw,igas,jgas) = interm
-                    elseif(jgas.eq.igas_He)then 
-                       interm = indi(nw,igas,jgas)
-                       call interpolateHeCH4(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)
-                       indi(nw,igas,jgas) = interm
-                    endif
-                    dtemp = dtemp + dtempc
-                 enddo
-
-              elseif(igas.eq.igas_H2O.and.T_cont.gt.100.0)then
-                 ! Compute self and foreign (with air) continuum of H2O 
-                 p_air = dble(PMID(k)*scalep) - p_cont ! note assumes background is air!
-                 interm = indi(nw,igas,igas)
-                 call interpolateH2O_self_foreign(wn_cont,T_cont,p_cont,p_air,dtemp,.false.,interm) ! MTCKD v3.3
-                 indi(nw,igas,igas) = interm
-
-              endif
-
-              DCONT = DCONT + dtemp
-
-           enddo
-
-           ! Oobleck test
-           !rho = PMID(k)*scalep / (TMID(k)*286.99)
-           !if(WNOI(nw).gt.300.0 .and. WNOI(nw).lt.500.0)then
-           !   DCONT = rho * 0.125 * 4.6e-4
-           !elseif(WNOI(nw).gt.500.0 .and. WNOI(nw).lt.700.0)then
-           !   DCONT = 1000*dpr(k) * 1.0 * 4.6e-4 / g
-           !   DCONT = rho * 1.0 * 4.6e-4
-           !elseif(WNOI(nw).gt.700.0 .and. WNOI(nw).lt.900.0)then
-           !   DCONT = rho * 0.125 * 4.6e-4
-           !endif
-
-           endif ! generic_continuum_database
            
            DCONT = DCONT*dz(k)