Changeset 3851 for LMDZ6/branches/LMDZ-tracers/libf/dyn3dmem

Timestamp:

Feb 22, 2021, 12:44:07 PM (3 years ago)

Author:

dcugnet

Message:

Update the branch to the current trunk.

Location:

LMDZ6/branches/LMDZ-tracers

Files:

• . (modified) (1 prop)
• libf/dyn3dmem/conf_gcm.F90 (modified) (1 diff)
• libf/dyn3dmem/dynredem_loc.F90 (modified) (2 diffs)
• libf/dyn3dmem/groupe_loc.F (modified) (1 diff)
• libf/dyn3dmem/guide_loc_mod.F90 (modified) (4 diffs)
• libf/dyn3dmem/qminimum_loc.F (modified) (12 diffs, 1 prop)
• libf/dyn3dmem/vlsplt_loc.F (modified) (13 diffs, 1 prop)
• libf/dyn3dmem/vlspltqs_loc.F (modified) (9 diffs, 1 prop)

LMDZ6/branches/LMDZ-tracers/libf/dyn3dmem/conf_gcm.F90

@@ -345 +345 @@
   ngroup=3
   CALL getin('ngroup',ngroup)
+  if (mod(iim, 2**ngroup) /= 0) &
+       call abort_gcm("conf_gcm", 'iim must be multiple of 2**ngroup', 1)
+  if (2**ngroup > jjm + 1) &
+       call abort_gcm("conf_gcm", '2**ngroup must be <= jjm + 1', 1)
 
   ! mode_top_bound : fields towards which sponge relaxation will be done:

LMDZ6/branches/LMDZ-tracers/libf/dyn3dmem/dynredem_loc.F90

@@ -11 +11 @@
   USE infotrac
   USE netcdf, ONLY: NF90_CREATE, NF90_DEF_DIM, NF90_INQ_VARID, NF90_GLOBAL,    &
-                    NF90_CLOSE,  NF90_PUT_ATT, NF90_UNLIMITED, NF90_CLOBBER
+                    NF90_CLOSE,  NF90_PUT_ATT, NF90_UNLIMITED, NF90_CLOBBER,   &
+                    NF90_64BIT_OFFSET
   USE dynredem_mod, ONLY: cre_var, put_var, err, modname, fil
   USE comvert_mod, ONLY: ap,bp,aps,bps,presnivs,pseudoalt,pa,preff, &
@@ -110 +111 @@
 
 !--- File creation
-  CALL err(NF90_CREATE(fichnom,NF90_CLOBBER,nid))
+  CALL err(NF90_CREATE(fichnom,IOR(NF90_CLOBBER,NF90_64BIT_OFFSET),nid))
 
 !--- Some global attributes

LMDZ6/branches/LMDZ-tracers/libf/dyn3dmem/groupe_loc.F

@@ -40 +40 @@
 c$OMP THREADPRIVATE(firstcall)
 
-      data firstcall/.true./
       integer ijb,ije,jjb,jje
       
-      if (firstcall) then
-         if(mod(iim,2**ngroup).ne.0) stop'probleme du nombre ede point'
-         firstcall=.false.
-      endif
-
 c   Champs 1D
 

LMDZ6/branches/LMDZ-tracers/libf/dyn3dmem/guide_loc_mod.F90

@@ -27 +27 @@
   LOGICAL, PRIVATE, SAVE  :: invert_p,invert_y,ini_anal
   LOGICAL, PRIVATE, SAVE  :: guide_2D,guide_sav,guide_modele
+!FC
+  LOGICAL, PRIVATE, SAVE  :: convert_Pa
   
   REAL, PRIVATE, SAVE     :: tau_min_u,tau_max_u
@@ -152 +154 @@
         guide_plevs=1
     ENDIF
+!FC
+    CALL getpar('convert_Pa',.true.,convert_Pa,'Convert Pressure levels in Pa')
     ! Fin raccord
     CALL getpar('ini_anal',.false.,ini_anal,'Etat initial = analyse')
@@ -1695 +1699 @@
              status=NF_GET_VARA_REAL(ncidpl,varidpl,1,nlevnc,apnc)
 #endif
-             apnc=apnc*100.! conversion en Pascals
+!FC Pour les corrections la pression est deja en Pascals on commente la ligne ci-dessous
+             IF(convert_Pa) apnc=apnc*100.! conversion en Pascals
              bpnc(:)=0.
          ENDIF
@@ -2193 +2198 @@
 ! ----------------------------------------------
 ! Ouverture du fichier
