source: LMDZ5/trunk/libf/dyn3dmem/qminimum_loc.F @ 2526

Last change on this file since 2526 was 2286, checked in by Ehouarn Millour, 9 years ago
  • fix in dyn3d the array out of bound issue (in qminimum) that was only corrected in dyn3dmem in rev 2285.
  • comment out many invasive debug writes.

EM

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
File size: 8.5 KB
Line 
1      SUBROUTINE qminimum_loc( q,nqtot,deltap )
2      USE parallel_lmdz
3      USE infotrac, ONLY: ok_isotopes,ntraciso,iqiso,ok_iso_verif
4      IMPLICIT none
5c
6c  -- Objet : Traiter les valeurs trop petites (meme negatives)
7c             pour l'eau vapeur et l'eau liquide
8c
9#include "dimensions.h"
10#include "paramet.h"
11#include "comvert.h"
12c
13      INTEGER nqtot ! CRisi: on remplace nq par nqtot
14      REAL q(ijb_u:ije_u,llm,nqtot), deltap(ijb_u:ije_u,llm)
15c
16      INTEGER iq_vap, iq_liq
17      PARAMETER ( iq_vap = 1 ) ! indice pour l'eau vapeur
18      PARAMETER ( iq_liq = 2 ) ! indice pour l'eau liquide
19      REAL seuil_vap, seuil_liq
20      PARAMETER ( seuil_vap = 1.0e-10 ) ! seuil pour l'eau vapeur
21      PARAMETER ( seuil_liq = 1.0e-11 ) ! seuil pour l'eau liquide
22c
23c  NB. ....( Il est souhaitable mais non obligatoire que les valeurs des
24c            parametres seuil_vap, seuil_liq soient pareilles a celles
25c            qui  sont utilisees dans la routine    ADDFI       )
26c     .................................................................
27c
28      INTEGER i, k, iq
29      REAL zx_defau, zx_abc, zx_pump(ijb_u:ije_u), pompe
30
31      real zx_defau_diag(ijb_u:ije_u,llm,2)
32      real q_follow(ijb_u:ije_u,llm,2)
33c
34      REAL SSUM
35      EXTERNAL SSUM
36c
37      INTEGER imprim
38      SAVE imprim
39      DATA imprim /0/
40c$OMP THREADPRIVATE(imprim)
41      INTEGER ijb,ije
42      INTEGER Index_pump(ij_end-ij_begin+1)
43      INTEGER nb_pump
44      INTEGER ixt
45      INTEGER iso_verif_noNaN_nostop
46c
47c Quand l'eau liquide est trop petite (ou negative), on prend
48c l'eau vapeur de la meme couche et la convertit en eau liquide
49c (sans changer la temperature !)
50c
51
52        !write(*,*) 'qminimum 52: entree'
53        if (ok_iso_verif) then
54           call check_isotopes(q,ij_begin,ij_end,'qminimum 52')   
55        endif !if (ok_iso_verif) then     
56
57      ijb=ij_begin
58      ije=ij_end
59
60      zx_defau_diag(ijb:ije,:,:)=0.0
61      q_follow(ijb:ije,:,1:2)=q(ijb:ije,:,1:2) 
62
63      !write(*,*) 'qminimum 57'
64c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)       
65      DO 1000 k = 1, llm
66      DO 1040 i = ijb, ije
67            if (seuil_liq - q(i,k,iq_liq) .gt. 0.d0 ) then
68
69              if (ok_isotopes) then
70                 zx_defau_diag(i,k,iq_liq)=AMAX1
71     :               ( seuil_liq - q(i,k,iq_liq), 0.0 )
72              endif !if (ok_isotopes) then
73
74               q(i,k,iq_vap) = q(i,k,iq_vap) + q(i,k,iq_liq) - seuil_liq
75               q(i,k,iq_liq) = seuil_liq
76            endif
77 1040 CONTINUE
78 1000 CONTINUE
79c$OMP END DO NOWAIT
80c$OMP BARRIER
81c --->  SYNCHRO OPENMP ICI
82
83
84c
85c Quand l'eau vapeur est trop faible (ou negative), on complete
86c le defaut en prennant de l'eau vapeur de la couche au-dessous.
87c
88      !write(*,*) 'qminimum 81'
89      iq = iq_vap
90c
91      DO k = llm, 2, -1
92ccc      zx_abc = dpres(k) / dpres(k-1)
93c$OMP DO SCHEDULE(STATIC)
94      DO i = ijb, ije
95
96         if ( seuil_vap - q(i,k,iq) .gt. 0.d0 ) then
97
98            if (ok_isotopes) then
99              zx_defau_diag(i,k,iq)=AMAX1( seuil_vap - q(i,k,iq), 0.0 )
100            endif !if (ok_isotopes) then
101
102            q(i,k-1,iq) =  q(i,k-1,iq) - ( seuil_vap - q(i,k,iq) ) *
103     &           deltap(i,k) / deltap(i,k-1)
104            q(i,k,iq)   =  seuil_vap 
105
106         endif
107      ENDDO
108c$OMP END DO NOWAIT
109      ENDDO
110c$OMP BARRIER
111
112c
113c Quand il s'agit de la premiere couche au-dessus du sol, on
114c doit imprimer un message d'avertissement (saturation possible).
115c
116      !write(*,*) 'qminimum 106'
117      nb_pump=0
118c$OMP DO SCHEDULE(STATIC)
119      DO i = ijb, ije
120         zx_pump(i) = AMAX1( 0.0, seuil_vap - q(i,1,iq) )
121         q(i,1,iq)  = AMAX1( q(i,1,iq), seuil_vap )
122         IF (zx_pump(i) > 0.0) THEN
123            nb_pump = nb_pump+1
124            Index_pump(nb_pump)=i
125         ENDIF
126      ENDDO
127c$OMP END DO 
128!      pompe = SSUM(ije-ijb+1,zx_pump(ijb),1)
129
130      IF (imprim.LE.100 .AND. nb_pump .GT. 0 ) THEN
131         PRINT *, 'ATT!:on pompe de l eau au sol'
132         DO i = 1, nb_pump
133               imprim = imprim + 1
134               PRINT*,'  en ',index_pump(i),zx_pump(index_pump(i))
135         ENDDO
136      ENDIF
137
138      !write(*,*) 'qminimum 128'
139      if (ok_isotopes) then
140      ! CRisi: traiter de même les traceurs d'eau
141      ! Mais il faut les prendre à l'envers pour essayer de conserver la
142      ! masse.
143      ! 1) pompage dans le sol 
144      ! On suppose que ce pompage se fait sans isotopes -> on ne modifie
145      ! rien ici et on croise les doigts pour que ça ne soit pas trop
146      ! génant
147      DO i = ijb, ije
148        if (zx_pump(i).gt.0.0) then
149          q_follow(i,1,iq_vap)=q_follow(i,1,iq_vap)+zx_pump(i)
150        endif !if (zx_pump(i).gt.0.0) then
151      enddo !DO i = ijb, ije 
152
153      ! 2) transfert de vap vers les couches plus hautes
154      !write(*,*) 'qminimum 139'
155      do k=2,llm
156        DO i = ijb, ije
157          if (zx_defau_diag(i,k,iq_vap).gt.0.0) then             
158              ! on ajoute la vapeur en k             
159              do ixt=1,ntraciso
160               q(i,k,iqiso(ixt,iq_vap))=q(i,k,iqiso(ixt,iq_vap))
161     :              +zx_defau_diag(i,k,iq_vap)
162     :              *q(i,k-1,iqiso(ixt,iq_vap))/q_follow(i,k-1,iq_vap)
163               
164              if (ok_iso_verif) then
165                if (iso_verif_noNaN_nostop(q(i,k,iqiso(ixt,iq_vap)),
166     :                   'qminimum 155').eq.1) then
167                   write(*,*) 'i,k,ixt=',i,k,ixt
168                   write(*,*) 'q_follow(i,k-1,iq_vap)=',
169     :                   q_follow(i,k-1,iq_vap)
170                   write(*,*) 'q(i,k,iqiso(ixt,iq_vap))=',
171     :                   q(i,k,iqiso(ixt,iq_vap))
172                   write(*,*) 'zx_defau_diag(i,k,iq_vap)=',
173     :                   zx_defau_diag(i,k,iq_vap)
174                   write(*,*) 'q(i,k-1,iqiso(ixt,iq_vap))=',
175     :                   q(i,k-1,iqiso(ixt,iq_vap))
176                   stop
177                endif
178              endif
179
180              ! et on la retranche en k-1
181               q(i,k-1,iqiso(ixt,iq_vap))=q(i,k-1,iqiso(ixt,iq_vap))
182     :              -zx_defau_diag(i,k,iq_vap)
183     :              *deltap(i,k)/deltap(i,k-1)
184     :              *q(i,k-1,iqiso(ixt,iq_vap))/q_follow(i,k-1,iq_vap)
185
186               if (ok_iso_verif) then
187                if (iso_verif_noNaN_nostop(q(i,k-1,iqiso(ixt,iq_vap)),
188     :                   'qminimum 175').eq.1) then
189                   write(*,*) 'k,i,ixt=',k,i,ixt
190                   write(*,*) 'q_follow(i,k-1,iq_vap)=',
191     :                   q_follow(i,k-1,iq_vap)
192                   write(*,*) 'q(i,k,iqiso(ixt,iq_vap))=',
193     :                   q(i,k,iqiso(ixt,iq_vap))
194                   write(*,*) 'zx_defau_diag(i,k,iq_vap)=',
195     :                   zx_defau_diag(i,k,iq_vap)
196                   write(*,*) 'q(i,k-1,iqiso(ixt,iq_vap))=',
197     :                   q(i,k-1,iqiso(ixt,iq_vap))
198                   stop
199                endif
200              endif
201
202              enddo !do ixt=1,niso
203              q_follow(i,k,iq_vap)=   q_follow(i,k,iq_vap)
204     :               +zx_defau_diag(i,k,iq_vap)
205              q_follow(i,k-1,iq_vap)=   q_follow(i,k-1,iq_vap)
206     :               -zx_defau_diag(i,k,iq_vap)
207     :              *deltap(i,k)/deltap(i,k-1)
208          endif !if (zx_defau_diag(i,k,iq_vap).gt.0.0) then
209        enddo !DO i = 1, ip1jmp1       
210       enddo !do k=2,llm
211
212        if (ok_iso_verif) then
213           call check_isotopes(q,ijb,ije,'qminimum 168')
214        endif !if (ok_iso_verif) then
215       
216     
217        ! 3) transfert d'eau de la vapeur au liquide
218        !write(*,*) 'qminimum 164'
219        do k=1,llm
220        DO i = ijb, ije
221          if (zx_defau_diag(i,k,iq_liq).gt.0.0) then
222
223              ! on ajoute eau liquide en k en k             
224              do ixt=1,ntraciso
225               q(i,k,iqiso(ixt,iq_liq))=q(i,k,iqiso(ixt,iq_liq))
226     :              +zx_defau_diag(i,k,iq_liq)
227     :              *q(i,k,iqiso(ixt,iq_vap))/q_follow(i,k,iq_vap)
228              ! et on la retranche à la vapeur en k
229               q(i,k,iqiso(ixt,iq_vap))=q(i,k,iqiso(ixt,iq_vap))
230     :              -zx_defau_diag(i,k,iq_liq)
231     :              *q(i,k,iqiso(ixt,iq_vap))/q_follow(i,k,iq_vap)   
232              enddo !do ixt=1,niso
233              q_follow(i,k,iq_liq)=   q_follow(i,k,iq_liq)
234     :               +zx_defau_diag(i,k,iq_liq)
235              q_follow(i,k,iq_vap)=   q_follow(i,k,iq_vap)
236     :               -zx_defau_diag(i,k,iq_liq)
237          endif !if (zx_defau_diag(i,k,iq_vap).gt.0.0) then
238        enddo !DO i = 1, ip1jmp1
239       enddo !do k=2,llm 
240
241        if (ok_iso_verif) then
242           call check_isotopes(q,ijb,ije,'qminimum 197')
243        endif !if (ok_iso_verif) then
244
245      endif !if (ok_isotopes) then
246      !write(*,*) 'qminimum 188'
247c
248      RETURN
249      END
Note: See TracBrowser for help on using the repository browser.