Changeset 3700 for trunk/LMDZ.TITAN

Timestamp:

Mar 28, 2025, 4:17:11 PM (3 months ago)

Author:

emoisan

Message:

Titan:

• Add safeguards not to enter microphisics when there are no clouds
• Removed a comment

EMo

Location:

trunk/LMDZ.TITAN/libf

Files:

• muphytitan/mm_haze.F90 (modified) (1 diff)
• phytitan/get_haze_and_cloud_opacity.F90 (modified) (1 diff)
• phytitan/optci.F90 (modified) (1 diff)
• phytitan/optcv.F90 (modified) (1 diff)

trunk/LMDZ.TITAN/libf/muphytitan/mm_haze.F90 ¶

@@ -792 +792 @@
       call EXIT(11)
 #endif
-#ifdef MESOSCALE
-      WRITE(*,*) "Aerosol production level out of model, no production of simple aerosols"
-#endif
     ENDIF
 

trunk/LMDZ.TITAN/libf/phytitan/get_haze_and_cloud_opacity.F90 ¶

@@ -200 +200 @@
   !-------------------------
 
-  rinit = 1.e-9                ! fin = rinit*step**(33-1) = 1.e-5
-  if(CTYPE.eq.0) rinit = 1.e-7 ! fin = rinit*step**(33-1) = 1.e-3
-
-  ! No alpha(k) !
-  ! Because it is taken into account in the optical tables ...
-  rc = (m3/m0)**(1./3.)
-  
-  step = 10.**(1./8.)
-  idx = int(log(rc/rinit) / log(step)) + 1
-
-  ! Spherical (cloud) :
-  !--------------------
-  IF(CTYPE.EQ.0) THEN
-
-    if(abs(rc_s(1)-rinit)/rc_s(1).gt.1.e-5 .or. &
-       abs(rc_s(nrad)-rinit*step**32)/rc_s(nrad).gt.1.e-5) then
-        write(*,*) 'rinit          = ',rinit
-        write(*,*) 'rc_s(1)        = ',rc_s(1)
-        write(*,*) 'rinit*step**32 = ',rinit*step**32
-        write(*,*) 'rc_s(nrad)     = ',rc_s(nrad)
-        stop
-    endif
-
-    if(rc.le.rc_s(1)) then
-      wbar   = wbar_s(ich,1)
-      gbar   = gbar_s(ich,1)
-      tauext = sext_s(ich,1)*(rc/rc_s(1))**3
-    endif
-
-    if(rc.ge.rc_s(nrad)) then
-      wbar   = wbar_s(ich,nrad)
-      gbar   = gbar_s(ich,nrad)
-      tauext = sext_s(ich,nrad)*(rc/rc_s(nrad))**3.
-    endif
-
-    if(rc.gt.rc_s(1) .and. rc.lt.rc_s(nrad)) then
-      xfact  = (log(rc)-log(rc_s(idx))) / (log(rc_s(idx+1))-log(rc_s(idx)))
-      tauext = log(sext_s(ich,idx))*(1.-xfact) + log(sext_s(ich,idx+1))*xfact 
-      tauext = exp(tauext)
-      wbar   = wbar_s(ich,idx)*(1.-xfact) + wbar_s(ich,idx+1)*xfact 
-      gbar   = gbar_s(ich,idx)*(1.-xfact) + gbar_s(ich,idx+1)*xfact 
-    endif
-
-  ENDIF ! (CTYPE.EQ.0)
-
-  ! Spherical (haze) :
-  !-------------------
-  IF(CTYPE.EQ.5) THEN
-
-    if(abs(rc_as(1)-rinit)/rc_as(1).gt.1.e-5 .or. &
-       abs(rc_as(nrad)-rinit*step**32)/rc_as(nrad).gt.1.e-5) then
-        write(*,*) 'rinit          = ',rinit
-        write(*,*) 'rc_as(1)       = ',rc_as(1)
-        write(*,*) 'rinit*step**32 = ',rinit*step**32
-        write(*,*) 'rc_as(nrad)    = ',rc_as(nrad)
-        stop
-    endif
-
-    if(rc.le.rc_as(1)) then
-      wbar   = wbar_as(ich,1)
-      gbar   = gbar_as(ich,1)
-      tauext = sext_as(ich,1)*(rc/rc_as(1))**3
-    endif
-
-    if(rc.ge.rc_as(nrad)) then
-      wbar   = wbar_as(ich,nrad)
-      gbar   = gbar_as(ich,nrad)
-      tauext = sext_as(ich,nrad)*(rc/rc_as(nrad))**3.
-    endif
-
-    if(rc.gt.rc_as(1) .and. rc.lt.rc_as(nrad)) then
-      xfact  = (log(rc)-log(rc_as(idx))) / (log(rc_as(idx+1))-log(rc_as(idx)))
-      tauext = log(sext_as(ich,idx))*(1.-xfact) + log(sext_as(ich,idx+1))*xfact
-      tauext = exp(tauext)
-      wbar   = wbar_as(ich,idx)*(1.-xfact) + wbar_as(ich,idx+1)*xfact
-      gbar   = gbar_as(ich,idx)*(1.-xfact) + gbar_as(ich,idx+1)*xfact
-    endif
-  
-  ENDIF ! (CTYPE.EQ.5)
-
-  ! Fractal (haze) :
-  !-----------------
-  IF(CTYPE.NE.0 .AND. CTYPE.NE.5) THEN
-
-    if(abs(rc_f(1)-rinit)/rc_f(1).gt.1.e-5 .or. &
-       abs(rc_f(nrad)-rinit*step**32)/rc_f(nrad).gt.1.e-5) then
-        write(*,*) 'rinit          = ',rinit
-        write(*,*) 'rc_f(1)        = ',rc_f(1)
-        write(*,*) 'rinit*step**32 = ',rinit*step**32
-        write(*,*) 'rc_f(nrad)     = ',rc_f(nrad)
-        stop
-    endif
-
-    if(rc.le.rc_f(1)) then
-      wbar   =  wbar_f(ich,1)
-      gbar   =  gbar_f(ich,1)
-      tauext = sext_f(ich,1)*(rc/rc_f(1))**3
-    endif
-
-    if(rc.ge.rc_f(nrad)) then
-      wbar   = wbar_f(ich,nrad)
-      gbar   = gbar_f(ich,nrad)
-      tauext = sext_f(ich,nrad)*(rc/rc_f(nrad))**3.
-    endif
-    
-    if(rc.gt.rc_f(1) .and. rc.lt.rc_f(nrad)) then
-      xfact  = (log(rc)-log(rc_f(idx))) / (log(rc_f(idx+1))-log(rc_f(idx)))
-      tauext = log(sext_f(ich,idx))*(1.-xfact) + log(sext_f(ich,idx+1))*xfact 
-      tauext = exp(tauext)
-      wbar   = wbar_f(ich,idx)*(1.-xfact) + wbar_f(ich,idx+1)*xfact 
-      gbar   = gbar_f(ich,idx)*(1.-xfact) + gbar_f(ich,idx+1)*xfact 
-    endif
-
-  ENDIF ! (CTYPE.NE.0 .AND. CTYPE.NE.5)
-
-  ! Sanity check :
-  !---------------
-  IF (CTYPE.LT.0 .AND. CTYPE.GT.5) THEN
-    write(*,*) 'CTYPE < 0 ET > 5  = ',CTYPE
-    call abort
-  ENDIF
-
-  tauext = tauext * m0
+  IF((m0 .gt. tiny(m0)) .and. (m3 .gt.tiny(m3))) THEN  ! if there are clouds !emoisan tests
+      rinit = 1.e-9                ! fin = rinit*step**(33-1) = 1.e-5
+      if(CTYPE.eq.0) rinit = 1.e-7 ! fin = rinit*step**(33-1) = 1.e-3
+    
+      ! No alpha(k) !
+      ! Because it is taken into account in the optical tables ...
+      rc = (m3/m0)**(1./3.)
+      
+      step = 10.**(1./8.)
+      idx = int(log(rc/rinit) / log(step)) + 1
+    
+      ! Spherical (cloud) :
+      !--------------------
+      IF(CTYPE.EQ.0) THEN
+    
+        if(abs(rc_s(1)-rinit)/rc_s(1).gt.1.e-5 .or. &
+           abs(rc_s(nrad)-rinit*step**32)/rc_s(nrad).gt.1.e-5) then
+            write(*,*) 'rinit          = ',rinit
+            write(*,*) 'rc_s(1)        = ',rc_s(1)
+            write(*,*) 'rinit*step**32 = ',rinit*step**32
+            write(*,*) 'rc_s(nrad)     = ',rc_s(nrad)
+            stop
+        endif
+    
+        if(rc.le.rc_s(1)) then
+          wbar   = wbar_s(ich,1)
+          gbar   = gbar_s(ich,1)
+          tauext = sext_s(ich,1)*(rc/rc_s(1))**3
+        endif
+    
+        if(rc.ge.rc_s(nrad)) then
+          wbar   = wbar_s(ich,nrad)
+          gbar   = gbar_s(ich,nrad)
+          tauext = sext_s(ich,nrad)*(rc/rc_s(nrad))**3.
+        endif
+    
+        if(rc.gt.rc_s(1) .and. rc.lt.rc_s(nrad)) then
+          xfact  = (log(rc)-log(rc_s(idx))) / (log(rc_s(idx+1))-log(rc_s(idx)))
+          tauext = log(sext_s(ich,idx))*(1.-xfact) + log(sext_s(ich,idx+1))*xfact 
+          tauext = exp(tauext)
+          wbar   = wbar_s(ich,idx)*(1.-xfact) + wbar_s(ich,idx+1)*xfact 
+          gbar   = gbar_s(ich,idx)*(1.-xfact) + gbar_s(ich,idx+1)*xfact 
+        endif
+    
+      ENDIF ! (CTYPE.EQ.0)
+    
+      ! Spherical (haze) :
+      !-------------------
+      IF(CTYPE.EQ.5) THEN
+    
+        if(abs(rc_as(1)-rinit)/rc_as(1).gt.1.e-5 .or. &
+           abs(rc_as(nrad)-rinit*step**32)/rc_as(nrad).gt.1.e-5) then
+            write(*,*) 'rinit          = ',rinit
+            write(*,*) 'rc_as(1)       = ',rc_as(1)
+            write(*,*) 'rinit*step**32 = ',rinit*step**32
+            write(*,*) 'rc_as(nrad)    = ',rc_as(nrad)
+            stop
+        endif
+    
+        if(rc.le.rc_as(1)) then
+          wbar   = wbar_as(ich,1)
+          gbar   = gbar_as(ich,1)
+          tauext = sext_as(ich,1)*(rc/rc_as(1))**3
+        endif
+    
+        if(rc.ge.rc_as(nrad)) then
+          wbar   = wbar_as(ich,nrad)
+          gbar   = gbar_as(ich,nrad)
+          tauext = sext_as(ich,nrad)*(rc/rc_as(nrad))**3.
+        endif
+    
+        if(rc.gt.rc_as(1) .and. rc.lt.rc_as(nrad)) then
+          xfact  = (log(rc)-log(rc_as(idx))) / (log(rc_as(idx+1))-log(rc_as(idx)))
+          tauext = log(sext_as(ich,idx))*(1.-xfact) + log(sext_as(ich,idx+1))*xfact
+          tauext = exp(tauext)
+          wbar   = wbar_as(ich,idx)*(1.-xfact) + wbar_as(ich,idx+1)*xfact
+          gbar   = gbar_as(ich,idx)*(1.-xfact) + gbar_as(ich,idx+1)*xfact
+        endif
+      
+      ENDIF ! (CTYPE.EQ.5)
+    
+      ! Fractal (haze) :
+      !-----------------
+      IF(CTYPE.NE.0 .AND. CTYPE.NE.5) THEN
+    
+        if(abs(rc_f(1)-rinit)/rc_f(1).gt.1.e-5 .or. &
+           abs(rc_f(nrad)-rinit*step**32)/rc_f(nrad).gt.1.e-5) then
+            write(*,*) 'rinit          = ',rinit
+            write(*,*) 'rc_f(1)        = ',rc_f(1)
+            write(*,*) 'rinit*step**32 = ',rinit*step**32
+            write(*,*) 'rc_f(nrad)     = ',rc_f(nrad)
+            stop
+        endif
+    
+        if(rc.le.rc_f(1)) then
+          wbar   =  wbar_f(ich,1)
+          gbar   =  gbar_f(ich,1)
+          tauext = sext_f(ich,1)*(rc/rc_f(1))**3
+        endif
+    
+        if(rc.ge.rc_f(nrad)) then
+          wbar   = wbar_f(ich,nrad)
+          gbar   = gbar_f(ich,nrad)
+          tauext = sext_f(ich,nrad)*(rc/rc_f(nrad))**3.
+        endif
+        
+        if(rc.gt.rc_f(1) .and. rc.lt.rc_f(nrad)) then
+          xfact  = (log(rc)-log(rc_f(idx))) / (log(rc_f(idx+1))-log(rc_f(idx)))
+          tauext = log(sext_f(ich,idx))*(1.-xfact) + log(sext_f(ich,idx+1))*xfact 
+          tauext = exp(tauext)
+          wbar   = wbar_f(ich,idx)*(1.-xfact) + wbar_f(ich,idx+1)*xfact 
+          gbar   = gbar_f(ich,idx)*(1.-xfact) + gbar_f(ich,idx+1)*xfact 
+        endif
+    
+      ENDIF ! (CTYPE.NE.0 .AND. CTYPE.NE.5)
+    
+      ! Sanity check :
+      !---------------
+      IF (CTYPE.LT.0 .AND. CTYPE.GT.5) THEN
+        write(*,*) 'CTYPE < 0 ET > 5  = ',CTYPE
+        call abort
+      ENDIF
+    
+      tauext = tauext * m0
+  ENDIF !if there are clouds !emoisan tests
+  !else just return the initialize values of tauext, wbar and gbar =0
 
   return