-        ierr=NF_CREATE("guide_ins.nc",NF_CLOBBER,nid)
+        ierr=NF_CREATE("guide_ins.nc",IOR(NF_CLOBBER,NF_64BIT_OFFSET),nid)
 ! Definition des dimensions
         ierr=NF_DEF_DIM(nid,"LONU",iip1,id_lonu) 

LMDZ6/branches/LMDZ-tracers/libf/dyn3dmem/qminimum_loc.F

@@ +1 @@
+!
+!     $Id$
+!
       SUBROUTINE qminimum_loc( q,nqtot,deltap )
       USE parallel_lmdz
-      USE infotrac, ONLY: ok_isotopes,ntraciso,iqiso,ok_iso_verif
+      USE infotrac, ONLY: ok_isotopes,ntraciso,iqiso,ok_iso_verif,             &
+     &   ratiomin,qperemin ! CRisi 23nov2020
       IMPLICIT none
 c
@@ -9 +13 @@
       include "dimensions.h"
       include "paramet.h"
+      include "iniprint.h"
 c
       INTEGER nqtot ! CRisi: on remplace nq par nqtot
@@ -49 +54 @@
 c
 
-        !write(*,*) 'qminimum 52: entree'
+        !write(lunout,*) 'qminimum 52: entree'
         if (ok_iso_verif) then
            call check_isotopes(q,ij_begin,ij_end,'qminimum 52')   
@@ -60 +65 @@
       q_follow(ijb:ije,:,1:2)=q(ijb:ije,:,1:2)  
 
-      !write(*,*) 'qminimum 57'
+      !write(lunout,*) 'qminimum 57'
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)       
       DO 1000 k = 1, llm
@@ -85 +90 @@
 c le defaut en prennant de l'eau vapeur de la couche au-dessous.
 c
-      !write(*,*) 'qminimum 81'
+      !write(lunout,*) 'qminimum 81'
       iq = iq_vap
 c
@@ -113 +118 @@
 c doit imprimer un message d'avertissement (saturation possible).
 c
-      !write(*,*) 'qminimum 106'
+      !write(lunout,*) 'qminimum 106'
       nb_pump=0
 c$OMP DO SCHEDULE(STATIC)
@@ -135 +140 @@
       ENDIF
 
-      !write(*,*) 'qminimum 128'
+      !write(lunout,*) 'qminimum 128'
       if (ok_isotopes) then
+              !write(lunout,*) 'qminimum 140'
       ! CRisi: traiter de même les traceurs d'eau
       ! Mais il faut les prendre à l'envers pour essayer de conserver la
@@ -144 +150 @@
       ! rien ici et on croise les doigts pour que ça ne soit pas trop
       ! génant
+      ! en fait, si, c'est genant quand les isotopes doivent eux même transporter des
+      ! traceurs -> apporter aussi un peu d'isotopes... Combien?
+      ! Essayer tnat/2 = -500 permil? C'est déjà mieux que -1000
+      ! permil...
+      ! pb: que faire pour les traceurs?
+c$OMP DO SCHEDULE(STATIC)      
       DO i = ijb, ije
         if (zx_pump(i).gt.0.0) then
@@ -149 +161 @@
         endif !if (zx_pump(i).gt.0.0) then
       enddo !DO i = ijb, ije  
+c$OMP END DO
 
       ! 2) transfert de vap vers les couches plus hautes
-      !write(*,*) 'qminimum 139'
+      !write(lunout,*) 'qminimum 158'
       do k=2,llm
+c$OMP DO SCHEDULE(STATIC)      
         DO i = ijb, ije
           if (zx_defau_diag(i,k,iq_vap).gt.0.0) then             
-              ! on ajoute la vapeur en k              
-              do ixt=1,ntraciso
+              ! on ajoute la vapeur en k     
+!              write(lunout,*) 'i,k,q_follow(i,k-1,iq_vap)=',
+!     :                 i,k,q_follow(i,k-1,iq_vap)          
+              if (q_follow(i,k-1,iq_vap).lt.qperemin) then
+                write(lunout,*) 'tmp qmin: on stoppe'
+                write(lunout,*) 'zx_pump(i)=',zx_pump(i)
+                write(lunout,*) 'q_follow(i,:,iq_vap)=',
+     :                   q_follow(i,:,iq_vap)
+                write(lunout,*) 'k=',k
+                call abort_gcm("qminimum","not enough vapor",1)
+              endif  
+            do ixt=1,ntraciso
+!                write(lunout,*) 'qmin 168: ixt=',ixt
+!                write(lunout,*) 'q(i,k,iqiso(ixt,iq_vap)=',
+!     :             q(i,k,iqiso(ixt,iq_vap))
+!                write(lunout,*) 'zx_defau_diag(i,k,iq_vap)=',
+!     :                  zx_defau_diag(i,k,iq_vap)
+!                write(lunout,*) 'q(i,k-1,iqiso(ixt,iq_vap)=',
+!     :                   q(i,k-1,iqiso(ixt,iq_vap))     
+
                q(i,k,iqiso(ixt,iq_vap))=q(i,k,iqiso(ixt,iq_vap))
      :              +zx_defau_diag(i,k,iq_vap)
