Changeset 3800 for LMDZ6/trunk/libf/dyn3dmem/vlsplt_loc.F

Timestamp:

Jan 15, 2021, 6:10:56 PM (3 years ago)

Author:

Laurent Fairhead

Message:

Modifications nécessaires pour les isotopes
CRisi

File:

• LMDZ6/trunk/libf/dyn3dmem/vlsplt_loc.F (modified) (13 diffs, 1 prop)

LMDZ6/trunk/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)