Ignore:
Timestamp:
Feb 12, 2013, 12:26:32 PM (12 years ago)
Author:
jleconte
Message:

12/02/2013 == JL

  • Follows previous commit by Aymeric about bilinear interpolations:
    • Extended to all existing continua
    • generalized bilinearbig to work for various size inputs
  • because N2 and H2O continua databases are smaller, improvement around 15% for

an earth case.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.GENERIC/libf/phystd/interpolateN2N2.F90

    r716 r878  
    1 subroutine interpolateN2N2(wn,temp,pres,abcoef,firstcall)
     1subroutine interpolateN2N2(wn,temp,pres,abcoef,firstcall,ind)
    22
    33  !==================================================================
     
    2121  double precision temp               ! temperature            (Kelvin)
    2222  double precision pres               ! pressure               (Pascals)
     23  integer :: ind
    2324
    2425  ! output
     
    104105     print*,'   pressure ',pres,' Pa'
    105106
    106      call bilinearN2N2(wn_arr,temp_arr,abs_arr,wn,temp,abcoef)
     107  endif
     108     call bilinearbig(nS,nT,wn_arr,temp_arr,abs_arr,wn,temp,abcoef,ind)
    107109
    108      print*,'the absorption is ',abcoef,' cm^5 molecule^-2'
    109      print*,'or ',abcoef*losch**2,' cm^-1 amagat^-2'
     110!     print*,'the absorption is ',abcoef,' cm^5 molecule^-2'
     111!     print*,'or ',abcoef*losch**2,' cm^-1 amagat^-2'
    110112
    111113     abcoef=abcoef*losch**2*100.0*amagat**2 ! convert to m^-1
     114!     abcoef=0.
    112115
    113      print*,'We have ',amagat,' amagats of N2'
    114      print*,'So the absorption is ',abcoef,' m^-1'
     116!     print*,'We have ',amagat,' amagats of N2'
     117!     print*,'So the absorption is ',abcoef,' m^-1'
    115118
    116   else
     119!     unlike for Rayleigh scattering, we do not currently weight by the BB function
     120!     however our bands are normally thin, so this is no big deal.
    117121
    118      call bilinearN2N2(wn_arr,temp_arr,abs_arr,wn,temp,abcoef)
    119      abcoef=abcoef*losch**2*100.0*amagat**2 ! convert to m^-1
    120 
    121      ! unlike for Rayleigh scattering, we do not currently weight by the BB function
    122      ! however our bands are normally thin, so this is no big deal.
    123 
    124   endif
    125122
    126123  return
    127124end subroutine interpolateN2N2
    128125
    129 
    130 !-------------------------------------------------------------------------
    131       subroutine bilinearN2N2(x_arr,y_arr,f2d_arr,x_in,y_in,f)
    132 
    133       implicit none
    134 
    135       integer nX,nY,i,j,a,b
    136       parameter(nX=582)
    137       parameter(nY=10)
    138 
    139       real*8 x_in,y_in,x,y,x1,x2,y1,y2
    140       real*8 f,f11,f12,f21,f22,fA,fB
    141       real*8 x_arr(nX)
    142       real*8 y_arr(nY)
    143       real*8 f2d_arr(nX,nY)
    144      
    145       integer strlen
    146       character*100 label
    147       label='subroutine bilinear'
    148 
    149       x=x_in
    150       y=y_in
    151 
    152 !     1st check we're within the wavenumber range
    153       if ((x.lt.x_arr(2)).or.(x.gt.x_arr(nX-2))) then
    154          f=0.0D+0
    155          return
    156       else
    157          
    158 !     in the x (wavenumber) direction 1st
    159          i=1
    160  10      if (i.lt.(nX+1)) then
    161             if (x_arr(i).gt.x) then
    162                x1=x_arr(i-1)
    163                x2=x_arr(i)
    164                a=i-1
    165                i=9999
    166             endif
    167             i=i+1
    168             goto 10
    169          endif
    170       endif
    171      
    172       if ((y.lt.y_arr(1)).or.(y.gt.y_arr(nY))) then
    173          write(*,*) 'Warning from bilinearN2N2:'
    174          write(*,*) 'Outside continuum temperature range!'
    175          if(y.lt.y_arr(1))then
    176             y=y_arr(1)+0.01
    177          endif
    178          if(y.gt.y_arr(nY))then
    179             y=y_arr(nY)-0.01
    180          endif
    181       else
    182 
    183 !     in the y (temperature) direction 2nd
    184          j=1
    185  20      if (j.lt.(nY+1)) then
    186             if (y_arr(j).gt.y) then
    187                y1=y_arr(j-1)
    188                y2=y_arr(j)
    189                b=j-1
    190                j=9999
    191             endif
    192             j=j+1
    193             goto 20
    194          endif
    195       endif
    196      
    197       f11=f2d_arr(a,b)
    198       f21=f2d_arr(a+1,b)
    199       f12=f2d_arr(a,b+1)
    200       f22=f2d_arr(a+1,b+1)
    201 
    202       call bilinear(f,f11,f21,f12,f22,x,x1,x2,y,y1,y2)
    203 
    204       return
    205     end subroutine bilinearN2N2
Note: See TracChangeset for help on using the changeset viewer.