Ignore:
Timestamp:
Sep 24, 2016, 7:14:59 PM (8 years ago)
Author:
oboucher
Message:

Quicker way to compute aerosol optical properties.
Does not converge with previous version but results are close.

Location:
LMDZ5/trunk/libf/phylmd
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • LMDZ5/trunk/libf/phylmd/aeropt_2bands.F90

    r2593 r2634  
    1212  USE phys_local_var_mod, only: absvisaer
    1313  USE pres2lev_mod
    14 
    1514
    1615  !    Yves Balkanski le 12 avril 2006
    … …  
    128127  REAL :: zrho
    129128  REAL :: fac
    130   REAL :: zdp1(klon,klev)
    131   REAL, PARAMETER ::  gravit = 9.80616    ! m2/s
     129  REAL :: zdh(klon,klev)
    132130  INTEGER, ALLOCATABLE, DIMENSION(:)   :: aerosol_name
    133131  INTEGER :: nb_aer
    134   REAL, DIMENSION(klon,klev,naero_tot) :: mass_temp
    135 !RAF
    136   REAL, DIMENSION(klon,klev,naero_tot) :: mass_temp_pi
    137132
    138133  !
    … …  
    587582  END IF ! firstcall
    588583
    589 
    590584  DO k=1, klev
    591585    DO i=1, klon
    592       zrho=pplay(i,k)/t_seri(i,k)/RD                  ! kg/m3
    593       mass_temp(i,k,:) = m_allaer(i,k,:) / zrho / 1.e+9
    594       mass_temp_pi(i,k,:) = m_allaer_pi(i,k,:) / zrho / 1.e+9
    595       zdp1(i,k)=pdel(i,k)/(gravit*delt)      ! air mass auxiliary  variable --> zdp1 [kg/(m^2 *s)]
     586      zrho=pplay(i,k)/t_seri(i,k)/RD    ! kg/m3
     587      zdh(i,k)=pdel(i,k)/(RG*zrho)      ! m
    596588    ENDDO
    597589  ENDDO
    … …  
    645637  !      compute optical_thickness_at_gridpoint_per_species
    646638
    647 
    648 
    649639!CDIR ON_ADB(fact_RH)
    650640!CDIR SHORTLOOP
    … …  
    660650      IF (rh(i,k).GT.85.) RH_num(i,k)=10
    661651      IF (rh(i,k).GT.90.) RH_num(i,k)=11
    662      
    663652      DELTA(i,k)=(rh(i,k)-RH_tab(RH_num(i,k)))*fact_RH(RH_num(i,k))
    664653    ENDDO
    … …  
    718707    used_aer(id)=.TRUE.
    719708
    720      
    721709    IF (soluble) THEN
    722710
    … …  
    739727!CDIR ON_ADB(C1_ASSSM_b2)
    740728!CDIR ON_ADB(C2_ASSSM_b2)
     729
    741730              DO i=1, KLON
     731
    742732                H=rh(i,k)/100
    743                 tmp_var=mass_temp(i,k,spsol)*1000.*zdp1(i,k)*delt*fac
    744                 tmp_var_pi=mass_temp_pi(i,k,spsol)*1000.*zdp1(i,k)*delt*fac
     733                tmp_var=m_allaer(i,k,spsol)/1.e6*zdh(i,k)*fac
     734                tmp_var_pi=m_allaer_pi(i,k,spsol)/1.e6*zdh(i,k)*fac
    745735
    746736                ! band 1
    … …  
    784774!CDIR ON_ADB(C1_CSSSM_b2)
    785775!CDIR ON_ADB(C2_CSSSM_b2)
     776
    786777              DO i=1, KLON
     778
    787779                H=rh(i,k)/100
    788                 tmp_var=mass_temp(i,k,spsol)*1000.*zdp1(i,k)*delt*fac
    789                 tmp_var_pi=mass_temp_pi(i,k,spsol)*1000.*zdp1(i,k)*delt*fac
     780                tmp_var=m_allaer(i,k,spsol)/1.e6*zdh(i,k)*fac
     781                tmp_var_pi=m_allaer_pi(i,k,spsol)/1.e6*zdh(i,k)*fac
    790782                ! band 1
    791783                tau_ae2b_int=A1_CSSSM_b1(k)+A2_CSSSM_b1(k)*H+A3_CSSSM_b1(k)/(H-1.05)
    … …  
    828820!CDIR ON_ADB(C1_SSSSM_b2)
    829821!CDIR ON_ADB(C2_SSSSM_b2)
     822
    830823              DO i=1, KLON
     824
    831825                H=rh(i,k)/100
    832                 tmp_var=mass_temp(i,k,spsol)*1000.*zdp1(i,k)*delt*fac
    833                 tmp_var_pi=mass_temp_pi(i,k,spsol)*1000.*zdp1(i,k)*delt*fac
     826                tmp_var=m_allaer(i,k,spsol)/1.e6*zdh(i,k)*fac
     827                tmp_var_pi=m_allaer_pi(i,k,spsol)/1.e6*zdh(i,k)*fac
    834828
    835829                ! band 1
    … …  
    864858          DO k=1, KLEV
    865859            DO i=1, KLON
    866               tmp_var=mass_temp(i,k,spsol)*1000.*zdp1(i,k)*delt*fac
    867               tmp_var_pi=mass_temp_pi(i,k,spsol)*1000.*zdp1(i,k)*delt*fac
     860              tmp_var=m_allaer(i,k,spsol)/1.e6*zdh(i,k)*fac
     861              tmp_var_pi=m_allaer_pi(i,k,spsol)/1.e6*zdh(i,k)*fac
    868862!CDIR UNROLL=nbands
    869863              DO inu=1,nbands
    … …  
    899893        DO k=1, KLEV
    900894          DO i=1, KLON
    901             tmp_var=mass_temp(i,k,naero_soluble+ spinsol)*1000.*zdp1(i,k)*delt*fac
    902             tmp_var_pi=mass_temp_pi(i,k,naero_soluble+spinsol)*1000.*zdp1(i,k)*delt*fac
     895            tmp_var=m_allaer(i,k,naero_soluble+spinsol)/1.e6*zdh(i,k)*fac
     896            tmp_var_pi=m_allaer_pi(i,k,naero_soluble+spinsol)/1.e6*zdh(i,k)*fac
    903897!CDIR UNROLL=nbands
    904898            DO inu=1,nbands
    … …  
    931925
    932926  DO inu=1, nbands
    933     DO mrfspecies=1,naero_grp
     927
     928    !!DO mrfspecies=1,naero_grp
     929    DO mrfspecies=2,3    !--only deal with total and natural aerosols
     930
    934931      IF (mrfspecies .EQ. 2) THEN             ! = total aerosol AER     
     932
    935933        DO k=1, KLEV
    936934          DO i=1, KLON
    … …  
    10261024                   
    10271025      ELSEIF (mrfspecies .EQ. 4) THEN             ! = BC
     1026
    10281027        DO k=1, KLEV
    10291028          DO i=1, KLON
    … …  
    11281127  ENDDO
    11291128
    1130   inu=1         ! visible wavaband
     1129  inu=1         ! visible waveband
    11311130  mrfspecies=2  ! total aerosol AER     
    11321131  DO i=1, KLON

  • LMDZ5/trunk/libf/phylmd/aeropt_5wv.F90

    r2550 r2634  
    1313  USE phys_local_var_mod, only: od550aer,od865aer,ec550aer,od550lt1aer
    1414  USE pres2lev_mod
    15 
    1615
    1716  !
    … …  
    6665  ! Output arguments:
    6766  !
    68   REAL, DIMENSION(klon), INTENT(out)          :: ai      ! POLDER aerosol index 
     67  REAL, DIMENSION(klon), INTENT(out)          :: ai      ! POLDER aerosol index
    6968  REAL, DIMENSION(klon,nwave,naero_tot), INTENT(out)      :: tausum
    7069  REAL, DIMENSION(klon,klev,nwave,naero_tot), INTENT(out) :: tau
    71 
    72 
    7370  !
    7471  ! Local
    … …  
    134131  REAL :: cg_ae5wv_int  ! Intermediate asymmetry parameter aerosol
    135132  REAL, PARAMETER :: RH_MAX=95.
    136   REAL :: taue670(KLON)       ! epaisseur optique aerosol absorption 550 nm
    137   REAL :: taue865(KLON)       ! epaisseur optique aerosol extinction 865 nm
     133  REAL :: taue670(klon)       ! epaisseur optique aerosol absorption 550 nm
     134  REAL :: taue865(klon)       ! epaisseur optique aerosol extinction 865 nm
    138135  REAL :: fac
    139   REAL :: zdp1(klon,klev)
    140   REAL, PARAMETER ::  gravit = 9.80616    ! m2/s
    141136  INTEGER, ALLOCATABLE, DIMENSION(:)  :: aerosol_name
    142137  INTEGER :: nb_aer
    143138 
    144   REAL :: tau3d(KLON,KLEV), piz3d(KLON,KLEV), cg3d(KLON,KLEV)
    145   REAL :: abs3d(KLON,KLEV)     ! epaisseur optique d'absorption
    146   REAL :: dh(KLON,KLEV)
     139  REAL :: tau3d(klon,klev), piz3d(klon,klev), cg3d(klon,klev)
     140  REAL :: abs3d(klon,klev)     ! epaisseur optique d'absorption
     141  REAL :: dh(klon,klev)
    147142 
    148143  REAL :: alpha_aers_5wv(nbre_RH,las,naero_soluble)   ! ext. coeff. Soluble comp. units *** m2/g
    … …  
    155150   !  1- BC soluble; 2- POM soluble; 3- SO4 acc.; 4- SO4 coarse; 5 seasalt super-C; 6 seasalt coarse; 7 seasalt acc.
    156151  REAL :: piz_aeri_5wv(las,naero_insoluble)           ! Insoluble comp. 1- Dust: 2- BC; 3- POM
    157 
    158   REAL, DIMENSION(klon,klev,naero_tot) :: mass_temp
    159  
    160152  !
    161153  ! Proprietes optiques
    162154  !
    163   REAL :: radry = 287.054
    164   REAL :: tau_tmp                     ! dry air mass constant
    165155  REAL :: fact_RH(nbre_RH)
    166156  LOGICAL :: used_tau(naero_spc)
    … …  
    303293!*********************************************************************
    304294!
    305 !
    306 !
    307 !
    308 ! 
    309295!
    310296! From here on we look at the optical parameters at 5 wavelengths: 
    … …  
    478464      0.737,0.750,0.765,0.775,0.787,0.803/
    479465 !
    480 
    481466  DATA cg_aeri_5wv/&
    482467     ! dust insoluble
    … …  
    616601  END IF ! firstcall
    617602
    618 
    619603  ! Initialisations
    620   ai(:) = 0.
    621   tausum(:,:,:) = 0.
    622 
     604  tau(:,:,:,:) =0.
     605  tausum(:,:,:)= 0.
     606  ai(:)=0.0
    623607
    624608  DO k=1, klev
    625609    DO i=1, klon
    626       zrho=pplay(i,k)/t_seri(i,k)/RD                  ! kg/m3
    627       dh(i,k)=pdel(i,k)/(gravit*zrho)
    628       mass_temp(i,k,:) = m_allaer(i,k,:) / zrho / 1.e+9
    629       zdp1(i,k)=pdel(i,k)/(gravit*delt)     ! air mass auxiliary  variable --> zdp1 [kg/(m^2 *s)]
    630 
     610      zrho=pplay(i,k)/t_seri(i,k)/RD     ! kg/m3
     611      dh(i,k)=pdel(i,k)/(RG*zrho)    ! m
    631612    ENDDO
    632613  ENDDO
    633 
    634614
    635615  IF (flag_aerosol .EQ. 1) THEN
    … …  
    680660  !      compute optical_thickness_at_gridpoint_per_species
    681661 
    682 
    683662  !
    684663  ! Calculations that need to be done since we are not in the subroutines INCA
    685664  !     
    686665
    687 !CDIR ON_ADB(RH_tab)
    688 !CDIR ON_ADB(fact_RH)
    689 !CDIR NOVECTOR
    690666  DO n=1,nbre_RH-1
    691667    fact_RH(n)=1./(RH_tab(n+1)-RH_tab(n))
    692668  ENDDO
    693669   
    694   DO k=1, KLEV
    695 !CDIR ON_ADB(RH_tab)
    696 !CDIR ON_ADB(fact_RH)
    697     DO i=1, KLON
     670  DO k=1, klev
     671    DO i=1, klon
    698672      rh(i,k)=MIN(RHcl(i,k)*100.,RH_MAX)
    699673      RH_num(i,k) = INT( rh(i,k)/10. + 1.)
     674!--test olivier pour pas de reindicage
     675!      RH_num(i,k) =1
    700676      IF (rh(i,k).GT.85.) RH_num(i,k)=10
    701677      IF (rh(i,k).GT.90.) RH_num(i,k)=11
    … …  
    704680  ENDDO
    705681
    706 !CDIR SHORTLOOP 
    707682  used_tau(:)=.FALSE.
    708683   
    … …  
    765740    DO la=1,las
    766741
     742    !--only 550 and 865 nm are used
     743    IF (la.NE.la550.AND.la.NE.la865) CYCLE
     744
    767745      IF (soluble) THEN
    768746
    769         IF((la.EQ.2).AND.(spss.NE.0)) THEN !la=2 corresponds to 550 nm
     747        IF ((la.EQ.2).AND.(spss.NE.0)) THEN !la=2 corresponds to 550 nm
    770748          IF (spss.EQ.1) THEN !accumulation mode
    771             DO k=1, KLEV
    772 !CDIR ON_ADB(A1_ASSSM)
    773 !CDIR ON_ADB(A2_ASSSM)
    774 !CDIR ON_ADB(A3_ASSSM)
    775               DO i=1, KLON
    776                 H=rh(i,k)/100
     749            DO k=1, klev
     750              DO i=1, klon
     751                H=rh(i,k)/100.
    777752                tau_ae5wv_int=A1_ASSSM(k)+A2_ASSSM(k)*H+A3_ASSSM(k)/(H-1.05)
    778                 tau(i,k,la,aerindex) = mass_temp(i,k,aerindex)*1000.*zdp1(i,k)   &
    779                                    *tau_ae5wv_int*delt*fac
     753                tau(i,k,la,aerindex) = m_allaer(i,k,aerindex)/1.e6*dh(i,k)*tau_ae5wv_int*fac
    780754                tausum(i,la,aerindex)=tausum(i,la,aerindex)+tau(i,k,la,aerindex)
    781755              ENDDO
    … …  
    784758 
    785759          IF (spss.EQ.2) THEN !coarse mode
    786             DO k=1, KLEV
    787 !CDIR ON_ADB(A1_CSSSM)
    788 !CDIR ON_ADB(A2_CSSSM)
    789 !CDIR ON_ADB(A3_CSSSM)
    790               DO i=1, KLON
    791                 H=rh(i,k)/100
     760            DO k=1, klev
     761              DO i=1, klon
     762                H=rh(i,k)/100.
    792763                tau_ae5wv_int=A1_CSSSM(k)+A2_CSSSM(k)*H+A3_CSSSM(k)/(H-1.05)
    793                 tau(i,k,la,aerindex) = mass_temp(i,k,aerindex)*1000.*zdp1(i,k)  &
    794                                    *tau_ae5wv_int*delt*fac
    795                 tausum(i,la,aerindex) = tausum(i,la,aerindex)+tau(i,k,la,aerindex)
    796               ENDDO
    797             ENDDO
    798           ENDIF
    799 
    800           IF (spss.EQ.3) THEN !super coarse mode
    801             DO k=1, KLEV
    802 !CDIR ON_ADB(A1_SSSSM)
    803 !CDIR ON_ADB(A2_SSSSM)
    804 !CDIR ON_ADB(A3_SSSSM)
    805               DO i=1, KLON
    806                 H=rh(i,k)/100
    807                 tau_ae5wv_int=A1_SSSSM(k)+A2_SSSSM(k)*H+A3_SSSSM(k)/(H-1.05)
    808                 tau(i,k,la,aerindex) = mass_temp(i,k,aerindex)*1000.*zdp1(i,k)  &
    809                                    *tau_ae5wv_int*delt*fac
     764                tau(i,k,la,aerindex) = m_allaer(i,k,aerindex)/1.e6*dh(i,k)*tau_ae5wv_int*fac
    810765                tausum(i,la,aerindex)=tausum(i,la,aerindex)+tau(i,k,la,aerindex)
    811766              ENDDO
    … …  
    813768          ENDIF
    814769
     770          IF (spss.EQ.3) THEN !super coarse mode
     771            DO k=1, klev
     772              DO i=1, klon
     773                H=rh(i,k)/100.
     774                tau_ae5wv_int=A1_SSSSM(k)+A2_SSSSM(k)*H+A3_SSSSM(k)/(H-1.05)
     775                tau(i,k,la,aerindex) = m_allaer(i,k,aerindex)/1.e6*dh(i,k)*tau_ae5wv_int*fac
     776                tausum(i,la,aerindex)=tausum(i,la,aerindex)+tau(i,k,la,aerindex)
     777              ENDDO
     778            ENDDO
     779          ENDIF
     780
    815781        ELSE
    816           DO k=1, KLEV
    817 !CDIR ON_ADB(alpha_aers_5wv)
    818             DO i=1, KLON
     782          DO k=1, klev
     783            DO i=1, klon
    819784              tau_ae5wv_int = alpha_aers_5wv(RH_num(i,k),la,spsol)+DELTA(i,k)* &
    820785                             (alpha_aers_5wv(RH_num(i,k)+1,la,spsol) - &
    821786                              alpha_aers_5wv(RH_num(i,k),la,spsol))
    822 
    823               tau(i,k,la,aerindex) = mass_temp(i,k,aerindex)*1000.*zdp1(i,k)   &
    824                                  *tau_ae5wv_int*delt*fac
     787              tau(i,k,la,aerindex) = m_allaer(i,k,aerindex)/1.e6*dh(i,k)*tau_ae5wv_int*fac
    825788              tausum(i,la,aerindex)=tausum(i,la,aerindex)+tau(i,k,la,aerindex)
    826789            ENDDO
    … …  
    829792
    830793      ELSE                                                  ! For insoluble aerosol
    831         DO k=1, KLEV
    832 !CDIR ON_ADB(alpha_aeri_5wv)
    833           DO i=1, KLON
     794
     795        DO k=1, klev
     796          DO i=1, klon
    834797            tau_ae5wv_int = alpha_aeri_5wv(la,spinsol)
    835             tau(i,k,la,aerindex) = mass_temp(i,k,aerindex)*1000.*zdp1(i,k)* &
    836                                                 tau_ae5wv_int*delt*fac
    837             tausum(i,la,aerindex)= tausum(i,la,aerindex)+tau(i,k,la,aerindex)
     798            tau(i,k,la,aerindex) = m_allaer(i,k,aerindex)/1.e6*dh(i,k)*tau_ae5wv_int*fac
     799            tausum(i,la,aerindex)=tausum(i,la,aerindex)+tau(i,k,la,aerindex)
    838800          ENDDO
    839801        ENDDO
     802
    840803      ENDIF
     804
    841805    ENDDO   ! boucle sur les longueurs d'onde
    842806  ENDDO     ! Boucle  sur les masses de traceurs
    … …  
    845809    IF (.NOT.used_tau(m)) tau(:,:,:,m)=0.
    846810  ENDDO 
    847 !
    848 !
    849 !  taue670(:) = SUM(tausum(:,la670,:),dim=2)
    850 !  taue865(:) = SUM(tausum(:,la865,:),dim=2)
    851 !
    852 !  DO i=1, klon
    853 !    ai(i)=-LOG(MAX(taue670(i),0.0001)/ &
    854 !       MAX(taue865(i),0.0001))/LOG(670./865.)
    855 !  ENDDO
    856 
    857   DO i=1, klon
    858      od550aer(i)=0.
    859      DO m=1,naero_spc
    860         od550aer(i)=od550aer(i)+tausum(i,2,m)
    861      ENDDO
    862   ENDDO
    863 
    864   DO i=1, klon
    865      od865aer(i)=0.
    866      DO m=1,naero_spc
    867         od865aer(i)=od865aer(i)+tausum(i,5,m)
    868      ENDDO
    869   ENDDO
    870 
    871   DO i=1, klon
    872      DO k=1, KLEV
    873         ec550aer(i,k)=0.
    874         DO m=1,naero_spc
    875            ec550aer(i,k)=ec550aer(i,k)+tau(i,k,2,m)/dh(i,k)
    876         ENDDO
    877      ENDDO
    878   ENDDO
     811
     812!--AOD calculations for diagnostics
     813  od550aer(:)=SUM(tausum(:,la550,:),dim=2)
     814  od865aer(:)=SUM(tausum(:,la865,:),dim=2)
     815
     816!--extinction coefficient for diagnostic
     817  ec550aer(:,:)=SUM(tau(:,:,la550,:),dim=3)/dh(:,:)
    879818 
    880    od550lt1aer(:)=tausum(:,2,id_ASSO4M_phy)+tausum(:,2,id_ASBCM_phy)+tausum(:,2,id_AIBCM_phy)+  &
    881                   tausum(:,2,id_ASPOMM_phy)+tausum(:,2,id_AIPOMM_phy)+tausum(:,2,id_ASSSM_phy)+ &
    882                   0.03*tausum(:,2,id_CSSSM_phy)+0.4*tausum(:,2,id_CIDUSTM_phy)
     819!--acc mode AOD calculation for diagnostic
     820  od550lt1aer(:)=tausum(:,la550,id_ASSO4M_phy)+tausum(:,la550,id_ASBCM_phy)+tausum(:,la550,id_AIBCM_phy)+  &
     821                 tausum(:,la550,id_ASPOMM_phy)+tausum(:,la550,id_AIPOMM_phy)+tausum(:,la550,id_ASSSM_phy)+ &
     822                 0.03*tausum(:,la550,id_CSSSM_phy)+0.4*tausum(:,la550,id_CIDUSTM_phy)
    883823
    884824  DEALLOCATE(aerosol_name)

  • LMDZ5/trunk/libf/phylmd/rrtm/aeropt_5wv_rrtm.F90

    r2550 r2634  
    9494  REAL :: od670aer(klon)  ! epaisseur optique aerosol extinction 670 nm
    9595  REAL :: fac
    96   REAL :: zdp1(klon,klev)
    9796  INTEGER, ALLOCATABLE, DIMENSION(:)  :: aerosol_name
    9897  INTEGER :: nb_aer, itau
    9998  LOGICAL :: ok_itau
    10099 
    101   REAL :: dh(KLON,KLEV)
     100  REAL :: zdh(klon,klev)
    102101 
    103102   ! Soluble components 1- BC soluble; 2- POM soluble; 3- SO4 acc.; 4- SO4 coarse; 5 seasalt super-coarse; 6 seasalt coarse; 7 seasalt acc.
    … …  
    106105  REAL :: alpha_aeri_5wv(las,naero_insoluble)         ! Ext. coeff. ** m2/g
    107106
    108   REAL, DIMENSION(klon,klev,naero_tot) :: mass_temp
    109  
    110107  !
    111108  ! Proprietes optiques
    112109  !
    113110  REAL :: fact_RH(nbre_RH)
    114   LOGICAL :: used_tau(naero_tot)
    115111  INTEGER :: n
    116112
    … …  
    209205  ai(:) = 0.
    210206  tausum(:,:,:) = 0.
     207  tau(:,:,:,:)=0.
    211208
    212209  DO k=1, klev
    213210    DO i=1, klon
    214211      zrho=pplay(i,k)/t_seri(i,k)/RD                  ! kg/m3
    215       dh(i,k)=pdel(i,k)/(RG*zrho)
    216 !CDIR UNROLL=naero_spc
    217       mass_temp(i,k,:) = m_allaer(i,k,:) / zrho / 1.e+9
    218       zdp1(i,k)=pdel(i,k)/(RG*delt)     ! air mass auxiliary  variable --> zdp1 [kg/(m^2 *s)]
     212      zdh(i,k)=pdel(i,k)/(RG*zrho)                    ! m
    219213    ENDDO
    220214  ENDDO
    … …  
    270264  !     
    271265
    272 !CDIR ON_ADB(RH_tab)
    273 !CDIR ON_ADB(fact_RH)
    274 !CDIR NOVECTOR
    275266  DO n=1,nbre_RH-1
    276267    fact_RH(n)=1./(RH_tab(n+1)-RH_tab(n))
    277268  ENDDO
    278269   
    279   DO k=1, KLEV
    280 !CDIR ON_ADB(RH_tab)
    281 !CDIR ON_ADB(fact_RH)
    282     DO i=1, KLON
     270  DO k=1, klev
     271    DO i=1, klon
    283272      rh(i,k)=MIN(RHcl(i,k)*100.,RH_MAX)
    284273      RH_num(i,k) = INT( rh(i,k)/10. + 1.)
    … …  
    289278  ENDDO
    290279
    291 !CDIR SHORTLOOP 
    292   used_tau(:)=.FALSE.
    293    
    294280  DO m=1,nb_aer   ! tau is only computed for each mass   
    295281    fac=1.0
    … …  
    330316    ENDIF
    331317
    332     IF (soluble) then
    333       used_tau(spsol)=.TRUE.
    334     ELSE
    335       used_tau(naero_soluble+spinsol)=.TRUE.
    336     ENDIF
    337 
    338318    aerindex=aerosol_name(m)
    339319
    340320    DO la=1,las
    341321
     322    !--only 550, 670 and 865 nm are used
     323    IF (la.NE.la550.AND.la.NE.la670.AND.la.NE.la865) CYCLE
     324
    342325      IF (soluble) THEN            ! For soluble aerosol
    343326
    344           DO k=1, KLEV
    345             DO i=1, KLON
     327          DO k=1, klev
     328            DO i=1, klon
    346329              tau_ae5wv_int = alpha_aers_5wv(RH_num(i,k),la,spsol)+DELTA(i,k)* &
    347330                             (alpha_aers_5wv(RH_num(i,k)+1,la,spsol) - &
    348331                              alpha_aers_5wv(RH_num(i,k),la,spsol))
    349               tau(i,k,la,aerindex) = mass_temp(i,k,aerindex)*1000.*zdp1(i,k)*   &
    350                                     tau_ae5wv_int*delt*fac
     332              tau(i,k,la,aerindex) = m_allaer(i,k,aerindex)/1.e6*zdh(i,k)*tau_ae5wv_int*fac
    351333              tausum(i,la,aerindex)=tausum(i,la,aerindex)+tau(i,k,la,aerindex)
    352334            ENDDO
    … …  
    355337      ELSE                         ! For insoluble aerosol
    356338
    357         DO k=1, KLEV
    358           DO i=1, KLON
     339        DO k=1, klev
     340          DO i=1, klon
    359341            tau_ae5wv_int = alpha_aeri_5wv(la,spinsol)
    360             tau(i,k,la,aerindex) = mass_temp(i,k,aerindex)*1000.*zdp1(i,k)* &
    361                                    tau_ae5wv_int*delt*fac
     342            tau(i,k,la,aerindex) = m_allaer(i,k,aerindex)/1.e6*zdh(i,k)*tau_ae5wv_int*fac
    362343            tausum(i,la,aerindex)= tausum(i,la,aerindex)+tau(i,k,la,aerindex)
    363344          ENDDO
    … …  
    369350  ENDDO     ! Boucle sur les masses de traceurs
    370351
    371   DO m=1,naero_tot
    372     IF (.NOT.used_tau(m)) tau(:,:,:,m)=0.
    373   ENDDO 
    374 
    375   DO i=1, klon
    376      od550aer(i)=0.
    377      DO m=1,naero_tot
    378         od550aer(i)=od550aer(i)+tausum(i,la550,m)
    379      END DO
    380   END DO
    381 
    382   DO i=1, klon
    383      od670aer(i)=0.
    384      DO m=1,naero_tot
    385         od670aer(i)=od670aer(i)+tausum(i,la670,m)
    386      END DO
    387   END DO
    388 
    389   DO i=1, klon
    390      od865aer(i)=0.
    391      DO m=1,naero_tot
    392         od865aer(i)=od865aer(i)+tausum(i,la865,m)
    393      END DO
    394   END DO
    395 
    396   DO i=1, klon
    397      DO k=1, KLEV
    398         ec550aer(i,k)=0.
    399         DO m=1,naero_tot
    400            ec550aer(i,k)=ec550aer(i,k)+tau(i,k,la550,m)/dh(i,k)
    401         END DO
    402      END DO
    403   END DO
    404  
    405   DO i=1, klon
    406     ai(i)=-LOG(MAX(od670aer(i),1.e-8)/MAX(od865aer(i),1.e-8))/LOG(670./865.)
    407   ENDDO
     352!--AOD calculations for diagnostics
     353  od550aer(:)=SUM(tausum(:,la550,:),dim=2)
     354  od670aer(:)=SUM(tausum(:,la670,:),dim=2)
     355  od865aer(:)=SUM(tausum(:,la865,:),dim=2)
     356
     357!--extinction coefficient for diagnostic
     358  ec550aer(:,:)=SUM(tau(:,:,la550,:),dim=3)/zdh(:,:)
     359
     360!--aerosol index
     361  ai(:)=-LOG(MAX(od670aer(:),1.e-8)/MAX(od865aer(:),1.e-8))/LOG(670./865.)
    408362
    409363  od550lt1aer(:)=tausum(:,la550,id_ASSO4M_phy)+tausum(:,la550,id_ASBCM_phy) +tausum(:,la550,id_AIBCM_phy)+ &

  • LMDZ5/trunk/libf/phylmd/rrtm/aeropt_6bands_rrtm.F90

    r2311 r2634  
    3636  !
    3737  ! Output arguments:
    38   ! 1= total aerosols
    39   ! 2= natural aerosols
     38  ! 2= total aerosols
     39  ! 1= natural aerosols
    4040  !
    4141  REAL, DIMENSION(klon,klev,2,nbands_sw_rrtm), INTENT(out) :: tau_allaer ! epaisseur optique aerosol
    … …  
    6666  REAL :: Fact_RH(nbre_RH)
    6767  REAL :: fac
    68   REAL :: zdp1(klon,klev)
     68  REAL :: zdh(klon,klev)
    6969  INTEGER, ALLOCATABLE, DIMENSION(:)   :: aerosol_name
    7070  INTEGER :: nb_aer
    7171
    72   REAL, DIMENSION(klon,klev,naero_tot) :: mass_temp
    73   REAL, DIMENSION(klon,klev,naero_tot) :: mass_temp_pi
    7472  REAL, DIMENSION(klon,klev,naero_tot,nbands_sw_rrtm) ::  tau_ae
    7573  REAL, DIMENSION(klon,klev,naero_tot,nbands_sw_rrtm) ::  tau_ae_pi
    7674  REAL, DIMENSION(klon,klev,naero_tot,nbands_sw_rrtm) ::  piz_ae
    7775  REAL, DIMENSION(klon,klev,naero_tot,nbands_sw_rrtm) ::  cg_ae
    78 
    79 
    8076  !
    8177  ! Proprietes optiques
    … …  
    8985
    9086  INTEGER :: id
    91   LOGICAL :: used_aer(naero_tot)
    9287  REAL :: tmp_var, tmp_var_pi
    9388
    … …  
    277272  spsol = 0
    278273  spinsol = 0
     274
    279275  IF (NSW.NE.nbands_sw_rrtm) THEN
    280276     print *,'Erreur NSW doit etre egal a 6 pour cette routine'
    … …  
    282278  ENDIF
    283279
    284   DO k=1, klev
    285     DO i=1, klon
    286 !CDIR UNROLL=naero_tot
    287       mass_temp(i,k,:) = m_allaer(i,k,:) / zrho(i,k) / 1.e+9  !--kg/kg
    288 !CDIR UNROLL=naero_tot
    289       mass_temp_pi(i,k,:) = m_allaer_pi(i,k,:) / zrho(i,k) / 1.e+9
    290       zdp1(i,k)=pdel(i,k)/(RG*delt)      ! air mass auxiliary  variable --> zdp1 [kg/(m^2 *s)]
    291     ENDDO
    292   ENDDO
     280  zdh(:,:)=pdel(:,:)/(RG*zrho(:,:))      ! m
    293281
    294282  IF (flag_aerosol .EQ. 1) THEN
    … …  
    339327  !      compute optical_thickness_at_gridpoint_per_species
    340328
    341 !!CDIR ON_ADB(RH_tab)
    342 !CDIR ON_ADB(fact_RH)
    343 !CDIR SHORTLOOP
    344329  DO n=1,nbre_RH-1
    345330    fact_RH(n)=1./(RH_tab(n+1)-RH_tab(n))
    346331  ENDDO
    347332   
    348   DO k=1, KLEV
    349 !CDIR ON_ADB(fact_RH)
    350     DO i=1, KLON
     333  DO k=1, klev
     334    DO i=1, klon
    351335      rh(i,k)=MIN(RHcl(i,k)*100.,RH_MAX)
    352336      RH_num(i,k) = INT(rh(i,k)/10. + 1.)
    … …  
    357341  ENDDO
    358342
    359   used_aer(:)=.FALSE.
     343  tau_ae(:,:,:,:)=0.
     344  tau_ae_pi(:,:,:,:)=0.
     345  piz_ae(:,:,:,:)=0.
     346  cg_ae(:,:,:,:)=0.
    360347   
    361348  DO m=1,nb_aer   ! tau is only computed for each mass
    … …  
    398385
    399386    id=aerosol_name(m)
    400     used_aer(id)=.TRUE.
    401387
    402388    IF (soluble) THEN
    403389
    404        DO k=1, KLEV
    405          DO i=1, KLON
    406            tmp_var=mass_temp(i,k,spsol)*1000.*zdp1(i,k)*delt*fac
    407            tmp_var_pi=mass_temp_pi(i,k,spsol)*1000.*zdp1(i,k)*delt*fac
     390       DO k=1, klev
     391         DO i=1, klon
     392           tmp_var=m_allaer(i,k,spsol)/1.e6*zdh(i,k)*fac
     393           tmp_var_pi=m_allaer_pi(i,k,spsol)/1.e6*zdh(i,k)*fac
    408394
    409395           DO inu=1,NSW
    … …  
    432418     ELSE    ! For all aerosol insoluble components
    433419
    434        DO k=1, KLEV
    435          DO i=1, KLON
    436            tmp_var=mass_temp(i,k,naero_soluble+ spinsol)*1000.*zdp1(i,k)*delt*fac
    437            tmp_var_pi=mass_temp_pi(i,k,naero_soluble+spinsol)*1000.*zdp1(i,k)*delt*fac
     420       DO k=1, klev
     421         DO i=1, klon
     422           tmp_var=m_allaer(i,k,naero_soluble+spinsol)/1.e6*zdh(i,k)*fac
     423           tmp_var_pi=m_allaer_pi(i,k,naero_soluble+spinsol)/1.e6*zdh(i,k)*fac
    438424
    439425           DO inu=1,NSW
    … …  
    454440  ENDDO  ! nb_aer 
    455441
    456   DO m=1,naero_tot
    457     IF (.NOT. used_aer(m)) THEN
    458       tau_ae(:,:,m,:)=0.
    459       tau_ae_pi(:,:,m,:)=0.
    460       piz_ae(:,:,m,:)=0.
    461       cg_ae(:,:,m,:)=0.
    462     ENDIF
    463   ENDDO
    464 
    465442  DO inu=1, NSW
    466      DO k=1, KLEV
    467        DO i=1, KLON
     443     DO k=1, klev
     444       DO i=1, klon
    468445!--anthropogenic aerosol
    469446         tau_allaer(i,k,2,inu)=tau_ae(i,k,id_ASSO4M_phy,inu)+tau_ae(i,k,id_CSSO4M_phy,inu)+ &
    … …  
    536513    ENDDO
    537514   
    538 !--???????
    539   inu=1
    540   DO i=1, KLON
    541      absvisaer(i)=SUM((1-piz_allaer(i,:,:,inu))*tau_allaer(i,:,:,inu))
    542   END DO       
     515!--waveband 2 and all aerosol
     516  inu=2
     517  DO i=1, klon
     518     absvisaer(i)=SUM((1-piz_allaer(i,:,2,inu))*tau_allaer(i,:,2,inu))
     519  ENDDO
    543520
    544521  DEALLOCATE(aerosol_name)
Note: See TracChangeset for help on using the changeset viewer.