Ignore:
Timestamp:
Mar 4, 2011, 6:39:18 PM (14 years ago)
Author:
aslmd
Message:

mars: nouvelle version amelioree de aeroptproperties.F dans phymars. permet de limiter les incertitudes et laisse la possibilite d'augmenter ngau. la derniere version de cette routine sur les sources reference de JB [dans /u/jbmlmd/LMDZ.VERSO/11W02_SP40_LMDZ.BETA.RECHERCHE] n'avait pas ete incluse dans la version de reference SVN [et LMDZ.MARS.BETA].

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/mars/libf/phymars/aeroptproperties.F

    r38 r89  
    3737
    3838c     Min. and max radius of the interpolation grid (in METERS)
    39       REAL, PARAMETER :: refftabmin = 5e-8 !2e-8
     39      REAL, PARAMETER :: refftabmin = 1e-8 !5e-8
    4040      REAL, PARAMETER :: refftabmax = 35e-6
    4141c     Log of the min and max variance of the interpolation grid
    … …  
    6060      REAL kint
    6161      REAL drad
    62       REAL weightgaus(5),radgaus(5)
     62      INTEGER, PARAMETER :: ngau = 10
     63      REAL weightgaus(ngau),radgaus(ngau)
    6364      SAVE weightgaus,radgaus
    64       DATA weightgaus/.2955242247,.2692667193,.2190863625,
    65      &                .1494513491,.0666713443/
    66       DATA radgaus/.1488743389,.4333953941,.6794095682,
    67      &             .8650633666,.9739065285/
     65c     DATA weightgaus/.2955242247,.2692667193,.2190863625,
     66c    &                .1494513491,.0666713443/
     67c     DATA radgaus/.1488743389,.4333953941,.6794095682,
     68c    &             .8650633666,.9739065285/
     69      DATA    radgaus/0.07652652113350,0.22778585114165,
     70     &                0.37370608871528,0.51086700195146,
     71     &                0.63605368072468,0.74633190646476,
     72     &                0.83911697181213,0.91223442826796,
     73     &                0.96397192726078,0.99312859919241/
     74
     75      DATA weightgaus/0.15275338723120,0.14917298659407,
     76     &                0.14209610937519,0.13168863843930,
     77     &                0.11819453196154,0.10193011980823,
     78     &                0.08327674160932,0.06267204829828,
     79     &                0.04060142982019,0.01761400714091/
    6880
    6981c     Indices
    … …  
    106118c     Variables used by the Gauss-Legendre integration:
    107119      REAL,SAVE :: normd(refftabsize,nuefftabsize,naerkind,2)
    108       REAL,SAVE :: dista(refftabsize,nuefftabsize,naerkind,2,5)
    109       REAL,SAVE :: distb(refftabsize,nuefftabsize,naerkind,2,5)
    110 
    111       REAL,SAVE :: radGAUSa(5,naerkind,2)
    112       REAL,SAVE :: radGAUSb(5,naerkind,2)
    113 
    114       REAL,SAVE :: qsqrefVISa(nsun,5,naerkind)
    115       REAL,SAVE :: qrefVISa(5,naerkind)
    116       REAL,SAVE :: qsqrefVISb(nsun,5,naerkind)
    117       REAL,SAVE :: qrefVISb(5,naerkind)
    118       REAL,SAVE :: omegVISa(nsun,5,naerkind)
    119       REAL,SAVE :: omegrefVISa(5,naerkind)
    120       REAL,SAVE :: omegVISb(nsun,5,naerkind)
    121       REAL,SAVE :: omegrefVISb(5,naerkind)
    122       REAL,SAVE :: gVISa(nsun,5,naerkind)
    123       REAL,SAVE :: gVISb(nsun,5,naerkind)
    124 
    125       REAL,SAVE :: qsqrefIRa(nir,5,naerkind)
    126       REAL,SAVE :: qrefIRa(5,naerkind)
    127       REAL,SAVE :: qsqrefIRb(nir,5,naerkind)
    128       REAL,SAVE :: qrefIRb(5,naerkind)
    129       REAL,SAVE :: omegIRa(nir,5,naerkind)
    130       REAL,SAVE :: omegrefIRa(5,naerkind)
    131       REAL,SAVE :: omegIRb(nir,5,naerkind)
    132       REAL,SAVE :: omegrefIRb(5,naerkind)
    133       REAL,SAVE :: gIRa(nir,5,naerkind)
    134       REAL,SAVE :: gIRb(nir,5,naerkind)
     120      REAL,SAVE :: dista(refftabsize,nuefftabsize,naerkind,2,ngau)
     121      REAL,SAVE :: distb(refftabsize,nuefftabsize,naerkind,2,ngau)
     122
     123      REAL,SAVE :: radGAUSa(ngau,naerkind,2)
     124      REAL,SAVE :: radGAUSb(ngau,naerkind,2)
     125
     126      REAL,SAVE :: qsqrefVISa(nsun,ngau,naerkind)
     127      REAL,SAVE :: qrefVISa(ngau,naerkind)
     128      REAL,SAVE :: qsqrefVISb(nsun,ngau,naerkind)
     129      REAL,SAVE :: qrefVISb(ngau,naerkind)
     130      REAL,SAVE :: omegVISa(nsun,ngau,naerkind)
     131      REAL,SAVE :: omegrefVISa(ngau,naerkind)
     132      REAL,SAVE :: omegVISb(nsun,ngau,naerkind)
     133      REAL,SAVE :: omegrefVISb(ngau,naerkind)
     134      REAL,SAVE :: gVISa(nsun,ngau,naerkind)
     135      REAL,SAVE :: gVISb(nsun,ngau,naerkind)
     136
     137      REAL,SAVE :: qsqrefIRa(nir,ngau,naerkind)
     138      REAL,SAVE :: qrefIRa(ngau,naerkind)
     139      REAL,SAVE :: qsqrefIRb(nir,ngau,naerkind)
     140      REAL,SAVE :: qrefIRb(ngau,naerkind)
     141      REAL,SAVE :: omegIRa(nir,ngau,naerkind)
     142      REAL,SAVE :: omegrefIRa(ngau,naerkind)
     143      REAL,SAVE :: omegIRb(nir,ngau,naerkind)
     144      REAL,SAVE :: omegrefIRb(ngau,naerkind)
     145      REAL,SAVE :: gIRa(nir,ngau,naerkind)
     146      REAL,SAVE :: gIRb(nir,ngau,naerkind)
    135147
    136148      REAL :: radiusm
    … …  
    162174      REAL :: omegaREFvis3d(ngridmx,nlayermx,naerkind)
    163175      REAL :: omegaREFir3d(ngridmx,nlayermx,naerkind)
     176
     177c     Tests
     178c     -----
     179
     180      LOGICAL,SAVE :: out_qwg = .false.
     181      INTEGER, PARAMETER :: out_iaer = 1
     182      REAL :: out_qext(ngridmx,nlayermx)
     183      REAL :: out_omeg(ngridmx,nlayermx)
     184      REAL :: out_g(ngridmx,nlayermx)
     185      INTEGER :: out_nchannel
     186      CHARACTER*1 :: out_str
    164187
    165188      DO iaer = 1, naerkind ! Loop on aerosol kind
    … …  
    225248
    226249c       1.5 Interpolating data at the Gauss quadrature points:
    227         DO gausind=1,5
     250        DO gausind=1,ngau
    228251          drad=radiusr*radgaus(gausind)
    229252          radGAUSa(gausind,iaer,idomain)=radiusm-drad
    … …  
    286309        ENDDO
    287310
    288         DO gausind=1,5
     311        DO gausind=1,ngau
    289312          drad=radiusr*radgaus(gausind)
    290313          radGAUSb(gausind,iaer,idomain)=radiusm+drad
    … …  
    391414
    392415                normd(j,1,iaer,idomain) = 1e-30
    393                 DO gausind=1,5
     416                DO gausind=1,ngau
    394417                  drad=radiusr*radgaus(gausind)
    395418                  dista(j,1,iaer,idomain,gausind) =
    … …  
    435458                  omegrefVISgrid(j,1,iaer)=0.
    436459
    437                   DO gausind=1,5
     460                  DO gausind=1,ngau
    438461                    DO m=1,nsun
    439462c                     Convolution:
    … …  
    551574                  omegrefIRgrid(j,1,iaer)=0.
    552575
    553                   DO gausind=1,5
     576                  DO gausind=1,ngau
    554577                    DO m=1,nir
    555578c                     Convolution:
    … …  
    772795
    773796                normd(j,k,iaer,idomain) = 1e-30
    774                 DO gausind=1,5
     797                DO gausind=1,ngau
    775798                  drad=radiusr*radgaus(gausind)
    776799
    … …  
    817840                  omegrefVISgrid(j,k,iaer)=0.
    818841
    819                   DO gausind=1,5
     842                  DO gausind=1,ngau
    820843                    DO m=1,nsun
    821844c                     Convolution:
    … …  
    932955                  omegrefIRgrid(j,k,iaer)=0.
    933956
    934                   DO gausind=1,5
     957                  DO gausind=1,ngau
    935958                    DO m=1,nir
    936959c                     Convolution:
    … …  
    11141137      ENDDO ! iaer (loop on aerosol kind)
    11151138
     1139c=====Radiative properties - TESTS=================================
     1140      IF (out_qwg) THEN
     1141      IF (ngrid.NE.1) THEN
     1142        DO out_nchannel = 1, 2
     1143c     -------------------------------------------------------------
     1144          DO lg = 1, nlayer
     1145            DO ig = 1, ngrid
     1146              out_qext(ig,lg) =
     1147     &           QVISsQREF3d(ig,lg,out_nchannel,out_iaer)*
     1148     &           QREFvis3d(ig,lg,out_iaer)
     1149              out_omeg(ig,lg) =
     1150     &           omegaVIS3d(ig,lg,out_nchannel,out_iaer)
     1151              out_g(ig,lg)    = gVIS3d(ig,lg,out_nchannel,out_iaer)
     1152            ENDDO ! ig
     1153          ENDDO ! lg
     1154          write(out_str(1:1),'(i1.1)') out_nchannel
     1155          call WRITEDIAGFI(ngrid,'qextvis'//out_str,"Ext.efficiency","",
     1156     &                     3,out_qext)
     1157          call WRITEDIAGFI(ngrid,'omegvis'//out_str,"Sing.Scat.Alb.","",
     1158     &                     3,out_omeg)
     1159          call WRITEDIAGFI(ngrid,'gvis'//out_str,"Asym.Factor","",
     1160     &                     3,out_g)
     1161c     -------------------------------------------------------------
     1162        ENDDO ! out_nchannel
     1163        DO out_nchannel = 2, 4
     1164c     -------------------------------------------------------------
     1165          DO lg = 1, nlayer
     1166            DO ig = 1, ngrid
     1167              out_qext(ig,lg) =
     1168     &          QIRsQREF3d(ig,lg,out_nchannel,out_iaer)*
     1169     &          QREFir3d(ig,lg,out_iaer)
     1170              out_omeg(ig,lg) =
     1171     &          omegaIR3d(ig,lg,out_nchannel,out_iaer)
     1172              out_g(ig,lg)    = gIR3d(ig,lg,out_nchannel,out_iaer)
     1173            ENDDO ! ig
     1174          ENDDO ! lg
     1175          write(out_str(1:1),'(i1.1)') out_nchannel
     1176          call WRITEDIAGFI(ngrid,'qextir'//out_str,"Ext.efficiency","",
     1177     &                     3,out_qext)
     1178          call WRITEDIAGFI(ngrid,'omegir'//out_str,"Sing.Scat.Alb.","",
     1179     &                     3,out_omeg)
     1180          call WRITEDIAGFI(ngrid,'gir'//out_str,"Asym.Factor","",
     1181     &                     3,out_g)
     1182c     -------------------------------------------------------------
     1183        ENDDO ! out_nchannel
     1184        call WRITEDIAGFI(ngrid,"omegvisref","Sing.Scat.Alb.","",
     1185     &                   3,omegaREFvis3d(1,1,out_iaer))
     1186        call WRITEDIAGFI(ngrid,"omegirref","Sing.Scat.Alb.","",
     1187     &                   3,omegaREFir3d(1,1,out_iaer))
     1188      ENDIF ! ngrid.EQ.1
     1189      ENDIF ! out_qwg
     1190c==================================================================
     1191
    11161192      RETURN
    11171193      END
    1118 
    1119 
    1120 
    1121 
Note: See TracChangeset for help on using the changeset viewer.