@@ -207 +239 @@
           endif !if (zx_defau_diag(i,k,iq_vap).gt.0.0) then
         enddo !DO i = 1, ip1jmp1        
-       enddo !do k=2,llm
+c$OMP END DO
+        enddo !do k=2,llm
 
         if (ok_iso_verif) then
@@ -217 +250 @@
         !write(*,*) 'qminimum 164'
         do k=1,llm
+c$OMP DO SCHEDULE(STATIC)
         DO i = ijb, ije
           if (zx_defau_diag(i,k,iq_liq).gt.0.0) then
@@ -235 +269 @@
      :               -zx_defau_diag(i,k,iq_liq)
           endif !if (zx_defau_diag(i,k,iq_vap).gt.0.0) then
-        enddo !DO i = 1, ip1jmp1
+        enddo !DO i = ijb, ije
+c$OMP END DO        
        enddo !do k=2,llm  
 

LMDZ6/branches/LMDZ-tracers/libf/dyn3dmem/vlsplt_loc.F

@@ -14 +14 @@
 c   --------------------------------------------------------------------
       USE parallel_lmdz
-      USE infotrac, ONLY : nqtot,nqfils,nqdesc,iqfils ! CRisi
+      USE infotrac, ONLY : nqtot,nqfils,nqdesc,iqfils, ! CRisi                 &
+     &                     qperemin,masseqmin,ratiomin ! MVals et CRisi
       IMPLICIT NONE
 c
@@ -329 +330 @@
 ! Il faut faire ça avant d'avoir mis à jour q et masse
 
-      !write(*,*) 'vlsplt 326: iq,ijb_x,nqfils(iq)=',iq,ijb_x,nqfils(iq) 
-
-      if (nqfils(iq).gt.0) then  
+       if (nqfils(iq).gt.0) then
        do ifils=1,nqdesc(iq)
+       !do ifils=1,nqfils(iq) ! modif C Risi 22nov2020
+        ! attention: comme Ratio est utilisé comme q dans l'appel
+        ! recursif, il doit contenir à lui seul tous les indices de tous
+        ! les descendants!
          iq2=iqfils(ifils,iq)
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
@@ -339 +342 @@
            ! On a besoin de q et masse seulement entre ijb et ije. On ne
            ! les calcule donc que de ijb à ije
-           masse(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq)
-           Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
+           !MVals: veiller a ce qu'on n'ait pas de denominateur nul
+           masse(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin)
+           if (q(ij,l,iq).gt.qperemin) then ! modif 13 nov 2020
+             Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
+           else
+             Ratio(ij,l,iq2)=ratiomin
+           endif
           enddo   
          enddo
@@ -352 +360 @@
 ! end CRisi
 
-      !write(*,*) 'vlsplt 360: iq,ijb_x=',iq,ijb_x
 
 c   calcul des tENDances
@@ -358 +365 @@
       DO l=1,llm
          DO ij=ijb+1,ije
-            new_m=masse(ij,l,iq)+u_m(ij-1,l)-u_m(ij,l)
+            !MVals: veiller a ce qu'on n'ait pas de denominateur nul
+            new_m=max(masse(ij,l,iq)+u_m(ij-1,l)-u_m(ij,l),masseqmin)
             q(ij,l,iq)=(q(ij,l,iq)*masse(ij,l,iq)+
      &        u_mq(ij-1,l)-u_mq(ij,l))
@@ -371 +379 @@
       ENDDO
 c$OMP END DO NOWAIT
-      !write(*,*) 'vlsplt 380: iq,ijb_x=',iq,ijb_x
 
 ! retablir les fils en rapport de melange par rapport a l'air:
@@ -378 +385 @@
       if (nqfils(iq).gt.0) then  
        do ifils=1,nqdesc(iq)
+       !do ifils=1,nqfils(iq) ! modif C Risi 22nov2020
          iq2=iqfils(ifils,iq)  
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)   
@@ -414 +422 @@
 c   --------------------------------------------------------------------
       USE parallel_lmdz
