Ignore:
Timestamp:
Feb 18, 2026, 10:28:23 AM (6 weeks ago)
Author:
mturbet
Message:

renaming and merging of radiative routines in phygeneric

Location:
trunk/LMDZ.GENERIC/libf/phygeneric
Files:
1 deleted
13 edited
29 moved

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.GENERIC/libf/phygeneric/aerosol_global_variables.F90

    r4076 r4077  
    11!==================================================================
    2 module aerosol_mod
     2module aerosol_global_variables
    33implicit none
    44
     
    77!  aerosol indexes: these are initialized to be 0 if the
    88!                 corresponding aerosol was not activated in callphys.def
    9 !                 -- otherwise a value is set via iniaerosol
     9!                 -- otherwise a value is set via aerosol_init
    1010      integer, save, protected :: iaero_co2 = 0
    1111      integer, save, protected :: iaero_h2o = 0
     
    4343contains
    4444
    45   SUBROUTINE iniaerosol
     45  SUBROUTINE aerosol_init
    4646
    4747  use mod_phys_lmdz_para, only : is_master
     
    181181
    182182! For the zero aerosol case, we currently make a dummy co2 aerosol which is zero everywhere.
    183 ! (See aeropacity.F90 for how this works). A better solution would be to turn off the
     183! (See aerosol_opacity.F90 for how this works). A better solution would be to turn off the
    184184! aerosol machinery in the no aerosol case, but this would be complicated. LK
    185185
     
    197197      print*, 'according to current options in callphys.def'
    198198      print*, 'or change/correct incompatible options there'
    199       print*, 'Abort in iniaerosol'
     199      print*, 'Abort in aerosol_init'
    200200    endif
    201201    call abort_physic("iniaerosl",'wrong number of aerosols',1)
    202202  endif ! of if (ia.ne.naerkind)
    203203
    204   END SUBROUTINE iniaerosol
    205 
    206 end module aerosol_mod
    207 !==================================================================
     204  END SUBROUTINE aerosol_init
     205
     206end module aerosol_global_variables
     207!==================================================================
  • trunk/LMDZ.GENERIC/libf/phygeneric/aerosol_opacity.F90

    r4076 r4077  
    1 module aeropacity_mod
     1module aerosol_opacity_mod
    22
    33implicit none
     
    55contains
    66
    7       Subroutine aeropacity(ngrid,nlayer,nq,pplay,pplev,pt,pq,zls, &
     7      Subroutine aerosol_opacity(ngrid,nlayer,nq,pplay,pplev,pt,pq,zls, &
    88         aerosol,reffrad,nueffrad, QREFvis3d,QREFir3d,tau_col, &
    99         cloudfrac,totcloudfrac,clearsky)
    1010
    1111       use radinc_h, only : L_TAUMAX,naerkind
    12        use aerosol_mod, only: iaero_nlay, iaero_generic, &
     12       use aerosol_global_variables , only: iaero_nlay, iaero_generic, &
    1313                              iaero_aurora, iaero_back2lay, iaero_co2, &
    1414                              iaero_dust, iaero_h2o, iaero_h2so4, &
     
    120120      IF (firstcall) THEN
    121121        ia =0
    122         write(*,*) "Tracers found in aeropacity:"
     122        write(*,*) "Tracers found in aerosol_opacity:"
    123123        do iq=1,nq
    124124          tracername=noms(iq)
     
    135135
    136136        if (noaero) then
    137           print*, "No active aerosols found in aeropacity"
     137          print*, "No active aerosols found in aerosol_opacity"
    138138        else
    139139          print*, "If you would like to use aerosols, make sure any old"
    140140          print*, "start files are updated in newstart using the option"
    141141          print*, "q=0"
    142           write(*,*) "Active aerosols found in aeropacity:"
     142          write(*,*) "Active aerosols found in aerosol_opacity:"
    143143        endif
    144144
     
    269269                  end do
    270270
    271                   call abort_physic("aeropacity", "Something wrong happened on water ice liquid opacity calculation",1)
     271                  call abort_physic("aerosol_opacity", "Something wrong happened on water ice liquid opacity calculation",1)
    272272               endif
    273273
     
    276276               do ig=1, ngrid
    277277                  !do l=1,nlayer-1 ! to stop the rad tran bug
    278                   do l=1,nlayer !JL18 if aerosols are present in the last layer we must account for them. Provides better upper boundary condition in the IR. They must however be put to zero in the sw (see optcv)
     278                  do l=1,nlayer !JL18 if aerosols are present in the last layer we must account for them. Provides better upper boundary condition in the IR. They must however be put to zero in the sw (see rad_correlatedk_opacities_stellar)
    279279                                ! same correction should b-probably be done for other aerosol types.
    280280                     aerosol(ig,l,iaer) =                                    & !modification by BC
     
    294294                  do ig=1, ngrid
    295295                     !do l=1,nlayer-1 ! to stop the rad tran bug
    296                      do l=1,nlayer !JL18 if aerosols are present in the last layer we must account for them. Provides better upper boundary condition in the IR. They must however be put to zero in the sw (see optcv)
     296                     do l=1,nlayer !JL18 if aerosols are present in the last layer we must account for them. Provides better upper boundary condition in the IR. They must however be put to zero in the sw (see rad_correlatedk_opacities_stellar)
    297297                        CLFtot  = max(totcloudfrac(ig),0.01)
    298298                        aerosol(ig,l,iaer)=aerosol(ig,l,iaer)/CLFtot
     
    303303                  do ig=1, ngrid
    304304                     !do l=1,nlayer-1 ! to stop the rad tran bug
    305                      do l=1,nlayer !JL18 if aerosols are present in the last layer we must account for them. Provides better upper boundary condition in the IR. They must however be put to zero in the sw (see optcv)
     305                     do l=1,nlayer !JL18 if aerosols are present in the last layer we must account for them. Provides better upper boundary condition in the IR. They must however be put to zero in the sw (see rad_correlatedk_opacities_stellar)
    306306                        CLFtot  = CLFfixval
    307307                        aerosol(ig,l,iaer)=aerosol(ig,l,iaer)/CLFtot
     
    10921092      ! end do
    10931093
    1094     end subroutine aeropacity
     1094    end subroutine aerosol_opacity
    10951095     
    1096 end module aeropacity_mod
     1096end module aerosol_opacity_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/aerosol_optical_properties.F90

    r4076 r4077  
    1 module aeroptproperties_mod
     1module aerosol_optical_properties_mod
    22
    33implicit none
     
    55contains
    66
    7       SUBROUTINE aeroptproperties(ngrid,nlayer,reffrad,nueffrad,   &
     7      SUBROUTINE aerosol_optical_properties(ngrid,nlayer,reffrad,  &
     8                                  nueffrad,                        &
    89                                  QVISsQREF3d,omegaVIS3d,gVIS3d,   &
    910                                  QIRsQREF3d,omegaIR3d,gIR3d,      &
    1011                                  QREFvis3d,QREFir3d)!,            &
    11 !                                  omegaREFvis3d,omegaREFir3d)
     12!                                 omegaREFvis3d,omegaREFir3d)
    1213
    1314      use radinc_h,    only: L_NSPECTI,L_NSPECTV,nsizemax,naerkind
     
    310311!       1.3 Effective variance
    311312        if(nuefftabsize.eq.1)then ! addded by RDW
    312            print*,'Warning: no variance range in aeroptproperties'
     313           print*,'Warning: no variance range in aerosol_optical_properties'
    313314           nuefftab(1)=0.2
    314315        else
     
    814815      ENDDO ! iaer (loop on aerosol kind)
    815816
    816     END SUBROUTINE aeroptproperties
    817 
    818 
    819 end module aeroptproperties_mod
     817    END SUBROUTINE aerosol_optical_properties
     818
     819
     820end module aerosol_optical_properties_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/aerosol_optical_properties_averaging.F

    r4076 r4077  
    1       SUBROUTINE aerave_new ( ndata,
     1      SUBROUTINE aerosol_optical_properties_averaging(ndata,
    22     & longdata,epdata,omegdata,gdata,         
    3      &            longref,epref,temp,nir,longir
    4      &            ,epir,omegir,gir,qref,omegaref        )
     3     & longref,epref,temp,nir,longir,
     4     & epir,omegir,gir,qref,omegaref)
    55
    66
     
    8787        IF (longdata(1).LT.longdata(ndata)) THEN
    8888          IF (.not.(longdata(idata).LT.longdata(idata+1))) THEN
    89            call abort_physic("aerave_new",
     89           call abort_physic("aerosol_optical_properties_averaging",
    9090     &     "Non descending order in longdata",1)
    9191          ENDIF
    9292        ELSEIF (longdata(1).GT.longdata(ndata)) THEN
    9393          IF (.not.(longdata(idata).GT.longdata(idata+1))) THEN
    94            call abort_physic("aerave_new",
     94           call abort_physic("aerosol_optical_properties_averaging",
    9595     &     "Non ascending order in longdata",1)
    9696          ENDIF
     
    174174c
    175175            long=longir(iir) + (ibande-0.5E+0) * deltalong
    176             CALL blackl(DBLE(long),DBLE(temp),tmp1)
     176            CALL rad_blackbody_planck_law_wavelength(DBLE(long),
     177     &      DBLE(temp),tmp1)
    177178            emit=REAL(tmp1)
    178179c
     
    240241c
    241242            long=longir(iir) + (ibande-0.5E+0) * deltalong
    242             CALL blackl(DBLE(long),DBLE(temp),tmp1)
     243            CALL rad_blackbody_planck_law_wavelength(DBLE(long)
     244     &      ,DBLE(temp),tmp1)
    243245            emit=REAL(tmp1)
    244246c
  • trunk/LMDZ.GENERIC/libf/phygeneric/aerosol_radius.F90

    r4076 r4077  
    11!==================================================================
    2 module radii_mod
     2module aerosol_radius
    33!==================================================================
    44!  module to centralize the radii calculations for aerosols
     
    1313!$OMP THREADPRIVATE(radfixed)
    1414
    15 !     water cloud optical properties (initialized in su_aer_radii below)
     15!     water cloud optical properties (initialized in aerosol_radius_init below)
    1616      real, save ::  rad_h2o
    1717      real, save ::  rad_h2o_ice
     
    2121
    2222      real,save :: nueff_iaero_h2o ! effective variance of H2O aerosol
    23                                    ! (initialized in su_aer_radii below)
     23                                   ! (initialized in aerosol_radius_init below)
    2424!$OMP THREADPRIVATE(nueff_iaero_h2o)
    2525! coefficients for a variable nueff() for h2o aerosol; disabled for now
     
    3232
    3333!==================================================================
    34    subroutine su_aer_radii(ngrid,nlayer,reffrad,nueffrad)
     34   subroutine aerosol_radius_init(ngrid,nlayer,reffrad,nueffrad)
    3535!==================================================================
    3636!     Purpose
     
    5050      use ioipsl_getin_p_mod, only: getin_p
    5151      use radinc_h, only: naerkind
    52       use aerosol_mod, only: iaero_back2lay, iaero_co2, iaero_dust, &
     52      use aerosol_global_variables , only: iaero_back2lay, iaero_co2, iaero_dust, &
    5353                             iaero_h2o, iaero_h2so4, iaero_nh3, iaero_nlay, &
    5454                             iaero_aurora, iaero_generic, i_rgcs_ice, &
     
    185185
    186186
    187    end subroutine su_aer_radii
     187   end subroutine aerosol_radius_init
    188188!==================================================================
    189189
     
    243243
    244244! For now only constant nueff is enabled (otherwise some specific handling
    245 ! of variable nueff is required in aeroptproperties)
     245! of variable nueff is required in aerosol_optical_properties)
    246246      nueffrad(1:ngrid,1:nlayer)=nueff_iaero_h2o
    247247
     
    296296
    297297!==================================================================
    298    subroutine co2_reffrad(ngrid,nlayer,nq,pq,reffrad)
     298   subroutine aerosol_radius_co2(ngrid,nlayer,nq,pq,reffrad)
    299299!==================================================================
    300300!     Purpose
     
    333333      end if
    334334
    335    end subroutine co2_reffrad
    336 !==================================================================
    337 
    338 
    339 
    340 !==================================================================
    341    subroutine dust_reffrad(ngrid,nlayer,reffrad)
     335   end subroutine aerosol_radius_co2
     336!==================================================================
     337
     338
     339
     340!==================================================================
     341   subroutine aerosol_radius_dust(ngrid,nlayer,reffrad)
    342342!==================================================================
    343343!     Purpose
     
    359359      reffrad(1:ngrid,1:nlayer) = 2.e-6 ! dust
    360360
    361    end subroutine dust_reffrad
    362 !==================================================================
    363 
    364 
    365 !==================================================================
    366    subroutine h2so4_reffrad(ngrid,nlayer,reffrad)
     361   end subroutine aerosol_radius_dust
     362!==================================================================
     363
     364
     365!==================================================================
     366   subroutine aerosol_radius_h2so4(ngrid,nlayer,reffrad)
    367367!==================================================================
    368368!     Purpose
     
    384384      reffrad(1:ngrid,1:nlayer) = 1.e-6 ! h2so4
    385385
    386    end subroutine h2so4_reffrad
    387 !==================================================================
    388 
    389 !==================================================================
    390    subroutine back2lay_reffrad(ngrid,reffrad,nlayer,pplev)
     386   end subroutine aerosol_radius_h2so4
     387!==================================================================
     388
     389!==================================================================
     390   subroutine aerosol_radius_back2lay(ngrid,reffrad,nlayer,pplev)
    391391!==================================================================
    392392!     Purpose
     
    426426          ENDDO
    427427
    428    end subroutine back2lay_reffrad
    429 !==================================================================
    430 
    431 end module radii_mod
    432 !==================================================================
     428   end subroutine aerosol_radius_back2lay
     429!==================================================================
     430
     431end module aerosol_radius
     432!==================================================================
  • trunk/LMDZ.GENERIC/libf/phygeneric/callsedim.F

    r3893 r4077  
    1010
    1111      use radinc_h, only : naerkind
    12       use radii_mod, only: h2o_reffrad
    13       use aerosol_mod, only : iaero_h2o
     12      use aerosol_radius, only: h2o_reffrad
     13      use aerosol_global_variables , only : iaero_h2o
    1414      USE tracer_h, only : igcm_co2_ice,igcm_h2o_ice,radius,rho_q
    1515      use comcstfi_mod, only: g
  • trunk/LMDZ.GENERIC/libf/phygeneric/condense_co2.F90

    r3893 r4077  
    88      use radinc_h, only : L_NSPECTV
    99      use gases_h, only: gfrac, igas_co2
    10       use radii_mod, only : co2_reffrad
    11       use aerosol_mod, only : iaero_co2
     10      use aerosol_radius, only : aerosol_radius_co2
     11      use aerosol_global_variables , only : iaero_co2
    1212      USE surfdat_h, only: emisice, emissiv
    1313      USE geometry_mod, only: latitude ! in radians
     
    273273         ! Gravitational sedimentation starts.
    274274           
    275          ! Sedimentation computed from radius computed from q in module radii_mod.
    276          call co2_reffrad(ngrid,nlayer,nq,zq,reffrad)
     275         ! Sedimentation computed from radius computed from q in module aerosol_radius.
     276         call aerosol_radius_co2(ngrid,nlayer,nq,zq,reffrad)
    277277         
    278278         DO  ilay=1,nlayer
  • trunk/LMDZ.GENERIC/libf/phygeneric/dyn1d/kcm1d.F90

    r3893 r4077  
    1111                          varspec, varspec_data, nvarlayer
    1212  use inifis_mod, only: inifis
    13   use aerosol_mod, only: iniaerosol
    14   use callcorrk_mod, only: callcorrk
     13  use aerosol_global_variables , only: aerosol_init
     14  use rad_correlatedk_mod, only: rad_correlatedk
    1515  use comcstfi_mod
    1616  use mod_grid_phy_lmdz, only : regular_lonlat
     
    410410         
    411411               
    412      call iniaerosol
     412     call aerosol_init
    413413     
    414414
     
    416416
    417417     !    Run radiative transfer
    418      call callcorrk(1,nlayer,q,nq,qsurf,                  &
     418     call rad_correlatedk(1,nlayer,q,nq,qsurf,                  &
    419419          albedo_wv,albedo_equivalent,                    &
    420420          emis,mu0,plev,play,temp,                        &
     
    470470  firstcall=.false.
    471471  lastcall=.true.
    472   call callcorrk(1,nlayer,q,nq,qsurf,                          &
     472  call rad_correlatedk(1,nlayer,q,nq,qsurf,                          &
    473473       albedo_wv,albedo_equivalent,emis,mu0,plev,play,temp,    &
    474474       tsurf,fract,dist_star,aerosol,muvar,                    &
  • trunk/LMDZ.GENERIC/libf/phygeneric/dyn1d/rcm1d.F

    r3995 r4077  
    11611161        IF (idt.eq.ndt) then       !test
    11621162         lastcall=.true.
    1163          call stellarlong(day*1.0,zls)
     1163         call ephemeris_stellar_longitude(day*1.0,zls)
    11641164!         write(103,*) 'Ls=',zls*180./pi
    11651165!         write(103,*) 'Lat=', latitude(1)*180./pi
  • trunk/LMDZ.GENERIC/libf/phygeneric/ephemeris_orbit.F

    r4076 r4077  
    1       subroutine orbite(pls,pdist_star,pdecli,pright_ascenc)
     1      subroutine ephemeris_orbit(pls,pdist_star,pdecli,pright_ascenc)
    22
    33      use planete_mod, only: p_elips, e_elips, timeperi, obliquit
  • trunk/LMDZ.GENERIC/libf/phygeneric/ephemeris_orbit_init.F

    r4076 r4077  
    1       SUBROUTINE iniorbit
     1      SUBROUTINE ephemeris_orbit_init
    22     $     (papoastr,pperiastr,pyear_day,pperi_day,pobliq)
    33     
     
    3232      peri_day=pperi_day
    3333
    34       PRINT*,'iniorbit: Periastron in AU  ',periastr
    35       PRINT*,'iniorbit: Apoastron in AU  ',apoastr
    36       PRINT*,'iniorbit: Obliquity in degrees  :',obliquit
     34      PRINT*,'ephemeris_orbit_init: Periastron in AU  ',periastr
     35      PRINT*,'ephemeris_orbit_init: Apoastron in AU  ',apoastr
     36      PRINT*,'ephemeris_orbit_init: Obliquity in degrees  :',obliquit
    3737
    3838
     
    4040      p_elips=0.5*(periastr+apoastr)*(1-e_elips*e_elips)
    4141
    42       print*,'iniorbit: e_elips',e_elips
    43       print*,'iniorbit: p_elips',p_elips
     42      print*,'ephemeris_orbit_init: e_elips',e_elips
     43      print*,'ephemeris_orbit_init: p_elips',p_elips
    4444
    4545!-----------------------------------------------------------------------
     
    5252      zanom=2.*pi*(zz-nint(zz))
    5353      zxref=abs(zanom)
    54       PRINT*,'iniorbit: zanom  ',zanom
     54      PRINT*,'ephemeris_orbit_init: zanom  ',zanom
    5555
    5656!  solve equation  zx0 - e * sin (zx0) = zxref for eccentric anomaly zx0
     
    6666      zx0=zx0+zdx
    6767      if(zanom.lt.0.) zx0=-zx0
    68       PRINT*,'iniorbit: zx0   ',zx0
     68      PRINT*,'ephemeris_orbit_init: zx0   ',zx0
    6969
    7070      timeperi=2.*atan(sqrt((1.+e_elips)/(1.-e_elips))*tan(zx0/2.))
    71       PRINT*,'iniorbit: Perihelion solar long. Ls (deg)=',
     71      PRINT*,'ephemeris_orbit_init: Perihelion solar long. Ls (deg)=',
    7272     &       360.-timeperi*180./pi
    7373
  • trunk/LMDZ.GENERIC/libf/phygeneric/ephemeris_stellar_angle.F

    r4076 r4077  
    1       subroutine stelang(kgrid,psilon,pcolon,psilat,pcolat,
     1      subroutine ephemeris_stellar_angle(kgrid,
     2     &                psilon,pcolon,psilat,pcolat,
    23     &                ptim1,ptim2,ptim3,pmu0,pfract, pflat)
    34      IMPLICIT NONE
     
    1213C**   INTERFACE.
    1314C     ----------
    14 C      SUBROUTINE STELANG ( KGRID )
     15C      SUBROUTINE ephemeris_stellar_angle ( KGRID )
    1516C
    1617C        EXPLICIT ARGUMENTS :
  • trunk/LMDZ.GENERIC/libf/phygeneric/ephemeris_stellar_longitude.F

    r4076 r4077  
    1       SUBROUTINE stellarlong(pday,pstellong)
     1      SUBROUTINE ephemeris_stellar_longitude(pday,pstellong)
    22     
    33      USE planete_mod, ONLY: year_day, peri_day, e_elips, timeperi
  • trunk/LMDZ.GENERIC/libf/phygeneric/inifis_mod.F90

    r4055 r4077  
    1212  use radinc_h, only: ini_radinc_h, naerkind
    1313  use radcommon_h, only: ini_radcommon_h
    14   use radii_mod, only: radfixed, Nmix_co2
     14  use aerosol_radius, only: radfixed, Nmix_co2
    1515  use datafile_mod, only: datadir
    1616  use comdiurn_h, only: sinlat, coslat, sinlon, coslon
     
    2727  use ioipsl_getin_p_mod, only : getin_p
    2828  use mod_phys_lmdz_para, only : is_parallel, is_master, bcast
    29   use newton_cooling_hotJ, only: planetary_suffix
     29  use rad_netwon_cooling_hot_jupiter, only: planetary_suffix
    3030
    3131!=======================================================================
     
    8080  CHARACTER(len=20) :: rname="inifis" ! routine name, for messages
    8181 
    82   EXTERNAL iniorbit,orbite
     82  EXTERNAL ephemeris_orbit_init,ephemeris_orbit
    8383  EXTERNAL SSUM
    8484  REAL SSUM
  • trunk/LMDZ.GENERIC/libf/phygeneric/initracer.F90

    r3893 r4077  
    44      USE tracer_h
    55      USE callkeys_mod, only: water
    6       USE recombin_corrk_mod, ONLY: ini_recombin
     6      USE rad_correlatedk_online_recombination_mod, ONLY: rad_correlatedk_recombination_init
    77      USE mod_phys_lmdz_para, only: is_master, bcast
    88      use generic_cloud_common_h
     
    503503!     Processing modern traceur options
    504504      if(moderntracdef) then
    505         call ini_recombin
     505        call rad_correlatedk_recombination_init
    506506      endif
    507507
  • trunk/LMDZ.GENERIC/libf/phygeneric/newsedim.F

    r3663 r4077  
    1515      use tracer_h, only : igcm_h2o_ice
    1616      use watercommon_h, only: T_h2O_ice_liq,T_h2O_ice_clouds     
    17       use radii_mod, only: h2o_cloudrad
     17      use aerosol_radius, only: h2o_cloudrad
    1818
    1919      IMPLICIT NONE
  • trunk/LMDZ.GENERIC/libf/phygeneric/physiq_mod.F90

    r4033 r4077  
    2323      use gases_h, only: gnom, gfrac, ngasmx
    2424      use radcommon_h, only: sigma, glat, grav, BWNV, WNOI, DWNI, DWNV, WNOV
    25       use suaer_corrk_mod, only: suaer_corrk
    26       use setspv_mod, only: setspv
    27       use radii_mod, only: h2o_reffrad, co2_reffrad
    28       use aerosol_mod, only: iniaerosol, iaero_co2, iaero_h2o
     25      use rad_correlatedk_ini_aerosol_mod, only: rad_correlatedk_ini_aerosol
     26      use rad_correlatedk_init_stellar_mod, only: rad_correlatedk_init_stellar
     27      use aerosol_radius, only: h2o_reffrad, aerosol_radius_co2
     28      use aerosol_global_variables , only: aerosol_init, iaero_co2, iaero_h2o
    2929      use surfdat_h, only: phisfi, zmea, zstd, zsig, zgam, zthe, &
    3030                           dryness
     
    8282      use conc_mod, only: rnew, cpnew, ini_conc_mod
    8383      use phys_state_var_mod
    84       use callcorrk_mod, only: callcorrk
     84      use rad_correlatedk_mod, only: rad_correlatedk
    8585      use conduction_mod, only: conduction
    8686      use molvis_mod, only: molvis
     
    9595      use condensation_generic_mod, only: condensation_generic
    9696      use datafile_mod, only: datadir
    97       use newton_cooling_hotJ, only: newtcool_MOCHA ! LT, adding for MOCHA protocol
     97      use rad_netwon_cooling_hot_jupiter, only: rad_newton_cooling_MOCHA_intercomparison ! LT, adding for MOCHA protocol
    9898
    9999#ifndef MESOSCALE
     
    337337      real dtmoist(ngrid,nlayer)                              ! Moistadj routine.
    338338      real dt_ekman(ngrid,nslay), dt_hdiff(ngrid,nslay), dt_gm(ngrid,nslay) ! Slab_ocean routine.
    339       real zdtsw1(ngrid,nlayer), zdtlw1(ngrid,nlayer)         ! Callcorrk routine.
     339      real zdtsw1(ngrid,nlayer), zdtlw1(ngrid,nlayer)         ! rad_correlatedk routine.
    340340      real zdtchim(ngrid,nlayer)                              ! Calchim routine.
    341341
     
    594594!        Initialize aerosol indexes.
    595595!        ~~~~~~~~~~~~~~~~~~~~~~~~~~~
    596          call iniaerosol
     596         call aerosol_init
    597597         ! allocate related local arrays
    598598         ! (need be allocated instead of automatic because of "naerkind")
     
    650650!        Initialize orbital calculation.
    651651!        ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    652          call iniorbit(apoastr,periastr,year_day,peri_day,obliquit)
     652         call ephemeris_orbit_init(apoastr,periastr,year_day,peri_day,obliquit)
    653653
    654654
     
    799799
    800800         call su_watercycle ! even if we don't have a water cycle, we might
    801                             ! need epsi for the wvp definitions in callcorrk.F
     801                            ! need epsi for the wvp definitions in rad_correlatedk.F
    802802                            ! or RETV, RLvCp for the thermal plume model
    803803
     
    812812         if (corrk) then
    813813            ! We initialise the spectral grid here instead of
    814             ! at firstcall of callcorrk so we can output XspecIR, XspecVI
     814            ! at firstcall of rad_correlatedk so we can output XspecIR, XspecVI
    815815            ! when using Dynamico
    816816            print*, "physiq_mod: Correlated-k data base folder:",trim(datadir)
     
    821821            banddir=trim(trim(adjustl(tmp1))//'x'//trim(adjustl(tmp2)))
    822822            banddir=trim(trim(adjustl(corrkdir))//'/'//trim(adjustl(banddir)))
    823             call setspi !Basic infrared properties.
    824             call setspv ! Basic visible properties.
    825             call sugas_corrk       ! Set up gaseous absorption properties.
    826             call suaer_corrk       ! Set up aerosol optical properties.
     823            call rad_correlatedk_init_thermal !Basic infrared properties.
     824            call rad_correlatedk_init_stellar ! Basic visible properties.
     825            call rad_correlatedk_read_opacity_tables        ! Set up gaseous absorption properties.
     826            call rad_correlatedk_ini_aerosol       ! Set up aerosol optical properties.
    827827         endif
    828828
     
    877877      ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    878878      if (season) then
    879          call stellarlong(zday,zls)
     879         call ephemeris_stellar_longitude(zday,zls)
    880880      else
    881          call stellarlong(noseason_day,zls)
     881         call ephemeris_stellar_longitude(noseason_day,zls)
    882882      end if
    883883
    884       call orbite(zls,dist_star,declin,right_ascen)
     884      call ephemeris_orbit(zls,dist_star,declin,right_ascen)
    885885
    886886      if (tlocked) then
     
    981981         ztim3=COS(declin)*SIN(zlss)
    982982
    983          call stelang(ngrid,sinlon,coslon,sinlat,coslat,    &
     983         call ephemeris_stellar_angle(ngrid,sinlon,coslon,sinlat,coslat,    &
    984984                        ztim1,ztim2,ztim3,mu0,fract, flatten)
    985985
     
    989989         ztim3=-COS(declin)*SIN(2.*pi*(zday-.5))
    990990
    991          call stelang(ngrid,sinlon,coslon,sinlat,coslat,    &
     991         call ephemeris_stellar_angle(ngrid,sinlon,coslon,sinlat,coslat,    &
    992992                        ztim1,ztim2,ztim3,mu0,fract, flatten)
    993993      else if(diurnal .eqv. .false.) then
     
    10081008            ! Eclipse incoming sunlight (e.g. Saturn ring shadowing).
    10091009            if(rings_shadow) then
    1010                 call call_rings(ngrid, ptime, pday, diurnal)
     1010                call rad_ring_shadowing(ngrid, ptime, pday, diurnal)
    10111011            endif
    10121012
     
    10631063               endif !(ok_slab_ocean)
    10641064
    1065                ! standard callcorrk
     1065               ! standard rad_correlatedk
    10661066               clearsky=.false.
    1067                call callcorrk(ngrid,nlayer,pq,nq,qsurf,zls,                        &
     1067               call rad_correlatedk(ngrid,nlayer,pq,nq,qsurf,zls,                 &
    10681068                              albedo,albedo_equivalent,emis,mu0,pplev,pplay,pt,   &
    10691069                              tsurf,fract,dist_star,aerosol,muvar,                &
     
    10941094               if(CLFvarying)then
    10951095
    1096                   ! ---> PROBLEMS WITH ALLOCATED ARRAYS : temporary solution in callcorrk: do not deallocate if CLFvarying ...
     1096                  ! ---> PROBLEMS WITH ALLOCATED ARRAYS : temporary solution in rad_correlatedk: do not deallocate if CLFvarying ...
    10971097                  clearsky=.true.
    1098                   call callcorrk(ngrid,nlayer,pq,nq,qsurf,zls,                       &
     1098                  call rad_correlatedk(ngrid,nlayer,pq,nq,qsurf,zls,                 &
    10991099                                 albedo,albedo_equivalent1,emis,mu0,pplev,pplay,pt,  &
    11001100                                 tsurf,fract,dist_star,aerosol,muvar,                &
     
    11621162! II.b Call Newtonian cooling scheme
    11631163! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    1164                ! call newtrelax(ngrid,nlayer,mu0,sinlat,zpopsk,pt,pplay,pplev,dtrad,firstcall)
    1165                call newtcool_MOCHA(ngrid,nlayer,coslon,coslat,pt,pplay,firstcall,lastcall,dtrad)
     1164               ! call rad_netwon_cooling(ngrid,nlayer,mu0,sinlat,zpopsk,pt,pplay,pplev,dtrad,firstcall)
     1165               call rad_newton_cooling_MOCHA_intercomparison(ngrid,nlayer,coslon,coslat,pt,pplay,firstcall,lastcall,dtrad)
    11661166
    11671167               zdtsurf(1:ngrid) = +(pt(1:ngrid,1)-tsurf(1:ngrid))/ptimestep
     
    22992299         reffcol(1:ngrid,1:naerkind)=0.0
    23002300         if(co2cond.and.(iaero_co2.ne.0))then
    2301             call co2_reffrad(ngrid,nlayer,nq,zq,reffrad(1,1,iaero_co2))
     2301            call aerosol_radius_co2(ngrid,nlayer,nq,zq,reffrad(1,1,iaero_co2))
    23022302            do ig=1,ngrid
    23032303               reffcol(ig,iaero_co2) = SUM(zq(ig,1:nlayer,igcm_co2_ice)*reffrad(ig,1:nlayer,iaero_co2)*mass(ig,1:nlayer))
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_blackbody.F

    r4076 r4077  
    1       subroutine blackl(blalong,blat,blae)
     1      subroutine rad_blackbody_planck_law_wavelength(blalong,blat,blae)
    22
    33      implicit double precision (a-h,o-z)
     
    2121      end
    2222
    23       subroutine blackn(blalong,blat,blae)
     23      subroutine rad_blackbody_planck_law_wavenumber(blalong,blat,blae)
    2424
    2525      implicit double precision (a-h,o-z)
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk.F90

    r4076 r4077  
    1 MODULE callcorrk_mod
     1MODULE rad_correlatedk_mod
    22
    33IMPLICIT NONE
     
    55CONTAINS
    66
    7       subroutine callcorrk(ngrid,nlayer,pq,nq,qsurf,zls,       &
     7      subroutine rad_correlatedk(ngrid,nlayer,pq,nq,qsurf,zls,       &
    88          albedo,albedo_equivalent,emis,mu0,pplev,pplay,pt,    &
    99          tsurf,fract,dist_star,aerosol,muvar,                 &
     
    2727      use ioipsl_getin_p_mod, only: getin_p
    2828      use gases_h, only: ngasmx
    29       use radii_mod, only : su_aer_radii,co2_reffrad,h2o_reffrad,dust_reffrad,h2so4_reffrad,back2lay_reffrad
    30       use aerosol_mod, only : iaero_co2,iaero_h2o,iaero_dust,iaero_h2so4, &
     29      use aerosol_radius, only : aerosol_radius_init,aerosol_radius_co2,h2o_reffrad, &
     30                    aerosol_radius_dust,aerosol_radius_h2so4,aerosol_radius_back2lay
     31      use aerosol_global_variables , only : iaero_co2,iaero_h2o,iaero_dust,iaero_h2so4, &
    3132                              iaero_back2lay, iaero_aurora,               &
    3233                              iaero_venus1,iaero_venus2,iaero_venus2p,    &
    3334                              iaero_venus3,iaero_venusUV
    34       use aeropacity_mod, only: aeropacity
    35       use aeroptproperties_mod, only: aeroptproperties
     35      use aerosol_opacity_mod, only: aerosol_opacity
     36      use aerosol_optical_properties_mod, only: aerosol_optical_properties
    3637      use tracer_h, only: igcm_h2o_ice, igcm_h2o_vap, igcm_co2_ice
    3738      use tracer_h, only: constants_epsi_generic
     
    4142                              CLFvarying,tplanckmin,tplanckmax,global1d, &
    4243                              generic_condensation, aerovenus, nvarlayer, varspec
    43       use optcv_mod, only: optcv
    44       use optci_mod, only: optci
    45       use sfluxi_mod, only: sfluxi
    46       use sfluxv_mod, only: sfluxv
    47       use recombin_corrk_mod, only: corrk_recombin, call_recombin
     44      use rad_correlatedk_opacities_stellar_mod, &
     45        only: rad_correlatedk_opacities_stellar
     46      use rad_correlatedk_opacities_thermal_mod, &
     47        only: rad_correlatedk_opacities_thermal
     48      use rad_correlatedk_fluxes_thermal_mod, &
     49        only: rad_correlatedk_fluxes_thermal
     50      use rad_correlatedk_fluxes_stellar_mod, &
     51        only: rad_correlatedk_fluxes_stellar
     52      use rad_correlatedk_online_recombination_mod, &
     53        only: corrk_recombin, &
     54        rad_correlatedk_recombination_main
    4855      use pindex_mod, only: pindex
    4956      use generic_cloud_common_h, only: Psat_generic, epsi_generic
     
    111118      REAL,INTENT(OUT) :: OSR_nu(ngrid,L_NSPECTV)        ! Outgoing SW radiation in each band (Normalized to the band width (W/m2/cm-1).
    112119      REAL,INTENT(OUT) :: GSR_nu(ngrid,L_NSPECTV)        ! Surface SW radiation in each band (Normalized to the band width (W/m2/cm-1).
    113       REAL,INTENT(OUT) :: tau_col(ngrid)                 ! Diagnostic from aeropacity.
     120      REAL,INTENT(OUT) :: tau_col(ngrid)                 ! Diagnostic from aerosol_opacity.
    114121      REAL,INTENT(OUT) :: albedo_equivalent(ngrid)       ! Spectrally Integrated Albedo. For Diagnostic. By MT2015
    115122      REAL,INTENT(OUT) :: totcloudfrac(ngrid)            ! Column Fraction of clouds (%).
     
    148155      REAL*8 tlevrad(L_LEVELS),plevrad(L_LEVELS)
    149156
    150       ! Optical values for the optci/cv subroutines
     157      ! Optical values for the rad_correlatedk_opacities_thermal/cv subroutines
    151158      REAL*8 stel(L_NSPECTV),stel_fract(L_NSPECTV)
    152159      ! NB: Arrays below are "save" to avoid reallocating them at every call
     
    200207      character(len=10) :: tmp2
    201208      character(len=100) :: message
    202       character(len=10),parameter :: subname="callcorrk"
     209      character(len=10),parameter :: subname="rad_correlatedk"
    203210
    204211      ! For fixed water vapour profiles.
     
    242249      if(firstcall) then
    243250
    244         ! test on allocated necessary because of CLFvarying (two calls to callcorrk in physiq)
     251        ! test on allocated necessary because of CLFvarying (two calls to rad_correlatedk in physiq)
    245252        if(.not.allocated(QVISsQREF3d)) then
    246253          allocate(QVISsQREF3d(ngrid,nlayer,L_NSPECTV,naerkind))
     
    308315        endif
    309316
    310          !!! ALLOCATED instances are necessary because of CLFvarying (strategy to call callcorrk twice in physiq...)
     317         !!! ALLOCATED instances are necessary because of CLFvarying (strategy to call rad_correlatedk twice in physiq...)
    311318         IF(.not.ALLOCATED(QREFvis3d))THEN
    312319           ALLOCATE(QREFvis3d(ngrid,nlayer,naerkind), stat=ok)
     
    343350#endif
    344351
    345          call su_aer_radii(ngrid,nlayer,reffrad,nueffrad)
     352         call aerosol_radius_init(ngrid,nlayer,reffrad,nueffrad)
    346353         
    347354         
     
    351358
    352359      !this block is now done at firstcall of physiq_mod
    353          ! print*, "callcorrk: Correlated-k data base folder:",trim(datadir)
     360         ! print*, "rad_correlatedk: Correlated-k data base folder:",trim(datadir)
    354361         ! call getin_p("corrkdir",corrkdir)
    355362         ! print*, "corrkdir = ",corrkdir
     
    359366         ! banddir=trim(adjustl(corrkdir))//'/'//trim(adjustl(banddir))
    360367
    361          ! call setspi            ! Basic infrared properties.
    362          ! call setspv            ! Basic visible properties.
    363          ! call sugas_corrk       ! Set up gaseous absorption properties.
    364          ! call suaer_corrk       ! Set up aerosol optical properties.
     368         ! call rad_correlatedk_init_thermal            ! Basic infrared properties.
     369         ! call rad_correlatedk_init_stellar            ! Basic visible properties.
     370         ! call rad_correlatedk_read_opacity_tables        ! Set up gaseous absorption properties.
     371         ! call rad_correlatedk_ini_aerosol       ! Set up aerosol optical properties.
    365372       
    366373
    367          ! now that L_NGAUSS has been initialized (by sugas_corrk)
     374         ! now that L_NGAUSS has been initialized (by rad_correlatedk_read_opacity_tables )
    368375         ! allocate related arrays
    369376         if(.not.allocated(dtaui)) then
     
    446453
    447454         if((igcm_h2o_vap.eq.0) .and. varactive .and. water)then
    448             message='varactive in callcorrk but no h2o_vap tracer.'
     455            message='varactive in rad_correlatedk but no h2o_vap tracer.'
    449456            call abort_physic(subname,message,1)
    450457         endif
     
    485492
    486493         if ((iaer.eq.iaero_co2).and.tracer.and.(igcm_co2_ice.gt.0)) then ! Treat condensed co2 particles.
    487             call co2_reffrad(ngrid,nlayer,nq,pq,reffrad(1,1,iaero_co2))
     494            call aerosol_radius_co2(ngrid,nlayer,nq,pq,reffrad(1,1,iaero_co2))
    488495
    489496            call planetwide_maxval(reffrad(:,:,iaero_co2),maxvalue)
     
    517524         
    518525         if(iaer.eq.iaero_dust)then
    519             call dust_reffrad(ngrid,nlayer,reffrad(1,1,iaero_dust))
     526            call aerosol_radius_dust(ngrid,nlayer,reffrad(1,1,iaero_dust))
    520527            if (is_master) then
    521528               print*,'Dust particle size = ',reffrad(1,1,iaer)/1.e-6,' um'
     
    524531         
    525532         if(iaer.eq.iaero_h2so4)then
    526             call h2so4_reffrad(ngrid,nlayer,reffrad(1,1,iaero_h2so4))
     533            call aerosol_radius_h2so4(ngrid,nlayer,reffrad(1,1,iaero_h2so4))
    527534            if (is_master) then
    528535               print*,'H2SO4 particle size =',reffrad(1,1,iaer)/1.e-6,' um'
     
    531538         
    532539          if(iaer.eq.iaero_back2lay)then
    533             call back2lay_reffrad(ngrid,reffrad(1,1,iaero_back2lay),nlayer,pplev)
    534          endif
    535 
    536          !  For n-layer aerosol size set once for all at firstcall in su_aer_radii
     540            call aerosol_radius_back2lay(ngrid,reffrad(1,1,iaero_back2lay),nlayer,pplev)
     541         endif
     542
     543         !  For n-layer aerosol size set once for all at firstcall in aerosol_radius_init
    537544
    538545!         if(iaer.eq.iaero_aurora)then
     
    549556
    550557      ! Get 3D aerosol optical properties.
    551       call aeroptproperties(ngrid,nlayer,reffrad,nueffrad,         &
     558      call aerosol_optical_properties(ngrid,nlayer,reffrad,nueffrad,         &
    552559           QVISsQREF3d,omegaVIS3d,gVIS3d,                          &
    553560           QIRsQREF3d,omegaIR3d,gIR3d,                             &
     
    555562
    556563      ! Get aerosol optical depths.
    557       call aeropacity(ngrid,nlayer,nq,pplay,pplev,pt,pq,zls,aerosol,      &
     564      call aerosol_opacity(ngrid,nlayer,nq,pplay,pplev,pt,pq,zls,aerosol,      &
    558565           reffrad,nueffrad,QREFvis3d,QREFir3d,                             &
    559566           tau_col,cloudfrac,totcloudfrac,clearsky)               
     
    692699                       (pplev(ig,L_NLAYRAD-k)-pplev(ig,L_NLAYRAD-k+1))
    693700               ! As 'aerosol' is at reference (visible) wavelenght we scale it as
    694                ! it will be multplied by qxi/v in optci/v
     701               ! it will be multplied by qxi/v in rad_correlatedk_opacities_thermal/v
    695702               temp=aerosol(ig,L_NLAYRAD-k,iaer)/QREFvis3d(ig,L_NLAYRAD-k,iaer)
    696703               tauaero(2*k+2,iaer)=max(temp*pweight,0.d0)
     
    918925            qvar(1)=qvar(2)
    919926
    920             write(*,*)trim(subname),' :Warning: reducing qvar in callcorrk.'
     927            write(*,*)trim(subname),' :Warning: reducing qvar in rad_correlatedk.'
    921928            write(*,*)trim(subname),' :Temperature profile no longer consistent ', &
    922929                   'with saturated H2O. qsat=',satval
     
    10181025      ! -- JVO 20 : Also add a sanity test checking that tlevrad is
    10191026      !             within Planck function temperature boundaries,
    1020       !             which would cause gfluxi/sfluxi to crash.
     1027      !             which would cause rad_correlatedk_fluxes_solver_thermal/rad_correlatedk_fluxes_thermal to crash.
    10211028      do k=1,L_LEVELS
    10221029
     
    11271134! Recombine reference corrk tables if needed - Added by JVO, 2020.
    11281135         if (corrk_recombin) then
    1129            call call_recombin(ig,nlayer,pq(ig,:,:),pplay(ig,:),pt(ig,:),qvar(:),tmid(:),pmid(:))
     1136           call rad_correlatedk_recombination_main(ig,nlayer,pq(ig,:,:),pplay(ig,:),pt(ig,:),qvar(:),tmid(:),pmid(:))
    11301137         endif
    11311138! ----------------------------------------------------------------
     
    11491156            endif
    11501157
    1151             call optcv(dtauv,tauv,taucumv,plevrad,                 &
    1152                  qxvaer,qsvaer,gvaer,wbarv,cosbv,tauaero,          &
     1158            call rad_correlatedk_opacities_stellar(dtauv,    &
     1159                 tauv,taucumv,plevrad,                       &
     1160                 qxvaer,qsvaer,gvaer,wbarv,cosbv,tauaero,    &
    11531161                 tmid,pmid,taugsurf,qvar,muvarrad,fracvari)
    11541162
    1155             call sfluxv(dtauv,tauv,taucumv,albv,dwnv,wbarv,cosbv,  &
    1156                  acosz,stel_fract,                                 &
    1157                  nfluxtopv,fluxtopvdn,nfluxoutv_nu,nfluxgndv_nu,   &
     1163            call rad_correlatedk_fluxes_stellar(dtauv,tauv,  &
     1164                 taucumv,albv,dwnv,wbarv,cosbv,              &
     1165                 acosz,stel_fract,                           &
     1166                 nfluxtopv,fluxtopvdn,                       &
     1167                 nfluxoutv_nu,nfluxgndv_nu,                  &
    11581168                 fmnetv,fluxupv,fluxdnv,fzerov,taugsurf)
    11591169
     
    11921202!-----------------------------------------------------------------------
    11931203
    1194          call optci(plevrad,tlevrad,dtaui,taucumi,                  &
    1195               qxiaer,qsiaer,giaer,cosbi,wbari,tauaero,tmid,pmid,    &
     1204         call rad_correlatedk_opacities_thermal(plevrad,          &
     1205              tlevrad,dtaui,taucumi,                              &
     1206              qxiaer,qsiaer,giaer,cosbi,wbari,tauaero,tmid,pmid,  &
    11961207              taugsurfi,qvar,muvarrad,fracvari)
    11971208
    1198          call sfluxi(plevrad,tlevrad,dtaui,taucumi,ubari,albi,      &
    1199               wnoi,dwni,cosbi,wbari,nfluxtopi,nfluxtopi_nu,         &
     1209         call rad_correlatedk_fluxes_thermal(plevrad,             &
     1210              tlevrad,dtaui,taucumi,ubari,albi,                   &
     1211              wnoi,dwni,cosbi,wbari,nfluxtopi,nfluxtopi_nu,       &
    12001212              fmneti,fluxupi,fluxdni,fluxupi_nu,fzeroi,taugsurfi)
    12011213
     
    13461358
    13471359
    1348     end subroutine callcorrk
    1349 
    1350 END MODULE callcorrk_mod
     1360    end subroutine rad_correlatedk
     1361
     1362END MODULE rad_correlatedk_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_continuum_interpolation.F90

    r4076 r4077  
    1 module interpolate_continuum_mod
     1module rad_correlatedk_continuum_interpolation_mod
    22
    33implicit none
     
    55contains
    66
    7      subroutine interpolate_continuum(filename,igas_X,igas_Y,c_WN,ind_WN,temp,pres_X,pres_Y,abs_coef,firstcall)
     7     subroutine rad_correlatedk_continuum_interpolation(filename,igas_X,igas_Y,c_WN,ind_WN,temp,pres_X,pres_Y,abs_coef,firstcall)
    88
    99!==================================================================
     
    5353
    5454      character(len=512) :: line
    55       character(len=21),parameter :: rname="interpolate_continuum"
     55      character(len=21),parameter :: rname="rad_correlatedk_continuum_interpolation"
    5656
    5757      integer i, pos, iT, iW, iB, count_norm, igas
     
    206206
    207207
    208       if(firstcall)then ! called by sugas_corrk only
     208      if(firstcall)then ! called by rad_correlatedk_read_opacity_tables only
    209209        if (is_master) print*,'----------------------------------------------------'
    210         if (is_master) print*,'Initialising continuum (interpolate_continuum routine) from ', trim(filename)
     210        if (is_master) print*,'Initialising continuum (rad_correlatedk_continuum_interpolation routine) from ', trim(filename)
    211211
    212212!$OMP MASTER
     
    216216        if (ios.ne.0) then        ! file not found
    217217          if (is_master) then
    218             write(*,*) 'Error from interpolate_continuum routine'
     218            write(*,*) 'Error from rad_correlatedk_continuum_interpolation routine'
    219219            write(*,*) 'Data file ',trim(filename),' not found.'
    220220            write(*,*) 'Check that your path to datagcm:',trim(datadir)
     
    650650      !print*,'So the absorption is ',abs_coef,' m^-1'
    651651     
    652     end subroutine interpolate_continuum
     652    end subroutine rad_correlatedk_continuum_interpolation
    653653   
    654654   
     
    836836    end subroutine interpolate_T_abs_coeff
    837837
    838 end module interpolate_continuum_mod
     838end module rad_correlatedk_continuum_interpolation_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_fluxes_solver_stellar.F

    r4076 r4077  
    1       module gfluxv_mod
     1      module rad_correlatedk_fluxes_solver_stellar_mod
    22     
    33      implicit none
     
    55      contains
    66     
    7       SUBROUTINE GFLUXV(DTDEL,TDEL,TAUCUMIN,WDEL,CDEL,UBAR0,F0PI,RSF,
    8      *                  BTOP,BSURF,FMIDP,FMIDM,DIFFV,FLUXUP,FLUXDN)
    9 
     7      SUBROUTINE rad_correlatedk_fluxes_solver_stellar(DTDEL,
     8     *    TDEL,TAUCUMIN,WDEL,CDEL,UBAR0,F0PI,RSF,BTOP,
     9     *    BSURF,FMIDP,FMIDM,DIFFV,FLUXUP,FLUXDN)
    1010
    1111C  THIS SUBROUTINE TAKES THE OPTICAL CONSTANTS AND BOUNDARY CONDITIONS
     
    1919C THE FLUXES DIRECTLY USING THE GENERALIZED NOTATION OF MEADOR AND WEAVOR
    2020C J.A.S., 37, 630-642, 1980.
    21 C THE TRI-DIAGONAL MATRIX SOLVER IS DSOLVER AND IS DOUBLE PRECISION SO MANY
    22 C VARIABLES ARE PASSED AS SINGLE THEN BECOME DOUBLE IN DSOLVER
     21C THE TRI-DIAGONAL MATRIX SOLVER IS rad_tridiagonal_matrix_solver AND IS DOUBLE PRECISION SO MANY
     22C VARIABLES ARE PASSED AS SINGLE THEN BECOME DOUBLE IN rad_tridiagonal_matrix_solver
    2323C
    2424C NLL           = NUMBER OF LEVELS (NAYER + 1) THAT WILL BE SOLVED
     
    200200      END DO
    201201
    202       CALL DSOLVER(NAYER,GAMA,CP,CM,CPM1,CMM1,E1,E2,E3,E4,BTOP,
    203      *             BSURF,RSF,XK1,XK2)
     202      CALL rad_tridiagonal_matrix_solver(NAYER,GAMA,CP,CM,CPM1,
     203     *             CMM1,E1,E2,E3,E4,BTOP,BSURF,RSF,XK1,XK2)
    204204
    205205C     NOW WE CALCULATE THE FLUXES AT THE MIDPOINTS OF THE LAYERS.
     
    335335
    336336
    337       END SUBROUTINE GFLUXV
    338 
    339       end module gfluxv_mod
     337      END SUBROUTINE rad_correlatedk_fluxes_solver_stellar
     338
     339      end module rad_correlatedk_fluxes_solver_stellar_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_fluxes_solver_thermal.F

    r4076 r4077  
    1       module gfluxi_mod
     1      module rad_correlatedk_fluxes_solver_thermal_mod
    22     
    33      implicit none
     
    55      contains
    66     
    7       SUBROUTINE GFLUXI(NLL,TLEV,NW,DW,DTAU,TAUCUM,W0,COSBAR,UBARI,
     7      SUBROUTINE rad_correlatedk_fluxes_solver_thermal(NLL,
     8     *                  TLEV,NW,DW,DTAU,TAUCUM,W0,COSBAR,UBARI,
    89     *                  RSF,BTOP,BSURF,FTOPUP,FMIDP,FMIDM)
    910     
     
    2324!  HAS LEVEL N ON TOP AND LEVEL N+1 ON BOTTOM. OPTICAL DEPTH INCREASES
    2425!  FROM TOP TO BOTTOM.  SEE C.P. MCKAY, TGM NOTES.
    25 !  THE TRI-DIAGONAL MATRIX SOLVER IS DSOLVER AND IS DOUBLE PRECISION SO MANY
    26 !  VARIABLES ARE PASSED AS SINGLE THEN BECOME DOUBLE IN DSOLVER
     26!  THE TRI-DIAGONAL MATRIX SOLVER IS rad_tridiagonal_matrix_solver AND IS DOUBLE PRECISION SO MANY
     27!  VARIABLES ARE PASSED AS SINGLE THEN BECOME DOUBLE IN rad_tridiagonal_matrix_solver
    2728!
    2829! NLL            = NUMBER OF LEVELS (NLAYERS + 1) MUST BE LESS THAT NL (101)
     
    8182!         open(888,file='W0')
    8283!           if ((W0(L).eq.0.).or.(W0(L).eq.1.)) then
    83 !             write(888,*) W0(L), L, 'gfluxi'
     84!             write(888,*) W0(L), L, 'rad_correlatedk_fluxes_solver_thermal'
    8485!           endif
    8586! Prevent this with an if statement:
     
    170171      END DO
    171172     
    172 !      B81=BTOP  ! RENAME BEFORE CALLING DSOLVER - used to be to set
     173!      B81=BTOP  ! RENAME BEFORE CALLING rad_tridiagonal_matrix_solver - used to be to set
    173174!      B82=BSURF ! them to real*8 - but now everything is real*8
    174175!      R81=RSF   ! so this may not be necessary
     
    176177! DOUBLE PRECISION TRIDIAGONAL SOLVER
    177178     
    178       CALL DSOLVER(NLAYER,GAMA,CP,CM,CPM1,CMM1,E1,E2,E3,E4,BTOP,
    179      *             BSURF,RSF,XK1,XK2)
     179      CALL rad_tridiagonal_matrix_solver(NLAYER,GAMA,CP,CM,CPM1,
     180     *             CMM1,E1,E2,E3,E4,BTOP,BSURF,RSF,XK1,XK2)
    180181     
    181182! NOW WE CALCULATE THE FLUXES AT THE MIDPOINTS OF THE LAYERS.
     
    242243     
    243244     
    244       END SUBROUTINE GFLUXI
    245 
    246       end module gfluxi_mod
     245      END SUBROUTINE rad_correlatedk_fluxes_solver_thermal
     246
     247      end module rad_correlatedk_fluxes_solver_thermal_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_fluxes_stellar.F

    r4076 r4077  
    1       module sfluxv_mod
     1      module rad_correlatedk_fluxes_stellar_mod
    22
    33      implicit none
     
    55      contains
    66     
    7       SUBROUTINE SFLUXV(DTAUV,TAUV,TAUCUMV,RSFV,DWNV,WBARV,COSBV,
     7      SUBROUTINE rad_correlatedk_fluxes_stellar(DTAUV,TAUV,
     8     *                  TAUCUMV,RSFV,DWNV,WBARV,COSBV,
    89     *                  UBAR0,STEL,NFLUXTOPV,FLUXTOPVDN,
    910     *                  NFLUXOUTV_nu,NFLUXGNDV_nu,
     
    1314      use radinc_h, only: L_NLAYRAD, L_NLEVRAD
    1415      use radcommon_h, only: tlimit, gweight
    15       use gfluxv_mod, only: gfluxv
     16      use rad_correlatedk_fluxes_solver_stellar_mod, only:
     17     *               rad_correlatedk_fluxes_solver_stellar
    1618
    1719      implicit none
     
    101103
    102104
    103           CALL GFLUXV(DTAUV(1,NW,NG),TAUV(1,NW,NG),TAUCUMV(1,NW,NG),
    104      *                WBARV(1,NW,NG),COSBV(1,NW,NG),UBAR0,F0PI,RSFV(NW),   
    105      *                BTOP,BSURF,FMUPV,FMDV,DIFFV,FLUXUP,FLUXDN)
     105          CALL rad_correlatedk_fluxes_solver_stellar(DTAUV(1,NW,NG),
     106     *           TAUV(1,NW,NG),TAUCUMV(1,NW,NG),
     107     *           WBARV(1,NW,NG),COSBV(1,NW,NG),UBAR0,F0PI,RSFV(NW),   
     108     *           BTOP,BSURF,FMUPV,FMDV,DIFFV,FLUXUP,FLUXDN)
    106109
    107110C         NOW CALCULATE THE CUMULATIVE VISIBLE NET FLUX
     
    159162C       RETURN FLUXES FOR A GIVEN NT
    160163
    161         CALL GFLUXV(DTAUV(1,NW,NG),TAUV(1,NW,NG),TAUCUMV(1,NW,NG),
    162      *              WBARV(1,NW,NG),COSBV(1,NW,NG),UBAR0,F0PI,RSFV(NW),
    163      *              BTOP,BSURF,FMUPV,FMDV,DIFFV,FLUXUP,FLUXDN)
     164        CALL rad_correlatedk_fluxes_solver_stellar(DTAUV(1,NW,NG),
     165     *          TAUV(1,NW,NG),TAUCUMV(1,NW,NG),
     166     *          WBARV(1,NW,NG),COSBV(1,NW,NG),UBAR0,F0PI,RSFV(NW),
     167     *          BTOP,BSURF,FMUPV,FMDV,DIFFV,FLUXUP,FLUXDN)
    164168
    165169
     
    193197
    194198
    195       END SUBROUTINE SFLUXV
    196 
    197       end module sfluxv_mod
    198      
     199      END SUBROUTINE rad_correlatedk_fluxes_stellar
     200
     201      end module rad_correlatedk_fluxes_stellar_mod
     202     
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_fluxes_thermal.F

    r4076 r4077  
    1       module sfluxi_mod
     1      module rad_correlatedk_fluxes_thermal_mod
    22     
    33      implicit none
     
    55      contains
    66     
    7       SUBROUTINE SFLUXI(PLEV,TLEV,DTAUI,TAUCUMI,UBARI,RSFI,WNOI,DWNI,
     7      SUBROUTINE rad_correlatedk_fluxes_thermal(PLEV,TLEV,
     8     *                  DTAUI,TAUCUMI,UBARI,RSFI,WNOI,DWNI,
    89     *                  COSBI,WBARI,NFLUXTOPI,NFLUXTOPI_nu,
    910     *                  FMNETI,fluxupi,fluxdni,fluxupi_nu,
     
    1415      use radcommon_h, only: planckir, tlimit,sigma, gweight
    1516      use comcstfi_mod, only: pi
    16       use gfluxi_mod, only: gfluxi
     17      use rad_correlatedk_fluxes_solver_thermal_mod,
     18     * only: rad_correlatedk_fluxes_solver_thermal
    1719     
    1820      implicit none
     
    132134! WITHIN EACH INTERVAL AT THE MIDPOINT WAVENUMBER
    133135           
    134             CALL GFLUXI(NLEVRAD,TLEV,NW,DWNI(NW),DTAUI(1,NW,NG),
     136            CALL rad_correlatedk_fluxes_solver_thermal(NLEVRAD,
     137     *                TLEV,NW,DWNI(NW),DTAUI(1,NW,NG),
    135138     *                TAUCUMI(1,NW,NG),
    136139     *                WBARI(1,NW,NG),COSBI(1,NW,NG),UBARI,RSFI,BTOP,
     
    175178! WITHIN EACH INTERVAL AT THE MIDPOINT WAVENUMBER
    176179         
    177          CALL GFLUXI(NLEVRAD,TLEV,NW,DWNI(NW),DTAUI(1,NW,NG),
     180         CALL rad_correlatedk_fluxes_solver_thermal(NLEVRAD,
     181     *                TLEV,NW,DWNI(NW),DTAUI(1,NW,NG),
    178182     *                TAUCUMI(1,NW,NG),
    179183     *                WBARI(1,NW,NG),COSBI(1,NW,NG),UBARI,RSFI,BTOP,
     
    200204! *** END OF MAJOR SPECTRAL INTERVAL LOOP IN THE INFRARED****
    201205     
    202       END SUBROUTINE SFLUXI
     206      END SUBROUTINE rad_correlatedk_fluxes_thermal
    203207
    204       end module sfluxi_mod
     208      end module rad_correlatedk_fluxes_thermal_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_ini_aerosol.F90

    r4076 r4077  
    1 module suaer_corrk_mod
     1module rad_correlatedk_ini_aerosol_mod
    22
    33implicit none
     
    55contains
    66
    7 subroutine suaer_corrk
     7subroutine rad_correlatedk_ini_aerosol
    88
    99      ! inputs
     
    1616      use radcommon_h, only: radiustab,nsize,tstellar
    1717      use radcommon_h, only: qrefvis,qrefir,omegarefir !,omegarefvis
    18       use aerosol_mod, only: noaero,iaero_co2,iaero_h2o,iaero_dust,iaero_h2so4
    19       use aerosol_mod, only: iaero_back2lay,iaero_nh3,iaero_nlay,iaero_aurora
    20       use aerosol_mod, only: iaero_venus1,iaero_venus2,iaero_venus2p
    21       use aerosol_mod, only: iaero_venus3,iaero_venusUV
    22       use aerosol_mod, only: iaero_generic,i_rgcs_ice
     18      use aerosol_global_variables , only: noaero,iaero_co2,iaero_h2o,iaero_dust,iaero_h2so4
     19      use aerosol_global_variables , only: iaero_back2lay,iaero_nh3,iaero_nlay,iaero_aurora
     20      use aerosol_global_variables , only: iaero_venus1,iaero_venus2,iaero_venus2p
     21      use aerosol_global_variables , only: iaero_venus3,iaero_venusUV
     22      use aerosol_global_variables , only: iaero_generic,i_rgcs_ice
    2323      use callkeys_mod, only: tplanet, optprop_back2lay_vis, optprop_back2lay_ir, &
    2424                              optprop_aeronlay_vis, optprop_aeronlay_ir,          &
     
    8989      CHARACTER(LEN=132) :: scanline ! ASCII file scanning line
    9090
    91 !     I/O  of "aerave" (subroutine that spectrally averages
     91!     I/O  of "aerosol_optical_properties_averaging" (subroutine that spectrally averages
    9292!     the single scattering parameters)
    9393
     
    342342            else
    343343! If you want to add another specie, copy,paste & adapt the elseif block up here to your new specie (LT 2022)
    344                call abort_physic("suaer_corrk", "Unknown specie in radiative generic condensable species",1)
     344               call abort_physic("rad_correlatedk_ini_aerosol", "Unknown specie in radiative generic condensable species",1)
    345345            endif
    346346         endif
     
    393393            ENDIF
    394394            IF(.NOT.file_ok) THEN
    395                write(*,*)'suaer_corrk: Problem opening ',&
     395               write(*,*)'rad_correlatedk_ini_aerosol: Problem opening ',&
    396396               TRIM(file_id(iaer,idomain))
    397397               write(*,*)'It should be in: ',TRIM(datadir)//'/'//TRIM(aerdir)
     
    405405               write(*,*)' http://www.lmd.jussieu.fr/',&
    406406               '~lmdz/planets/generic/datagcm/'
    407                CALL abort_physic("suaer_corrk", "Unable to read file",1)
     407               CALL abort_physic("rad_correlatedk_ini_aerosol", "Unable to read file",1)
    408408            ENDIF
    409409
     
    425425               endwhile = .true.
    426426            CASE DEFAULT reading1_seq ! ----------------------------
    427                CALL abort_physic("suaer_corrk","Error while loading optical properties",1)
     427               CALL abort_physic("rad_correlatedk_ini_aerosol","Error while loading optical properties",1)
    428428            END SELECT reading1_seq ! ==============================
    429429         ENDIF
     
    505505         endwhile = .true.
    506506      CASE DEFAULT reading2_seq ! ----------------------------
    507          CALL abort_physic("suaer_corrk","Error while loading optical properties",1)
     507         CALL abort_physic("rad_correlatedk_ini_aerosol","Error while loading optical properties",1)
    508508      END SELECT reading2_seq   ! ==============================
    509509      ENDIF
     
    551551!     Averaged optical properties (GCM channels)
    552552
    553             CALL aerave_new ( nwvl,&
     553            CALL aerosol_optical_properties_averaging ( nwvl,&
    554554            wvl(:),ep(:,isize),omeg(:,isize),gfactor(:,isize),&
    555555            lamref,epref,tstellar,&
     
    580580
    581581!     epav is <QIR>/Qext(lamrefir) since epref=1
    582 !     Note: aerave also computes the extinction coefficient at
     582!     Note: aerosol_optical_properties_averaging also computes the extinction coefficient at
    583583!     the reference wavelength. This is called QREFvis or QREFir
    584584!     (not epref, which is a different parameter).
     
    586586!     radcommon_h.F90
    587587
    588             CALL aerave_new ( nwvl,&
     588            CALL aerosol_optical_properties_averaging ( nwvl,&
    589589            wvl(:),ep(:,isize),omeg(:,isize),gfactor(:,isize),&
    590590            lamref,epref,tplanet,&
     
    622622  deallocate(file_id)
    623623
    624 end subroutine suaer_corrk
     624end subroutine rad_correlatedk_ini_aerosol
    625625     
    626 end module suaer_corrk_mod
     626end module rad_correlatedk_ini_aerosol_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_init_stellar.F90

    r4076 r4077  
    1       module setspv_mod
     1      module rad_correlatedk_init_stellar_mod
    22     
    33      implicit none
     
    55      contains
    66     
    7       subroutine setspv
     7      subroutine rad_correlatedk_init_stellar
    88
    99!==================================================================
     
    1515!     Authors
    1616!     -------
    17 !     Adapted from setspv in the NASA Ames radiative code by
     17!     Adapted from rad_correlatedk_init_stellar in the NASA Ames radiative code by
    1818!     Robin Wordsworth (2009).
    1919!
    2020!     Called by
    2121!     ---------
    22 !     callcorrk.F
     22!     rad_correlatedk.F
    2323!     
    2424!     Calls
    2525!     -----
    26 !     ave_stelspec.F
     26!     rad_correlatedk_stellar_spectrum .F
    2727!     
    2828!==================================================================
     
    3333      use datafile_mod, only: datadir
    3434      use callkeys_mod, only: Fat1AU
    35       use ave_stelspec_mod, only: ave_stelspec
     35      use rad_correlatedk_stellar_spectrum_mod, only: rad_correlatedk_stellar_spectrum
    3636
    3737      implicit none
     
    6666      if(.not.file_ok) then
    6767         write(*,*)'The file ',TRIM(file_path)
    68          write(*,*)'was not found by setspv.F90, exiting.'
     68         write(*,*)'was not found by rad_correlatedk_init_stellar.F90, exiting.'
    6969         write(*,*)'Check that your path to datagcm:',trim(datadir)
    7070         write(*,*)' is correct. You can change it in callphys.def with:'
    7171         write(*,*)' datadir = /absolute/path/to/datagcm'
    7272         write(*,*)'Also check that the corrkdir you chose in callphys.def exists.'
    73          call abort_physic("setspv", "Unable to read file",1)
     73         call abort_physic("rad_correlatedk_init_stellar", "Unable to read file",1)
    7474      endif
    7575       
     
    8686      close(131)
    8787
    88       write(*,*) 'setspv: L_NSPECTV = ',L_NSPECTV, 'in the model '
     88      write(*,*) 'rad_correlatedk_init_stellar: L_NSPECTV = ',L_NSPECTV, 'in the model '
    8989      write(*,*) '        there are   ',nb, 'entries in ',TRIM(file_path)
    9090      if(nb.ne.L_NSPECTV) then
    9191         write(*,*) 'MISMATCH !! I stop here'
    92          call abort_physic("setspv","The number of entries in narrowbands_VI.in does not match with L_NSPECTV",1)
     92         call abort_physic("rad_correlatedk_init_stellar","The number of entries in narrowbands_VI.in does not match with L_NSPECTV",1)
    9393      endif
    9494
     
    106106!$OMP BARRIER
    107107
    108       print*,'setspv: VI band limits:'
     108      print*,'rad_correlatedk_init_stellar: VI band limits:'
    109109      do M=1,L_NSPECTV+1
    110110         print*,m,'-->',BWNV(M),' cm^-1'
     
    127127!     Set up stellar spectrum
    128128
    129       write(*,*)'setspv: Interpolating stellar spectrum from the hires data...'
    130       call ave_stelspec(STELLAR)
     129      write(*,*)'rad_correlatedk_init_stellar: Interpolating stellar spectrum from the hires data...'
     130      call rad_correlatedk_stellar_spectrum (STELLAR)
    131131
    132132!     Sum the stellar flux, and write out the result. 
     
    136136         sum         = sum+STELLARF(N)
    137137      end do
    138       write(6,'("setspv: Stellar flux at 1 AU = ",f9.2," W m-2")') sum
     138      write(6,'("rad_correlatedk_init_stellar: Stellar flux at 1 AU = ",f9.2," W m-2")') sum
    139139      print*,' '
    140140
    141     END subroutine setspv
     141    END subroutine rad_correlatedk_init_stellar
    142142   
    143     end module setspv_mod
     143    end module rad_correlatedk_init_stellar_mod
    144144   
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_init_thermal.F90

    r4076 r4077  
    1       subroutine setspi
     1      subroutine rad_correlatedk_init_thermal
    22
    33!==================================================================
     
    99!     Authors
    1010!     -------
    11 !     Adapted from setspi in the NASA Ames radiative code by
     11!     Adapted from rad_correlatedk_init_thermal in the NASA Ames radiative code by
    1212!     Robin Wordsworth (2009).
    1313!     
    1414!     Called by
    1515!     ---------
    16 !     callcorrk.F
     16!     rad_correlatedk.F
    1717!     
    1818!     Calls
     
    8282      if(.not.file_ok) then
    8383         write(*,*)'The file ',TRIM(file_path)
    84          write(*,*)'was not found by setspi.F90, exiting.'
     84         write(*,*)'was not found by rad_correlatedk_init_thermal.F90, exiting.'
    8585         write(*,*)'Check that your path to datagcm:',trim(datadir)
    8686         write(*,*)' is correct. You can change it in callphys.def with:'
    8787         write(*,*)' datadir = /absolute/path/to/datagcm'
    8888         write(*,*)'Also check that the corrkdir you chose in callphys.def exists.'
    89          call abort_physic("setspi","Unable to read file",1)
     89         call abort_physic("rad_correlatedk_init_thermal","Unable to read file",1)
    9090      endif
    9191   
     
    9898      do while (ierr==0)
    9999        read(131,*,iostat=ierr) dummy
    100 !        write(*,*) 'setspi: file_entries:',dummy,'ierr=',ierr
     100!        write(*,*) 'rad_correlatedk_init_thermal: file_entries:',dummy,'ierr=',ierr
    101101        if (ierr==0) nb=nb+1
    102102      enddo
    103103      close(131)
    104104
    105       write(*,*) 'setspi: L_NSPECTI = ',L_NSPECTI, 'in the model '
     105      write(*,*) 'rad_correlatedk_init_thermal: L_NSPECTI = ',L_NSPECTI, 'in the model '
    106106      write(*,*) '        there are   ',nb, 'entries in ',TRIM(file_path)
    107107      if(nb.ne.L_NSPECTI) then
    108108         write(*,*) 'MISMATCH !! I stop here'
    109          call abort_physic("setspi","The number of entries in narrowbands_IR.in does not match with L_NSPECTI",1)
     109         call abort_physic("rad_correlatedk_init_thermal","The number of entries in narrowbands_IR.in does not match with L_NSPECTI",1)
    110110      endif
    111111
     
    124124
    125125      print*,''
    126       print*,'setspi: IR band limits:'
     126      print*,'rad_correlatedk_init_thermal: IR band limits:'
    127127      do M=1,L_NSPECTI+1
    128128         print*,m,'-->',BWNI(M),' cm^-1'
     
    148148
    149149      print*,''
    150       print*,'setspi: Current Planck integration range:'
     150      print*,'rad_correlatedk_init_thermal: Current Planck integration range:'
    151151      print*,'T = ',dble(NTstart)/NTfac, ' to ',dble(NTstop)/NTfac,' K.'
    152152
     
    182182      ! force planck=sigma*eps*T^4 for each temperature in array
    183183      if(forceEC)then
    184          print*,'setspi: Force F=sigma*eps*T^4 for all values of T!'
     184         print*,'rad_correlatedk_init_thermal: Force F=sigma*eps*T^4 for all values of T!'
    185185         do nt=NTstart,NTstop
    186186            plancksum=0.0D0
     
    207207            plancksum=plancksum+planckir(NW,nt-NTstart+1)*DWNI(NW)*pi
    208208         end do
    209          print*,'setspi: At lower limit:'
     209         print*,'rad_correlatedk_init_thermal: At lower limit:'
    210210         print*,'in model sig*T^4 = ',plancksum,' W m^-2'
    211211         print*,'actual sig*T^4   = ',sigma*(dble(nt)/NTfac)**4,' W m^-2'
     
    217217            plancksum=plancksum+planckir(NW,nt-NTstart+1)*DWNI(NW)*pi
    218218         end do
    219          print*,'setspi: At upper limit:'
     219         print*,'rad_correlatedk_init_thermal: At upper limit:'
    220220         print*,'in model sig*T^4 = ',plancksum,' W m^-2'
    221221         print*,'actual sig*T^4   = ',sigma*(dble(nt)/NTfac)**4,' W m^-2'
     
    224224
    225225      return
    226     end subroutine setspi
     226    end subroutine rad_correlatedk_init_thermal
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_online_recombination.F90

    r4076 r4077  
    1 MODULE recombin_corrk_mod
     1MODULE rad_correlatedk_online_recombination_mod
    22
    33    !
     
    55    !
    66    ! This module contains the following subroutines :
    7     ! - ini_recombin    : From modern traceur.def options check if we will use recombining and for which species. Called by initracer.
    8     ! - su_recombin     : Initialise tables. Called by sugas_corrk
    9     ! - call_recombin   : Test profile of species and decide whether to call recombin_corrk. Called by callcork
    10     ! - recombin_corrk  : The core algorithm properly recombining corrk tables. Called by callrecombin_corrk.
     7    ! - rad_correlatedk_recombination_init    : From modern traceur.def options check if we will use recombining and for which species. Called by initracer.
     8    ! - rad_correlatedk_recombination_setup     : Initialise tables. Called by rad_correlatedk_read_opacity_tables
     9    ! - rad_correlatedk_recombination_main   : Test profile of species and decide whether to call rad_correlatedk_online_recombination. Called by callcork
     10    ! - rad_correlatedk_online_recombination  : The core algorithm properly recombining corrk tables.
    1111    !
    1212 
     
    3636  !$OMP THREADPRIVATE(all_otf)
    3737 
    38     ! Following arrays are allocated in su_recombin (excepted otf2tot_idx, in ini_recombin) and deallocated in callcork lastcall
     38    ! Following arrays are allocated in rad_correlatedk_recombination_setup (excepted otf2tot_idx, in rad_correlatedk_recombination_init) and deallocated in callcork lastcall
    3939    REAL*8, save,  DIMENSION(:,:,:,:,:), ALLOCATABLE :: gasi_recomb, gasv_recomb
    4040    REAL*8, save,  DIMENSION(:,:,:,:,:), ALLOCATABLE :: gasi_otf, gasv_otf
     
    5252 
    5353 
    54     SUBROUTINE ini_recombin
     54    SUBROUTINE rad_correlatedk_recombination_init
    5555 
    5656      USE tracer_h
     
    101101           
    102102          write(*,*)
    103           write(*,*) 'ini_recombin: found specie : Name = ',trim(noms(iq)), &
     103          write(*,*) 'rad_correlatedk_recombination_init: found specie : Name = ',trim(noms(iq)), &
    104104                     ' ; Predefined vmr=', is_recomb_qset(iq),               &
    105105                     ' ; On-the-fly=',     is_recomb_qotf(iq)
     
    142142      ! Summary
    143143      write(*,*)
    144       write(*,*) 'ini_recombin: Total species found for corrk recombination =', nrecomb_tot
     144      write(*,*) 'rad_correlatedk_recombination_init: Total species found for corrk recombination =', nrecomb_tot
    145145     
    146146      if (corrk_recombin) then
    147147        if (use_premix) then
    148           write(*,*) 'ini_recombin: .. Total radiative species matching total species for recombination...'
    149           write(*,*) 'ini_recombin: .. Any pre-mixed set of opacities will be ignored.'
     148          write(*,*) 'rad_correlatedk_recombination_init: .. Total radiative species matching total species for recombination...'
     149          write(*,*) 'rad_correlatedk_recombination_init: .. Any pre-mixed set of opacities will be ignored.'
    150150        else
    151           write(*,*) 'ini_recombin: .. Found less species to recombine than total radiative species..'
    152           write(*,*) 'ini_recombin: .. Recombination will occur ontop of premix set of opacities'
     151          write(*,*) 'rad_correlatedk_recombination_init: .. Found less species to recombine than total radiative species..'
     152          write(*,*) 'rad_correlatedk_recombination_init: .. Recombination will occur ontop of premix set of opacities'
    153153        endif
    154154      else
    155         write(*,*) 'ini_recombin: .. No species found for recombination, I will use premix set only.'
     155        write(*,*) 'rad_correlatedk_recombination_init: .. No species found for recombination, I will use premix set only.'
    156156      endif
    157157      write(*,*)
    158158 
    159     END SUBROUTINE ini_recombin
     159    END SUBROUTINE rad_correlatedk_recombination_init
    160160     
    161161   
    162162   
    163     SUBROUTINE su_recombin
     163    SUBROUTINE rad_correlatedk_recombination_setup
    164164      USE radinc_h
    165165      USE radcommon_h, only: gweight, gasi, gasv
     
    178178      IF(.NOT. ALLOCATED(wtwospec))    ALLOCATE(wtwospec(L_NGAUSS*L_NGAUSS))   
    179179     
    180       ! Init. for recombin_corrk firstcall
     180      ! Init. for rad_correlatedk_online_recombination firstcall
    181181      permut_idx = (/(i, i=1,L_NGAUSS*L_NGAUSS)/)
    182182     
     
    221221      gasv_recomb(:,:,:,:,:) = gasv(:,:,1:L_REFVAR,:,:) ! non-zero init (=kappa_vi)
    222222 
    223     END SUBROUTINE su_recombin
    224  
    225  
    226  
    227     SUBROUTINE call_recombin(igrid,nlayer,pq,pplay,pt,qvar,tmid,pmid)
     223    END SUBROUTINE rad_correlatedk_recombination_setup
     224 
     225 
     226 
     227    SUBROUTINE rad_correlatedk_recombination_main(igrid,nlayer,pq,pplay,pt,qvar,tmid,pmid)
    228228 
    229229      USE comcstfi_mod, only: mugaz
     
    283283            if (is_recomb_qset(rcb2tot_idx(iq))==0 .and. is_recomb_qotf(rcb2tot_idx(iq))==0) then
    284284              print*, 'Recombining tracer ', noms(rcb2tot_idx(iq)),' requires an input profile, this is not implemented yet !!'
    285               call abort_physic('call_recombin','Missing implementation',1)
     285              call abort_physic('rad_correlatedk_recombination_main','Missing implementation',1)
    286286              ! Read pqr(:,iq)
    287287            endif
     
    289289 
    290290          ! Recombine for all T-P-Q as we do it only once for all.
    291           call recombin_corrk_ini(pqr)
     291          call rad_correlatedk_online_recombination_init(pqr)
    292292        ELSE
    293293          all_otf=.true.
     
    304304      ! calculate a gasi/v_recomb variable on the reference corrk-k T,P,X grid (only for T,P,X values
    305305      ! who match the atmospheric conditions) which is then processed as a standard pre-mix in
    306       ! optci/v routines, but updated every time tracers on the ref P grid have varied > 1%.
     306      ! rad_correlatedk_opacities_thermal/v routines, but updated every time tracers on the ref P grid have varied > 1%.
    307307 
    308308      ! Extract tracers for species which are recombined on-the-fly
     
    341341      ! The following useptx is a trick to call recombine only for the reference T-P-X
    342342      ! reference grid points that are useful given the temperature range (and variable specie amount) at a given altitude
    343       ! (cf optci/optcv routines where we interpolate corrk calling tpindex)
     343      ! (cf rad_correlatedk_opacities_thermal/rad_correlatedk_opacities_stellar routines where we interpolate corrk calling tpindex)
    344344      ! It saves a looot of time - JVO 18
    345345 
     
    350350
    351351      if (.not.all_otf) then
    352          call recombin_corrk_mix(igrid,pqr,useptx)
     352         call rad_correlatedk_online_recombination_mix(igrid,pqr,useptx)
    353353       else
    354354         if (nrecomb_qotf.gt.1) then
    355            call recombin_corrk_mix_allotf(igrid,pqr,useptx)
     355           call rad_correlatedk_online_recombination_mix_allotf(igrid,pqr,useptx)
    356356         else
    357357           do ix=1,L_REFVAR
     
    368368       endif
    369369 
    370     END SUBROUTINE call_recombin
    371 
    372     SUBROUTINE recombin_corrk_ini(pqr)
     370    END SUBROUTINE rad_correlatedk_recombination_main
     371
     372    SUBROUTINE rad_correlatedk_online_recombination_init(pqr)
    373373
    374374      USE radinc_h
     
    587587      enddo ! ip=1,L_REFVAR
    588588 
    589     END SUBROUTINE recombin_corrk_ini
    590 
    591     SUBROUTINE recombin_corrk_mix(igrid,pqr,useptx)
     589    END SUBROUTINE rad_correlatedk_online_recombination_init
     590
     591    SUBROUTINE rad_correlatedk_online_recombination_mix(igrid,pqr,useptx)
    592592
    593593      USE radinc_h
     
    808808      enddo ! ip=1,L_REFVAR
    809809 
    810     END SUBROUTINE recombin_corrk_mix
    811 
    812     SUBROUTINE recombin_corrk_mix_allotf(igrid,pqr,useptx)
     810    END SUBROUTINE rad_correlatedk_online_recombination_mix
     811
     812    SUBROUTINE rad_correlatedk_online_recombination_mix_allotf(igrid,pqr,useptx)
    813813
    814814      USE radinc_h
     
    10271027      enddo ! ip=1,L_REFVAR
    10281028 
    1029     END SUBROUTINE recombin_corrk_mix_allotf
    1030 
    1031   END MODULE recombin_corrk_mod
    1032  
     1029    END SUBROUTINE rad_correlatedk_online_recombination_mix_allotf
     1030
     1031  END MODULE rad_correlatedk_online_recombination_mod
     1032 
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_opacities_stellar.F90

    r4076 r4077  
    1 MODULE optcv_mod
     1MODULE rad_correlatedk_opacities_stellar_mod
    22
    33IMPLICIT NONE
     
    55CONTAINS
    66
    7 SUBROUTINE OPTCV(DTAUV,TAUV,TAUCUMV,PLEV,  &
     7SUBROUTINE rad_correlatedk_opacities_stellar(DTAUV,TAUV,TAUCUMV,PLEV,  &
    88     QXVAER,QSVAER,GVAER,WBARV,COSBV,       &
    99     TAUAERO,TMID,PMID,TAUGSURF,QVAR,MUVAR,FRACVAR)
     
    1616  use callkeys_mod, only: kastprof,continuum,graybody,callgasvis,varspec, &
    1717                          rayleigh
    18   use recombin_corrk_mod, only: corrk_recombin, gasv_recomb
     18  use rad_correlatedk_online_recombination_mod, only: corrk_recombin, gasv_recomb
    1919  use tpindex_mod, only: tpindex
    20   use interpolate_continuum_mod, only: interpolate_continuum
    21   use calc_rayleigh_mod, only: calc_rayleigh
     20  use rad_correlatedk_continuum_interpolation_mod, only: rad_correlatedk_continuum_interpolation
     21  use rad_correlatedk_rayleigh_scattering_opacity_mod, only: rad_correlatedk_rayleigh_scattering_opacity
    2222
    2323  implicit none
     
    138138     else
    139139        dz(k) = dpr(k)*R*TMID(K)/(glat_ig*PMID(K))*mugaz/muvar(k)
    140         U(k)  = Cmk*DPR(k)*mugaz/muvar(k)     ! only Cmk line in optci.F 
     140        U(k)  = Cmk*DPR(k)*mugaz/muvar(k)     ! only Cmk line in rad_correlatedk_opacities_thermal.F 
    141141            !JL13 the mugaz/muvar factor takes into account water meanmolecular weight if water is present
    142142     endif
     
    154154            !   but visible does not handle very well diffusion in first layer.
    155155            !   The tauaero and tauray are thus set to 0 (a small value for rayleigh because the code crashes otherwise)
    156             !   in the 4 first semilayers in optcv, but not optci.
     156            !   in the 4 first semilayers in rad_correlatedk_opacities_stellar, but not rad_correlatedk_opacities_thermal.
    157157            !   This solves random variations of the sw heating at the model top.
    158158  do iaer=1,naerkind
     
    170170
    171171      if(rayleigh) then
    172          call calc_rayleigh(QVAR,MUVAR,PMID,TMID,TAURAY)
     172         call rad_correlatedk_rayleigh_scattering_opacity(QVAR,MUVAR,PMID,TMID,TAURAY)
    173173      else
    174          print*,'setspv: No Rayleigh scattering, check for NaN in output!'
     174         print*,'rad_correlatedk_init_stellar: No Rayleigh scattering, check for NaN in output!'
    175175         do NW=1,L_NSPECTV
    176176            TAURAY(:,NW) = 1E-16
     
    243243                     ((igas .eq. igas_CO2) .and. (jgas .eq. igas_CH4)) ) then
    244244
    245                   call interpolate_continuum('',igas,jgas,'VI',nw,T_cont,p_cont,p_cross,dtemp,.false.)
     245                  call rad_correlatedk_continuum_interpolation('',igas,jgas,'VI',nw,T_cont,p_cont,p_cross,dtemp,.false.)
    246246
    247247                endif
     
    268268           
    269269              ! JVO 2017 : added tmpk because the repeated calls to gasi/v increased dramatically
    270               ! the execution time of optci/v -> ~ factor 2 on the whole radiative
     270              ! the execution time of rad_correlatedk_opacities_thermal/v -> ~ factor 2 on the whole radiative
    271271              ! transfer on the tested simulations !
    272272
     
    334334            !JL18 It seems to be good to have aerosols in the first "radiative layer" of the gcm in the IR
    335335            !   but not in the visible
    336             !   The tauaero is thus set to 0 in the 4 first semilayers in optcv, but not optci.
     336            !   The tauaero is thus set to 0 in the 4 first semilayers in rad_correlatedk_opacities_stellar, but not rad_correlatedk_opacities_thermal.
    337337            !   This solves random variations of the sw heating at the model top.
    338338  do iaer=1,naerkind
     
    407407
    408408
    409 end subroutine optcv
    410 
    411 END MODULE optcv_mod
     409end subroutine rad_correlatedk_opacities_stellar
     410
     411END MODULE rad_correlatedk_opacities_stellar_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_opacities_thermal.F90

    r4076 r4077  
    1 MODULE optci_mod
     1MODULE rad_correlatedk_opacities_thermal_mod
    22
    33IMPLICIT NONE
     
    55CONTAINS
    66
    7 subroutine optci(PLEV,TLEV,DTAUI,TAUCUMI,      &
     7subroutine rad_correlatedk_opacities_thermal(PLEV,TLEV,DTAUI,TAUCUMI,      &
    88     QXIAER,QSIAER,GIAER,COSBI,WBARI,TAUAERO,  &
    99     TMID,PMID,TAUGSURF,QVAR,MUVAR,FRACVAR)
     
    1616  use comcstfi_mod, only: g, r, mugaz
    1717  use callkeys_mod, only: kastprof,continuum,graybody,varspec
    18   use recombin_corrk_mod, only: corrk_recombin, gasi_recomb
     18  use rad_correlatedk_online_recombination_mod, only: corrk_recombin, gasi_recomb
    1919  use tpindex_mod, only: tpindex
    20   use interpolate_continuum_mod, only: interpolate_continuum
     20  use rad_correlatedk_continuum_interpolation_mod, only: rad_correlatedk_continuum_interpolation
    2121
    2222  implicit none
     
    150150     else
    151151        dz(k) = dpr(k)*R*TMID(K)/(glat_ig*PMID(K))*mugaz/muvar(k)
    152         U(k)  = Cmk*DPR(k)*mugaz/muvar(k)     ! only Cmk line in optci.F 
     152        U(k)  = Cmk*DPR(k)*mugaz/muvar(k)     ! only Cmk line in rad_correlatedk_opacities_thermal.F 
    153153            !JL13 the mugaz/muvar factor takes into account water meanmolecular weight if water is present
    154154     endif
     
    222222                     ((igas .eq. igas_CO2) .and. (jgas .eq. igas_CH4))  ) then
    223223
    224                    call interpolate_continuum('',igas,jgas,'IR',nw,T_cont,p_cont,p_cross,dtemp,.false.)
     224                   call rad_correlatedk_continuum_interpolation('',igas,jgas,'IR',nw,T_cont,p_cont,p_cross,dtemp,.false.)
    225225
    226226                 endif
     
    247247           
    248248              ! JVO 2017 : added tmpk because the repeated calls to gasi/v increased dramatically
    249               ! the execution time of optci/v -> ~ factor 2 on the whole radiative
     249              ! the execution time of rad_correlatedk_opacities_thermal/v -> ~ factor 2 on the whole radiative
    250250              ! transfer on the tested simulations !
    251251
     
    457457!  call abort
    458458
    459 end subroutine optci
    460 
    461 END MODULE optci_mod
    462 
     459end subroutine rad_correlatedk_opacities_thermal
     460
     461END MODULE rad_correlatedk_opacities_thermal_mod
     462
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_rayleigh_scattering_opacity.F90

    r4076 r4077  
    1 module calc_rayleigh_mod
     1module rad_correlatedk_rayleigh_scattering_opacity_mod
    22
    33implicit none
     
    55contains
    66
    7       subroutine calc_rayleigh(qvar,muvar,PMID,TMID,tauray)
     7      subroutine rad_correlatedk_rayleigh_scattering_opacity(qvar,muvar,PMID,TMID,tauray)
    88
    99!==================================================================
     
    2424!     Called by
    2525!     ---------
    26 !     setspv.F
     26!     rad_correlatedk_init_stellar.F
    2727!     
    2828!     Calls
     
    9797      ! mu is the mean molecular weight
    9898      ! x_i is the mass fraction of the ith gas
    99       ! The pressure P dependence is calculated in optcv.F90
     99      ! The pressure P dependence is calculated in rad_correlatedk_opacities_stellar.F90
    100100     
    101101      if(firstcall) then
     
    325325            enddo !ngasmx
    326326
    327             call blackl(dble(wl*1e-6),dble(tstellar),df)
     327            call rad_blackbody_planck_law_wavelength(dble(wl*1e-6),dble(tstellar),df)
    328328            df=df*bwidth/Nfine
    329329            tauwei=tauwei+df
     
    337337
    338338
    339    end subroutine calc_rayleigh
    340 
    341 end module calc_rayleigh_mod
     339   end subroutine rad_correlatedk_rayleigh_scattering_opacity
     340
     341end module rad_correlatedk_rayleigh_scattering_opacity_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_read_opacity_tables.F90

    r4076 r4077  
    1       subroutine sugas_corrk
     1      subroutine rad_correlatedk_read_opacity_tables
    22
    33!==================================================================
     
    3636                continuum
    3737      use tracer_h, only : nqtot, moderntracdef, is_recomb, noms
    38       use recombin_corrk_mod, only: su_recombin,        &
     38      use rad_correlatedk_online_recombination_mod, only: rad_correlatedk_recombination_setup,        &
    3939                corrk_recombin, use_premix, nrecomb_tot, rcb2tot_idx
    40       use interpolate_continuum_mod, only: interpolate_continuum
     40      use rad_correlatedk_continuum_interpolation_mod, only: rad_correlatedk_continuum_interpolation
    4141      implicit none
    4242
     
    7777      if(.not.file_ok) then
    7878         write(*,*)'The file ',TRIM(file_path)
    79          write(*,*)'was not found by sugas_corrk.F90, exiting.'
     79         write(*,*)'was not found by rad_correlatedk_read_opacity_tables .F90, exiting.'
    8080         write(*,*)'Check that your path to datagcm:',trim(datadir)
    8181         write(*,*)' is correct. You can change it in callphys.def with:'
    8282         write(*,*)' datadir = /absolute/path/to/datagcm'
    8383         write(*,*)'Also check that the corrkdir you chose in callphys.def exists.'
    84          call abort_physic("sugas_corrk", "Unable to read file", 1)
     84         call abort_physic("rad_correlatedk_read_opacity_tables ", "Unable to read file", 1)
    8585      endif
    8686
     
    9696           print*,'Number of gases in radiative transfer data (',ngas,') does not', &
    9797                  'match that in gases.def (',ngasmx,'), exiting.'
    98            call abort_physic("sugas_corrk", "Number of gases in radiative transfer data does not match that in gases.def", 1)
     98           call abort_physic("rad_correlatedk_read_opacity_tables ", "Number of gases in radiative transfer data does not match that in gases.def", 1)
    9999        endif
    100100      endif
     
    104104                     corrkdir(1:LEN_TRIM(corrkdir)),           &
    105105                '] has no variable species, exiting.'
    106          call abort_physic("sugas_corrk", "You have varactive/fixed=.true. but the database has no variable species",1)
     106         call abort_physic("rad_correlatedk_read_opacity_tables ", "You have varactive/fixed=.true. but the database has no variable species",1)
    107107      elseif(ngas.gt.5 .or. ngas.lt.1)then
    108108         print*,ngas,' species in database [',               &
    109109                     corrkdir(1:LEN_TRIM(corrkdir)),           &
    110110                '], radiative code cannot handle this.'
    111          call abort_physic("sugas_corrk", "No gas or too many gases for radiative code", 1)
     111         call abort_physic("rad_correlatedk_read_opacity_tables ", "No gas or too many gases for radiative code", 1)
    112112      endif
    113113
     
    135135                     corrkdir(1:LEN_TRIM(corrkdir)),           &
    136136                '] has a variable species.'
    137          call abort_physic("sugas_corrk", "You have varactive and varfixed=.false. and the database has a variable species",1)
     137         call abort_physic("rad_correlatedk_read_opacity_tables ", "You have varactive and varfixed=.false. and the database has a variable species",1)
    138138      endif
    139139
     
    149149                   'match that in gases.def (',trim(gnom(igas)),'), exiting. You should compare ', &
    150150                   'gases.def with Q.dat in your radiative transfer directory.'
    151               call abort_physic("sugas_corrk", "Name of a gas in radiative transfer data does not match that in gases.def",1)
     151              call abort_physic("rad_correlatedk_read_opacity_tables ", "Name of a gas in radiative transfer data does not match that in gases.def",1)
    152152           endif
    153153        enddo
     
    177177      if(.not.file_ok) then
    178178         write(*,*)'The file ',TRIM(file_path)
    179          write(*,*)'was not found by sugas_corrk.F90, exiting.'
     179         write(*,*)'was not found by rad_correlatedk_read_opacity_tables .F90, exiting.'
    180180         write(*,*)'Check that your path to datagcm:',trim(datadir)
    181181         write(*,*)' is correct. You can change it in callphys.def with:'
    182182         write(*,*)' datadir = /absolute/path/to/datagcm'
    183183         write(*,*)'Also check that the corrkdir you chose in callphys.def exists.'
    184          call abort_physic("sugas_corrk", "Unable to read file", 1)
     184         call abort_physic("rad_correlatedk_read_opacity_tables ", "Unable to read file", 1)
    185185      endif
    186186     
     
    213213      if(.not.file_ok) then
    214214         write(*,*)'The file ',TRIM(file_path)
    215          write(*,*)'was not found by sugas_corrk.F90, exiting.'
     215         write(*,*)'was not found by rad_correlatedk_read_opacity_tables .F90, exiting.'
    216216         write(*,*)'Check that your path to datagcm:',trim(datadir)
    217217         write(*,*)' is correct. You can change it in callphys.def with:'
    218218         write(*,*)' datadir = /absolute/path/to/datagcm'
    219219         write(*,*)'Also check that the corrkdir you chose in callphys.def exists.'
    220          call abort_physic("sugas_corrk", "Unable to read file", 1)
     220         call abort_physic("rad_correlatedk_read_opacity_tables ", "Unable to read file", 1)
    221221      endif
    222222     
     
    259259      if(.not.file_ok) then
    260260         write(*,*)'The file ',TRIM(file_path)
    261          write(*,*)'was not found by sugas_corrk.F90, exiting.'
     261         write(*,*)'was not found by rad_correlatedk_read_opacity_tables .F90, exiting.'
    262262         write(*,*)'Check that your path to datagcm:',trim(datadir)
    263263         write(*,*)' is correct. You can change it in callphys.def with:'
    264264         write(*,*)' datadir = /absolute/path/to/datagcm'
    265265         write(*,*)'Also check that the corrkdir you chose in callphys.def exists.'
    266          call abort_physic("sugas_corrk", "Unable to read file",1)
     266         call abort_physic("rad_correlatedk_read_opacity_tables ", "Unable to read file",1)
    267267      endif
    268268
     
    379379              if(.not.file_ok) then
    380380                 write(*,*)'The file ',TRIM(file_path)
    381                  write(*,*)'was not found by sugas_corrk.F90.'
     381                 write(*,*)'was not found by rad_correlatedk_read_opacity_tables .F90.'
    382382                 write(*,*)'Are you sure you have absorption data for these bands?'
    383                  call abort_physic("sugas_corrk", "No absorption data found", 1)
     383                 call abort_physic("rad_correlatedk_read_opacity_tables ", "No absorption data found", 1)
    384384              endif
    385385           
     
    403403            if(.not.file_ok) then
    404404               write(*,*)'The file ',TRIM(file_path)
    405                write(*,*)'was not found by sugas_corrk.F90.'
     405               write(*,*)'was not found by rad_correlatedk_read_opacity_tables .F90.'
    406406               write(*,*)'Are you sure you have absorption data for this specie at these bands?'
    407                call abort_physic("sugas_corrk", "No absorption data found", 1)
     407               call abort_physic("rad_correlatedk_read_opacity_tables ", "No absorption data found", 1)
    408408            endif
    409409         
     
    449449           if(.not.file_ok) then
    450450              write(*,*)'The file ',TRIM(file_path)
    451               write(*,*)'was not found by sugas_corrk.F90.'
     451              write(*,*)'was not found by rad_correlatedk_read_opacity_tables .F90.'
    452452              write(*,*)'Are you sure you have absorption data for these bands?'
    453               call abort_physic("sugas_corrk", "No absorption data found", 1)
     453              call abort_physic("rad_correlatedk_read_opacity_tables ", "No absorption data found", 1)
    454454           endif
    455455         
     
    516516          if(.not.file_ok) then
    517517             write(*,*)'The file ',TRIM(file_path)
    518              write(*,*)'was not found by sugas_corrk.F90.'
     518             write(*,*)'was not found by rad_correlatedk_read_opacity_tables .F90.'
    519519             write(*,*)'Are you sure you have absorption data for this specie at these bands?'
    520              call abort_physic("sugas_corrk", "No absorption data found", 1)
     520             call abort_physic("rad_correlatedk_read_opacity_tables ", "No absorption data found", 1)
    521521          endif
    522522!$OMP MASTER           
     
    710710
    711711! Allocate and initialise arrays for corrk_recombin
    712       if (corrk_recombin) call su_recombin
     712      if (corrk_recombin) call rad_correlatedk_recombination_setup
    713713
    714714!=======================================================================
     
    721721         file_id='/continuum_data/' // 'N2-N2_continuum_70-500K_2025.dat'
    722722         file_path=TRIM(datadir)//TRIM(file_id)
    723          call interpolate_continuum(file_path,igas_N2,igas_N2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     723         call rad_correlatedk_continuum_interpolation(file_path,igas_N2,igas_N2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    724724         do jgas=1,ngasmx
    725725          if (jgas .eq. igas_H2) then
    726726           file_id='/continuum_data/' // 'H2-N2_continuum_40-600K_2025.dat'
    727727           file_path=TRIM(datadir)//TRIM(file_id)
    728            call interpolate_continuum(file_path,igas_N2,igas_H2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     728           call rad_correlatedk_continuum_interpolation(file_path,igas_N2,igas_H2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    729729          elseif (jgas .eq. igas_O2) then
    730730           file_id='/continuum_data/' // 'O2-N2_continuum_100-500K_2025.dat'
    731731           file_path=TRIM(datadir)//TRIM(file_id)
    732            call interpolate_continuum(file_path,igas_N2,igas_O2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     732           call rad_correlatedk_continuum_interpolation(file_path,igas_N2,igas_O2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    733733          elseif (jgas .eq. igas_CH4) then
    734734           file_id='/continuum_data/' // 'CH4-N2_continuum_40-600K_2025.dat'
    735735           file_path=TRIM(datadir)//TRIM(file_id)
    736            call interpolate_continuum(file_path,igas_N2,igas_CH4,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     736           call rad_correlatedk_continuum_interpolation(file_path,igas_N2,igas_CH4,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    737737          endif
    738738         enddo
     
    740740         file_id='/continuum_data/' // 'O2-O2_continuum_100-400K_2025.dat'
    741741         file_path=TRIM(datadir)//TRIM(file_id)
    742          call interpolate_continuum(file_path,igas_O2,igas_O2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     742         call rad_correlatedk_continuum_interpolation(file_path,igas_O2,igas_O2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    743743         do jgas=1,ngasmx
    744744          if (jgas .eq. igas_CO2) then
    745745           file_id='/continuum_data/' // 'CO2-O2_continuum_150-600K_2025.dat'
    746746           file_path=TRIM(datadir)//TRIM(file_id)
    747            call interpolate_continuum(file_path,igas_CO2,igas_O2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     747           call rad_correlatedk_continuum_interpolation(file_path,igas_CO2,igas_O2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    748748          endif
    749749         enddo
     
    751751         file_id='/continuum_data/' // 'H2-H2_continuum_40-7000K_2025.dat'
    752752         file_path=TRIM(datadir)//TRIM(file_id)
    753          call interpolate_continuum(file_path,igas_H2,igas_H2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     753         call rad_correlatedk_continuum_interpolation(file_path,igas_H2,igas_H2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    754754         do jgas=1,ngasmx
    755755          if (jgas .eq. igas_CH4) then
    756756           file_id='/continuum_data/' // 'H2-CH4_continuum_40-600K_2025.dat'
    757757           file_path=TRIM(datadir)//TRIM(file_id)
    758            call interpolate_continuum(file_path,igas_H2,igas_CH4,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     758           call rad_correlatedk_continuum_interpolation(file_path,igas_H2,igas_CH4,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    759759          elseif (jgas .eq. igas_He) then
    760760           file_id='/continuum_data/' // 'H2-He_continuum_40-5500K_2025.dat'
    761761           file_path=TRIM(datadir)//TRIM(file_id)
    762            call interpolate_continuum(file_path,igas_H2,igas_He,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     762           call rad_correlatedk_continuum_interpolation(file_path,igas_H2,igas_He,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    763763          endif
    764764         enddo   
     
    766766         file_id='/continuum_data/' // 'CH4-CH4_continuum_40-500K_2025.dat'
    767767         file_path=TRIM(datadir)//TRIM(file_id)
    768          call interpolate_continuum(file_path,igas_CH4,igas_CH4,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     768         call rad_correlatedk_continuum_interpolation(file_path,igas_CH4,igas_CH4,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    769769        elseif (igas .eq. igas_H2O) then
    770770         file_id='/continuum_data/' // 'H2O-H2O_continuum_100-2000K_2025.dat'
    771771         file_path=TRIM(datadir)//TRIM(file_id)
    772          call interpolate_continuum(file_path,igas_H2O,igas_H2O,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     772         call rad_correlatedk_continuum_interpolation(file_path,igas_H2O,igas_H2O,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    773773         do jgas=1,ngasmx
    774774          if (jgas .eq. igas_N2) then
    775775           file_id='/continuum_data/' // 'H2O-N2_continuum_100-2000K_2025.dat'
    776776           file_path=TRIM(datadir)//TRIM(file_id)
    777            call interpolate_continuum(file_path,igas_H2O,igas_N2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     777           call rad_correlatedk_continuum_interpolation(file_path,igas_H2O,igas_N2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    778778          elseif (jgas .eq. igas_O2) then
    779779           file_id='/continuum_data/' // 'H2O-O2_continuum_100-2000K_2025.dat'
    780780           file_path=TRIM(datadir)//TRIM(file_id)
    781            call interpolate_continuum(file_path,igas_H2O,igas_O2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     781           call rad_correlatedk_continuum_interpolation(file_path,igas_H2O,igas_O2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    782782          elseif (jgas .eq. igas_CO2) then
    783783           file_id='/continuum_data/' // 'H2O-CO2_continuum_100-800K_2025.dat'
    784784           file_path=TRIM(datadir)//TRIM(file_id)
    785            call interpolate_continuum(file_path,igas_H2O,igas_CO2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     785           call rad_correlatedk_continuum_interpolation(file_path,igas_H2O,igas_CO2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    786786          endif
    787787         enddo   
     
    789789         file_id='/continuum_data/' // 'CO2-CO2_continuum_100-800K_2025.dat'
    790790         file_path=TRIM(datadir)//TRIM(file_id)
    791          call interpolate_continuum(file_path,igas_CO2,igas_CO2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     791         call rad_correlatedk_continuum_interpolation(file_path,igas_CO2,igas_CO2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    792792         do jgas=1,ngasmx
    793793          if (jgas .eq. igas_H2) then
    794794           file_id='/continuum_data/' // 'CO2-H2_continuum_100-800K_2025.dat'
    795795           file_path=TRIM(datadir)//TRIM(file_id)
    796            call interpolate_continuum(file_path,igas_CO2,igas_H2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     796           call rad_correlatedk_continuum_interpolation(file_path,igas_CO2,igas_H2,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    797797          elseif (jgas .eq. igas_CH4) then
    798798           file_id='/continuum_data/' // 'CO2-CH4_continuum_100-800K_2025.dat'
    799799           file_path=TRIM(datadir)//TRIM(file_id)
    800            call interpolate_continuum(file_path,igas_CO2,igas_CH4,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
     800           call rad_correlatedk_continuum_interpolation(file_path,igas_CO2,igas_CH4,'IR',1,300.D+0,10000.D+0,20000.D+0,testcont,.true.)
    801801          endif
    802802         enddo
     
    822822!$OMP BARRIER
    823823
    824     end subroutine sugas_corrk
     824    end subroutine rad_correlatedk_read_opacity_tables
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_correlatedk_stellar_spectrum.F90

    r4076 r4077  
    1       module ave_stelspec_mod
     1      module rad_correlatedk_stellar_spectrum_mod
    22     
    33      implicit none
     
    55      contains
    66     
    7       subroutine ave_stelspec(STELLAR)
     7      subroutine rad_correlatedk_stellar_spectrum (STELLAR)
    88
    99!==================================================================
     
    2222!     Called by
    2323!     ---------
    24 !     setspv.F
     24!     rad_correlatedk_init_stellar.F
    2525!     
    2626!     Calls
     
    7171            do ifine=1,Nfineband
    7272               lam_temp=lamm+(lamp-lamm)*(ifine-1.)/(Nfineband)
    73                call blackl(dble(lam_temp*1e-6),dble(tstellar),stel_temp)
     73               call rad_blackbody_planck_law_wavelength(dble(lam_temp*1e-6),dble(tstellar),stel_temp)
    7474               STELLAR(band)=STELLAR(band)+stel_temp*dl
    7575            enddo           
     
    8686           write(*,*)'Error: tstellar (effective stellar temperature) needs to be specified'
    8787           write(*,*)'in callphys.def: tstellar=...'
    88            call abort_physic("ave_stelspec", "tstellar needs to be specified",1)
     88           call abort_physic("rad_correlatedk_stellar_spectrum", "tstellar needs to be specified",1)
    8989         end if
    9090         
     
    121121           write(*,*)'available stellar spectra here : '
    122122           write(*,*)'https://web.lmd.jussieu.fr/~lmdz/planets/generic/datagcm/stellar_spectra/'
    123            call abort_physic("ave_stelspec", "Unable to read stellar flux file", 1)
     123           call abort_physic("rad_correlatedk_stellar_spectrum", "Unable to read stellar flux file", 1)
    124124         end if
    125125
     
    174174      endif !stelbbody
    175175
    176       end subroutine ave_stelspec
     176      end subroutine rad_correlatedk_stellar_spectrum
    177177     
    178       end module ave_stelspec_mod
     178      end module rad_correlatedk_stellar_spectrum_mod
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_newton_cooling.F90

    r4076 r4077  
    1 subroutine newtrelax(ngrid,nlayer,mu0,sinlat,popsk,temp,pplay,pplev,dtrad,firstcall)
     1subroutine rad_netwon_cooling(ngrid,nlayer,mu0,sinlat,popsk,temp,pplay,pplev,dtrad,firstcall)
    22       
    33  use comcstfi_mod, only: rcp, pi
     
    126126  !call writediagfi(ngrid,'ThetaZ','stellar zenith angle','deg',2,mu0)
    127127
    128 end subroutine newtrelax
     128end subroutine rad_netwon_cooling
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_newton_cooling_hot_jupiter.F90

    r4076 r4077  
    1 module newton_cooling_hotJ
     1module rad_netwon_cooling_hot_jupiter
    22   
    33    !==========================================================================================
     
    1010    !
    1111    ! We aim at having a generic code but you never know, it might need improving at some point.
    12     ! The current (at time of writing) newtrelax.F90 routine is hardcoded for telluric temperate planets and untested.
     12    ! The current (at time of writing) rad_netwon_cooling.F90 routine is hardcoded for telluric temperate planets and untested.
    1313    ! Thus, we don't use it and use this one instead.
    1414    !
     
    3030    contains
    3131
    32     subroutine newtcool_MOCHA(ngrid,nlayer,coslon,coslat,temp,pplay,firstcall,lastcall,dtrad)
     32    subroutine rad_newton_cooling_MOCHA_intercomparison(ngrid,nlayer,coslon,coslat,temp,pplay,firstcall,lastcall,dtrad)
    3333
    3434        ! use callkeys_mod, only: planetary_suffix ! this is to know which profiles to load for the T0, the delta Teq and the tau_rad
     
    131131        endif
    132132
    133     end subroutine newtcool_MOCHA
     133    end subroutine rad_newton_cooling_MOCHA_intercomparison
    134134
    135135    subroutine read_input(nlayer,filename, field)
     
    163163        open(401,form='formatted',status='old',file=trim(filename) ,iostat=ierr)
    164164            if (ierr /=0) then
    165                 print*,"Problem in newton_cooling_hotJ.F90"
     165                print*,"Problem in rad_netwon_cooling_hot_jupiter.F90"
    166166                print*,"I have an issue opening file ",trim(filename)
    167167                call abort_physic("newton_cooling_hot_J", "Unable to read input file", 1)
     
    180180    end subroutine read_input
    181181
    182 end module newton_cooling_hotJ
     182end module rad_netwon_cooling_hot_jupiter
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_ring_shadowing.F90

    r4076 r4077  
    1 subroutine call_rings(ngrid, ptime, pday, diurnal)
     1subroutine rad_ring_shadowing(ngrid, ptime, pday, diurnal)
    22    ! A subroutine to compute the day fraction in case of rings shadowing.
    33
     
    4141        do m=1, nb_hours
    4242            ptime_day = m*pas
    43             call stellarlong(pday+ptime_day,tmp_zls)
    44             call orbite(tmp_zls,tmp_dist_star,tmp_declin,tmp_right_ascen)
     43            call ephemeris_stellar_longitude(pday+ptime_day,tmp_zls)
     44            call ephemeris_orbit(tmp_zls,tmp_dist_star,tmp_declin,tmp_right_ascen)
    4545           
    4646            ztim1=SIN(tmp_declin)
     
    4848            ztim3=-COS(tmp_declin)*SIN(2.*pi*(pday+ptime_day-.5))
    4949
    50             call stelang(ngrid,sinlon,coslon,sinlat,coslat,    &
     50            call ephemeris_stellar_angle(ngrid,sinlon,coslon,sinlat,coslat,    &
    5151                        ztim1,ztim2,ztim3,tmp_mu0,tmp_fract, flatten)       
    52             call rings(ngrid, tmp_declin, ptime_day, rad, flatten, eclipse)
     52            call saturn_rings(ngrid, tmp_declin, ptime_day, rad, flatten, eclipse)
    5353            fract(:) = fract(:) + (1.-eclipse(:))*tmp_fract(:) !! fract takes into account the rings shadow and the day/night alternation
    5454
     
    6161                 
    6262     else   ! instant insolation is weighted by the rings shadow
    63             call rings(ngrid, declin, ptime, rad, 0., eclipse)
     63            call saturn_rings(ngrid, declin, ptime, rad, 0., eclipse)
    6464            fract(:) = fract(:) * (1.-eclipse)
    6565    endif
     
    6767    IF (ALLOCATED(eclipse)) DEALLOCATE(eclipse)
    6868
    69 end subroutine call_rings
     69end subroutine rad_ring_shadowing
     70
     71
     72SUBROUTINE saturn_rings(ngrid, declin, ptime, rad, flat, eclipse)
     73! Calculates Saturn's rings shadowing
     74! Includes rings opacities measured by Cassini/UVIS
     75! Authors: M. Sylvestre, M. Capderou, S. Guerlet, A. Spiga
     76
     77    use comdiurn_h, only: sinlat, sinlon, coslat, coslon
     78    use geometry_mod, only: latitude ! (rad)
     79 
     80    implicit none   
     81
     82    INTEGER, INTENT(IN) :: ngrid  ! horizontal grid dimension
     83    REAL, INTENT(IN) :: declin    ! latitude of the subsolar point
     84    REAL, INTENT(IN) :: ptime     ! UTC time in sol fraction : ptime=0.5 at noon
     85    REAL, INTENT(IN) :: rad       ! equatorial radius of the planet
     86    REAL, INTENT(IN) :: flat      ! flattening of the planet
     87    REAL, DIMENSION(ngrid), INTENT(OUT) :: eclipse ! absorption of the light by the rings   
     88   
     89    REAL :: rpol   ! polar radius of the planet
     90    REAL :: e      ! shape excentricity of the planet : (1-e*e) = (1-f)*(1-f)   
     91    INTEGER, PARAMETER :: nb_a = 4 ! number of subdivisions of the A ring
     92    INTEGER, PARAMETER :: nb_b = 3 ! number of subdivisions of the B ring
     93    INTEGER, PARAMETER :: nb_c = 3 ! number of subdivisions of the C ring
     94    INTEGER, PARAMETER :: nb_ca = 2 ! number of subdivisions in the Cassini division
     95    INTEGER :: i
     96
     97    ! arrays for the rings. TBD: dynamical?
     98    REAL, DIMENSION(nb_a) :: A_Rint ! internal radii of the subdivisions of the A ring
     99    REAL, DIMENSION(nb_a) :: A_Rext ! external radii of the subdivisions of the A ring
     100    REAL, DIMENSION(nb_b) :: B_Rint ! internal radii of the subdivisions of the B ring
     101    REAL, DIMENSION(nb_b) :: B_Rext ! external radii of the subdivisions of the B ring
     102    REAL, DIMENSION(nb_c) :: C_Rint ! internal radii of the subdivisions of the C ring
     103    REAL, DIMENSION(nb_c) :: C_Rext ! external radii of the subdivisions of the C ring
     104    REAL, DIMENSION(nb_ca) :: Ca_Rint ! internal radii of the subdivisions of the Cassini Division
     105    REAL, DIMENSION(nb_ca) :: Ca_Rext ! external radii of the subdivisions of the Cassini Division
     106
     107    ! Opacities of the rings : for each one we can give different opacities for each part
     108    REAL, DIMENSION(nb_a) :: tau_A ! opacity of the A ring
     109    REAL, DIMENSION(nb_b) :: tau_B ! opacity of the B ring
     110    REAL, DIMENSION(nb_c) :: tau_C ! opacity of the C ring
     111    REAL, DIMENSION(nb_ca) :: tau_Ca ! opacity of the Cassini Division
     112
     113    ! Parameters used to calculate if a point is under a ring subdivision's shadow
     114    REAL :: phi_S                             ! subsolar point longitude
     115    REAL, PARAMETER :: pi=acos(-1.0)   
     116    REAL, DIMENSION(:), ALLOCATABLE:: x, y, z ! cartesian coordinates of the points on the planet
     117    REAL :: xs, ys, zs                        ! cartesian coordinates of the points of the subsolar point
     118    REAL, DIMENSION(:), ALLOCATABLE :: k
     119    REAL, DIMENSION(:), ALLOCATABLE :: N      ! parameter to compute cartesian coordinates on a ellipsoidal planet
     120    REAL, DIMENSION(:), ALLOCATABLE :: r      ! distance at which the incident ray of sun crosses the equatorial plane
     121                                              ! measured from the center of the planet   
     122    REAL :: Ns                                ! (same for the subsolar point)
     123   
     124    ! equinox --> no shadow (AS: why is this needed?)
     125    if(declin .eq. 0.) then
     126        eclipse(:) = 0.
     127        return
     128    endif
     129
     130! 1) INITIALIZATION
     131
     132    ! Generic
     133    rpol = (1.- flat)*rad
     134    e = sqrt(2*flat - flat**2)
     135    ALLOCATE(x(ngrid))
     136    ALLOCATE(y(ngrid))
     137    ALLOCATE(z(ngrid))
     138    ALLOCATE(k(ngrid))
     139    ALLOCATE(N(ngrid))
     140    ALLOCATE(r(ngrid))
     141    eclipse(:) = 2000.
     142
     143! Model of the rings with Cassini/UVIS opacities
     144
     145    ! Size of the rings
     146    A_Rint(1) = 2.03*rad
     147    A_Rext(1) = 2.06*rad
     148    A_Rint(2) = 2.06*rad
     149    A_Rext(2) = 2.09*rad
     150    A_Rint(3) = 2.09*rad
     151    A_Rext(3) = 2.12*rad
     152    A_Rint(4) = 2.12*rad
     153    A_Rext(4) = 2.27*rad
     154
     155    B_Rint(1) = 1.53*rad
     156    B_Rext(1) = 1.64*rad
     157    B_Rint(2) = 1.64*rad
     158    B_Rext(2) = 1.83*rad
     159    B_Rint(3) = 1.83*rad
     160    B_Rext(3) = 1.95*rad
     161   
     162    C_Rint(1) = 1.24*rad
     163    C_Rext(1) = 1.29*rad
     164    C_Rint(2) = 1.29*rad
     165    C_Rext(2) = 1.43*rad
     166    C_Rint(3) = 1.43*rad
     167    C_Rext(3) = 1.53*rad
     168
     169    Ca_Rint(1) = 1.95*rad
     170    Ca_Rext(1) = 1.99*rad
     171    Ca_Rint(2) = 1.99*rad
     172    Ca_Rext(2) = 2.03*rad
     173
     174
     175    ! Opacities of the rings
     176    tau_A(1) = 1.24
     177    tau_A(2) = 0.81
     178    tau_A(3) = 0.67
     179    tau_A(4) = 0.58
     180               
     181    tau_B(1) = 1.29
     182    tau_B(2) = 5.13
     183    tau_B(3) = 2.84
     184   
     185    tau_C(1) = 0.06
     186    tau_C(2) = 0.10
     187    tau_C(3) = 0.14
     188
     189    tau_Ca(1) = 0.06
     190    tau_Ca(2) = 0.24
     191
     192    ! Convert to cartesian coordinates
     193    N(:) = rad / sqrt(1-(e**2)*sinlat(:)**2)
     194    x(:) = N(:)*coslat(:)*coslon(:)
     195    y(:) = N(:)*coslat(:)*sinlon(:)
     196    z(:) = N(:)*(1-e**2)*sinlat(:)
     197
     198! 2) LOCATION OF THE SUBSOLAR POINT
     199 
     200    ! subsolar longitude is deduced from time fraction ptime
     201    ! SG: the minus sign is important! ... otherwise subsolar point adopts a reverse rotation
     202    phi_S = -(ptime - 0.5)*2.*pi
     203!    write(*,*) 'subsol point coords : ', declin*180./pi, phi_S*180./pi
     204
     205    ! subsolar latitude is declin (declination of the sun)
     206    ! now convert in cartesian coordinates :
     207    Ns = rad/sqrt(1-(e**2)*sin(declin)**2)
     208    xs = Ns*cos(declin)*cos(phi_S)
     209    ys = Ns*cos(declin)*sin(phi_S)
     210    zs = Ns*(1-e**2)*sin(declin)
     211
     212! 3) WHERE DOES THE INCIDENT RAY OF SUN CROSS THE EQUATORIAL PLAN ?
     213
     214    k(:) = -z(:)/zs
     215    r(:) = (k(:)*xs + x(:))**2 + (k(:)*ys + y(:))**2
     216    r(:) = sqrt(r(:))
     217
     218! 4) SO WHERE ARE THE SHADOW OF THESE RINGS ?
     219
     220    ! Summer hemisphere is not under the shadow of the rings
     221    where(latitude(:)*declin .gt. 0.)
     222       eclipse(:) = 1000.
     223    end where
     224
     225    ! No shadow of the rings by night
     226    where(x(:)*xs + y(:)*ys + z(:)*zs .lt. 0.)
     227       eclipse(:) = 1000.
     228    end where
     229
     230    ! if the incident rays of sun cross a ring, there is a shadow
     231    do i=1, nb_A
     232        where(r(:) .ge. A_Rint(i) .and. r(:) .le. A_Rext(i) .and. eclipse(:) .ne. 1000.)
     233            eclipse(:) = 1. - exp(-tau_A(i)/abs(sin(declin)))
     234        end where
     235    end do
     236
     237    do i=1, nb_B
     238        where(r(:) .ge. B_Rint(i) .and. r(:) .le. B_Rext(i) .and. eclipse(:) .ne. 1000.)
     239            eclipse(:) = 1. - exp(-tau_B(i)/abs(sin(declin)))
     240        end where
     241    enddo
     242   
     243    do i=1, nb_C
     244        where(r(:) .ge. C_Rint(i) .and. r(:) .le. C_Rext(i) .and. eclipse(:) .ne. 1000.)
     245            eclipse(:) = 1. - exp(-tau_C(i)/abs(sin(declin)))
     246        end where
     247    enddo
     248
     249    do i=1, nb_ca
     250        where(r(:) .ge. Ca_Rint(i) .and. r(:) .le. Ca_Rext(i) .and. eclipse(:) .ne. 1000.)
     251            eclipse(:) = 1. - exp(-tau_Ca(i)/abs(sin(declin)))
     252        end where
     253    enddo
     254
     255    ! At the other places and the excluded ones, eclipse is 0.
     256    where(eclipse(:) .eq. 2000. .or. eclipse(:) .eq. 1000.)
     257        eclipse(:) = 0.
     258    end where
     259
     260! 5) CLEAN THE PLACE
     261    DEALLOCATE(x)
     262    DEALLOCATE(y)
     263    DEALLOCATE(z)
     264    DEALLOCATE(k)
     265    DEALLOCATE(N)
     266    DEALLOCATE(r)
     267
     268END SUBROUTINE saturn_rings
  • trunk/LMDZ.GENERIC/libf/phygeneric/rad_tridiagonal_matrix_solver.F

    r4076 r4077  
    1       SUBROUTINE DSOLVER(NL,GAMA,CP,CM,CPM1,CMM1,E1,E2,E3,E4,BTOP,
    2      *                   BSURF,RSF,XK1,XK2)
     1      SUBROUTINE rad_tridiagonal_matrix_solver(NL,
     2     *           GAMA,CP,CM,CPM1,CMM1,E1,E2,E3,E4,
     3     *           BTOP,BSURF,RSF,XK1,XK2)
    34
    45C  GCM2.0  Feb 2003
  • trunk/LMDZ.GENERIC/libf/phygeneric/radcommon_h.F90

    r4055 r4077  
    6363!
    6464
    65       REAL*8 BWNI(L_NSPECTI+1), WNOI(L_NSPECTI), DWNI(L_NSPECTI), WAVEI(L_NSPECTI) !BWNI read by master in setspi
    66       REAL*8 BWNV(L_NSPECTV+1), WNOV(L_NSPECTV), DWNV(L_NSPECTV), WAVEV(L_NSPECTV) !BWNV read by master in setspv
     65      REAL*8 BWNI(L_NSPECTI+1), WNOI(L_NSPECTI), DWNI(L_NSPECTI), WAVEI(L_NSPECTI) !BWNI read by master in rad_correlatedk_init_thermal
     66      REAL*8 BWNV(L_NSPECTV+1), WNOV(L_NSPECTV), DWNV(L_NSPECTV), WAVEV(L_NSPECTV) !BWNV read by master in rad_correlatedk_init_stellar
    6767      REAL*8 STELLARF(L_NSPECTV)
    6868!$OMP THREADPRIVATE(WNOI,DWNI,WAVEI,&
     
    8181      real*8 pgasmin, pgasmax
    8282      real*8 tgasmin, tgasmax
    83 !$OMP THREADPRIVATE(gasi,gasv,&  !wrefvar,pgasref,tgasref,pfgasref read by master in sugas_corrk
    84         !$OMP FZEROI,FZEROV)     !pgasmin,pgasmax,tgasmin,tgasmax read by master in sugas_corrk
     83!$OMP THREADPRIVATE(gasi,gasv,&  !wrefvar,pgasref,tgasref,pfgasref read by master in rad_correlatedk_read_opacity_tables
     84        !$OMP FZEROI,FZEROV)     !pgasmin,pgasmax,tgasmin,tgasmax read by master in rad_correlatedk_read_opacity_tables
    8585
    8686      real,allocatable,save :: QVISsQREF(:,:,:)
     
    102102
    103103      integer,allocatable,save :: nsize(:,:)
    104 !$OMP THREADPRIVATE(nsize) ! nsize filled by suaer_corrk
     104!$OMP THREADPRIVATE(nsize) ! nsize filled by rad_correlatedk_ini_aerosol
    105105
    106106! Particle size axis (depend on the kind of aerosol and the
     
    111111
    112112! Extinction coefficient at reference wavelengths;
    113 !   These wavelengths are defined in aeroptproperties, and called
     113!   These wavelengths are defined in aerosol_optical_properties, and called
    114114!   longrefvis and longrefir.
    115115
     
    127127      real*8,parameter :: UBARI = 0.5D0
    128128
    129 !$OMP THREADPRIVATE(QREFvis,QREFir,omegaREFir,&         ! gweight read by master in sugas_corrk
     129!$OMP THREADPRIVATE(QREFvis,QREFir,omegaREFir,&         ! gweight read by master in rad_correlatedk_read_opacity_tables
    130130                !$OMP tstellar,planckir,PTOP)
    131131
  • trunk/LMDZ.GENERIC/libf/phygeneric/radinc_h.F90

    r2972 r4077  
    6464!$OMP THREADPRIVATE(L_NLAYRAD,L_LEVELS,L_NLEVRAD)
    6565
    66       ! These are set in sugas_corrk
     66      ! These are set in rad_correlatedk_read_opacity_tables
    6767      ! [uses allocatable arrays] -- AS 12/2011
    68       integer :: L_NPREF, L_NTREF, L_REFVAR, L_PINT, L_NGAUSS  !L_NPREF, L_NTREF, L_REFVAR, L_PINT, L_NGAUSS read by master in sugas_corrk
     68      integer :: L_NPREF, L_NTREF, L_REFVAR, L_PINT, L_NGAUSS  !L_NPREF, L_NTREF, L_REFVAR, L_PINT, L_NGAUSS read by master in rad_correlatedk_read_opacity_tables
    6969
    7070      integer, parameter :: L_NSPECTI = NBinfrared
  • trunk/LMDZ.GENERIC/libf/phygeneric/rain.F90

    r3893 r4077  
    44  use ioipsl_getin_p_mod, only: getin_p
    55  use watercommon_h, only: T_h2O_ice_liq,T_h2O_ice_clouds, RLVTT, RCPD, RCPV, RW, RVTMP2,Psat_water,Tsat_water,rhowater
    6   use radii_mod, only: h2o_cloudrad
     6  use aerosol_radius, only: h2o_cloudrad
    77  USE tracer_h, only: igcm_h2o_vap, igcm_h2o_ice
    88  use comcstfi_mod, only: g, r
  • trunk/LMDZ.GENERIC/libf/phygeneric/rain_generic.F90

    r3893 r4077  
    66   use watercommon_h, only: T_h2O_ice_liq,T_h2O_ice_clouds,rhowater
    77   ! T_h2O_ice_clouds,rhowater  are only used for precip_scheme_generic >=2
    8    use radii_mod, only: h2o_cloudrad ! only used for precip_scheme_generic >=2
     8   use aerosol_radius, only: h2o_cloudrad ! only used for precip_scheme_generic >=2
    99   use tracer_h
    1010   use comcstfi_mod, only: g, r, cpp
  • trunk/LMDZ.GENERIC/libf/phygeneric/su_gases.F90

    r3893 r4077  
    145145
    146146     if(count.ne.ngasmx)then
    147         print*,'Mismatch between ngas and number of recognised gases in sugas_corrk.F90.'
     147        print*,'Mismatch between ngas and number of recognised gases in rad_correlatedk_read_opacity_tables .F90.'
    148148        print*,'Either we haven`t managed to assign all the gases, or there are duplicates.'
    149149        print*,'Please try again.'
Note: See TracChangeset for help on using the changeset viewer.