Changeset 3793


Ignore:
Timestamp:
Jun 5, 2025, 1:10:40 PM (3 days ago)
Author:
mturbet
Message:

fix a bug and update new continuum interpolation routine in optci and optcv

Location:
trunk/LMDZ.GENERIC
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.GENERIC/changelog.txt

    r3746 r3793  
    20802080Minor fix in interpolate_continuum: handle the extreme case when temperature
    20812081is below minimum value in the table.
     2082
     2083== 05/06/2025 == MT
     2084Update + fix a bug for the new continuum interpolate routine (in optci and optcv routines)
  • trunk/LMDZ.GENERIC/libf/phystd/optci.F90

    r3693 r3793  
    203203               endif
    204204             
    205                dtemp=0.0
     205               do jgas=1,ngasmx
     206                 if(gfrac(jgas).eq.-1)then ! variable
     207                   p_cross  = dble(PMID(k)*scalep*QVAR(k)) ! qvar = mol/mol
     208                 elseif(varspec) then
     209                   p_cross  = dble(PMID(k)*scalep*FRACVAR(jgas,k)*(1.-QVAR(k)))
     210                 else
     211                   p_cross  = dble(PMID(k)*scalep*gfrac(jgas)*(1.-QVAR(k)))
     212                 endif
    206213               
    207                if (igas .eq. igas_N2) then
    208                 call interpolate_continuum('',igas_N2,igas_N2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    209                 do jgas=1,ngasmx
    210                  if (jgas .eq. igas_H2) then
    211                   call interpolate_continuum('',igas_N2,igas_H2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    212                  elseif (jgas .eq. igas_O2) then
    213                   call interpolate_continuum('',igas_N2,igas_O2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    214                  elseif (jgas .eq. igas_CH4) then
    215                   call interpolate_continuum('',igas_N2,igas_CH4,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    216                  endif
    217                 enddo
    218                elseif (igas .eq. igas_O2) then
    219                 call interpolate_continuum('',igas_O2,igas_O2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    220                 do jgas=1,ngasmx
    221                  if (jgas .eq. igas_CO2) then
    222                   call interpolate_continuum('',igas_CO2,igas_O2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
     214                 dtemp=0.0
     215
     216                 if ( ((igas .eq. igas_N2) .and. (jgas .eq. igas_N2)) .or.   &
     217                     ((igas .eq. igas_N2) .and. (jgas .eq. igas_H2)) .or.    &
     218                     ((igas .eq. igas_N2) .and. (jgas .eq. igas_O2)) .or.    &
     219                     ((igas .eq. igas_N2) .and. (jgas .eq. igas_CH4)) .or.   &
     220                     ((igas .eq. igas_O2) .and. (jgas .eq. igas_O2)) .or.    &
     221                     ((igas .eq. igas_CO2) .and. (jgas .eq. igas_O2)) .or.   &
     222                     ((igas .eq. igas_H2) .and. (jgas .eq. igas_H2)) .or.    &
     223                     ((igas .eq. igas_H2) .and. (jgas .eq. igas_CH4)) .or.   &
     224                     ((igas .eq. igas_H2) .and. (jgas .eq. igas_He)) .or.    &
     225                     ((igas .eq. igas_CH4) .and. (jgas .eq. igas_CH4)) .or.  &
     226                     ((igas .eq. igas_H2O) .and. (jgas .eq. igas_H2O)) .or.  &
     227                     ((igas .eq. igas_H2O) .and. (jgas .eq. igas_N2)) .or.   &
     228                     ((igas .eq. igas_H2O) .and. (jgas .eq. igas_O2)) .or.   &
     229                     ((igas .eq. igas_H2O) .and. (jgas .eq. igas_CO2)) .or.  &
     230                     ((igas .eq. igas_CO2) .and. (jgas .eq. igas_CO2)) .or.  &
     231                     ((igas .eq. igas_CO2) .and. (jgas .eq. igas_H2)) .or.   &
     232                     ((igas .eq. igas_CO2) .and. (jgas .eq. igas_CH4))  ) then
     233
     234                   call interpolate_continuum('',igas,jgas,'IR',nw,T_cont,p_cont,p_cross,dtemp,.false.)
     235
    223236                 endif
    224                 enddo
    225                elseif (igas .eq. igas_H2) then
    226                 call interpolate_continuum('',igas_H2,igas_H2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    227                 do jgas=1,ngasmx
    228                  if (jgas .eq. igas_CH4) then
    229                   call interpolate_continuum('',igas_H2,igas_CH4,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    230                  elseif (jgas .eq. igas_He) then
    231                   call interpolate_continuum('',igas_H2,igas_He,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    232                  endif
    233                 enddo   
    234                elseif (igas .eq. igas_CH4) then
    235                 call interpolate_continuum('',igas_CH4,igas_CH4,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    236                elseif (igas .eq. igas_H2O) then
    237                 call interpolate_continuum('',igas_H2O,igas_H2O,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    238                 do jgas=1,ngasmx
    239                  if (jgas .eq. igas_N2) then
    240                   call interpolate_continuum('',igas_H2O,igas_N2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    241                  elseif (jgas .eq. igas_O2) then
    242                   call interpolate_continuum('',igas_H2O,igas_O2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    243                  elseif (jgas .eq. igas_CO2) then
    244                   call interpolate_continuum('',igas_H2O,igas_CO2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    245                  endif
    246                 enddo   
    247                elseif (igas .eq. igas_CO2) then
    248                 call interpolate_continuum('',igas_CO2,igas_CO2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    249                 do jgas=1,ngasmx
    250                  if (jgas .eq. igas_H2) then
    251                   call interpolate_continuum('',igas_CO2,igas_H2,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    252                  elseif (jgas .eq. igas_CH4) then
    253                   call interpolate_continuum('',igas_CO2,igas_CH4,'IR',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    254                  endif
    255                 enddo
    256                endif
    257237               
    258                DCONT = DCONT + dtemp
     238                 DCONT = DCONT + dtemp
     239                 
     240               enddo ! jgas=1,ngasmx
    259241               
    260242             enddo ! igas=1,ngasmx
     243             
    261244           else ! generic_continuum_database
     245           
    262246            wn_cont = dble(wnoi(nw))
    263247            T_cont  = dble(TMID(k))
  • trunk/LMDZ.GENERIC/libf/phystd/optcv.F90

    r3696 r3793  
    195195 
    196196  !     we ignore K=1...
     197 
    197198  do K=2,L_LEVELS
    198199
    199200     do NW=1,L_NSPECTV
    200 
     201     
    201202        DRAYAER = TRAY(K,NW)
    202203        !     DRAYAER is Tau RAYleigh scattering, plus AERosol opacity
     
    222223              endif
    223224             
    224               dtemp=0.0
    225              
    226               if (igas .eq. igas_N2) then
    227                call interpolate_continuum('',igas_N2,igas_N2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    228                do jgas=1,ngasmx
    229                 if (jgas .eq. igas_H2) then
    230                  call interpolate_continuum('',igas_N2,igas_H2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    231                 elseif (jgas .eq. igas_O2) then
    232                  call interpolate_continuum('',igas_N2,igas_O2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    233                 elseif (jgas .eq. igas_CH4) then
    234                  call interpolate_continuum('',igas_N2,igas_CH4,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
     225              do jgas=1,ngasmx
     226                if(gfrac(jgas).eq.-1)then ! variable
     227                  p_cross  = dble(PMID(k)*scalep*QVAR(k)) ! qvar = mol/mol
     228                elseif(varspec) then
     229                  p_cross  = dble(PMID(k)*scalep*FRACVAR(jgas,k)*(1.-QVAR(k)))
     230                else
     231                  p_cross  = dble(PMID(k)*scalep*gfrac(jgas)*(1.-QVAR(k)))
    235232                endif
    236                enddo
    237               elseif (igas .eq. igas_O2) then
    238                call interpolate_continuum('',igas_O2,igas_O2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    239                do jgas=1,ngasmx
    240                 if (jgas .eq. igas_CO2) then
    241                  call interpolate_continuum('',igas_CO2,igas_O2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
     233               
     234                dtemp=0.0
     235
     236                if ( ((igas .eq. igas_N2) .and. (jgas .eq. igas_N2)) .or.    &
     237                     ((igas .eq. igas_N2) .and. (jgas .eq. igas_H2)) .or.    &
     238                     ((igas .eq. igas_N2) .and. (jgas .eq. igas_O2)) .or.    &
     239                     ((igas .eq. igas_N2) .and. (jgas .eq. igas_CH4)) .or.   &
     240                     ((igas .eq. igas_O2) .and. (jgas .eq. igas_O2)) .or.    &
     241                     ((igas .eq. igas_CO2) .and. (jgas .eq. igas_O2)) .or.   &
     242                     ((igas .eq. igas_H2) .and. (jgas .eq. igas_H2)) .or.    &
     243                     ((igas .eq. igas_H2) .and. (jgas .eq. igas_CH4)) .or.   &
     244                     ((igas .eq. igas_H2) .and. (jgas .eq. igas_He)) .or.    &
     245                     ((igas .eq. igas_CH4) .and. (jgas .eq. igas_CH4)) .or.  &
     246                     ((igas .eq. igas_H2O) .and. (jgas .eq. igas_H2O)) .or.  &
     247                     ((igas .eq. igas_H2O) .and. (jgas .eq. igas_N2)) .or.   &
     248                     ((igas .eq. igas_H2O) .and. (jgas .eq. igas_O2)) .or.   &
     249                     ((igas .eq. igas_H2O) .and. (jgas .eq. igas_CO2)) .or.  &
     250                     ((igas .eq. igas_CO2) .and. (jgas .eq. igas_CO2)) .or.  &
     251                     ((igas .eq. igas_CO2) .and. (jgas .eq. igas_H2)) .or.   &
     252                     ((igas .eq. igas_CO2) .and. (jgas .eq. igas_CH4)) ) then
     253
     254                  call interpolate_continuum('',igas,jgas,'VI',nw,T_cont,p_cont,p_cross,dtemp,.false.)
     255
    242256                endif
    243                enddo
    244               elseif (igas .eq. igas_H2) then
    245                call interpolate_continuum('',igas_H2,igas_H2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    246                do jgas=1,ngasmx
    247                 if (jgas .eq. igas_CH4) then
    248                  call interpolate_continuum('',igas_H2,igas_CH4,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    249                 elseif (jgas .eq. igas_He) then
    250                  call interpolate_continuum('',igas_H2,igas_He,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    251                 endif
    252                enddo     
    253               elseif (igas .eq. igas_CH4) then
    254                call interpolate_continuum('',igas_CH4,igas_CH4,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    255               elseif (igas .eq. igas_H2O) then
    256                call interpolate_continuum('',igas_H2O,igas_H2O,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    257                do jgas=1,ngasmx
    258                 if (jgas .eq. igas_N2) then
    259                  call interpolate_continuum('',igas_H2O,igas_N2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    260                 elseif (jgas .eq. igas_O2) then
    261                  call interpolate_continuum('',igas_H2O,igas_O2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    262                 elseif (jgas .eq. igas_CO2) then
    263                  call interpolate_continuum('',igas_H2O,igas_CO2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    264                 endif
    265                enddo     
    266               elseif (igas .eq. igas_CO2) then
    267                call interpolate_continuum('',igas_CO2,igas_CO2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    268                do jgas=1,ngasmx
    269                 if (jgas .eq. igas_H2) then
    270                  call interpolate_continuum('',igas_CO2,igas_H2,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    271                 elseif (jgas .eq. igas_CH4) then
    272                  call interpolate_continuum('',igas_CO2,igas_CH4,'VI',nw,T_cont,p_cont,p_cont,dtemp,.false.)
    273                 endif
    274                enddo
    275               endif
    276                
    277               DCONT = DCONT + dtemp
     257               
     258                DCONT = DCONT + dtemp
     259               
     260              enddo ! jgas=1,ngasmx
    278261               
    279262            enddo ! igas=1,ngasmx
     
    456439  end do
    457440
    458 
    459441  !=======================================================================
    460442  !     Now the full treatment for the layers, where besides the opacity
Note: See TracChangeset for help on using the changeset viewer.