-      USE infotrac, ONLY : nqtot,nqfils,nqdesc,iqfils ! CRisi
+      USE infotrac, ONLY : nqtot,nqfils,nqdesc,iqfils, ! CRisi                 &
+     &                     qperemin,masseqmin,ratiomin ! MVals et CRisi    
       USE comconst_mod, ONLY: pi
       IMPLICIT NONE
@@ -468 +477 @@
       DATA temps0,temps1,temps2,temps3,temps4,temps5/0.,0.,0.,0.,0.,0./
       INTEGER ijb,ije
+      INTEGER ijbm,ijem
 
       ijb=ij_begin-2*iip1
@@ -726 +736 @@
       ijb=ij_begin-2*iip1
       ije=ij_end+2*iip1
+      ijbm=ij_begin-iip1
+      ijem=ij_end+iip1
       if (pole_nord) ijb=ij_begin
-      if (pole_sud)  ije=ij_end
-   
+      if (pole_sud)  ije=ij_end  
+      if (pole_nord) ijbm=ij_begin
+      if (pole_sud)  ijem=ij_end
+
       if (nqfils(iq).gt.0) then  
        do ifils=1,nqdesc(iq)
+       !do ifils=1,nqfils(iq) ! modif C Risi 22nov2020
          iq2=iqfils(ifils,iq)
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK) 
          DO l=1,llm
+          ! modif des bornes: CRisi 16 nov 2020
+          ! d'abord masse avec bornes corrigées
+          DO ij=ijbm,ijem 
+           !MVals: veiller a ce qu'on n'ait pas de denominateur nul
+           masse(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin)
+          enddo !DO ij=ijbm,ijem 
+
+          ! ensuite Ratio avec anciennes bornes
          DO ij=ijb,ije
-           masse(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq)
-           Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)     
-          enddo   
-         enddo
+           !MVals: veiller a ce qu'on n'ait pas de denominateur nul
+           if (q(ij,l,iq).gt.qperemin) then ! modif 13 nov 2020
+             Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
+           else
+             Ratio(ij,l,iq2)=ratiomin  
+           endif     
+          enddo !DO ij=ijbm,ijem  
+         enddo !DO l=1,llm
 c$OMP END DO NOWAIT
         enddo !do ifils=1,nqdesc(iq)
@@ -868 +895 @@
       USE parallel_lmdz
       USE vlz_mod
-      USE infotrac, ONLY : nqtot,nqfils,nqdesc,iqfils ! CRisi  
+      USE infotrac, ONLY : nqtot,nqfils,nqdesc,iqfils, ! CRisi                 &
+     &                     qperemin,masseqmin,ratiomin ! MVals et CRisi
+      
       IMPLICIT NONE
 c
@@ -1084 +1113 @@
                 lorig(ij,l)=lorig(ij,l)-1
              ENDIF
+             ! CRisi 24nov2020: ajout d'un message d'erreur clair au lieu d'un plantage
+             ! pour seg fault
+             if (lorig(ij,l).eq.0) then
+                call abort_gcm("vlz in vlsplt_loc",
+     :           "unfixable violation of CFL",1)
+             endif
              morig(ij,l)=masse(ij,lorig(ij,l),iq)
              qorig(ij,l)=q(ij,lorig(ij,l),iq)
@@ -1127 +1162 @@
       if (nqfils(iq).gt.0) then  
        do ifils=1,nqdesc(iq)
+       !do ifils=1,nqfils(iq) ! modif C Risi 22 nov 2020
          iq2=iqfils(ifils,iq)
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
          DO l=1,llm
           DO ij=ijb,ije
-           masse(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq)
-           Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
+           !MVals: veiller a ce qu'on n'ait pas de denominateur nul
+           masse(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin)
+           if (q(ij,l,iq).gt.qperemin) then
+             Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
+           else
+             Ratio(ij,l,iq2)=ratiomin
+           endif
            !wq(ij,l,iq2)=wq(ij,l,iq) ! correction bug le 15mai2015
            w(ij,l,iq2)=wq(ij,l,iq)

LMDZ6/branches/LMDZ-tracers/libf/dyn3dmem/vlspltqs_loc.F

@@ +1 @@
+!
+!     $Id$
+!
       SUBROUTINE vlxqs_loc(q,pente_max,masse,u_m,qsat,ijb_x,ije_x,iq)
 c
@@ -9 +12 @@
 c   --------------------------------------------------------------------
       USE parallel_lmdz
