Changeset 1237 for LMDZ4


Ignore:
Timestamp:
Sep 3, 2009, 2:03:33 PM (15 years ago)
Author:
Laurent Fairhead
Message:

Des modifications sur la lecture des aerosols par Michael
Correction du test sur le jour de lecture des aerosols qui ne marchait
pas avec le nouveau calendrier (a revoir?)
Menage sur quelques prints
SD/MAF

Location:
LMDZ4/branches/LMDZ4-dev/libf/phylmd
Files:
8 edited

Legend:

Unmodified
Added
Removed

  • LMDZ4/branches/LMDZ4-dev/libf/phylmd/aeropt_2bands.F90

    r1221 r1237  
     1!
     2! $Id$
     3!
    14SUBROUTINE AEROPT_2BANDS( &
    25     pdel, m_allaer, delt, RHcl, &
    … …  
    4447  REAL, DIMENSION(klon,klev,naero_tot,nbands) ::  cg_ae
    4548  LOGICAL ::  soluble
    46   INTEGER :: i, k, inu, m, mrfspecies
     49  INTEGER :: i, k, ierr, inu, m, mrfspecies
    4750  INTEGER :: spsol, spinsol, spss
    4851  INTEGER :: RH_num
    … …  
    5053
    5154  INTEGER, PARAMETER :: nbre_RH=12
    52   INTEGER, PARAMETER :: nbsol_compaer=3    ! 1- BC soluble; 2- POM soluble; 3- SO4.
    53   INTEGER, PARAMETER :: nbinsol_compaer=3  ! 1- Dust; 2- BC insoluble; 3- POM insoluble
     55  INTEGER, PARAMETER :: naero_soluble=7    ! 1- BC soluble; 2- POM soluble; 3- SO4.
     56  INTEGER, PARAMETER :: naero_insoluble=3  ! 1- Dust; 2- BC insoluble; 3- POM insoluble
    5457  LOGICAL, SAVE :: firstcall=.TRUE.
    5558
    … …  
    7073
    7174! Coefficient optiques interpole sur le nombre de niveau du modele
    72   REAL, DIMENSION(klev) :: A1_ASSSM_b1, A2_ASSSM_b1, A3_ASSSM_b1,&
     75  REAL, ALLOCATABLE, DIMENSION(:), SAVE :: &
     76          A1_ASSSM_b1, A2_ASSSM_b1, A3_ASSSM_b1,&
    7377          B1_ASSSM_b1, B2_ASSSM_b1, C1_ASSSM_b1, C2_ASSSM_b1,&
    7478          A1_CSSSM_b1, A2_CSSSM_b1, A3_CSSSM_b1,&
    … …  
    101105  ! Proprietes optiques
    102106  !
    103   REAL:: alpha_aers_2bands(nbre_RH,nbands,nbsol_compaer)   !--unit m2/g SO4
    104   REAL:: alpha_aeri_2bands(nbands,nbinsol_compaer)
    105   REAL:: cg_aers_2bands(nbre_RH,nbands,nbsol_compaer)      !--unit
    106   REAL:: cg_aeri_2bands(nbands,nbinsol_compaer)
    107   REAL:: piz_aers_2bands(nbre_RH,nbands,nbsol_compaer)     !-- unit
    108   REAL:: piz_aeri_2bands(nbands,nbinsol_compaer)           !-- unit
     107  REAL:: alpha_aers_2bands(nbre_RH,nbands,naero_soluble)   !--unit m2/g SO4
     108  REAL:: alpha_aeri_2bands(nbands,naero_insoluble)
     109  REAL:: cg_aers_2bands(nbre_RH,nbands,naero_soluble)      !--unit
     110  REAL:: cg_aeri_2bands(nbands,naero_insoluble)
     111  REAL:: piz_aers_2bands(nbre_RH,nbands,naero_soluble)     !-- unit
     112  REAL:: piz_aeri_2bands(nbands,naero_insoluble)           !-- unit
    109113
    110114
    … …  
    354358       7.556,8.613,10.687,12.265,16.32,21.692, &
    355359       1.107,1.239,1.381,1.490,1.635,1.8030,   &
    356        2.071,2.407,3.126,3.940,5.539,7.921/
    357 
     360       2.071,2.407,3.126,3.940,5.539,7.921,    &
     361                                ! sulfate coarse
     362       4.681,5.062,5.460,5.798,6.224,6.733,    &
     363       7.556,8.613,10.687,12.265,16.32,21.692, &
     364       1.107,1.239,1.381,1.490,1.635,1.8030,   &
     365       2.071,2.407,3.126,3.940,5.539,7.921,    &
     366                                ! seasalt Super Coarse Soluble (SS)
     367       0.5090,0.6554,0.7129,0.7767,0.8529,1.2728, &
     368       1.3820,1.5792,1.9173,2.2002,2.7173,4.1487, &
     369       0.5167,0.6613,0.7221,0.7868,0.8622,1.3027, &
     370       1.4227,1.6317,1.9887,2.2883,2.8356,4.3453, &
     371                                ! seasalt  Coarse Soluble (CS)
     372       0.5090,0.6554,0.7129,0.7767,0.8529,1.2728, &
     373       1.3820,1.5792,1.9173,2.2002,2.7173,4.1487, &
     374       0.5167,0.6613,0.7221,0.7868,0.8622,1.3027, &
     375       1.4227,1.6317,1.9887,2.2883,2.8356,4.3453, &
     376                                ! seasalt  Accumulation Soluble (AS)
     377       4.125, 4.674, 5.005, 5.434, 5.985, 10.006, &
     378       11.175,13.376,17.264,20.540,26.604, 42.349,&
     379       4.187, 3.939, 3.919, 3.937, 3.995,  5.078, &
     380       5.511, 6.434, 8.317,10.152,14.024, 26.537/
    358381
    359382  DATA alpha_aeri_2bands/  &
    … …  
    364387       ! pom insoluble
    365388       3.741,0.606/
    366 
    367389
    368390  DATA cg_aers_2bands/ &
    … …  
    381403       .719, .733, .752, .760, .773, .786, &
    382404       .544, .555, .565, .573, .583, .593, &
    383        .610, .628, .655, .666, .692, .719/
    384 
     405       .610, .628, .655, .666, .692, .719, &
     406                                ! sulfate coarse
     407       .658, .669, .680, .688, .698, .707, &
     408       .719, .733, .752, .760, .773, .786, &
     409       .544, .555, .565, .573, .583, .593, &
     410       .610, .628, .655, .666, .692, .719, &
     411                                ! seasalt Super Coarse soluble (SS)
     412       .727, .747, .755, .761, .770, .788, &
     413       .792, .799, .805, .809, .815, .826, &
     414       .717, .738, .745, .752, .761, .779, &
     415       .781, .786, .793, .797, .803, .813, &
     416                                ! seasalt Coarse soluble (CS)
     417       .727, .747, .755, .761, .770, .788, &
     418       .792, .799, .805, .809, .815, .826, &
     419       .717, .738, .745, .752, .761, .779, &
     420       .781, .786, .793, .797, .803, .813, &
     421                                ! Sesalt Accumulation Soluble (AS)
     422       .727, .741, .748, .754, .761, .782, &
     423       .787, .792, .797, .799, .801, .799, &
     424       .606, .645, .658, .669, .681, .726, &
     425       .734, .746, .761, .770, .782, .798/
    385426
    386427  DATA cg_aeri_2bands/ &
    … …  
    407448       1.000,1.000,1.000,1.000,1.000,1.000, &
    408449       .992, .988, .988, .987, .986, .985,  &
    409        .985, .985, .984, .984, .984, .984 /
    410 
     450       .985, .985, .984, .984, .984, .984,  &
     451                                ! sulfate coarse
     452       1.000,1.000,1.000,1.000,1.000,1.000, &
     453       1.000,1.000,1.000,1.000,1.000,1.000, &
     454       .992, .988, .988, .987, .986, .985,  &
     455       .985, .985, .984, .984, .984, .984,  &
     456                                ! seasalt Super Coarse Soluble (SS)
     457       1.000,1.000,1.000,1.000,1.000,1.000, &
     458       1.000,1.000,1.000,1.000,1.000,1.000, &
     459       0.992,0.989,0.987,0.986,0.986,0.980, &
     460       0.980,0.978,0.976,0.976,0.974,0.971, &
     461                                ! seasalt Coarse soluble (CS)
     462       1.000,1.000,1.000,1.000,1.000,1.000, &
     463       1.000,1.000,1.000,1.000,1.000,1.000, &
     464       0.992,0.989,0.987,0.986,0.986,0.980, &
     465       0.980,0.978,0.976,0.976,0.974,0.971, &
     466                                ! seasalt Accumulation Soluble (AS)
     467       1.000, 1.000, 1.000, 1.000, 1.000, 1.000, &
     468       1.000, 1.000, 1.000, 1.000, 1.000, 1.000, &
     469       0.970, 0.975, 0.976, 0.977, 0.978, 0.982, &
     470       0.982, 0.983, 0.984, 0.984, 0.985, 0.985/
    411471
    412472  DATA piz_aeri_2bands/ &
    … …  
    418478       .966, .859/
    419479
    420 
    421480! Interpolation des coefficients optiques de 19 niveaux vers le nombre des niveaux du model
    422481  IF (firstcall) THEN
    423482     firstcall=.FALSE.
     483     
     484     IF (.NOT. ALLOCATED(A1_ASSSM_b1)) THEN
     485        ALLOCATE(A1_ASSSM_b1(klev),A2_ASSSM_b1(klev), A3_ASSSM_b1(klev),&
     486          B1_ASSSM_b1(klev), B2_ASSSM_b1(klev), C1_ASSSM_b1(klev), C2_ASSSM_b1(klev),&
     487          A1_CSSSM_b1(klev), A2_CSSSM_b1(klev), A3_CSSSM_b1(klev),&
     488          B1_CSSSM_b1(klev), B2_CSSSM_b1(klev), C1_CSSSM_b1(klev), C2_CSSSM_b1(klev),&
     489          A1_SSSSM_b1(klev), A2_SSSSM_b1(klev), A3_SSSSM_b1(klev),&
     490          B1_SSSSM_b1(klev), B2_SSSSM_b1(klev), C1_SSSSM_b1(klev), C2_SSSSM_b1(klev),&
     491          A1_ASSSM_b2(klev), A2_ASSSM_b2(klev), A3_ASSSM_b2(klev),&
     492          B1_ASSSM_b2(klev), B2_ASSSM_b2(klev), C1_ASSSM_b2(klev), C2_ASSSM_b2(klev),&
     493          A1_CSSSM_b2(klev), A2_CSSSM_b2(klev), A3_CSSSM_b2(klev),&
     494          B1_CSSSM_b2(klev), B2_CSSSM_b2(klev), C1_CSSSM_b2(klev), C2_CSSSM_b2(klev),&
     495          A1_SSSSM_b2(klev), A2_SSSSM_b2(klev), A3_SSSSM_b2(klev),&
     496          B1_SSSSM_b2(klev), B2_SSSSM_b2(klev), C1_SSSSM_b2(klev), C2_SSSSM_b2(klev), stat=ierr)
     497        IF (ierr /= 0) CALL abort_gcm('aeropt_2bands', 'pb in allocation 1',1)
     498     END IF
    424499     
    425500! bande 1
    … …  
    478553  DO k=1, klev
    479554     DO i=1, klon
    480         IF (t_seri(i,k).EQ.0.) THEN
    481            WRITE(lunout,*) 't_seri(i,k)=0 for i=',i,'k=',k
    482            CALL abort_gcm('aeropt_2bands','t_seri=0',1)
    483         END IF
    484         IF (pplay(i,k).EQ.0.) THEN
    485            WRITE(lunout,*) 'pplay(i,k)=0 for i=',i,'k=',k
    486            CALL abort_gcm('aeropt_2bands','pplay=0',1)
    487         END IF
    488 
    489         zrho=pplay(i,k)/t_seri(i,k)/RD                  ! kg/m3
     555!        IF (t_seri(i,k).EQ.0.) THEN
     556!           WRITE(lunout,*) 't_seri(i,k)=0 for i=',i,'k=',k
     557!           CALL abort_gcm('aeropt_2bands','t_seri=0',1)
     558!        END IF
     559!        IF (pplay(i,k).EQ.0.) THEN
     560!           WRITE(lunout,*) 'pplay(i,k)=0 for i=',i,'k=',k
     561!           CALL abort_gcm('aeropt_2bands','pplay=0',1)
     562!        END IF
     563        zrho=pplay(i,k)/t_seri(i,k)/RD                    ! kg/m3
    490564        mass_temp(i,k,:) = m_allaer(i,k,:) / zrho / 1.e+9
    491 
    492565     ENDDO
    493566  ENDDO
    … …  
    563636         spss=0
    564637     ELSEIF (aerosol_name(m).EQ.id_ASPOMM) THEN
    565          soluble=.TRUE.
    566          spsol=2
    567          spss=0
    568      ELSEIF ((aerosol_name(m).EQ.id_ASSO4M) .OR.  (aerosol_name(m).EQ.id_CSSO4M)) THEN
    569          soluble=.TRUE.
    570          spsol=3
    571          spss=0
    572          fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for OD
     638        soluble=.TRUE.
     639        spsol=2
     640        spss=0
     641     ELSEIF (aerosol_name(m).EQ.id_ASSO4M) THEN
     642        soluble=.TRUE.
     643        spsol=3
     644        spss=0
     645        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for OD
     646     ELSEIF  (aerosol_name(m).EQ.id_CSSO4M) THEN
     647        soluble=.TRUE.
     648        spsol=4
     649        spss=0
     650        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for OD
    573651     ELSEIF (aerosol_name(m).EQ.id_SSSSM) THEN
    574652         soluble=.TRUE.
    … …  
    604682     cg_ae2b_int(:,:,:)=0.
    605683
    606      DO k=1, KLEV
    607         DO i=1, KLON
    608 
    609            rh=MIN(RHcl(i,k)*100.,RH_MAX)
    610            RH_num = INT( rh/10. + 1.)
    611 
    612            IF (rh.GT.85.) RH_num=10
    613            IF (rh.GT.90.) RH_num=11
    614            DELTA=(rh-RH_tab(RH_num))/(RH_tab(RH_num+1)-RH_tab(RH_num))
    615 
    616            DO inu=1,nbands
     684     DO inu=1,nbands
     685        DO k=1, KLEV
     686           DO i=1, KLON
     687
     688              rh=MIN(RHcl(i,k)*100.,RH_MAX)
     689              RH_num = INT( rh/10. + 1.)
     690
     691              IF (rh.GT.85.) RH_num=10
     692              IF (rh.GT.90.) RH_num=11
     693              DELTA=(rh-RH_tab(RH_num))/(RH_tab(RH_num+1)-RH_tab(RH_num))
     694
     695!           DO inu=1,nbands
    617696              IF (soluble) THEN
    618697
     698              ! First optical parameters are computed for seasalt
    619699                  IF (spss.NE.0) THEN
    620700                      H=rh/100
    … …  
    664744                         DELTA* (cg_aers_2bands(RH_num+1,inu,spsol) - &
    665745                         cg_aers_2bands(RH_num,inu,spsol))
    666                      
    667746                  ENDIF
    668747
    669                   tau_ae(i,k,aerosol_name(m),inu) = &
    670                      mass_temp(i,k,spsol)*1000.*zdp1(i,k)*delt*tau_ae2b_int(i,k,inu)*fac
    671 
    672               ELSE
     748                 tau_ae(i,k,aerosol_name(m),inu) = &
     749                      mass_temp(i,k,spsol)*1000.*zdp1(i,k)*delt* &
     750                                         tau_ae2b_int(i,k,inu)*fac
     751
     752              ELSE                                                    ! For all aerosol insoluble components
    673753                 tau_ae2b_int(i,k,inu) = alpha_aeri_2bands(inu,spinsol)
    674754                 piz_ae2b_int(i,k,inu) = piz_aeri_2bands(inu,spinsol)
    … …  
    676756
    677757                 tau_ae(i,k,aerosol_name(m),inu) = &
    678                     mass_temp(i,k,7+ spinsol)*1000.*zdp1(i,k)*delt*tau_ae2b_int(i,k,inu)*fac
     758                      mass_temp(i,k,naero_soluble+ spinsol)*1000.*zdp1(i,k)* &
     759                                               delt*tau_ae2b_int(i,k,inu)*fac
    679760              ENDIF
    680761
    … …  
    682763
    683764              cg_ae(i,k,aerosol_name(m),inu)= cg_ae2b_int(i,k,inu)
    684 
    685765
    686766           ENDDO    ! nbands : boucle sur les bandes spectrale
    … …  
    726806                      (tau_allaer(i,k,mrfspecies,inu)*piz_allaer(i,k,mrfspecies,inu))
    727807
    728 
    729808              ELSEIF (mrfspecies .EQ. 3) THEN             ! = natural aerosol NAT
    730809                 tau_allaer(i,k,mrfspecies,inu)=tau_ae(i,k,id_ASSO4M,inu)*fractnat_allaer(i,id_ASSO4M)+ &
    … …  
    767846                      (tau_allaer(i,k,mrfspecies,inu)*piz_allaer(i,k,mrfspecies,inu))
    768847
    769 
    770848              ELSEIF (mrfspecies .EQ. 4) THEN             ! = BC
    771849                 tau_allaer(i,k,mrfspecies,inu)=tau_ae(i,k,id_ASBCM,inu)+tau_ae(i,k,id_AIBCM,inu)
    … …  
    778856                      +tau_ae(i,k,id_AIBCM,inu)*piz_ae(i,k,id_AIBCM,inu)*cg_ae(i,k,id_AIBCM,inu))/ &
    779857                      (tau_allaer(i,k,mrfspecies,inu)*piz_allaer(i,k,mrfspecies,inu))
     858
    780859              ELSEIF (mrfspecies .EQ. 5) THEN             ! = SO4
    781860                 tau_allaer(i,k,mrfspecies,inu)=tau_ae(i,k,id_ASSO4M,inu)+tau_ae(i,k,id_CSSO4M,inu)
    … …  
    799878                      +tau_ae(i,k,id_AIPOMM,inu)*piz_ae(i,k,id_AIPOMM,inu)*cg_ae(i,k,id_AIPOMM,inu))/ &
    800879                      (tau_allaer(i,k,mrfspecies,inu)*piz_allaer(i,k,mrfspecies,inu))
     880
    801881              ELSEIF (mrfspecies .EQ. 7) THEN             ! = DUST
    802882                 tau_allaer(i,k,mrfspecies,inu)=tau_ae(i,k,id_CIDUSTM,inu)

  • LMDZ4/branches/LMDZ4-dev/libf/phylmd/aeropt_5wv.F90

    r1224 r1237  
    1 
     1!
     2! $Id$
     3!
    24
    35SUBROUTINE AEROPT_5WV(&
    … …  
    911  USE DIMPHY
    1012  USE aero_mod
    11 
     13  USE mod_phys_lmdz_para, ONLY : mpi_rank
    1214  !
    1315  !    Yves Balkanski le 12 avril 2006
    … …  
    7274  LOGICAL :: soluble
    7375 
    74   INTEGER :: i, k, m
     76  INTEGER :: i, k, ierr, m
    7577  INTEGER :: spsol, spinsol, spss, la
    7678  INTEGER :: RH_num
    … …  
    8183  INTEGER, PARAMETER :: la865 = 5
    8284  INTEGER, PARAMETER :: nbre_RH=12
    83   INTEGER, PARAMETER :: nbsol_compaer=3   !  1- BC soluble; 2- POM soluble; 3- SO4.
    84   INTEGER, PARAMETER :: nbinsol_compaer=3 !  1- Dust; 2- BC insoluble; 3- POM insoluble
     85  INTEGER, PARAMETER :: naero_soluble=7   !  1- BC soluble; 2- POM soluble; 3- SO4 acc.
     86                                          ! 4- SO4 coarse; 5 seasalt super-C; 6 seasalt coarse; 7 seasalt acc.
     87  INTEGER, PARAMETER :: naero_insoluble=3 !  1- Dust; 2- BC insoluble; 3- POM insoluble
    8588  INTEGER, PARAMETER :: nb_level = 19 ! number of vertical levels
    8689  LOGICAL, SAVE :: firstcall=.TRUE.
    … …  
    98101
    99102  ! Coefficient optiques interpole sur le nombre de niveau du modele
    100   REAL, DIMENSION(klev) :: A1_ASSSM, A2_ASSSM, A3_ASSSM,&
     103  REAL, ALLOCATABLE,  DIMENSION(:), SAVE :: &
     104          A1_ASSSM, A2_ASSSM, A3_ASSSM,&
    101105          B1_ASSSM, B2_ASSSM, C1_ASSSM, C2_ASSSM,&
    102106          A1_CSSSM, A2_CSSSM, A3_CSSSM,&
    … …  
    124128
    125129 
    126   REAL :: alpha_aers_5wv(nbre_RH,las,nbsol_compaer)   ! ext. coeff. Soluble comp. units *** m2/g
    127                                                       ! 1- Seasalt AS: 2- Sesalt CS; 3- BC; 4- POM; 5- SO4.
    128   REAL :: alpha_aeri_5wv(las,nbinsol_compaer)         ! ext. coeff. Insoluble comp. 1- Dust: 2- BC; 3- POM
    129   REAL :: cg_aers_5wv(nbre_RH,las,nbsol_compaer)      ! Asym. param. soluble comp.
    130                                                       ! 1- Seasalt AS: 2- Sesalt CS; 3- BC; 4- POM; 5- SO4.
    131   REAL :: cg_aeri_5wv(las,nbinsol_compaer)            ! Asym. param. insoluble comp. 1- Dust: 2- BC; 3- POM
    132   REAL :: piz_aers_5wv(nbre_RH,las,nbsol_compaer)   
    133                                                       ! 1- Seasalt AS: 2- Sesalt CS; 3- BC; 4- POM; 5- SO4.
    134   REAL :: piz_aeri_5wv(las,nbinsol_compaer)           ! Insoluble comp. 1- Dust: 2- BC; 3- POM
     130  REAL :: alpha_aers_5wv(nbre_RH,las,naero_soluble)   ! ext. coeff. Soluble comp. units *** m2/g
     131   !  1- BC soluble; 2- POM soluble; 3- SO4 acc.; 4- SO4 coarse; 5 seasalt super-C; 6 seasalt coarse; 7 seasalt acc.
     132  REAL :: alpha_aeri_5wv(las,naero_insoluble)         ! ext. coeff. Insoluble comp. 1- Dust: 2- BC; 3- POM
     133  REAL :: cg_aers_5wv(nbre_RH,las,naero_soluble)      ! Asym. param. soluble comp.
     134   !  1- BC soluble; 2- POM soluble; 3- SO4 acc.; 4- SO4 coarse; 5 seasalt super-C; 6 seasalt coarse; 7 seasalt acc.
     135  REAL :: cg_aeri_5wv(las,naero_insoluble)            ! Asym. param. insoluble comp. 1- Dust: 2- BC; 3- POM
     136  REAL :: piz_aers_5wv(nbre_RH,las,naero_soluble)   
     137   !  1- BC soluble; 2- POM soluble; 3- SO4 acc.; 4- SO4 coarse; 5 seasalt super-C; 6 seasalt coarse; 7 seasalt acc.
     138  REAL :: piz_aeri_5wv(las,naero_insoluble)           ! Insoluble comp. 1- Dust: 2- BC; 3- POM
    135139
    136140  REAL, DIMENSION(klon,klev,naero_spc) :: mass_temp
    … …  
    284288  !                                   le 12 AVRIL 2006
    285289  !
    286   DATA alpha_aers_5wv/ &
    287      ! bc soluble
    288      7.930,7.930,7.930,7.930,7.930,7.930,     &
    289      7.930,7.930,10.893,12.618,14.550,16.613, &
    290      7.658,7.658,7.658,7.658,7.658,7.658,     &
    291      7.658,7.658,10.351,11.879,13.642,15.510, &
    292      7.195,7.195,7.195,7.195,7.195,7.195,     &
    293      7.195,7.195,9.551,10.847,12.381,13.994,  &
    294      6.736,6.736,6.736,6.736,6.736,6.736,     &
    295      6.736,6.736,8.818,9.938,11.283,12.687,   &
    296      6.277,6.277,6.277,6.277,6.277,6.277,     &
    297      6.277,6.277,8.123,9.094,10.275,11.501,   &
    298      ! pom soluble
    299      6.676,6.676,6.676,6.676,6.710,6.934,   &
    300      7.141,7.569,8.034,8.529,9.456,10.511,  &
    301      5.109,5.109,5.109,5.109,5.189,5.535,   &
    302      5.960,6.852,8.008,9.712,12.897,19.676, &
    303      3.718,3.718,3.718,3.718,3.779,4.042,   &
    304      4.364,5.052,5.956,7.314,9.896,15.688,  &
    305      2.849,2.849,2.849,2.849,2.897,3.107,   &
    306      3.365,3.916,4.649,5.760,7.900,12.863,  &
    307      2.229,2.229,2.229,2.229,2.268,2.437,   &
    308      2.645,3.095,3.692,4.608,6.391,10.633,  &
    309      ! Sulfate
    310      5.751,6.215,6.690,7.024,7.599,8.195,      &
    311      9.156,10.355,12.660,14.823,18.908,24.508, &
    312      4.320,4.675,5.052,5.375,5.787,6.274,      &
    313      7.066,8.083,10.088,12.003,15.697,21.133,  &
    314      3.079,3.351,3.639,3.886,4.205,4.584,      &
    315      5.206,6.019,7.648,9.234,12.391,17.220,    &
    316      2.336,2.552,2.781,2.979,3.236,3.540,      &
    317      4.046,4.711,6.056,7.388,10.093,14.313,    &
    318      1.777,1.949,2.134,2.292,2.503,2.751,      &
    319      3.166,3.712,4.828,5.949,8.264,11.922/
     290 DATA alpha_aers_5wv/ &
     291                                ! bc soluble
     292       7.930,7.930,7.930,7.930,7.930,7.930,     &
     293       7.930,7.930,10.893,12.618,14.550,16.613, &
     294       7.658,7.658,7.658,7.658,7.658,7.658,     &
     295       7.658,7.658,10.351,11.879,13.642,15.510, &
     296       7.195,7.195,7.195,7.195,7.195,7.195,     &
     297       7.195,7.195,9.551,10.847,12.381,13.994,  &
     298       6.736,6.736,6.736,6.736,6.736,6.736,     &
     299       6.736,6.736,8.818,9.938,11.283,12.687,   &
     300       6.277,6.277,6.277,6.277,6.277,6.277,     &
     301       6.277,6.277,8.123,9.094,10.275,11.501,   &
     302                                ! pom soluble
     303       6.676,6.676,6.676,6.676,6.710,6.934,   &
     304       7.141,7.569,8.034,8.529,9.456,10.511,  &
     305       5.109,5.109,5.109,5.109,5.189,5.535,   &
     306       5.960,6.852,8.008,9.712,12.897,19.676, &
     307       3.718,3.718,3.718,3.718,3.779,4.042,   &
     308       4.364,5.052,5.956,7.314,9.896,15.688,  &
     309       2.849,2.849,2.849,2.849,2.897,3.107,   &
     310       3.365,3.916,4.649,5.760,7.900,12.863,  &
     311       2.229,2.229,2.229,2.229,2.268,2.437,   &
     312       2.645,3.095,3.692,4.608,6.391,10.633,  &
     313                                ! Sulfate (Accumulation)
     314       5.751,6.215,6.690,7.024,7.599,8.195,      &
     315       9.156,10.355,12.660,14.823,18.908,24.508, &
     316       4.320,4.675,5.052,5.375,5.787,6.274,      &
     317       7.066,8.083,10.088,12.003,15.697,21.133,  &
     318       3.079,3.351,3.639,3.886,4.205,4.584,      &
     319       5.206,6.019,7.648,9.234,12.391,17.220,    &
     320       2.336,2.552,2.781,2.979,3.236,3.540,      &
     321       4.046,4.711,6.056,7.388,10.093,14.313,    &
     322       1.777,1.949,2.134,2.292,2.503,2.751,      &
     323       3.166,3.712,4.828,5.949,8.264,11.922,     &
     324                                ! Sulfate (Coarse)
     325       5.751,6.215,6.690,7.024,7.599,8.195,      &
     326       9.156,10.355,12.660,14.823,18.908,24.508, &
     327       4.320,4.675,5.052,5.375,5.787,6.274,      &
     328       7.066,8.083,10.088,12.003,15.697,21.133,  &
     329       3.079,3.351,3.639,3.886,4.205,4.584,      &
     330       5.206,6.019,7.648,9.234,12.391,17.220,    &
     331       2.336,2.552,2.781,2.979,3.236,3.540,      &
     332       4.046,4.711,6.056,7.388,10.093,14.313,    &
     333       1.777,1.949,2.134,2.292,2.503,2.751,      &
     334       3.166,3.712,4.828,5.949,8.264,11.922,     &
     335                                ! Seasalt soluble super_coarse (computed below for 550nm)
     336       0.50,0.90,1.05,1.21,1.40,2.41, &
     337       2.66,3.11,3.88,4.52,5.69,8.84, &
     338       0.000,0.000,0.000,0.000,0.000,0.000, &
     339       0.000,0.000,0.000,0.000,0.000,0.000, &
     340     0.52,0.93,1.08,1.24,1.43,2.47, &
     341     2.73,3.20,3.99,4.64,5.84,9.04, &
     342     0.52,0.93,1.09,1.25,1.44,2.50, &
     343     2.76,3.23,4.03,4.68,5.89,9.14, &
     344     0.52,0.94,1.09,1.26,1.45,2.51, &
     345     2.78,3.25,4.06,4.72,5.94,9.22, &
     346                                ! seasalt soluble coarse (computed below for 550nm)
     347       0.50,0.90,1.05,1.21,1.40,2.41, &
     348       2.66,3.11,3.88,4.52,5.69,8.84, &
     349       0.000,0.000,0.000,0.000,0.000,0.000, &
     350       0.000,0.000,0.000,0.000,0.000,0.000, &
     351     0.52,0.93,1.08,1.24,1.43,2.47, &
     352     2.73,3.20,3.99,4.64,5.84,9.04, &
     353     0.52,0.93,1.09,1.25,1.44,2.50, &
     354     2.76,3.23,4.03,4.68,5.89,9.14, &
     355     0.52,0.94,1.09,1.26,1.45,2.51, &
     356     2.78,3.25,4.06,4.72,5.94,9.22, &
     357                                ! seasalt soluble accumulation (computed below for 550nm)
     358     4.28, 7.17, 8.44, 9.85,11.60,22.44,  &
     359     25.34,30.54,39.38,46.52,59.33,91.77, &
     360       0.000,0.000,0.000,0.000,0.000,0.000, &
     361       0.000,0.000,0.000,0.000,0.000,0.000, &
     362     2.48, 4.22, 5.02, 5.94, 7.11,15.29,  &
     363     17.70,22.31,30.73,38.06,52.15,90.59, &
     364     1.90, 3.29, 3.94, 4.69, 5.65, 12.58, &
     365     14.68,18.77,26.41,33.25,46.77,85.50, &
     366     1.47, 2.59, 3.12, 3.74, 4.54, 10.42, &
     367     12.24,15.82,22.66,28.91,41.54,79.33/
    320368
    321369  DATA alpha_aeri_5wv/ &
    322      ! dust insoluble
    323      0.759, 0.770, 0.775, 0.775, 0.772, &
    324      !!jb bc insoluble
    325      11.536,10.033, 8.422, 7.234, 6.270, &
    326      ! pom insoluble
    327      5.042, 3.101, 1.890, 1.294, 0.934/
    328 
    329   DATA cg_aers_5wv/ &
    330      ! bc soluble
    331      .651, .651, .651, .651, .651, .651, &
    332      .651, .651, .738, .764, .785, .800, &
    333      .597, .597, .597, .597, .597, .597, &
    334      .597, .597, .695, .725, .751, .770, &
    335      .543, .543, .543, .543, .543, .543, &
    336      .543, .543, .650, .684, .714, .736, &
    337      .504, .504, .504, .504, .504, .504, &
    338      .504, .504, .614, .651, .683, .708, &
    339      .469, .469, .469, .469, .469, .469, &
    340      .469, .469, .582, .620, .655, .681, &
    341      ! pom soluble
    342      .679, .679, .679, .679, .683, .691, &
    343      .703, .720, .736, .751, .766, .784, &
    344      .656, .656, .656, .656, .659, .669, &
    345      .681, .699, .717, .735, .750, .779, &
    346      .623, .623, .623, .623, .627, .637, &
    347      .649, .668, .688, .709, .734, .762, &
    348      .592, .592, .592, .592, .595, .605, &
    349      .618, .639, .660, .682, .711, .743, &
    350      .561, .561, .561, .561, .565, .575, &
    351      .588, .609, .632, .656, .688, .724, &
    352      ! sulfate
    353      .671, .684, .697, .704, .714, .723, &
    354      .734, .746, .762, .771, .781, .789, &
    355      .653, .666, .678, .687, .697, .707, &
    356      .719, .732, .751, .762, .775, .789, &
    357      .622, .635, .648, .657, .667, .678, &
    358      .691, .705, .728, .741, .758, .777, &
    359      .591, .604, .617, .627, .638, .650, &
    360      .664, .679, .704, .719, .739, .761, &
    361      .560, .574, .587, .597, .609, .621, &
    362      .637, .653, .680, .697, .719, .745/
    363   !
     370                                ! dust insoluble
     371       0.759, 0.770, 0.775, 0.775, 0.772, &
     372                                !!jb bc insoluble
     373       11.536,10.033, 8.422, 7.234, 6.270, &
     374                                ! pom insoluble
     375       5.042, 3.101, 1.890, 1.294, 0.934/
     376  !
     377 DATA cg_aers_5wv/ &
     378                                ! bc soluble
     379       .651, .651, .651, .651, .651, .651, &
     380       .651, .651, .738, .764, .785, .800, &
     381       .597, .597, .597, .597, .597, .597, &
     382       .597, .597, .695, .725, .751, .770, &
     383       .543, .543, .543, .543, .543, .543, &
     384       .543, .543, .650, .684, .714, .736, &
     385       .504, .504, .504, .504, .504, .504, &
     386       .504, .504, .614, .651, .683, .708, &
     387       .469, .469, .469, .469, .469, .469, &
     388       .469, .469, .582, .620, .655, .681, &
     389                                ! pom soluble
     390       .679, .679, .679, .679, .683, .691, &
     391       .703, .720, .736, .751, .766, .784, &
     392       .656, .656, .656, .656, .659, .669, &
     393       .681, .699, .717, .735, .750, .779, &
     394       .623, .623, .623, .623, .627, .637, &
     395       .649, .668, .688, .709, .734, .762, &
     396       .592, .592, .592, .592, .595, .605, &
     397       .618, .639, .660, .682, .711, .743, &
     398       .561, .561, .561, .561, .565, .575, &
     399       .588, .609, .632, .656, .688, .724, &
     400                                ! Accumulation sulfate
     401       .671, .684, .697, .704, .714, .723, &
     402       .734, .746, .762, .771, .781, .789, &
     403       .653, .666, .678, .687, .697, .707, &
     404       .719, .732, .751, .762, .775, .789, &
     405       .622, .635, .648, .657, .667, .678, &
     406       .691, .705, .728, .741, .758, .777, &
     407       .591, .604, .617, .627, .638, .650, &
     408       .664, .679, .704, .719, .739, .761, &
     409       .560, .574, .587, .597, .609, .621, &
     410       .637, .653, .680, .697, .719, .745, &
     411                                ! Coarse sulfate
     412       .671, .684, .697, .704, .714, .723, &
     413       .734, .746, .762, .771, .781, .789, &
     414       .653, .666, .678, .687, .697, .707, &
     415       .719, .732, .751, .762, .775, .789, &
     416       .622, .635, .648, .657, .667, .678, &
     417       .691, .705, .728, .741, .758, .777, &
     418       .591, .604, .617, .627, .638, .650, &
     419       .664, .679, .704, .719, .739, .761, &
     420       .560, .574, .587, .597, .609, .621, &
     421       .637, .653, .680, .697, .719, .745, &
     422                                ! For super coarse seasalt (computed below for 550nm!)
     423      0.730,0.753,0.760,0.766,0.772,0.793, &
     424      0.797,0.802,0.809,0.813,0.820,0.830, &
     425       0.000,0.000,0.000,0.000,0.000,0.000, &
     426       0.000,0.000,0.000,0.000,0.000,0.000, &
     427     0.721,0.744,0.750,0.756,0.762,0.784, &
     428     0.787,0.793,0.800,0.804,0.811,0.822, &
     429     0.717,0.741,0.747,0.753,0.759,0.780, &
     430     0.784,0.789,0.795,0.800,0.806,0.817, &
     431     0.715,0.739,0.745,0.751,0.757,0.777, & 
     432     0.781,0.786,0.793,0.797,0.803,0.814, &
     433                                ! For coarse-soluble seasalt (computed below for 550nm!)
     434     0.730,0.753,0.760,0.766,0.772,0.793, &
     435     0.797,0.802,0.809,0.813,0.820,0.830, &
     436       0.000,0.000,0.000,0.000,0.000,0.000, &
     437       0.000,0.000,0.000,0.000,0.000,0.000, &
     438     0.721,0.744,0.750,0.756,0.762,0.784, &
     439     0.787,0.793,0.800,0.804,0.811,0.822, &
     440     0.717,0.741,0.747,0.753,0.759,0.780, &
     441     0.784,0.789,0.795,0.800,0.806,0.817, &
     442     0.715,0.739,0.745,0.751,0.757,0.777, & 
     443     0.781,0.786,0.793,0.797,0.803,0.814, &
     444                                ! accumulation-seasalt soluble (computed below for 550nm!)
     445     0.698,0.722,0.729,0.736,0.743,0.765, &
     446     0.768,0.773,0.777,0.779,0.781,0.779, &
     447       0.000,0.000,0.000,0.000,0.000,0.000, &
     448       0.000,0.000,0.000,0.000,0.000,0.000, &
     449     0.658,0.691,0.701,0.710,0.720,0.756, &
     450     0.763,0.771,0.782,0.788,0.795,0.801, &
     451     0.632,0.668,0.679,0.690,0.701,0.743, &
     452     0.750,0.762,0.775,0.783,0.792,0.804, &
     453     0.605,0.644,0.656,0.669,0.681,0.729, &
     454     0.737,0.750,0.765,0.775,0.787,0.803/
    364455
    365456  DATA cg_aeri_5wv/&
    … …  
    372463  !
    373464  DATA piz_aers_5wv/&
    374      ! bc soluble
    375      .445, .445, .445, .445, .445, .445, &
    376      .445, .445, .470, .487, .508, .531, &
    377      .442, .442, .442, .442, .442, .442, &
    378      .442, .442, .462, .481, .506, .533, &
    379      .427, .427, .427, .427, .427, .427, &
    380      .427, .427, .449, .470, .497, .526, &
    381      .413, .413, .413, .413, .413, .413, &
    382      .413, .413, .437, .458, .486, .516, &
    383      .399, .399, .399, .399, .399, .399, &
    384      .399, .399, .423, .445, .473, .506, &
    385      ! pom soluble
    386      .975, .975, .975, .975, .975, .977, &
    387      .979, .982, .984, .987, .990, .994, &
    388      .972, .972, .972, .972, .973, .974, &
    389      .977, .980, .983, .986, .989, .993, &
    390      .963, .963, .963, .963, .964, .966, &
    391      .969, .974, .977, .982, .986, .991, &
    392      .955, .955, .955, .955, .955, .958, &
    393      .962, .967, .972, .977, .983, .989, &
    394      .944, .944, .944, .944, .944, .948, &
    395      .952, .959, .962, .972, .979, .987, &
    396      ! sulfate
    397      1.000,1.000,1.000,1.000,1.000,1.000, &
    398      1.000,1.000,1.000,1.000,1.000,1.000, &
    399      1.000,1.000,1.000,1.000,1.000,1.000, &
    400      1.000,1.000,1.000,1.000,1.000,1.000, &
    401      1.000,1.000,1.000,1.000,1.000,1.000, &
    402      1.000,1.000,1.000,1.000,1.000,1.000, &
    403      1.000,1.000,1.000,1.000,1.000,1.000, &
    404      1.000,1.000,1.000,1.000,1.000,1.000, &
    405      1.000,1.000,1.000,1.000,1.000,1.000, &
    406      1.000,1.000,1.000,1.000,1.000,1.000/
     465                                ! bc soluble
     466       .445, .445, .445, .445, .445, .445, &
     467       .445, .445, .470, .487, .508, .531, &
     468       .442, .442, .442, .442, .442, .442, &
     469       .442, .442, .462, .481, .506, .533, &
     470       .427, .427, .427, .427, .427, .427, &
     471       .427, .427, .449, .470, .497, .526, &
     472       .413, .413, .413, .413, .413, .413, &
     473       .413, .413, .437, .458, .486, .516, &
     474       .399, .399, .399, .399, .399, .399, &
     475       .399, .399, .423, .445, .473, .506, &
     476                                ! pom soluble
     477       .975, .975, .975, .975, .975, .977, &
     478       .979, .982, .984, .987, .990, .994, &
     479       .972, .972, .972, .972, .973, .974, &
     480       .977, .980, .983, .986, .989, .993, &
     481       .963, .963, .963, .963, .964, .966, &
     482       .969, .974, .977, .982, .986, .991, &
     483       .955, .955, .955, .955, .955, .958, &
     484       .962, .967, .972, .977, .983, .989, &
     485       .944, .944, .944, .944, .944, .948, &
     486       .952, .959, .962, .972, .979, .987, &
     487                                ! sulfate soluble accumulation
     488       1.000,1.000,1.000,1.000,1.000,1.000, &
     489       1.000,1.000,1.000,1.000,1.000,1.000, &
     490       1.000,1.000,1.000,1.000,1.000,1.000, &
     491       1.000,1.000,1.000,1.000,1.000,1.000, &
     492       1.000,1.000,1.000,1.000,1.000,1.000, &
     493       1.000,1.000,1.000,1.000,1.000,1.000, &
     494       1.000,1.000,1.000,1.000,1.000,1.000, &
     495       1.000,1.000,1.000,1.000,1.000,1.000, &
     496       1.000,1.000,1.000,1.000,1.000,1.000, &
     497       1.000,1.000,1.000,1.000,1.000,1.000, &
     498                                ! sulfate soluble coarse
     499       1.000,1.000,1.000,1.000,1.000,1.000, &
     500       1.000,1.000,1.000,1.000,1.000,1.000, &
     501       1.000,1.000,1.000,1.000,1.000,1.000, &
     502       1.000,1.000,1.000,1.000,1.000,1.000, &
     503       1.000,1.000,1.000,1.000,1.000,1.000, &
     504       1.000,1.000,1.000,1.000,1.000,1.000, &
     505       1.000,1.000,1.000,1.000,1.000,1.000, &
     506       1.000,1.000,1.000,1.000,1.000,1.000, &
     507       1.000,1.000,1.000,1.000,1.000,1.000, &
     508       1.000,1.000,1.000,1.000,1.000,1.000, &
     509                                ! seasalt super coarse (computed below for 550nm)
     510       1.000,1.000,1.000,1.000,1.000,1.000, &
     511       1.000,1.000,1.000,1.000,1.000,1.000, &
     512       1.000,1.000,1.000,1.000,1.000,1.000, &
     513       1.000,1.000,1.000,1.000,1.000,1.000, &
     514       1.000,1.000,1.000,1.000,1.000,1.000, &
     515       1.000,1.000,1.000,1.000,1.000,1.000, &
     516       1.000,1.000,1.000,1.000,1.000,1.000, &
     517       1.000,1.000,1.000,1.000,1.000,1.000, &
     518       1.000,1.000,1.000,1.000,1.000,1.000, &
     519       1.000,1.000,1.000,1.000,1.000,1.000, &
     520                                ! seasalt coarse (computed below for 550nm)
     521       1.000,1.000,1.000,1.000,1.000,1.000, &
     522       1.000,1.000,1.000,1.000,1.000,1.000, &
     523       1.000,1.000,1.000,1.000,1.000,1.000, &
     524       1.000,1.000,1.000,1.000,1.000,1.000, &
     525       1.000,1.000,1.000,1.000,1.000,1.000, &
     526       1.000,1.000,1.000,1.000,1.000,1.000, &
     527       1.000,1.000,1.000,1.000,1.000,1.000, &
     528       1.000,1.000,1.000,1.000,1.000,1.000, &
     529       1.000,1.000,1.000,1.000,1.000,1.000, &
     530       1.000,1.000,1.000,1.000,1.000,1.000, &
     531                                ! seasalt soluble accumulation (computed below for 550nm)
     532       1.000,1.000,1.000,1.000,1.000,1.000, &
     533       1.000,1.000,1.000,1.000,1.000,1.000, &
     534       1.000,1.000,1.000,1.000,1.000,1.000, &
     535       1.000,1.000,1.000,1.000,1.000,1.000, &
     536       1.000,1.000,1.000,1.000,1.000,1.000, &
     537       1.000,1.000,1.000,1.000,1.000,1.000, &
     538       1.000,1.000,1.000,1.000,1.000,1.000, &
     539       1.000,1.000,1.000,1.000,1.000,1.000, &
     540       1.000,1.000,1.000,1.000,1.000,1.000, &
     541       1.000,1.000,1.000,1.000,1.000,1.000/
    407542  !
    408543  DATA piz_aeri_5wv/&
    … …  
    417552  IF (firstcall) THEN
    418553     firstcall=.FALSE.
     554! Allocation
     555    IF (.NOT. ALLOCATED(A1_ASSSM)) THEN
     556        ALLOCATE(A1_ASSSM(klev),A2_ASSSM(klev), A3_ASSSM(klev),&
     557          B1_ASSSM(klev), B2_ASSSM(klev), C1_ASSSM(klev), C2_ASSSM(klev),&
     558          A1_CSSSM(klev), A2_CSSSM(klev), A3_CSSSM(klev),&
     559          B1_CSSSM(klev), B2_CSSSM(klev), C1_CSSSM(klev), C2_CSSSM(klev),&
     560          A1_SSSSM(klev), A2_SSSSM(klev), A3_SSSSM(klev),&
     561          B1_SSSSM(klev), B2_SSSSM(klev), C1_SSSSM(klev), C2_SSSSM(klev), stat=ierr)
     562        IF (ierr /= 0) CALL abort_gcm('aeropt_5mw', 'pb in allocation 1',1)
     563     END IF
     564
    419565!Accumulation mode
    420566     CALL pres2lev(A1_ASSSM_19, A1_ASSSM, nb_level, klev, presnivs_19, presnivs, 1, 1, .FALSE.)
    … …  
    456602  DO k=1, klev
    457603    DO i=1, klon
    458       IF (t_seri(i,k).EQ.0) stop 'stop aeropt_5wv T '
    459       IF (pplay(i,k).EQ.0) stop  'stop aeropt_5wv p '
    460 
     604!      IF (t_seri(i,k).EQ.0) stop 'stop aeropt_5wv T '
     605!      IF (pplay(i,k).EQ.0) stop  'stop aeropt_5wv p '
    461606      zrho=pplay(i,k)/t_seri(i,k)/RD                  ! kg/m3
    462607      mass_temp(i,k,:) = m_allaer(i,k,:) / zrho / 1.e+9
    463 
    464608    ENDDO
    465609  ENDDO
    … …  
    519663  zdp1=pdel(:,:)/(gravit*delt)
    520664 
    521   DO m=1,nb_aer   ! tau is only computed for each mass
    522    
     665  DO m=1,nb_aer   ! tau is only computed for each mass   
    523666    fac=1.0
    524667    IF (aerosol_name(m).EQ.id_ASBCM) THEN
    … …  
    530673        spsol=2
    531674        spss=0
    532     ELSEIF ((aerosol_name(m).EQ.id_ASSO4M) .OR.  (aerosol_name(m).EQ.id_CSSO4M)) THEN
     675    ELSEIF (aerosol_name(m).EQ.id_ASSO4M) THEN
    533676        soluble=.TRUE.
    534677        spsol=3
     678        spss=0
     679        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for OD
     680    ELSEIF (aerosol_name(m).EQ.id_CSSO4M) THEN
     681        soluble=.TRUE.
     682        spsol=4
    535683        spss=0
    536684        fac=1.375    ! (NH4)2-SO4/SO4 132/96 mass conversion factor for OD
    … …  
    581729          IF (soluble) THEN
    582730
    583 
    584731              IF((la.EQ.2).AND.(spss.NE.0)) THEN !la=2 corresponds to 550 nm
    585732                  H=rh/100
    … …  
    620767                 mass_temp(i,k,spsol)*1000.*zdp1(i,k)*tau_ae5wv_int(i,k,la)*delt*fac
    621768
    622           ELSE
     769               tausum(i,la,spsol)=tausum(i,la,spsol)+tau3d(i,k)
     770           tau(i,k,la,spsol)=tau3d(i,k)
     771
     772          ELSE                                           ! Case insoluble aerosol
    623773              tau_ae5wv_int(i,k,la) = alpha_aeri_5wv(la,spinsol)
    624774              piz_ae5wv_int(i,k,la) = piz_aeri_5wv(la,spinsol)
    … …  
    626776
    627777              tau3d(i,k) = &
    628                  mass_temp(i,k,7+spinsol)*1000.*zdp1(i,k)*tau_ae5wv_int(i,k,la)*delt*fac
     778                 mass_temp(i,k,naero_soluble+spinsol)*1000.*zdp1(i,k)* &
     779                      tau_ae5wv_int(i,k,la)*delt*fac
     780
     781               tausum(i,la,naero_soluble+spinsol)=tausum(i,la,spsol)+tau3d(i,k)
     782           tau(i,k,la,naero_soluble+spinsol)=tau3d(i,k)
     783
    629784          ENDIF
    630          
    631          
     785
    632786        ENDDO     ! Boucle sur les points géographiques (grille horizontale)
    633787      ENDDO     ! Boucle sur les niveaux verticaux
    634 
    635       IF (soluble) THEN
    636 
    637           tau(:,:,la,spsol)=tau3d(:,:)
    638      
    639           DO k=1, KLEV
    640             DO i=1,KLON
    641               tausum(i,la,spsol)=tausum(i,la,spsol)+tau3d(i,k)
    642             ENDDO
    643           ENDDO
    644       ELSE
    645           tau(:,:,la,spsol)=tau3d(:,:)
    646      
    647           DO k=1, KLEV
    648             DO i=1,KLON
    649               tausum(i,la,5+spinsol)=tausum(i,la,5+spinsol)+tau3d(i,k)
    650             ENDDO
    651           ENDDO
    652       ENDIF
    653 
    654 
    655 
    656788    ENDDO   ! boucle sur les longueurs d'onde
    657789  ENDDO     ! Boucle  sur les masses de traceurs

  • LMDZ4/branches/LMDZ4-dev/libf/phylmd/hines_gwd.F

    r1001 r1237  
     1!
     2! $Id$
     3!
    14      SUBROUTINE HINES_GWD(NLON,NLEV,DTIME,paphm1x, papm1x,
    25     I      rlat,tx,ux,vx,
    … …  
    16661669C  the variances.
    16671670C
    1668       DO 80 N = 1,NAZ
    1669         DO 70 I = IL1,IL2
    1670           IF (I_ALPHA(I,N).LT.0.)  THEN
    1671             WRITE (6,*)
    1672             WRITE (6,*) '******************************'
    1673             WRITE (6,*) 'Hines integral I_ALPHA < 0 '
    1674             WRITE (6,*) '  longitude I=',I
    1675             WRITE (6,*) '  azimuth   N=',N
    1676             WRITE (6,*) '  level   LEV=',LEV
    1677             WRITE (6,*) '  I_ALPHA =',I_ALPHA(I,N)
    1678             WRITE (6,*) '  V_ALPHA =',V_ALPHA(I,LEV,N)
    1679             WRITE (6,*) '  M_ALPHA =',M_ALPHA(I,LEV,N)
    1680             WRITE (6,*) '  Q_ALPHA =',V_ALPHA(I,LEV,N) / BVFB(I)
    1681             WRITE (6,*) '  QM      =',V_ALPHA(I,LEV,N) / BVFB(I)
    1682      ^                                * M_ALPHA(I,LEV,N)
    1683             WRITE (6,*) '******************************'
    1684           END IF
    1685  70     CONTINUE
    1686  80   CONTINUE
     1671c      DO 80 N = 1,NAZ
     1672c        DO 70 I = IL1,IL2
     1673c          IF (I_ALPHA(I,N).LT.0.)  THEN
     1674c            WRITE (6,*)
     1675c            WRITE (6,*) '******************************'
     1676c            WRITE (6,*) 'Hines integral I_ALPHA < 0 '
     1677c            WRITE (6,*) '  longitude I=',I
     1678c            WRITE (6,*) '  azimuth   N=',N
     1679c            WRITE (6,*) '  level   LEV=',LEV
     1680c            WRITE (6,*) '  I_ALPHA =',I_ALPHA(I,N)
     1681c            WRITE (6,*) '  V_ALPHA =',V_ALPHA(I,LEV,N)
     1682c            WRITE (6,*) '  M_ALPHA =',M_ALPHA(I,LEV,N)
     1683c            WRITE (6,*) '  Q_ALPHA =',V_ALPHA(I,LEV,N) / BVFB(I)
     1684c            WRITE (6,*) '  QM      =',V_ALPHA(I,LEV,N) / BVFB(I)
     1685c     ^                                * M_ALPHA(I,LEV,N)
     1686c            WRITE (6,*) '******************************'
     1687c          END IF
     1688c 70     CONTINUE
     1689c 80   CONTINUE
    16871690C
    16881691      RETURN

  • LMDZ4/branches/LMDZ4-dev/libf/phylmd/newmicro.F

    r1220 r1237  
    164164                                !
    165165                  cdnc_pi(i,k) = 10.**(bl95_b0+bl95_b1*
    166      &               log(MAX(mass_solu_aero_pi(i,k),1.e-4))/log(10.))*1.e6 !-m-3
     166     &             log(MAX(mass_solu_aero_pi(i,k),1.e-4))/log(10.))*1.e6 !-m-3
    167167                  cdnc_pi(i,k)=MIN(1000.e6,MAX(20.e6,cdnc_pi(i,k)))
    168168               ENDDO
    … …  
    178178     &                 *((pqlwp(i,k)*pplay(i,k)/(RD * T(i,k))) 
    179179     &                 /(4./3*RPI*1000.*cdnc(i,k)) )**(1./3.)
    180                   rad_chaud_tab(i,k) = MAX(rad_chaud_tab(i,k) * 1e6, 5.)
     180                  rad_chaud_tab(i,k) = MAX(rad_chaud_tab(i,k) *1.e6, 5.)
    181181               ENDDO           
    182182            ENDDO

  • LMDZ4/branches/LMDZ4-dev/libf/phylmd/physiq.F

    r1229 r1237  
    26712671         IF (.NOT. aerosol_couple)
    26722672     &        CALL readaerosol_optic(
    2673      &        debut, new_aod, flag_aerosol, jD_cur-jD_ref, pdtphys,
    2674      &        pplay, paprs, t_seri, rhcl, presnivs,
     2673     &        debut, new_aod, flag_aerosol, itap, jD_cur-jD_ref,
     2674     &        pdtphys, pplay, paprs, t_seri, rhcl, presnivs,
    26752675     &        mass_solu_aero, mass_solu_aero_pi,
    26762676     &        tau_aero, piz_aero, cg_aero,

  • LMDZ4/branches/LMDZ4-dev/libf/phylmd/readaerosol_interp.F90

    r1216 r1237  
    11! $Id$
    22!
    3 SUBROUTINE readaerosol_interp(id_aero, r_day, first, pplay, paprs, t_seri, mass_out, pi_mass_out)
     3SUBROUTINE readaerosol_interp(id_aero, itap, pdtphys, r_day, first, pplay, paprs, t_seri, mass_out, pi_mass_out)
    44!
    55! This routine will return the mass concentration at actual day(mass_out) and
    … …  
    1818  USE aero_mod, ONLY : naero_spc, name_aero
    1919  USE write_field_phy
     20  USE phys_cal_mod
    2021
    2122  IMPLICIT NONE
    … …  
    3233!****************************************************************************************
    3334  INTEGER, INTENT(IN)                    :: id_aero! Identity number for the aerosol to treat
     35  INTEGER, INTENT(IN)                    :: itap   ! Physic step count
     36  REAL, INTENT(IN)                       :: pdtphys! Physic day step
    3437  REAL, INTENT(IN)                       :: r_day  ! Day of integration
    3538  LOGICAL, INTENT(IN)                    :: first  ! First model timestep
    … …  
    4649!****************************************************************************************
    4750  INTEGER                         :: i, k, ierr
    48   INTEGER                         :: iday, iyr
     51  INTEGER                         :: iday, iyr, lmt_pas
    4952  INTEGER                         :: im, day1, day2, im2
    5053  INTEGER                         :: pi_klev_src ! Only for testing purpose
    … …  
    7780
    7881  LOGICAL            :: lnewday      ! Indicates if first time step at a new day
     82  LOGICAL            :: OLDNEWDAY
    7983  LOGICAL,SAVE       :: vert_interp  ! Indicates if vertical interpolation will be done
    8084  LOGICAL,SAVE       :: debug=.FALSE.! Debugging in this subroutine
    … …  
    8892
    8993! Calculation to find if it is a new day
    90   iday = INT(r_day)
     94
     95  IF(mpi_rank == 0)then
     96     PRINT*,'CONTROL PANEL REGARDING TIME STEPING'
     97  ENDIF
     98
     99  ! Use phys_cal_mod
     100  !iday= day_cur
     101  !iyr = year_cur
     102  !im  = mth_cur
     103
     104  iday = INT(r_day)
    91105  iyr  = iday/360
    92106  iday = iday-iyr*360         ! day of the actual year
    … …  
    94108  im   = iday/30 +1           ! the actual month
    95109
    96   ! 0.02 is about 0.5/24, namly less than half an hour
    97   lnewday = (r_day-FLOAT(iday) < 0.02)
     110  IF(MOD(itap-1,NINT(86400./pdtphys)) == 0)THEN
     111     lnewday=.TRUE.
     112  ENDIF
     113
     114  IF(mpi_rank == 0)then
     115     ! 0.02 is about 0.5/24, namly less than half an hour
     116     OLDNEWDAY = (r_day-FLOAT(iday) < 0.02)
     117     ! Once per day, update aerosol fields
     118     lmt_pas = NINT(86400./pdtphys)
     119     PRINT*,'r_day-FLOAT(iday) =',r_day-FLOAT(iday)
     120     PRINT*,'itap =',itap
     121     PRINT*,'pdtphys =',pdtphys
     122     PRINT*,'lmt_pas =',lmt_pas
     123     PRINT*,'iday =',iday
     124     PRINT*,'r_day =',r_day
     125     PRINT*,'NINT(86400./pdtphys) =',NINT(86400./pdtphys)
     126     PRINT*,'MOD(0,1) =',MOD(0,1)
     127     PRINT*,'lnewday =',lnewday
     128     PRINT*,'OLDNEWDAY =',OLDNEWDAY
     129  ENDIF
    98130
    99131  IF (.NOT. ALLOCATED(var_day)) THEN

  • LMDZ4/branches/LMDZ4-dev/libf/phylmd/readaerosol_optic.F90

    r1221 r1237  
    11! $Id$
    22!
    3 SUBROUTINE readaerosol_optic(debut, new_aod, flag_aerosol, rjourvrai, pdtphys, &
    4      pplay, paprs, t_seri, rhcl, presnivs, &
     3SUBROUTINE readaerosol_optic(debut, new_aod, flag_aerosol, itap, rjourvrai, &
     4     pdtphys, pplay, paprs, t_seri, rhcl, presnivs, &
    55     mass_solu_aero, mass_solu_aero_pi, &
    66     tau_aero, piz_aero, cg_aero, &
    … …  
    2121  LOGICAL, INTENT(IN)                      :: new_aod
    2222  INTEGER, INTENT(IN)                      :: flag_aerosol
     23  INTEGER, INTENT(IN)                      :: itap
    2324  REAL, INTENT(IN)                         :: rjourvrai
    2425  REAL, INTENT(IN)                         :: pdtphys
    … …  
    7475       flag_aerosol .EQ. 6 ) THEN
    7576
    76      CALL readaerosol_interp(id_ASSO4M, rjourvrai, debut, pplay, paprs, t_seri, sulfate, sulfate_pi)
     77     CALL readaerosol_interp(id_ASSO4M, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, sulfate, sulfate_pi)
    7778  ELSE
    7879     sulfate(:,:) = 0. ; sulfate_pi(:,:) = 0.
    … …  
    8485
    8586     ! Get bc aerosol distribution
    86      CALL readaerosol_interp(id_ASBCM, rjourvrai, debut, pplay, paprs, t_seri, bcsol, bcsol_pi )
    87      CALL readaerosol_interp(id_AIBCM, rjourvrai, debut, pplay, paprs, t_seri, bcins, bcins_pi )
     87     CALL readaerosol_interp(id_ASBCM, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, bcsol, bcsol_pi )
     88     CALL readaerosol_interp(id_AIBCM, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, bcins, bcins_pi )
    8889  ELSE
    8990     bcsol(:,:) = 0. ; bcsol_pi(:,:) = 0.
    … …  
    9697       flag_aerosol .EQ. 6 ) THEN
    9798
    98      CALL readaerosol_interp(id_ASPOMM, rjourvrai, debut, pplay, paprs, t_seri, pomsol, pomsol_pi)
    99      CALL readaerosol_interp(id_AIPOMM, rjourvrai, debut, pplay, paprs, t_seri, pomins, pomins_pi)
     99     CALL readaerosol_interp(id_ASPOMM, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, pomsol, pomsol_pi)
     100     CALL readaerosol_interp(id_AIPOMM, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, pomins, pomins_pi)
    100101  ELSE
    101102     pomsol(:,:) = 0. ; pomsol_pi(:,:) = 0.
    … …  
    108109      flag_aerosol .EQ. 6 ) THEN
    109110
    110       CALL readaerosol_interp(id_SSSSM ,rjourvrai, debut, pplay, paprs, t_seri, sssupco, sssupco_pi)
    111       CALL readaerosol_interp(id_CSSSM ,rjourvrai, debut, pplay, paprs, t_seri, sscoarse,sscoarse_pi)
    112       CALL readaerosol_interp(id_ASSSM ,rjourvrai, debut, pplay, paprs, t_seri, ssacu, ssacu_pi)
     111      CALL readaerosol_interp(id_SSSSM ,itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, sssupco, sssupco_pi)
     112      CALL readaerosol_interp(id_CSSSM ,itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, sscoarse,sscoarse_pi)
     113      CALL readaerosol_interp(id_ASSSM ,itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, ssacu, ssacu_pi)
    113114
    114115  ELSE
    … …  
    122123      flag_aerosol .EQ. 6 ) THEN
    123124
    124       CALL readaerosol_interp(id_CIDUSTM, rjourvrai, debut, pplay, paprs, t_seri, cidust, cidust_pi)
     125      CALL readaerosol_interp(id_CIDUSTM, itap, pdtphys, rjourvrai, debut, pplay, paprs, t_seri, cidust, cidust_pi)
    125126
    126127  ELSE
    … …  
    163164
    164165     fractnat_allaer(:,:) = 0.
    165      CALL aeropt_2bands( &
     166     CALL aeropt_2bands(                 &
    166167          pdel, m_allaer, pdtphys, rhcl, &
    167168          tau_aero, piz_aero, cg_aero,   &
    … …  
    170171     
    171172     ! aeropt_5wv only for validation and diagnostics.
    172      CALL aeropt_5wv( &
    173           pdel, m_allaer, &
    174           pdtphys, rhcl, aerindex, &
    175           flag_aerosol, pplay, t_seri, &
     173     CALL aeropt_5wv(                    &
     174          pdel, m_allaer,                &
     175          pdtphys, rhcl, aerindex,       &
     176          flag_aerosol, pplay, t_seri,   &
    176177          tausum_aero, tau3d_aero, presnivs)
    177178  ELSE

  • LMDZ4/branches/LMDZ4-dev/libf/phylmd/sw_aeroAR4.F90

    r1231 r1237  
     1!
     2! $Id$
     3!
    14SUBROUTINE SW_AEROAR4(PSCT, PRMU0, PFRAC, &
    25     PPMB, PDP, &
    … …  
    169172  !$OMP THREADPRIVATE(ZFSDN0_AERO)
    170173
     174!
     175  LOGICAL :: AEROSOLFEEDBACK_ACTIVE=.true.
    171176
    172177  IF(.NOT.initialized) THEN
    … …  
    207212        DO JL = 1, KDLON
    208213           ZCLDSW0(JL,JK) = 0.0
    209            ZOZ(JL,JK) = POZON(JL,JK) / dobson_u / 1e3 / RG * PDP(JL,JK)
     214           ZOZ(JL,JK) = POZON(JL,JK)*46.6968/RG &
     215                *PDP(JL,JK)*(101325.0/PPSOL(JL))
    210216        ENDDO
    211217     ENDDO
    … …  
    232238     DO JK = 1 , KFLEV+1
    233239        DO JL = 1, KDLON
    234            ZFSUP0(JL,JK) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    235            ZFSDN0(JL,JK) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    236            ZFSUP0_AERO(JL,JK,1) = ZFSUP0(JL,JK)
    237            ZFSDN0_AERO(JL,JK,1) = ZFSDN0(JL,JK)
     240           ZFSUP0_AERO(JL,JK,1) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
     241           ZFSDN0_AERO(JL,JK,1) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    238242        ENDDO
    239243     ENDDO
    … …  
    261265     DO JK = 1 , KFLEV+1
    262266        DO JL = 1, KDLON
    263            ZFSUP(JL,JK) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    264            ZFSDN(JL,JK) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    265            ZFSUP_AERO(JL,JK,1) = ZFSUP(JL,JK)
    266            ZFSDN_AERO(JL,JK,1) = ZFSDN(JL,JK)
     267           ZFSUP_AERO(JL,JK,1) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
     268           ZFSDN_AERO(JL,JK,1) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    267269        ENDDO
    268270     ENDDO
    269271
    270      ZFSUP0_AERO(:,:,2) = ZFSUP0_AERO(:,:,1)
    271      ZFSDN0_AERO(:,:,2) = ZFSDN0_AERO(:,:,1)
    272      ZFSUP_AERO(:,:,2) = ZFSUP_AERO(:,:,1)
    273      ZFSDN_AERO(:,:,2) = ZFSDN_AERO(:,:,1)
    274272
    275273
    … …  
    298296        DO JK = 1 , KFLEV+1
    299297           DO JL = 1, KDLON
    300               ZFSUPAD0_AERO(JL,JK) = ZFSUP0(JL,JK)
    301               ZFSDNAD0_AERO(JL,JK) = ZFSDN0(JL,JK)
    302               ZFSUP0(JL,JK) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    303               ZFSDN0(JL,JK) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    304               ZFSUP0_AERO(JL,JK,2) = ZFSUP0(JL,JK)
    305               ZFSDN0_AERO(JL,JK,2) = ZFSDN0(JL,JK)
     298              ZFSUP0_AERO(JL,JK,2) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
     299              ZFSDN0_AERO(JL,JK,2) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    306300           ENDDO
    307301        ENDDO
    … …  
    329323        DO JK = 1 , KFLEV+1
    330324           DO JL = 1, KDLON
    331               ZFSUPAD_AERO(JL,JK) = ZFSUP(JL,JK)
    332               ZFSDNAD_AERO(JL,JK) = ZFSDN(JL,JK)
    333               ZFSUP(JL,JK) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    334               ZFSDN(JL,JK) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    335               ZFSUP_AERO(JL,JK,2) = ZFSUP(JL,JK)
    336               ZFSDN_AERO(JL,JK,2) = ZFSDN(JL,JK)
     325              ZFSUP_AERO(JL,JK,2) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
     326              ZFSDN_AERO(JL,JK,2) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    337327           ENDDO
    338328        ENDDO
     329
     330
    339331
    340332        !CAS NAT
    … …  
    358350             ZFDOWN, ZFUP)
    359351
     352! Natural aerosol fluxes
    360353        DO JK = 1 , KFLEV+1
    361354           DO JL = 1, KDLON
    … …  
    392385        ENDDO
    393386
    394         ! clear sky (Yves 01/12/2007)
    395         ! CAS BC (4)
    396         flag_aer=1.0
    397         CALL SWU_LMDAR4(PSCT,ZCLDSW0,PPMB,PPSOL,&
    398              PRMU0,PFRAC,PTAVE,PWV,&
    399              ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
    400         INU = 1
    401         CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
    402              tauaero(:,:,4,:), pizaero(:,:,4,:), cgaero(:,:,4,:),&
    403              PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    404              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    405              ZFD, ZFU)
    406         INU = 2
    407         CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
    408              tauaero(:,:,4,:), pizaero(:,:,4,:), cgaero(:,:,4,:),&
    409              ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    410              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    411              PWV, PQS,&
    412              ZFDOWN, ZFUP)
    413 
    414         DO JK = 1 , KFLEV+1
    415            DO JL = 1, KDLON
    416               ZFSUP0_AERO(JL,JK,4) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    417               ZFSDN0_AERO(JL,JK,4) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    418            ENDDO
    419         ENDDO
    420 
    421         ! cloudy-sky + aerosol dir OB
    422         ! CAS BC (4)
    423         flag_aer=1.0
    424         CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
    425              PRMU0,PFRAC,PTAVE,PWV,&
    426              ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
    427         INU = 1
    428         CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
    429              tauaero(:,:,4,:), pizaero(:,:,4,:), cgaero(:,:,4,:),&
    430              PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    431              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    432              ZFD, ZFU)
    433         INU = 2
    434         CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
    435              tauaero(:,:,4,:), pizaero(:,:,4,:), cgaero(:,:,4,:),&
    436              ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    437              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    438              PWV, PQS,&
    439              ZFDOWN, ZFUP)
    440 
    441         DO JK = 1 , KFLEV+1
    442            DO JL = 1, KDLON
    443               ZFSUP_AERO(JL,JK,4) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    444               ZFSDN_AERO(JL,JK,4) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    445            ENDDO
    446         ENDDO
    447 
    448         ! clear sky (Yves 13/12/2007)
    449         ! CAS SO4 (5)
    450         flag_aer=1.0
    451         CALL SWU_LMDAR4(PSCT,ZCLDSW0,PPMB,PPSOL,&
    452              PRMU0,PFRAC,PTAVE,PWV,&
    453              ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
    454         INU = 1
    455         CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
    456              tauaero(:,:,5,:), pizaero(:,:,5,:), cgaero(:,:,5,:),&
    457              PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    458              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    459              ZFD, ZFU)
    460         INU = 2
    461         CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
    462              tauaero(:,:,5,:), pizaero(:,:,5,:), cgaero(:,:,5,:),&
    463              ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    464              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    465              PWV, PQS,&
    466              ZFDOWN, ZFUP)
    467 
    468         DO JK = 1 , KFLEV+1
    469            DO JL = 1, KDLON
    470               ZFSUP0_AERO(JL,JK,5) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    471               ZFSDN0_AERO(JL,JK,5) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    472            ENDDO
    473         ENDDO
    474 
    475         ! cloudy-sky + aerosol dir OB
    476         ! CAS SO4 (5)
    477         flag_aer=1.0
    478         CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
    479              PRMU0,PFRAC,PTAVE,PWV,&
    480              ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
    481         INU = 1
    482         CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
    483              tauaero(:,:,5,:), pizaero(:,:,5,:), cgaero(:,:,5,:),&
    484              PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    485              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    486              ZFD, ZFU)
    487         INU = 2
    488         CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
    489              tauaero(:,:,5,:), pizaero(:,:,5,:), cgaero(:,:,5,:),&
    490              ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    491              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    492              PWV, PQS,&
    493              ZFDOWN, ZFUP)
    494 
    495         DO JK = 1 , KFLEV+1
    496            DO JL = 1, KDLON
    497               ZFSUP_AERO(JL,JK,5) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    498               ZFSDN_AERO(JL,JK,5) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    499            ENDDO
    500         ENDDO
    501 
    502         ! clear sky (Yves 13/12/2007)
    503         ! CAS POM (6)
    504         flag_aer=1.0
    505         CALL SWU_LMDAR4(PSCT,ZCLDSW0,PPMB,PPSOL,&
    506              PRMU0,PFRAC,PTAVE,PWV,&
    507              ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
    508         INU = 1
    509         CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
    510              tauaero(:,:,6,:), pizaero(:,:,6,:), cgaero(:,:,6,:),&
    511              PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    512              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    513              ZFD, ZFU)
    514         INU = 2
    515         CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
    516              tauaero(:,:,6,:), pizaero(:,:,6,:), cgaero(:,:,6,:),&
    517              ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    518              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    519              PWV, PQS,&
    520              ZFDOWN, ZFUP)
    521 
    522         DO JK = 1 , KFLEV+1
    523            DO JL = 1, KDLON
    524               ZFSUP0_AERO(JL,JK,6) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    525               ZFSDN0_AERO(JL,JK,6) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    526            ENDDO
    527         ENDDO
    528 
    529         ! cloudy-sky + aerosol dir OB
    530         ! CAS POM (6)
    531         flag_aer=1.0
    532         CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
    533              PRMU0,PFRAC,PTAVE,PWV,&
    534              ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
    535         INU = 1
    536         CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
    537              tauaero(:,:,6,:), pizaero(:,:,6,:), cgaero(:,:,6,:),&
    538              PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    539              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    540              ZFD, ZFU)
    541         INU = 2
    542         CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
    543              tauaero(:,:,6,:), pizaero(:,:,6,:), cgaero(:,:,6,:),&
    544              ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    545              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    546              PWV, PQS,&
    547              ZFDOWN, ZFUP)
    548 
    549         DO JK = 1 , KFLEV+1
    550            DO JL = 1, KDLON
    551               ZFSUP_AERO(JL,JK,6) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    552               ZFSDN_AERO(JL,JK,6) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    553            ENDDO
    554         ENDDO
    555 
    556         ! clear sky (Yves 13/12/2007)
    557         ! CAS DUST (7)
    558         flag_aer=1.0
    559         CALL SWU_LMDAR4(PSCT,ZCLDSW0,PPMB,PPSOL,&
    560              PRMU0,PFRAC,PTAVE,PWV,&
    561              ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
    562         INU = 1
    563         CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
    564              tauaero(:,:,7,:), pizaero(:,:,7,:), cgaero(:,:,7,:),&
    565              PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    566              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    567              ZFD, ZFU)
    568         INU = 2
    569         CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
    570              tauaero(:,:,7,:), pizaero(:,:,7,:), cgaero(:,:,7,:),&
    571              ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    572              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    573              PWV, PQS,&
    574              ZFDOWN, ZFUP)
    575 
    576         DO JK = 1 , KFLEV+1
    577            DO JL = 1, KDLON
    578               ZFSUP0_AERO(JL,JK,7) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    579               ZFSDN0_AERO(JL,JK,7) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    580            ENDDO
    581         ENDDO
    582 
    583         ! cloudy-sky + aerosol dir OB
    584         ! CAS DUST (7)
    585         flag_aer=1.0
    586         CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
    587              PRMU0,PFRAC,PTAVE,PWV,&
    588              ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
    589         INU = 1
    590         CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
    591              tauaero(:,:,7,:), pizaero(:,:,7,:), cgaero(:,:,7,:),&
    592              PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    593              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    594              ZFD, ZFU)
    595         INU = 2
    596         CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
    597              tauaero(:,:,7,:), pizaero(:,:,7,:), cgaero(:,:,7,:),&
    598              ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    599              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    600              PWV, PQS,&
    601              ZFDOWN, ZFUP)
    602 
    603         DO JK = 1 , KFLEV+1
    604            DO JL = 1, KDLON
    605               ZFSUP_AERO(JL,JK,7) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    606               ZFSDN_AERO(JL,JK,7) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    607            ENDDO
    608         ENDDO
    609 
    610         ! clear sky (Yves 13/12/2007)
    611         ! CAS Seasalt SS (8)
    612         flag_aer=1.0
    613         CALL SWU_LMDAR4(PSCT,ZCLDSW0,PPMB,PPSOL,&
    614              PRMU0,PFRAC,PTAVE,PWV,&
    615              ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
    616         INU = 1
    617         CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
    618              tauaero(:,:,8,:), pizaero(:,:,8,:), cgaero(:,:,8,:),&
    619              PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    620              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    621              ZFD, ZFU)
    622         INU = 2
    623         CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
    624              tauaero(:,:,8,:), pizaero(:,:,8,:), cgaero(:,:,8,:),&
    625              ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    626              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    627              PWV, PQS,&
    628              ZFDOWN, ZFUP)
    629 
    630         DO JK = 1 , KFLEV+1
    631            DO JL = 1, KDLON
    632               ZFSUP0_AERO(JL,JK,8) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    633               ZFSDN0_AERO(JL,JK,8) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    634            ENDDO
    635         ENDDO
    636 
    637         ! cloudy-sky + aerosol dir OB
    638         ! CAS Seasalt SS (8)
    639         flag_aer=1.0
    640         CALL SWU_LMDAR4(PSCT,PCLDSW,PPMB,PPSOL,&
    641              PRMU0,PFRAC,PTAVE,PWV,&
    642              ZAKI,ZCLD,ZCLEAR,ZDSIG,ZFACT,ZRMU,ZSEC,ZUD)
    643         INU = 1
    644         CALL SW1S_LMDAR4(INU, PAER, flag_aer,&
    645              tauaero(:,:,8,:), pizaero(:,:,8,:), cgaero(:,:,8,:),&
    646              PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    647              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    648              ZFD, ZFU)
    649         INU = 2
    650         CALL SW2S_LMDAR4(INU, PAER, flag_aer,&
    651              tauaero(:,:,8,:), pizaero(:,:,8,:), cgaero(:,:,8,:),&
    652              ZAKI, PALBD, PALBP, PCG, ZCLD, ZCLEAR, PCLDSW,&
    653              ZDSIG, POMEGA, ZOZ, ZRMU, ZSEC, PTAU, ZUD,&
    654              PWV, PQS,&
    655              ZFDOWN, ZFUP)
    656 
    657         DO JK = 1 , KFLEV+1
    658            DO JL = 1, KDLON
    659               ZFSUP_AERO(JL,JK,8) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    660               ZFSDN_AERO(JL,JK,8) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    661            ENDDO
    662         ENDDO
    663387
    664388     ENDIF ! ok_ade
    … …  
    686410        DO JK = 1 , KFLEV+1
    687411           DO JL = 1, KDLON
    688               ZFSUPAI_AERO(JL,JK) = ZFSUP(JL,JK)
    689               ZFSDNAI_AERO(JL,JK) = ZFSDN(JL,JK)         
    690               ZFSUP(JL,JK) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
    691               ZFSDN(JL,JK) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    692               ZFSUP_AERO(JL,JK,2) = ZFSUP(JL,JK)
    693               ZFSDN_AERO(JL,JK,2) = ZFSDN(JL,JK)
     412              ZFSUP_AERO(JL,JK,4) = (ZFUP(JL,JK)   + ZFU(JL,JK)) * ZFACT(JL)
     413              ZFSDN_AERO(JL,JK,4) = (ZFDOWN(JL,JK) + ZFD(JL,JK)) * ZFACT(JL)
    694414           ENDDO
    695415        ENDDO
    … …  
    699419  ENDIF
    700420  itapsw = itapsw + 1
     421
     422
     423  IF ( ok_ade .and. ok_aie  ) THEN
     424    ZFSUP(:,:) =    ZFSUP_AERO(:,:,4)
     425    ZFSDN(:,:) =    ZFSDN_AERO(:,:,4)
     426    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,2)
     427    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,2)
     428  ENDIF
     429  IF ( ok_ade .and. (.not. ok_aie) )  THEN
     430    ZFSUP(:,:) =    ZFSUP_AERO(:,:,2)
     431    ZFSDN(:,:) =    ZFSDN_AERO(:,:,2)
     432    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,2)
     433    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,2)
     434  ENDIF
     435!MS the following combination would include the direct aerosol effect in cloud regions
     436!   because it takes the total aerosol effect
     437  IF ( (.not. ok_ade) .and. ok_aie  )  THEN
     438    ZFSUP(:,:) =    ZFSUP_AERO(:,:,4)
     439    ZFSDN(:,:) =    ZFSDN_AERO(:,:,4)
     440    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,1)
     441    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,1)
     442  ENDIF
     443  IF ((.not. ok_ade) .and. (.not. ok_aie)) THEN
     444    ZFSUP(:,:) =    ZFSUP_AERO(:,:,1)
     445    ZFSDN(:,:) =    ZFSDN_AERO(:,:,1)
     446    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,1)
     447    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,1)
     448  ENDIF
     449
     450! MS the following allows to compute the forcing diagostics without
     451! letting the aerosol forcing act on the meteorology
     452! assuming that the no-aerosol case creates the reference meteorological state
     453! for the natural aerosol state use: *_AERO(:,:3)
     454  IF  (.not. AEROSOLFEEDBACK_ACTIVE) THEN
     455    ZFSUP(:,:) =    ZFSUP_AERO(:,:,1)
     456    ZFSDN(:,:) =    ZFSDN_AERO(:,:,1)
     457    ZFSUP0(:,:) =   ZFSUP0_AERO(:,:,1)
     458    ZFSDN0(:,:) =   ZFSDN0_AERO(:,:,1)
     459  ENDIF
     460 
    701461
    702462  DO k = 1, KFLEV
    … …  
    704464     DO i = 1, KDLON
    705465
    706         PHEAT(i,k) = -(ZFSUP_AERO(i,kpl1,2)-ZFSUP_AERO(i,k,2)) &
    707              -(ZFSDN_AERO(i,k,2)-ZFSDN_AERO(i,kpl1,2))
     466        PHEAT(i,k) = -(ZFSUP(i,kpl1)-ZFSUP(i,k))-(ZFSDN(i,k)-ZFSDN(i,kpl1))
    708467        PHEAT(i,k) = PHEAT(i,k) * RDAY*RG/RCPD / PDP(i,k)
    709         PHEAT0(i,k) = -(ZFSUP0_AERO(i,kpl1,2)-ZFSUP0_AERO(i,k,2)) &
    710              -(ZFSDN0_AERO(i,k,2)-ZFSDN0_AERO(i,kpl1,2))
     468        PHEAT0(i,k) = -(ZFSUP0(i,kpl1)-ZFSUP0(i,k))-(ZFSDN0(i,k)-ZFSDN0(i,kpl1))
    711469        PHEAT0(i,k) = PHEAT0(i,k) * RDAY*RG/RCPD / PDP(i,k)
    712470
    713471     ENDDO
    714472  ENDDO
     473
    715474  DO i = 1, KDLON
    716475     PALBPLA(i) = ZFSUP(i,KFLEV+1)/(ZFSDN(i,KFLEV+1)+1.0e-20)
    … …  
    722481     PTOPSW(i) = ZFSDN(i,KFLEV+1) - ZFSUP(i,KFLEV+1)
    723482
    724      PSOLSW0AERO(i,:) = ZFSDN0_AERO(i,1,:) - ZFSUP0_AERO(i,1,:)
    725      PTOPSW0AERO(i,:) = &
    726           ZFSDN0_AERO(i,KFLEV+1,:) - ZFSUP0_AERO(i,KFLEV+1,:)
    727 
    728      PSOLSWAERO(i,:) = ZFSDN_AERO(i,1,:) - ZFSUP_AERO(i,1,:)
    729      PTOPSWAERO(i,:) = &
    730           ZFSDN_AERO(i,KFLEV+1,:) - ZFSUP_AERO(i,KFLEV+1,:)
    731 
    732      PSOLSWADAERO(i) = ZFSDNAD_AERO(i,1) - ZFSUPAD_AERO(i,1)
    733      PTOPSWADAERO(i) = &
    734           ZFSDNAD_AERO(i,KFLEV+1) - ZFSUPAD_AERO(i,KFLEV+1)
    735 
    736      PSOLSWAD0AERO(i) = ZFSDNAD0_AERO(i,1) - ZFSUPAD0_AERO(i,1)
    737      PTOPSWAD0AERO(i) = &
    738           ZFSDNAD0_AERO(i,KFLEV+1) - ZFSUPAD0_AERO(i,KFLEV+1)
    739 
    740      PSOLSWAIAERO(i) = ZFSDNAI_AERO(i,1) - ZFSUPAI_AERO(i,1)
    741      PTOPSWAIAERO(i) = &
    742           ZFSDNAI_AERO(i,KFLEV+1) - ZFSUPAI_AERO(i,KFLEV+1)
     483! MS the following is not output, so commented
     484!     PSOLSW0AERO(i,:) = ZFSDN0_AERO(i,1,:) - ZFSUP0_AERO(i,1,:)
     485!     PTOPSW0AERO(i,:) = &
     486!          ZFSDN0_AERO(i,KFLEV+1,:) - ZFSUP0_AERO(i,KFLEV+1,:)
     487
     488!     PSOLSWAERO(i,:) = ZFSDN_AERO(i,1,:) - ZFSUP_AERO(i,1,:)
     489!     PTOPSWAERO(i,:) = &
     490!          ZFSDN_AERO(i,KFLEV+1,:) - ZFSUP_AERO(i,KFLEV+1,:)
     491
     492
     493if (ok_ade) then
     494     PSOLSWADAERO(i) = (ZFSDN_AERO(i,1,2) - ZFSUP_AERO(i,1,2))-(ZFSDN_AERO(i,1,3) - ZFSUP_AERO(i,1,3))
     495     PTOPSWADAERO(i) = (ZFSDN_AERO(i,KFLEV+1,2) - ZFSUP_AERO(i,KFLEV+1,2))- (ZFSDN_AERO(i,KFLEV+1,3) - ZFSUP_AERO(i,KFLEV+1,3))
     496
     497     PSOLSWAD0AERO(i) = (ZFSDN0_AERO(i,1,2) - ZFSUP0_AERO(i,1,2))-(ZFSDN0_AERO(i,1,3) - ZFSUP0_AERO(i,1,3))
     498     PTOPSWAD0AERO(i) = (ZFSDN0_AERO(i,KFLEV+1,2) - ZFSUP0_AERO(i,KFLEV+1,2))-(ZFSDN0_AERO(i,KFLEV+1,3) - ZFSUP0_AERO(i,KFLEV+1,3))
     499endif
     500
     501if (ok_aie) then
     502     PSOLSWAIAERO(i) = (ZFSDN_AERO(i,1,4) - ZFSUP_AERO(i,1,4))-(ZFSDN_AERO(i,1,2) - ZFSUP_AERO(i,1,2))
     503     PTOPSWAIAERO(i) = (ZFSDN_AERO(i,KFLEV+1,4) - ZFSUP_AERO(i,KFLEV+1,4))-(ZFSDN_AERO(i,KFLEV+1,2) - ZFSUP_AERO(i,KFLEV+1,2))
     504endif
    743505
    744506  ENDDO
Note: See TracChangeset for help on using the changeset viewer.