Changeset 541

Timestamp:

Feb 23, 2012, 2:20:06 AM (13 years ago)

Author:

aslmd

Message:

LMDZ.MARS. Added small optimizations in watercycle. 10% improvement.

Location:

trunk/LMDZ.MARS

Files:

• README (modified) (1 diff)
• libf/phymars/improvedclouds.F (modified) (5 diffs)
• libf/phymars/watercloud.F (modified) (2 diffs)

trunk/LMDZ.MARS/README

@@ -1407 +1407 @@
 >> Code 35% faster with a Valles Marineris mesoscale run in springtime (set imicro=25)
    ... probably even better at cloudier seasons
+
+== 23/02/12 == AS
+>> In improvedclouds and waterclouds, added optimizations
+   ... moved calculations to firstcall
+   ... moved 'if' statement so that only useful calculations are performed
+>> Code about 10% faster [tests carried out with mesoscale summer Tharsis run]
+>> Check with a GCM run that results similar

trunk/LMDZ.MARS/libf/phymars/improvedclouds.F

@@ -192 +192 @@
         write(*,'(i2,3x,e12.6)') nbin_cld+1,rb_cld(nbin_cld+1)
         print*,'-----------------------------------'
+
+        do i=1,nbin_cld+1
+            rb_cld(i) = dlog(rb_cld(i))  !! we save that so that it is not computed
+                                         !! at each timestep and gridpoint
+        enddo
 
 c       Contact parameter of water ice on dust ( m=cos(theta) )
@@ -265 +270 @@
 
       Cste = ptimestep * 4. * pi * rho_ice
+      dev2 = 1. / ( sqrt(2.) * sigma_ice )
 
       call watersat(ngridmx*nlayermx,zt,pplay,zqsat)
@@ -286 +292 @@
         Mo = zq(ig,l,igcm_dust_mass)* tauscaling(ig)
         No = zq(ig,l,igcm_dust_number)* tauscaling(ig)+ 1.e-30
-        dev2 = 1. / ( sqrt(2.) * sigma_ice )
         Rn = rdust(ig,l)
         Rn = -dlog(Rn) 
         Rm = Rn - 3. * sigma_ice*sigma_ice  
-        yeah = dlog(rb_cld(1))
-        yeahn = derf( (yeah+Rn) *dev2)
-        yeahm = derf( (yeah+Rm) *dev2)
+        yeahn = derf( (rb_cld(1)+Rn) *dev2)
+        yeahm = derf( (rb_cld(1)+Rm) *dev2)
         do i = 1, nbin_cld
           n_aer(i) = -0.5 * No * yeahn !! this ith previously computed
           m_aer(i) = -0.5 * Mo * yeahm !! this ith previously computed
-          yeah = dlog(rb_cld(i+1))     !! this (i+1)th now computed
-          yeahn = derf( (yeah+Rn) *dev2)
-          yeahm = derf( (yeah+Rm) *dev2)
+          yeahn = derf( (rb_cld(i+1)+Rn) *dev2)
+          yeahm = derf( (rb_cld(i+1)+Rm) *dev2)
           n_aer(i) = n_aer(i) + 0.5 * No * yeahn
           m_aer(i) = m_aer(i) + 0.5 * Mo * yeahm
@@ -363 +366 @@
      &                      / Mo * rho_dust
           rhocloud(ig,l) = min(max(rhocloud(ig,l),rho_ice),rho_dust)
-          rice(ig,l) =
-     &      CBRT( real(Mo)/real(No) * 0.75 / pi / rhocloud(ig,l) ) !**(1./3.)
+          if ((Mo.lt.1.e-20) .or. (No.le.1)) then
+              rice(ig,l) = 1.e-8
+          else
+              rice(ig,l) =
+     &          CBRT( real(Mo)/real(No) * 0.75 / pi / rhocloud(ig,l) ) !**(1./3.)
+          endif
 c          nuice(ig,l)=nuice_ref ! used for rad. transfer calculations
-          if ((Mo.lt.1.e-20) .or. (No.le.1)) rice(ig,l) = 1.e-8
           seq  = exp( 2.*sig(zt(ig,l))*mh2o / 
      &           (rho_ice*rgp*zt(ig,l)*rice(ig,l)) )
@@ -483 +489 @@
           
           rice_out(ig,l)=rice(ig,l)
-          rice(ig,l) =
-     &      CBRT( real(Mo)/real(No) * 0.75 / pi / rhocloud(ig,l) ) !**(1./3.)
-          if ((Mo.lt.1.e-20) .or. (No.le.1)) rice(ig,l) = 1.e-8
+          if ((Mo.lt.1.e-20) .or. (No.le.1)) then
+              rice(ig,l) = 1.e-8
+          else
+              rice(ig,l) =
+     &          CBRT( real(Mo)/real(No) * 0.75 / pi / rhocloud(ig,l) ) !**(1./3.)
+          endif
           rice_out(ig,l)=rice(ig,l)-rice_out(ig,l)
           

trunk/LMDZ.MARS/libf/phymars/watercloud.F

@@ -402 +402 @@
      &                   *tauscaling(ig)/ Mo * rho_dust
           rhocloud(ig,l) = min(max(rhocloud(ig,l),rho_ice),rho_dust)
-          rice(ig,l)=(Mo / No * 0.75 / pi / rhocloud(ig,l))**(1./3.)
-           if ((Mo.lt.1.e-15) .or. (No.le.50)) rice(ig,l) = 1.e-8
+          if ((Mo.lt.1.e-15) .or. (No.le.50)) then
+              rice(ig,l) = 1.e-8
+          else
+              !! AS: only perform computations if conditions not met
+              rice(ig,l)=(Mo / No * 0.75 / pi / rhocloud(ig,l))**(1./3.)
+          endif
          ENDDO
         ENDDO
@@ -414 +418 @@
             rice(ig,l)=max( CBRT ( ((pq(ig,l,igcm_h2o_ice)
      &       + pdqcloud(ig,l,igcm_h2o_ice)*ptimestep)/rho_ice
-     &      +ccntyp*(4./3.)*pi*rdust(ig,l)**3.)
+     &      +ccntyp*(4./3.)*pi*rdust(ig,l)*rdust(ig,l)*rdust(ig,l))
      &      /(ccntyp*4./3.*pi) ), rdust(ig,l))
           ENDDO