-      USE infotrac, ONLY : nqtot,nqfils,nqdesc,iqfils ! CRisi  
+      USE infotrac, ONLY : nqtot,nqfils,nqdesc,iqfils, ! CRisi                 &
+     &                     qperemin,masseqmin,ratiomin ! MVals et CRisi
       IMPLICIT NONE
 c
@@ -342 +346 @@
          DO l=1,llm
           DO ij=ijb,ije
-           masse(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq)
-           Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)      
+           !MVals: veiller a ce qu'on n'ait pas de denominateur nul
+           masse(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin)
+           if (q(ij,l,iq).gt.qperemin) then ! modif 13 nov 2020
+             Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
+           else
+             Ratio(ij,l,iq2)=ratiomin
+           endif 
           enddo   
          enddo
@@ -362 +371 @@
       DO l=1,llm
          DO ij=ijb+1,ije
-            new_m=masse(ij,l,iq)+u_m(ij-1,l)-u_m(ij,l)
+            !MVals: veiller a ce qu'on n'ait pas de denominateur nul
+            new_m=max(masse(ij,l,iq)+u_m(ij-1,l)-u_m(ij,l),masseqmin)
             q(ij,l,iq)=(q(ij,l,iq)*masse(ij,l,iq)+
      &      u_mq(ij-1,l)-u_mq(ij,l))
@@ -416 +426 @@
 c   --------------------------------------------------------------------
       USE parallel_lmdz
-      USE infotrac, ONLY : nqtot,nqfils,nqdesc,iqfils ! CRisi  
+      USE infotrac, ONLY : nqtot,nqfils,nqdesc,iqfils, ! CRisi                 &
+     &                     qperemin,masseqmin,ratiomin ! MVals et CRisi
       USE comconst_mod, ONLY: pi
       IMPLICIT NONE
@@ -423 +434 @@
       include "paramet.h"
       include "comgeom.h"
+      include "iniprint.h"  
 c
 c
@@ -464 +476 @@
       DATA first/.true./
       INTEGER ijb,ije
+      INTEGER ijbm,ijem
 
       ijb=ij_begin-2*iip1
@@ -724 +737 @@
       ijb=ij_begin-2*iip1
       ije=ij_end+2*iip1
+      ijbm=ij_begin-iip1
+      ijem=ij_end+iip1
       if (pole_nord) ijb=ij_begin
       if (pole_sud)  ije=ij_end  
-
+      if (pole_nord) ijbm=ij_begin
+      if (pole_sud)  ijem=ij_end
+
+      !write(lunout,*) 'vlspltqs 737: iq,ijb,ije=',iq,ijb,ije
+      !write(lunout,*) 'ij_begin,ij_end=',ij_begin,ij_end
+      !write(lunout,*) 'pole_nord,pole_sud=',pole_nord,pole_sud
       if (nqfils(iq).gt.0) then  
        do ifils=1,nqdesc(iq)
@@ -732 +752 @@
 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)   
          DO l=1,llm
+          ! modif des bornes: CRisi 16 nov 2020
+          ! d'abord masse avec bornes corrigées
+          DO ij=ijbm,ijem 
+           !MVals: veiller a ce qu'on n'ait pas de denominateur nul
+           masse(ij,l,iq2)=max(masse(ij,l,iq)*q(ij,l,iq),masseqmin)
+          enddo !DO ij=ijbm,ijem 
+
+          ! ensuite Ratio avec anciennes bornes
           DO ij=ijb,ije
-           masse(ij,l,iq2)=masse(ij,l,iq)*q(ij,l,iq)
-           Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)   
-          enddo   
-         enddo
+           !MVals: veiller a ce qu'on n'ait pas de denominateur nul
+           !write(lunout,*) 'ij,l,q(ij,l,iq)=',ij,l,q(ij,l,iq)
+           if (q(ij,l,iq).gt.qperemin) then ! modif 13 nov 2020
+             Ratio(ij,l,iq2)=q(ij,l,iq2)/q(ij,l,iq)
+           else
+             Ratio(ij,l,iq2)=ratiomin   
+           endif
+          enddo !DO ij=ijbm,ijem
+         enddo !DO l=1,llm
 c$OMP END DO NOWAIT
         enddo !do ifils=1,nqdesc(iq)
         do ifils=1,nqfils(iq)
          iq2=iqfils(ifils,iq)
+         !write(lunout,*) 'vly: appel recursiv vly iq2=',iq2
          call vly_loc(Ratio,pente_max,masse,qbyv,iq2)
         enddo !do ifils=1,nqfils(iq)