trunk/LMDZ.TITAN/libf/phytitan/optci.F90 ¶

@@ -313 +313 @@
                
                ! For small dropplets, opacity of nucleus dominates...
-               dtau_cld = (dtau_cld*m3ccn + dtau_cld*m3cld) / (m3ccn + m3cld)
-               ssa_cld(nw) = (ssa_ccn(nw)*m3ccn + ssa_cld(nw)*m3cld) / (m3ccn + m3cld)
+               IF ((m3ccn + m3cld) .le. tiny(m3ccn)) THEN !no cloud !emoisan tests
+                  dtau_cld = 0.
+                  ssa_cld(nw) = 0.
+               ELSE
+                  dtau_cld = (dtau_cld*m3ccn + dtau_cld*m3cld) / (m3ccn + m3cld)
+                  ssa_cld(nw) = (ssa_ccn(nw)*m3ccn + ssa_cld(nw)*m3cld) / (m3ccn + m3cld)
+               ENDIF
                
                ! Tuning of optical properties for clouds :

trunk/LMDZ.TITAN/libf/phytitan/optcv.F90 ¶

@@ -337 +337 @@
 
                ! For small dropplets, opacity of nucleus dominates
-               dtau_cld = (dtau_cld*m3ccn + dtau_cld*m3cld) / (m3ccn + m3cld)
-               ssa_cld(nw) = (ssa_ccn(nw)*m3ccn + ssa_cld(nw)*m3cld) / (m3ccn + m3cld)
+               IF ((m3ccn + m3cld) .le. tiny(m3ccn)) THEN !no cloud !emoisan tests
+                  dtau_cld = 0.
+                  ssa_cld(nw) = 0.
+               ELSE
+                  dtau_cld = (dtau_cld*m3ccn + dtau_cld*m3cld) / (m3ccn + m3cld)
+                  ssa_cld(nw) = (ssa_ccn(nw)*m3ccn + ssa_cld(nw)*m3cld) / (m3ccn + m3cld)
+               ENDIF
 
                ! Tuning of optical properties for clouds :