Changeset 4055
- Timestamp:
- Feb 6, 2026, 10:20:21 PM (8 weeks ago)
- Location:
- trunk/LMDZ.GENERIC
- Files:
-
- 12 deleted
- 11 edited
-
changelog.txt (modified) (1 diff)
-
deftank/callphys.GJ581d (modified) (1 diff)
-
deftank/callphys.earlymars (modified) (1 diff)
-
deftank/callphys.earthslab (modified) (1 diff)
-
deftank/callphys.kcm1d (modified) (1 diff)
-
libf/phygeneric/bilinear.F90 (deleted)
-
libf/phygeneric/bilinearbig.F90 (deleted)
-
libf/phygeneric/callkeys_mod.F90 (modified) (1 diff)
-
libf/phygeneric/inifis_mod.F90 (modified) (1 diff)
-
libf/phygeneric/interpolateCH4CH4.F90 (deleted)
-
libf/phygeneric/interpolateCO2CH4.F90 (deleted)
-
libf/phygeneric/interpolateCO2H2.F90 (deleted)
-
libf/phygeneric/interpolateH2CH4.F90 (deleted)
-
libf/phygeneric/interpolateH2H2.F90 (deleted)
-
libf/phygeneric/interpolateH2He.F90 (deleted)
-
libf/phygeneric/interpolateH2O_self_foreign.F90 (deleted)
-
libf/phygeneric/interpolateHeCH4.F90 (deleted)
-
libf/phygeneric/interpolateN2H2.F90 (deleted)
-
libf/phygeneric/interpolateN2N2.F90 (deleted)
-
libf/phygeneric/optci.F90 (modified) (5 diffs)
-
libf/phygeneric/optcv.F90 (modified) (5 diffs)
-
libf/phygeneric/radcommon_h.F90 (modified) (1 diff)
-
libf/phygeneric/sugas_corrk.F90 (modified) (4 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/LMDZ.GENERIC/changelog.txt
r4037 r4055 2168 2168 the black body function. 2169 2169 2170 == 06/02/2026 == MT 2171 Make the 'new generic continuum database' the standard and only way to handle continuum absorption. 2172 Remove all the interpolateXY.F90 routines. Remove the bilinear and bilinearbig routines. RIP. 2173 More details here: https://lmdz-forge.lmd.jussieu.fr/mediawiki/Planets/index.php/Continuum_Database -
trunk/LMDZ.GENERIC/deftank/callphys.GJ581d
r3976 r4055 32 32 # Include continuum absorption in radiative transfer (note CO2 is treated separately) ? 33 33 continuum = .true. 34 generic_continuum_database = .true.35 34 # folder in which correlated-k data is stored ? 36 35 corrkdir = CO2_H2Ovar -
trunk/LMDZ.GENERIC/deftank/callphys.earlymars
r3976 r4055 32 32 # Include continuum absorption in radiative transfer (note CO2 is treated separately) ? 33 33 continuum = .true. 34 generic_continuum_database = .true.35 34 # folder in which correlated-k data is stored ? 36 35 corrkdir = earlymars -
trunk/LMDZ.GENERIC/deftank/callphys.earthslab
r3976 r4055 51 51 # call continuum in radiative transfer ? 52 52 continuum = .true. 53 generic_continuum_database = .true.54 53 # Include Rayleigh scattering in the visible ? 55 54 rayleigh = .true. -
trunk/LMDZ.GENERIC/deftank/callphys.kcm1d
r3976 r4055 16 16 # Include continuum absorption in radiative transfer (note CO2 is treated separately) ? 17 17 continuum = .true. 18 generic_continuum_database = .true.19 18 # folder in which correlated-k data is stored ? 20 19 corrkdir = N2_CO2step_H2Ovar/1E3ppmCO2 -
trunk/LMDZ.GENERIC/libf/phygeneric/callkeys_mod.F90
r4033 r4055 13 13 logical,save :: season,diurnal,tlocked,rings_shadow,lwrite 14 14 !$OMP THREADPRIVATE(season,diurnal,tlocked,rings_shadow,lwrite) 15 logical,save :: callgasvis,continuum,g eneric_continuum_database,graybody16 !$OMP THREADPRIVATE(callgasvis,continuum,g eneric_continuum_database,graybody)15 logical,save :: callgasvis,continuum,graybody 16 !$OMP THREADPRIVATE(callgasvis,continuum,graybody) 17 17 logical,save :: strictboundcorrk 18 18 !$OMP THREADPRIVATE(strictboundcorrk) -
trunk/LMDZ.GENERIC/libf/phygeneric/inifis_mod.F90
r4033 r4055 381 381 call getin_p("continuum",continuum) 382 382 if (is_master) write(*,*) trim(rname)//": continuum = ",continuum 383 384 if (is_master) write(*,*) trim(rname)//&385 ": use generic continuum opacity database ?"//&386 " (matters only if callrad=T)"387 generic_continuum_database=.true. ! default value388 call getin_p("generic_continuum_database",generic_continuum_database)389 if (is_master) write(*,*) trim(rname)//": generic_continuum_database = ", &390 generic_continuum_database391 383 392 384 if (is_master) write(*,*) trim(rname)//": version for H2H2 CIA file ?" -
trunk/LMDZ.GENERIC/libf/phygeneric/optci.F90
r3797 r4055 11 11 use radinc_h, only: L_LEVELS, L_NLAYRAD, L_NSPECTI, L_NGAUSS, & 12 12 L_NLEVRAD, L_REFVAR, naerkind 13 use radcommon_h, only: gasi,tlimit,wrefVAR,Cmk,tgasref,pfgasref,wnoi,scalep, indi,glat_ig13 use radcommon_h, only: gasi,tlimit,wrefVAR,Cmk,tgasref,pfgasref,wnoi,scalep,glat_ig 14 14 use gases_h, only: gfrac, ngasmx, igas_N2, igas_He, igas_H2O, igas_H2, & 15 15 igas_CH4, igas_CO2, igas_O2 16 16 use comcstfi_mod, only: g, r, mugaz 17 use callkeys_mod, only: kastprof,continuum,graybody,varspec, & 18 generic_continuum_database 17 use callkeys_mod, only: kastprof,continuum,graybody,varspec 19 18 use recombin_corrk_mod, only: corrk_recombin, gasi_recomb 20 19 use tpindex_mod, only: tpindex … … 98 97 99 98 integer igas, jgas 100 101 integer interm102 99 103 100 logical :: firstcall=.true. … … 126 123 endif ! of if (firstcall) 127 124 128 !! AS: to save time in computing continuum (see bilinearbig)129 IF (.not.ALLOCATED(indi)) THEN130 ALLOCATE(indi(L_NSPECTI,ngasmx,ngasmx))131 indi = -9999 ! this initial value means "to be calculated"132 ENDIF133 134 125 !======================================================================= 135 126 ! Determine the total gas opacity throughout the column, for each … … 191 182 ! include continua if necessary 192 183 193 if(generic_continuum_database)then194 184 T_cont = dble(TMID(k)) 195 185 do igas=1,ngasmx … … 241 231 242 232 enddo ! igas=1,ngasmx 243 244 else ! generic_continuum_database245 246 wn_cont = dble(wnoi(nw))247 T_cont = dble(TMID(k))248 do igas=1,ngasmx249 250 if(gfrac(igas).eq.-1)then ! variable251 p_cont = dble(PMID(k)*scalep*QVAR(k)) ! qvar = mol/mol252 elseif(varspec) then253 p_cont = dble(PMID(k)*scalep*FRACVAR(igas,k)*(1.-QVAR(k)))254 else255 p_cont = dble(PMID(k)*scalep*gfrac(igas)*(1.-QVAR(k)))256 endif257 258 dtemp=0.0d0259 if(igas.eq.igas_N2)then260 261 interm = indi(nw,igas,igas)262 call interpolateN2N2(wn_cont,T_cont,p_cont,dtemp,.false.,interm)263 indi(nw,igas,igas) = interm264 265 elseif(igas.eq.igas_H2)then266 267 ! first do self-induced absorption268 interm = indi(nw,igas,igas)269 call interpolateH2H2(wn_cont,T_cont,p_cont,dtemp,.false.,interm)270 indi(nw,igas,igas) = interm271 272 ! then cross-interactions with other gases273 do jgas=1,ngasmx274 if(varspec) then275 p_cross = dble(PMID(k)*scalep*FRACVAR(jgas,k)*(1.-QVAR(k)))276 else277 p_cross = dble(PMID(k)*scalep*gfrac(jgas)*(1.-QVAR(k)))278 endif279 dtempc = 0.0d0280 if(jgas.eq.igas_N2)then281 interm = indi(nw,igas,jgas)282 call interpolateN2H2(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)283 indi(nw,igas,jgas) = interm284 elseif(jgas.eq.igas_CO2)then285 interm = indi(nw,igas,jgas)286 call interpolateCO2H2(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)287 indi(nw,igas,jgas) = interm288 elseif(jgas.eq.igas_He)then289 interm = indi(nw,igas,jgas)290 call interpolateH2He(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)291 indi(nw,igas,jgas) = interm292 endif293 dtemp = dtemp + dtempc294 enddo295 296 elseif(igas.eq.igas_CH4)then297 298 ! first do self-induced absorption299 interm = indi(nw,igas,igas)300 call interpolateCH4CH4(wn_cont,T_cont,p_cont,dtemp,.false.,interm)301 indi(nw,igas,igas) = interm302 303 ! then cross-interactions with other gases304 do jgas=1,ngasmx305 if(varspec) then306 p_cross = dble(PMID(k)*scalep*FRACVAR(jgas,k)*(1.-QVAR(k)))307 else308 p_cross = dble(PMID(k)*scalep*gfrac(jgas)*(1.-QVAR(k)))309 endif310 dtempc = 0.0d0311 if(jgas.eq.igas_H2)then312 interm = indi(nw,igas,jgas)313 call interpolateH2CH4(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)314 indi(nw,igas,jgas) = interm315 elseif(jgas.eq.igas_CO2)then316 interm = indi(nw,igas,jgas)317 call interpolateCO2CH4(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)318 indi(nw,igas,jgas) = interm319 elseif(jgas.eq.igas_He)then320 interm = indi(nw,igas,jgas)321 call interpolateHeCH4(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)322 indi(nw,igas,jgas) = interm323 endif324 dtemp = dtemp + dtempc325 enddo326 327 elseif(igas.eq.igas_H2O.and.T_cont.gt.100.0)then328 ! Compute self and foreign (with air) continuum of H2O329 p_air = dble(PMID(k)*scalep) - p_cont ! note assumes background is air!330 interm = indi(nw,igas,igas)331 call interpolateH2O_self_foreign(wn_cont,T_cont,p_cont,p_air,dtemp,.false.,interm) ! MTCKD v3.3332 indi(nw,igas,igas) = interm333 334 endif335 336 DCONT = DCONT + dtemp337 338 enddo339 340 ! Oobleck test341 !rho = PMID(k)*scalep / (TMID(k)*286.99)342 !if(WNOI(nw).gt.300.0 .and. WNOI(nw).lt.500.0)then343 ! DCONT = rho * 0.125 * 4.6e-4344 !elseif(WNOI(nw).gt.500.0 .and. WNOI(nw).lt.700.0)then345 ! DCONT = 1000*dpr(k) * 1.0 * 4.6e-4 / g346 ! DCONT = rho * 1.0 * 4.6e-4347 !elseif(WNOI(nw).gt.700.0 .and. WNOI(nw).lt.900.0)then348 ! DCONT = rho * 0.125 * 4.6e-4349 !endif350 351 endif ! generic_continuum_database352 233 353 234 DCONT = DCONT*dz(k) -
trunk/LMDZ.GENERIC/libf/phygeneric/optcv.F90
r3797 r4055 10 10 11 11 use radinc_h, only: L_NLAYRAD, L_NLEVRAD, L_LEVELS, L_NSPECTV, L_NGAUSS, L_REFVAR, NAERKIND 12 use radcommon_h, only: gasv, tlimit, wrefVAR, Cmk, tgasref, pfgasref,wnov,scalep, indv,glat_ig12 use radcommon_h, only: gasv, tlimit, wrefVAR, Cmk, tgasref, pfgasref,wnov,scalep,glat_ig 13 13 use gases_h, only: gfrac, ngasmx, igas_H2, igas_H2O, igas_He, igas_N2, & 14 14 igas_CH4, igas_CO2, igas_O2 15 15 use comcstfi_mod, only: g, r, mugaz 16 16 use callkeys_mod, only: kastprof,continuum,graybody,callgasvis,varspec, & 17 generic_continuum_database,rayleigh17 rayleigh 18 18 use recombin_corrk_mod, only: corrk_recombin, gasv_recomb 19 19 use tpindex_mod, only: tpindex … … 107 107 integer igas, jgas 108 108 109 integer interm110 111 109 logical :: firstcall=.true. 112 110 !$OMP THREADPRIVATE(firstcall) … … 120 118 firstcall=.false. 121 119 endif ! of if (firstcall) 122 123 !! AS: to save time in computing continuum (see bilinearbig)124 IF (.not.ALLOCATED(indv)) THEN125 ALLOCATE(indv(L_NSPECTV,ngasmx,ngasmx))126 indv = -9999 ! this initial value means "to be calculated"127 ENDIF128 120 129 121 !======================================================================= … … 211 203 ! include continua if necessary 212 204 213 if(generic_continuum_database)then214 205 T_cont = dble(TMID(k)) 215 206 do igas=1,ngasmx … … 261 252 262 253 enddo ! igas=1,ngasmx 263 264 else ! generic_continuum_database265 266 267 wn_cont = dble(wnov(nw))268 T_cont = dble(TMID(k))269 do igas=1,ngasmx270 271 if(gfrac(igas).eq.-1)then ! variable272 p_cont = dble(PMID(k)*scalep*QVAR(k)) ! qvar = mol/mol273 elseif(varspec) then274 p_cont = dble(PMID(k)*scalep*FRACVAR(igas,k)*(1.-QVAR(k)))275 else276 p_cont = dble(PMID(k)*scalep*gfrac(igas)*(1.-QVAR(k)))277 endif278 279 dtemp=0.0280 if(igas.eq.igas_N2)then281 282 interm = indv(nw,igas,igas)283 ! call interpolateN2N2(wn_cont,T_cont,p_cont,dtemp,.false.,interm)284 indv(nw,igas,igas) = interm285 ! only goes to 500 cm^-1, so unless we're around a cold brown dwarf, this is irrelevant in the visible286 287 elseif(igas.eq.igas_H2)then288 289 ! first do self-induced absorption290 interm = indv(nw,igas,igas)291 call interpolateH2H2(wn_cont,T_cont,p_cont,dtemp,.false.,interm)292 indv(nw,igas,igas) = interm293 294 ! then cross-interactions with other gases295 do jgas=1,ngasmx296 if(varspec) then297 p_cross = dble(PMID(k)*scalep*FRACVAR(jgas,k)*(1.-QVAR(k)))298 else299 p_cross = dble(PMID(k)*scalep*gfrac(jgas)*(1.-QVAR(k)))300 endif301 dtempc = 0.0302 if(jgas.eq.igas_N2)then303 interm = indv(nw,igas,jgas)304 call interpolateN2H2(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)305 indv(nw,igas,jgas) = interm306 ! should be irrelevant in the visible307 elseif(jgas.eq.igas_CO2)then308 interm = indv(nw,igas,jgas)309 call interpolateCO2H2(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)310 indv(nw,igas,jgas) = interm311 ! might not be relevant in the visible312 elseif(jgas.eq.igas_He)then313 interm = indv(nw,igas,jgas)314 call interpolateH2He(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)315 indv(nw,igas,jgas) = interm316 endif317 dtemp = dtemp + dtempc318 enddo319 320 elseif(igas.eq.igas_CH4)then321 322 ! first do self-induced absorption323 interm = indv(nw,igas,igas)324 call interpolateCH4CH4(wn_cont,T_cont,p_cont,dtemp,.false.,interm)325 indv(nw,igas,igas) = interm326 327 ! then cross-interactions with other gases328 do jgas=1,ngasmx329 if(varspec) then330 p_cross = dble(PMID(k)*scalep*FRACVAR(jgas,k)*(1.-QVAR(k)))331 else332 p_cross = dble(PMID(k)*scalep*gfrac(jgas)*(1.-QVAR(k)))333 endif334 dtempc = 0.0d0335 if(jgas.eq.igas_H2)then336 interm = indv(nw,igas,jgas)337 call interpolateH2CH4(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)338 indv(nw,igas,jgas) = interm339 elseif(jgas.eq.igas_CO2)then340 interm = indv(nw,igas,jgas)341 call interpolateCO2CH4(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)342 indv(nw,igas,jgas) = interm343 ! might not be relevant in the visible344 elseif(jgas.eq.igas_He)then345 interm = indv(nw,igas,jgas)346 call interpolateHeCH4(wn_cont,T_cont,p_cross,p_cont,dtempc,.false.,interm)347 indv(nw,igas,jgas) = interm348 endif349 dtemp = dtemp + dtempc350 enddo351 352 elseif(igas.eq.igas_H2O.and.T_cont.gt.100.0)then353 ! Compute self and foreign (with air) continuum of H2O354 p_air = dble(PMID(k)*scalep) - p_cont ! note assumes background is air!355 interm = indv(nw,igas,igas)356 call interpolateH2O_self_foreign(wn_cont,T_cont,p_cont,p_air,dtemp,.false.,interm) ! MTCKD v3.3357 indv(nw,igas,igas) = interm358 359 endif360 361 DCONT = DCONT + dtemp362 363 enddo364 365 endif ! generic_continuum_database366 254 367 255 DCONT = DCONT*dz(k) -
trunk/LMDZ.GENERIC/libf/phygeneric/radcommon_h.F90
r3693 r4055 73 73 REAL*8 blamv(L_NSPECTV+1) ! these are needed by suaer.F90 74 74 !$OMP THREADPRIVATE(blami,blamv) 75 76 !! AS: introduced to avoid doing same computations again for continuum77 INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: indi78 INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: indv79 !$OMP THREADPRIVATE(indi,indv)80 75 81 76 !!! ALLOCATABLE STUFF SO THAT DIMENSIONS ARE READ in *.dat FILES -- AS 12/2011 -
trunk/LMDZ.GENERIC/libf/phygeneric/sugas_corrk.F90
r3893 r4055 34 34 use ioipsl_getin_p_mod, only: getin_p 35 35 use callkeys_mod, only: varactive,varfixed,graybody,callgasvis,& 36 continuum ,generic_continuum_database36 continuum 37 37 use tracer_h, only : nqtot, moderntracdef, is_recomb, noms 38 38 use recombin_corrk_mod, only: su_recombin, & … … 62 62 63 63 double precision testcont ! for continuum absorption initialisation 64 65 integer :: dummy66 64 67 65 if (.not. moderntracdef) use_premix=.true. ! Added by JVO for compatibility with 'old' traceur.def … … 717 715 ! Initialise the continuum absorption data 718 716 if(continuum)then 719 if(generic_continuum_database)then717 720 718 do igas=1,ngasmx ! we loop on all pairs of molecules that have data available 721 719 ! data can be downloaded from https://web.lmd.jussieu.fr/~lmdz/planets/generic/datagcm/continuum_data/ … … 806 804 enddo ! igas=1,ngasmx 807 805 808 else ! generic_continuum_database tag809 do igas=1,ngasmx810 if (igas .eq. igas_N2) then811 812 dummy = -9999813 call interpolateN2N2(100.D+0,250.D+0,17500.D+0,testcont,.true.,dummy)814 815 elseif (igas .eq. igas_H2) then816 817 ! first do self-induced absorption818 dummy = -9999819 call interpolateH2H2(500.D+0,250.D+0,17500.D+0,testcont,.true.,dummy)820 ! then cross-interactions with other gases821 do jgas=1,ngasmx822 if (jgas .eq. igas_N2) then823 dummy = -9999824 call interpolateN2H2(592.D+0,278.15D+0,200000.D+0,10000.D+0,testcont,.true.,dummy)825 elseif (jgas .eq. igas_CO2) then826 dummy = -9999827 call interpolateCO2H2(592.D+0,278.15D+0,200000.D+0,10000.D+0,testcont,.true.,dummy)828 elseif (jgas .eq. igas_He) then829 dummy = -9999830 call interpolateH2He(500.D+0,250.D+0,200000.D+0,10000.D+0,testcont,.true.,dummy)831 endif832 enddo833 834 835 elseif (igas .eq. igas_CH4) then836 837 ! first do self-induced absorption838 dummy = -9999839 call interpolateCH4CH4(600.0,66.7,400.0,testcont,.true.,dummy)840 ! then cross-interactions with other gases841 do jgas=1,ngasmx842 if (jgas .eq. igas_H2) then843 dummy = -9999844 call interpolateH2CH4(500.D+0,250.D+0,200000.D+0,10000.D+0,testcont,.true.,dummy)845 elseif (jgas .eq. igas_CO2) then846 dummy = -9999847 call interpolateCO2CH4(592.D+0,278.15D+0,200000.D+0,10000.D+0,testcont,.true.,dummy)848 elseif (jgas .eq. igas_He) then849 dummy = -9999850 call interpolateHeCH4(500.D+0,250.D+0,200000.D+0,10000.D+0,testcont,.true.,dummy)851 endif852 enddo853 854 855 elseif (igas .eq. igas_H2O) then856 857 ! Compute self and foreign (with air) continuum of H2O858 dummy = -9999859 call interpolateH2O_self_foreign(990.D+0,296.D+0,683.2D+0*2,0.D+0,testcont,.true.,dummy)860 861 endif862 863 enddo ! igas=1,ngasmx864 endif ! generic_continuum_database tag865 806 endif ! continuum flag 866 807
Note: See TracChangeset
for help on using the changeset viewer.
