Changeset 1410 for trunk/LMDZ.MARS/libf/phymars/aeropacity.F

Timestamp:

Apr 7, 2015, 3:41:32 PM (10 years ago)

Author:

aslmd

Message:

cosmetic changes for DUSTSTORM mode (better comments and location of code changes)

File:

• trunk/LMDZ.MARS/libf/phymars/aeropacity.F (modified) (10 diffs)

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

@@ -123 +123 @@
       real l_top ! layer of the storm's top
       REAL zalt(ngrid, nlayer) ! useful factor to compute l_top 
-      REAL zdqnorm(ngrid, nlayer,2)
 #endif
 
@@ -459 +458 @@
 #ifdef DUSTSTORM
 c -----------------------------------------------------------------
-c the quantity of dust to add at the first time step is calculated
-c to match a tunable opacity perturbation.
+! Calculate reference opacity without perturbation
 c -----------------------------------------------------------------
       IF (firstcall) THEN
@@ -472 +470 @@
         ENDDO
         tauref(:) = tauref(:) * odpref / pplev(:,1)
+
 c--------------------------------------------------
-c  Parameters of the opacity perturbation
+c Get parameters of the opacity perturbation
 c--------------------------------------------------
-
-      iaer=1  !!!! PROVISOIRE !!!!
+        iaer=1  ! just change dust
 
         write(*,*) "Add a local storm ?"
@@ -518 +516 @@
               write(*,*) " reffstorm = ",reffstorm
 
+!! LOOP: modify opacity
       DO ig=1,ngrid
 
-
-c---------------------------------------
-c        distance to the center:
-c-----------------------------------------
-
+      !! distance to the center:
       ray(ig)=SQRT((lati(ig)*180./pi-latloc)**2 +
      &          (long(ig)*180./pi -lonloc)**2)
 
-
       !! transition factor for storm
-      !! increase factor ray diff for steepness
+      !! factor is hardcoded -- increase it to steepen
       yeah = (TANH(2.+(radloc-ray(ig))*10.)+1.)/2.
 
-c-------------------------------------------------
-c           Tau's new map:
-c------------------------------------------
-
-      tauuser(ig)=max(tauref(ig) * pplev(ig,1) /odpref  , 
-     &          taulocref * yeah)
-
-c---------------------------------------------------------
-c           compute l_top
-c----------------------------------------------------------
-
+      !! new opacity field
+      !! -- add an opacity set to taulocref
+      !! -- the additional reference opacity will
+      !!      thus be taulocref*odpref/pplev
+      tauuser(ig)=max( tauref(ig) * pplev(ig,1) /odpref , 
+     &          taulocref * yeah )
+
+      !! compute l_top
           DO l=1,nlayer
             zalt(ig,l) = LOG( pplev(ig,1)/pplev(ig,l) )
@@ -552 +543 @@
           ENDDO
 
-c---------------------------------------------------------
-c           change reffrad if ever needed
-c----------------------------------------------------------
-
+     !! change reffrad if ever needed
       IF (reffstorm .gt. 0.) THEN
           DO l=1,nlayer
@@ -565 +553 @@
       ENDIF
 
-
-c---------------------------------------------------------------------------
-      ENDDO !! boucle sur ngridmx
-c---------------------------------------------------------------------------
-c----------------------------------------------------------------------------
-c    compute int_factor
-c---------------------------------------------------------------------------
-
+      ENDDO
+!! END LOOP
+
+      !! compute perturbation in each layer (equation 8 in Spiga et al. JGR 2013)
       DO ig=1,ngrid
           int_factor(ig)=0.
@@ -585 +569 @@
           ENDDO
           DO l=1, nlayer
-
-c-------------------------------------------------------------------------
-c Mass mixing ratio perturbation due to the local dust storm in each
-c layer
-c-------------------------------------------------------------------------
+          !! Mass mixing ratio perturbation due to local dust storm in each layer
           more_dust(ig,l,1)=
      &                     (tauuser(ig)-(tauref(ig)
@@ -603 +583 @@
       ENDDO
 
-c--------------------------------------------------------------------------------------
-c   quantity of dust which is added at the first time step 
-c     in dynamical core.
-c--------------------------------------------------------------------------------------
-
-        DO l=1, l_top
-          zdqnorm(:,l,1) = more_dust(:,l,1)
-          zdqnorm(:,l,2) = more_dust(:,l,2)
+      !! quantity of dust for each layer with the addition of the perturbation
+      DO l=1, l_top
           pq(:,l,igcm_dust_mass)= pq(:,l,igcm_dust_mass)
-     .            + zdqnorm(:,l,1)
+     .            + more_dust(:,l,1)
           pq(:,l,igcm_dust_number)= pq(:,l,igcm_dust_number)
-     .            + zdqnorm(:,l,2)
-        ENDDO
-        ENDIF
-        tauref(:)=0.
-      ENDIF  !firstcall
+     .            + more_dust(:,l,2)
+      ENDDO
+      ENDIF !! IF (localstorm)
+      tauref(:)=0.
+      ENDIF !! IF (firstcall)
 #endif
 
@@ -643 +617 @@
         ENDDO
 
-      ENDIF
+      ENDIF ! IF (freedust)
 
 c     Opacity computation
@@ -666 +640 @@
           DO l=1,nlayer
             DO ig=1,ngrid
+#ifdef DUSTSTORM
+      !! recalculate opacity because storm perturbation has been added
+      IF (firstcall) THEN
+              aerosol(ig,l,iaer) =
+     &          (  0.75 * QREFvis3d(ig,l,iaer) /
+     &          ( rho_dust * reffrad(ig,l,iaer) )  ) *
+     &          pq(ig,l,igcm_dust_mass) *
+     &          ( pplev(ig,l) - pplev(ig,l+1) ) / g
+      ENDIF
+#endif
               tauref(ig) = tauref(ig) + 
      &                    aerosol(ig,l,iaerdust(iaer))