Changeset 117 for trunk/mars/libf/phymars

Timestamp:

May 10, 2011, 6:33:29 PM (14 years ago)

Author:

aslmd

Message:

M 116 mesoscale/LMD_MM_MARS/SRC/WRFV2/phys/module_lmd_driver.F
* CORRECTION MAJEURE DE BUG : PSFC et TRACEURS PAS BIEN MIS A JOUR *

• pdq au lieu de pdq*dt (dt = pas de temps dynamique)
• idem pour pdpsrf

IMPACT SURTOUT SUR LES SIMULATIONS LONGUES ET LES SIMULATIONS TACHE DE POUSSIERE (stage Julien)
TEST A EFFECTUER POUR LES SIMULATIONS CYCLE DE L'EAU

M 116 mars/libf/phymars/newsedim.F
M 116 mars/README
Correction de deux bugs dans newsedim.F: 1-e(-x) trop faible et endif mal place

M 116 mesoscale/LMDZ.MARS.new/myGCM/start.nc
A 0 mesoscale/LMDZ.MARS.new/myGCM/DEFS_JB/SP51_HR_dq_r3.0n0.5a0.5_MY26_TM
A 0 mesoscale/LMDZ.MARS.new/myGCM/DEFS_JB/SP51_HR_dq_r3.0n0.5a0.5_MY26_TM/run.def
A 0 mesoscale/LMDZ.MARS.new/myGCM/DEFS_JB/SP51_HR_dq_r3.0n0.5a0.5_MY26_TM/traceur.def
A 0 mesoscale/LMDZ.MARS.new/myGCM/DEFS_JB/SP51_HR_dq_r3.0n0.5a0.5_MY26_TM/startfi72.nc
A 0 mesoscale/LMDZ.MARS.new/myGCM/DEFS_JB/SP51_HR_dq_r3.0n0.5a0.5_MY26_TM/start72.nc
A 0 mesoscale/LMDZ.MARS.new/myGCM/DEFS_JB/SP51_HR_dq_r3.0n0.5a0.5_MY26_TM/callphys.def
M 116 mesoscale/LMDZ.MARS.new/myGCM/traceur.def
M 116 mesoscale/LMDZ.MARS.new/myGCM/startfi.nc
M 116 mesoscale/LMDZ.MARS.new/myGCM/callphys.def
M 116 mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/compile
M 116 mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/launch_gcm
M 116 mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/run_mcd_3days
Nouvelle base de donnees d'etats initiaux sans les nuages radiativement actifs

A 0 mesoscale/PLOT/SPEC/LES/turb_period_psfc.pro
A 0 mesoscale/PLOT/SPEC/LES/turb_inc.pro.new
M 116 mesoscale/PLOT/SPEC/GW/gravitwave2.pro
M 116 mesoscale/PLOT/SPEC/GW/gravitwaveprof.pro
Petites MAJ routines graphiques

File:

• trunk/mars/libf/phymars/newsedim.F (modified) (4 diffs)

trunk/mars/libf/phymars/newsedim.F

@@ -67 +67 @@
 
 
-
 c    ** un petit test de coherence
 c       --------------------------
@@ -138 +137 @@
              k=0
 
+            w(ig,l) = 0. !! JF+AS ajout initialisation
 c **************************************************************
 c            Simple Method
-             w(ig,l) =
-     &       (1- exp(-dztop*g/(r*pt(ig,l))))*pplev(ig,l) / g
-cc           write(*,*) 'OK simple method l,w =', l, w(ig,l)
-cc           write(*,*) 'OK simple method dztop =', dztop
+
+cc             w(ig,l) =
+cc     &       (1.- exp(-dztop*g/(r*pt(ig,l))))*pplev(ig,l) / g
+cccc           write(*,*) 'OK simple method l,w =', l, w(ig,l)
+cccc           write(*,*) 'OK simple method dztop =', dztop
+
+             w(ig,l) = 1. - exp(-dztop*g/(r*pt(ig,l)))
+             !!! Diagnostic: JF. Fix: AS. Date: 05/11 
+             !!! Probleme arrondi avec la quantite ci-dessus
+             !!! ---> vaut 0 pour -dztop*g/(r*pt(ig,l)) trop petit
+             !!! ---> dans ce cas on utilise le developpement limite !
+             !!! ---> exp(-x) = 1 - x lorsque x --> 0 avec une erreur de x^2 / 2            
+             IF ( w(ig,l) .eq. 0. ) THEN
+                w(ig,l) = ( dztop*g/(r*pt(ig,l)) ) * pplev(ig,l) / g
+             ELSE
+                w(ig,l) = w(ig,l) * pplev(ig,l) / g  
+             ENDIF
+
+
 c **************************************************************
 cccc         Complex method :
-            if (dztop.gt.epaisseur(ig,l)) then
+            if (dztop.gt.epaisseur(ig,l)) then                !!!if on traverse plus d'une couche
 cccc            Cas ou on "epuise" la couche l : On calcule le flux
 cccc            Venant de dessus en tenant compte de la variation de Vstokes
-
+c **************************************************************
                Ep= epaisseur(ig,l)
                Stra= traversee(ig,l)
@@ -159 +174 @@
                enddo 
                Ep = Ep - epaisseur(ig,l+k)
-             end if
-             ptop=pplev(ig,l+k)*exp(-(dztop-Ep)*g/(r*pt(ig,l+k)))
-             w(ig,l) = (pplev(ig,l) -Ptop)/g
-c
+             !ptop=pplev(ig,l+k)*exp(-(dztop-Ep)*g/(r*pt(ig,l+k)))
+
+             !!! JF+AS 05/11 Probleme arrondi potentiel, meme solution que ci-dessus
+             ptop=exp(-(dztop-Ep)*g/(r*pt(ig,l+k)))
+             IF ( ptop .eq. 1. ) THEN
+                PRINT*, 'exposant trop petit ', ig, l
+                ptop=pplev(ig,l+k) * ( 1. - (dztop-Ep)*g/(r*pt(ig,l+k)))
+             ELSE
+                ptop=pplev(ig,l+k) * ptop 
+             ENDIF
+             w(ig,l) = (pplev(ig,l) - Ptop)/g
+
+            endif                !!!!!if complex method
+
+
 cc           write(*,*) 'OK new    method l,w =', l, w(ig,l)
 cc           write(*,*) 'OK new    method dztop =', dztop
@@ -171 +197 @@
 cc            if(l.eq.8)write(*,*)'l=8,k, w',k, w(1,l)
 c **************************************************************
+
+
         end do
       end do