Changeset 2634 for LMDZ5/trunk/libf/phylmd/rrtm/aeropt_5wv_rrtm.F90

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.

File:

• LMDZ5/trunk/libf/phylmd/rrtm/aeropt_5wv_rrtm.F90 (modified) (8 diffs)

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

@@ -94 +94 @@
   REAL :: od670aer(klon)  ! epaisseur optique aerosol extinction 670 nm
   REAL :: fac
-  REAL :: zdp1(klon,klev) 
   INTEGER, ALLOCATABLE, DIMENSION(:)  :: aerosol_name
   INTEGER :: nb_aer, itau
   LOGICAL :: ok_itau
   
-  REAL :: dh(KLON,KLEV)
+  REAL :: zdh(klon,klev)
   
    ! Soluble components 1- BC soluble; 2- POM soluble; 3- SO4 acc.; 4- SO4 coarse; 5 seasalt super-coarse; 6 seasalt coarse; 7 seasalt acc.
@@ -106 +105 @@
   REAL :: alpha_aeri_5wv(las,naero_insoluble)         ! Ext. coeff. ** m2/g 
 
-  REAL, DIMENSION(klon,klev,naero_tot) :: mass_temp
-  
   !
   ! Proprietes optiques
   !
   REAL :: fact_RH(nbre_RH)
-  LOGICAL :: used_tau(naero_tot)
   INTEGER :: n
 
@@ -209 +205 @@
   ai(:) = 0.
   tausum(:,:,:) = 0.
+  tau(:,:,:,:)=0.
 
   DO k=1, klev
     DO i=1, klon
       zrho=pplay(i,k)/t_seri(i,k)/RD                  ! kg/m3
-      dh(i,k)=pdel(i,k)/(RG*zrho)
-!CDIR UNROLL=naero_spc
-      mass_temp(i,k,:) = m_allaer(i,k,:) / zrho / 1.e+9
-      zdp1(i,k)=pdel(i,k)/(RG*delt)     ! air mass auxiliary  variable --> zdp1 [kg/(m^2 *s)]
+      zdh(i,k)=pdel(i,k)/(RG*zrho)                    ! m
     ENDDO
   ENDDO
@@ -270 +264 @@
   !      
 
-!CDIR ON_ADB(RH_tab)
-!CDIR ON_ADB(fact_RH)
-!CDIR NOVECTOR
   DO n=1,nbre_RH-1
     fact_RH(n)=1./(RH_tab(n+1)-RH_tab(n))
   ENDDO
    
-  DO k=1, KLEV
-!CDIR ON_ADB(RH_tab)
-!CDIR ON_ADB(fact_RH)
-    DO i=1, KLON
+  DO k=1, klev
+    DO i=1, klon
       rh(i,k)=MIN(RHcl(i,k)*100.,RH_MAX)
       RH_num(i,k) = INT( rh(i,k)/10. + 1.)
@@ -289 +278 @@
   ENDDO
 
-!CDIR SHORTLOOP  
-  used_tau(:)=.FALSE.
-    
   DO m=1,nb_aer   ! tau is only computed for each mass    
     fac=1.0
@@ -330 +316 @@
     ENDIF
 
-    IF (soluble) then
-      used_tau(spsol)=.TRUE.
-    ELSE
-      used_tau(naero_soluble+spinsol)=.TRUE.
-    ENDIF
-
     aerindex=aerosol_name(m)
 
     DO la=1,las
 
+    !--only 550, 670 and 865 nm are used
+    IF (la.NE.la550.AND.la.NE.la670.AND.la.NE.la865) CYCLE
+
       IF (soluble) THEN            ! For soluble aerosol
 
-          DO k=1, KLEV
-            DO i=1, KLON
+          DO k=1, klev
+            DO i=1, klon
               tau_ae5wv_int = alpha_aers_5wv(RH_num(i,k),la,spsol)+DELTA(i,k)* &
                              (alpha_aers_5wv(RH_num(i,k)+1,la,spsol) - & 
                               alpha_aers_5wv(RH_num(i,k),la,spsol))
-              tau(i,k,la,aerindex) = mass_temp(i,k,aerindex)*1000.*zdp1(i,k)*   &
-                                    tau_ae5wv_int*delt*fac
+              tau(i,k,la,aerindex) = m_allaer(i,k,aerindex)/1.e6*zdh(i,k)*tau_ae5wv_int*fac
               tausum(i,la,aerindex)=tausum(i,la,aerindex)+tau(i,k,la,aerindex)
             ENDDO
@@ -355 +337 @@
       ELSE                         ! For insoluble aerosol
 
-        DO k=1, KLEV
-          DO i=1, KLON
+        DO k=1, klev
+          DO i=1, klon
             tau_ae5wv_int = alpha_aeri_5wv(la,spinsol)
-            tau(i,k,la,aerindex) = mass_temp(i,k,aerindex)*1000.*zdp1(i,k)* &
-                                   tau_ae5wv_int*delt*fac
+            tau(i,k,la,aerindex) = m_allaer(i,k,aerindex)/1.e6*zdh(i,k)*tau_ae5wv_int*fac
             tausum(i,la,aerindex)= tausum(i,la,aerindex)+tau(i,k,la,aerindex)
           ENDDO
@@ -369 +350 @@
   ENDDO     ! Boucle sur les masses de traceurs
 
-  DO m=1,naero_tot
-    IF (.NOT.used_tau(m)) tau(:,:,:,m)=0.
-  ENDDO  
-
-  DO i=1, klon
-     od550aer(i)=0.
-     DO m=1,naero_tot
-        od550aer(i)=od550aer(i)+tausum(i,la550,m)
-     END DO
-  END DO
-
-  DO i=1, klon
-     od670aer(i)=0.
-     DO m=1,naero_tot
-        od670aer(i)=od670aer(i)+tausum(i,la670,m)
-     END DO
-  END DO
-
-  DO i=1, klon
-     od865aer(i)=0.
-     DO m=1,naero_tot
-        od865aer(i)=od865aer(i)+tausum(i,la865,m)
-     END DO
-  END DO
-
-  DO i=1, klon
-     DO k=1, KLEV
-        ec550aer(i,k)=0.
-        DO m=1,naero_tot
-           ec550aer(i,k)=ec550aer(i,k)+tau(i,k,la550,m)/dh(i,k)
-        END DO
-     END DO
-  END DO
-  
-  DO i=1, klon
-    ai(i)=-LOG(MAX(od670aer(i),1.e-8)/MAX(od865aer(i),1.e-8))/LOG(670./865.)
-  ENDDO
+!--AOD calculations for diagnostics
+  od550aer(:)=SUM(tausum(:,la550,:),dim=2)
+  od670aer(:)=SUM(tausum(:,la670,:),dim=2)
+  od865aer(:)=SUM(tausum(:,la865,:),dim=2)
+
+!--extinction coefficient for diagnostic
+  ec550aer(:,:)=SUM(tau(:,:,la550,:),dim=3)/zdh(:,:)
+
+!--aerosol index
+  ai(:)=-LOG(MAX(od670aer(:),1.e-8)/MAX(od865aer(:),1.e-8))/LOG(670./865.)
 
   od550lt1aer(:)=tausum(:,la550,id_ASSO4M_phy)+tausum(:,la550,id_ASBCM_phy) +tausum(:,la550,id_AIBCM_phy)+ &