source: LMDZ6/branches/contrails/libf/phylmdiso/isotopes_verif_mod.F90 @ 5456

Last change on this file since 5456 was 4982, checked in by crisi, 6 months ago

suppress isotope_params.def + update physiq_mod + proof of concept of 3rd dimension with reevap routine

  • Property svn:executable set to *
File size: 79.7 KB
Line 
1
2#ifdef ISOVERIF
3! $Id: $
4
5MODULE isotopes_verif_mod
6!use isotopes_mod, ONLY:
7!#ifdef ISOTRAC
8!   USE isotrac_mod, ONLY: nzone
9!#endif
10USE infotrac_phy, ONLY: ntraciso=>ntiso, niso, itZonIso, nzone
11implicit none
12save
13
14! variables de vérifications
15real errmax ! erreur maximale en absolu.
16parameter (errmax=1e-8)
17real errmax_sol ! erreur maximale en absolu.
18parameter (errmax_sol=5e-7)
19      real errmaxrel ! erreur maximale en relatif autorisée
20!      parameter (errmaxrel=1e10)
21      parameter (errmaxrel=1e-3)
22      real borne ! valeur maximale que n'importe quelle variable peut
23                 ! atteindre, en val abs; utile pour vérif que pas NAN
24      parameter (borne=1e19)
25      real, save ::  deltalim ! deltalim est le maximum de deltaD qu'on
26                             ! autorise dans la vapeur, au-dela, en suppose que c'est abérrant.
27                             ! dans le liquide, on autorise deltalim*faccond.
28!$OMP THREADPRIVATE(deltalim)
29!      parameter (deltalim=1e10)
30!      parameter (deltalim=300.0)
31       ! maintenant défini dans iso.def
32
33       real, save :: deltalimtrac ! max de deltaD dans les traceurs, défini dans iso.def
34!$OMP THREADPRIVATE(deltalimtrac)
35
36      real, save ::  deltalim_snow ! deltalim est le maximum de deltaD qu'on
37                             ! autorise dans la neige, au-dela, en suppose que c'est abérrant.
38!$OMP THREADPRIVATE(deltalim_snow)
39!      parameter (deltalim_snow=500.0) ! initalisé dans iso_init
40
41        real, save ::  deltaDmin
42!$OMP THREADPRIVATE(deltaDmin)
43!       parameter (deltaDmin=-900.0)
44    ! maintentant, défini dans iso.def
45
46      real, save ::  o17excess_bas,o17excess_haut ! bornes inf et sup de l'O17-excess
47!      parameter(o17excess_bas=-200.0,o17excess_haut=120)
48!$OMP THREADPRIVATE(o17excess_bas,o17excess_haut)
49      integer nlevmaxO17
50!$OMP THREADPRIVATE(nlevmaxO17)
51
52      logical, save ::  O17_verif
53!$OMP THREADPRIVATE(O17_verif)
54!      parameter (O17_verif=.true.)
55
56        real, save ::  dexcess_min,dexcess_max
57!$OMP THREADPRIVATE(dexcess_min,dexcess_max)
58
59      real faccond  ! dans le liquide, on autorise R(deltalim)*faccond.
60      parameter (faccond=1.1)
61     
62!      logical bidouille_anti_divergence ! si true, alors on fait un
63!                                        ! rappel régulier des xt4 vers les q pour
64!                                        !éviter accumulation d'érreurs et  divergences
65!      parameter (bidouille_anti_divergence=.true.)
66!      parameter (bidouille_anti_divergence=.false.)
67    ! bidouille_anti_divergence a été déplacé dans wateriso2.h et est lue à l'execution
68
69       
70        real deltaDfaible ! deltaD en dessous duquel la vapeur est tellement faible
71        !que on peut accepter que la remise à l'équilibre du sol avec
72        ! cette vapeur donne des deltaDevap aberrants.
73        parameter (deltaDfaible=-300)
74        real deltaDfaible_lax ! deltaD en dessous duquel la vapeur est tellement faible
75        !que on peut accepter que la remise à l'équilibre du sol avec
76        ! cette vapeur donne des deltaDevap aberrants.
77        parameter (deltaDfaible_lax=-180)
78
79        real faible_evap ! faible évaporation: on est plus laxiste
80        !pour les deltaD aberrant dans le cas de l'évap venant d'orchidee
81        parameter (faible_evap=3.0)
82
83        real Tmin_verif
84        parameter (Tmin_verif=5.0) ! en K
85        real Tmax_verif
86        parameter (Tmax_verif=1000.0) ! en K
87
88
89! subroutines de vérifications génériques, à ne pas modifier
90
91 
92CONTAINS
93
94        SUBROUTINE iso_verif_init()
95        use ioipsl_getin_p_mod, ONLY : getin_p
96        use isotopes_mod, ONLY: iso_O17, iso_O18, iso_HDO
97        implicit none
98
99        write(*,*) 'iso_verif_init 99: entree'
100        deltalim=300.0
101        deltalim_snow=500.0
102        deltaDmin=-900.0
103        deltalimtrac=2000.0
104        O17_verif=.false.
105        o17excess_bas=-200.0
106        o17excess_haut=120.0
107        dexcess_min=-100.0
108        dexcess_max=1000.0
109
110call getin_p('deltalim',deltalim)
111deltalim_snow=deltalim+200.0
112call getin_p('deltaDmin',deltaDmin)
113call getin_p('deltalimtrac',deltalimtrac)
114
115write(*,*) 'iso_O17,iso_O18,iso_HDO=',iso_O17,iso_O18,iso_HDO
116if ((iso_O17.gt.0).and.(iso_O18.gt.0)) then
117  call getin_p('O17_verif',O17_verif) 
118  call getin_p('o17excess_bas',o17excess_bas) 
119  call getin_p('o17excess_haut',o17excess_haut) 
120  call getin_p('nlevmaxO17',nlevmaxO17) 
121endif
122if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then
123  call getin_p('dexcess_min',dexcess_min) 
124  call getin_p('dexcess_max',dexcess_max) 
125endif
126
127write(*,*) 'deltalim=',deltalim
128write(*,*) 'deltaDmin=',deltaDmin     
129write(*,*) 'deltalimtrac=',deltalimtrac 
130write(*,*) 'O17_verif=',O17_verif
131write(*,*) 'o17excess_bas=',o17excess_bas
132write(*,*) 'o17excess_haut=',o17excess_haut 
133write(*,*) 'dexcess_min=',dexcess_min
134write(*,*) 'dexcess_max=',dexcess_max
135
136        end SUBROUTINE iso_verif_init
137
138      subroutine iso_verif_egalite(a,b,err_msg)
139        implicit none
140        ! compare a et b. Si pas egal à errmax près, on affiche message
141        ! d'erreur et stoppe
142
143        ! input:
144        real a, b
145        character*(*) err_msg ! message d''erreur à afficher
146
147        ! local
148        !integer iso_verif_egalite_choix_nostop
149
150
151        if (iso_verif_egalite_choix_nostop(a,b,err_msg, &
152     &           errmax,errmaxrel).eq.1) then
153                stop
154        endif
155       
156#ifdef ISOVERIF
157#else
158        write(*,*) 'err_msg=',err_msg,': pourquoi verif?'
159        stop
160#endif           
161
162        return
163        end subroutine iso_verif_egalite
164
165        !*****************
166
167        function iso_verif_egalite_nostop(a,b,err_msg)
168        implicit none
169        ! compare a et b. Si pas egal à errmax près, on affiche message
170        ! d'erreur et stoppe
171
172        ! input:
173        real a, b
174        character*(*) err_msg ! message d''erreur à afficher
175        !ouptut
176        integer iso_verif_egalite_nostop
177        ! local
178        !integer iso_verif_egalite_choix_nostop
179
180
181        iso_verif_egalite_nostop=iso_verif_egalite_choix_nostop &
182     &          (a,b,err_msg,errmax,errmaxrel)
183
184#ifdef ISOVERIF
185#else
186        write(*,*) 'err_msg=',err_msg,': pourquoi verif?'
187        stop
188#endif                 
189
190        return
191        end function iso_verif_egalite_nostop
192
193
194
195        !************************************
196
197        subroutine iso_verif_aberrant(R,err_msg)
198        use isotopes_mod, ONLY: ridicule, iso_HDO
199        implicit none
200        ! si le rapprot iso R est plus grand que deltaDlim, on affiche message
201        ! d'erreur et stoppe
202
203        ! input:
204        real R
205        character*(*) err_msg ! message d''erreur à afficher
206        !real deltaD
207        !integer iso_verif_aberrant_choix_nostop
208
209
210        if (iso_verif_aberrant_choix_nostop(R,1.0,ridicule, &
211     &           deltalim,err_msg).eq.1) then
212             stop
213        endif
214
215#ifdef ISOVERIF
216        if (.not.(iso_HDO.gt.0)) then   
217         write(*,*) 'iso86: err_msg=',err_msg,': pourquoi verif?'
218          stop
219        endif
220#else       
221        write(*,*) 'iso 90: err_msg=',err_msg,': pourquoi verif?'
222        stop
223#endif                   
224               
225        end subroutine iso_verif_aberrant
226
227        subroutine iso_verif_aberrant_encadre(R,err_msg)
228        use isotopes_mod, ONLY: ridicule, iso_HDO
229        implicit none
230        ! si le rapprot iso R est plus grand que deltaDlim, on affiche message
231        ! d'erreur et stoppe
232
233        ! input:
234        real R
235        character*(*) err_msg ! message d''erreur à afficher
236        !real deltaD
237
238        !integer iso_verif_aberrant_enc_choix_nostop
239
240
241        if (iso_verif_aberrant_enc_choix_nostop(R,1.0,ridicule, &
242     &           deltalim,err_msg).eq.1) then
243             write(*,*) 'deltaD=',deltaD(R)
244             call abort_physic('isotopes_verif_mod > iso_verif_aberrant_encadre',err_msg,1)
245             !stop             
246        endif
247
248#ifdef ISOVERIF
249        if (.not.(iso_HDO.gt.0)) then   
250         write(*,*) 'iso86: err_msg=',err_msg,': pourquoi verif?'
251          stop
252        endif
253#else       
254        write(*,*) 'iso 90: err_msg=',err_msg,': pourquoi verif?'
255        stop
256#endif                   
257       
258        end subroutine iso_verif_aberrant_encadre
259
260        !************************************
261
262        subroutine iso_verif_aberrant_choix(xt,q,qmin,deltaDmax,err_msg)
263        use isotopes_mod, ONLY: iso_HDO
264        implicit none
265        ! si le rapprot iso R est plus grand que deltaDlim, on affiche message
266        ! d'erreur et stoppe
267
268        ! input:
269        real xt,q,qmin,deltaDmax
270        character*(*) err_msg ! message d''erreur à afficher
271        !real deltaD
272
273        ! locals
274        !integer iso_verif_aberrant_choix_nostop
275
276
277        if (iso_verif_aberrant_choix_nostop(xt,q,qmin, &
278     &           deltaDmax,err_msg).eq.1) then
279             stop
280        endif
281
282#ifdef ISOVERIF
283        if (.not.(iso_HDO.gt.0)) then   
284         write(*,*) 'iso122: err_msg=',err_msg,': pourquoi verif?'
285          stop
286        endif
287#else
288        write(*,*) 'iso126: err_msg=',err_msg,': pourquoi verif?'
289        stop
290#endif                   
291       
292        end subroutine iso_verif_aberrant_choix
293
294        !************************************
295
296        function iso_verif_aberrant_nostop(R,err_msg)
297        use isotopes_mod, ONLY: ridicule,iso_HDO
298        implicit none
299        ! si le rapprot iso R est plus grand que deltaDlim, on affiche message
300        ! d'erreur et stoppe
301
302        ! input:
303        real R
304        character*(*) err_msg ! message d''erreur à afficher
305        integer iso_verif_aberrant_nostop ! output: 1 si erreur, 0 sinon
306        !real deltaD
307
308        ! locals
309        !integer iso_verif_aberrant_choix_nostop
310
311        iso_verif_aberrant_nostop=iso_verif_aberrant_choix_nostop &
312     &           (R,1.0,ridicule,deltalim,err_msg)
313
314#ifdef ISOVERIF
315        if (.not.(iso_HDO.gt.0)) then   
316         write(*,*) 'iso156: err_msg=',err_msg,': pourquoi verif?'
317          stop
318        endif
319#else
320        write(*,*) 'iso160: err_msg=',err_msg,': pourquoi verif?'
321        stop
322#endif           
323       
324        return
325        end function iso_verif_aberrant_nostop
326
327        function iso_verif_aberrant_enc_nostop(R,err_msg)
328        use isotopes_mod, ONLY: ridicule,iso_HDO
329        implicit none
330        ! si le rapprot iso R est plus grand que deltaDlim, on affiche message
331        ! d'erreur et stoppe
332
333        ! input:
334        real R
335        character*(*) err_msg ! message d''erreur à afficher
336        integer iso_verif_aberrant_enc_nostop ! output: 1 si erreur, 0 sinon
337        !real deltaD
338
339        ! locals
340        !integer iso_verif_aberrant_enc_choix_nostop
341
342        iso_verif_aberrant_enc_nostop= &
343     &           iso_verif_aberrant_enc_choix_nostop &
344     &           (R,1.0,ridicule,deltalim,err_msg)
345
346#ifdef ISOVERIF
347        if (.not.(iso_HDO.gt.0)) then   
348         write(*,*) 'iso156: err_msg=',err_msg,': pourquoi verif?'
349          stop
350        endif
351#else
352        write(*,*) 'iso160: err_msg=',err_msg,': pourquoi verif?'
353        stop
354#endif                   
355        return
356        end function iso_verif_aberrant_enc_nostop
357
358        !************************************
359
360        function iso_verif_aberrant_choix_nostop(xt,q,   &
361     &            qmin,deltaDmax,err_msg)
362
363        use isotopes_mod, ONLY: iso_HDO
364        implicit none
365        ! si le rapprot iso R est plus grand que deltaDlim, on affiche message
366        ! d'erreur et stoppe
367
368        ! input:
369        real xt,q,qmin,deltaDmax
370        character*(*) err_msg ! message d''erreur à afficher
371        ! output
372        integer iso_verif_aberrant_choix_nostop
373        ! locals
374        !real deltaD
375        !integer iso_verif_noNaN_nostop       
376
377
378        iso_verif_aberrant_choix_nostop=0
379
380#ifdef ISOVERIF       
381        if (iso_verif_noNaN_nostop(q,err_msg).eq.1) then
382            write(*,*) 'q=',q
383            iso_verif_aberrant_choix_nostop=1
384        endif     
385        if (iso_verif_noNaN_nostop(xt,err_msg).eq.1) then
386            write(*,*) 'xt=',xt
387            iso_verif_aberrant_choix_nostop=1
388        endif
389#endif
390
391        if (q.gt.qmin) then
392        if ((deltaD(xt/q).gt.deltaDmax).or. &
393     &          (deltaD(xt/q).lt.-borne).or. &
394     &          (xt.lt.-borne).or. &
395     &          (xt.gt.borne)) then                 
396            write(*,*) 'erreur detectee par '// &
397     &             'iso_verif_aberrant_choix_nostop:'
398            write(*,*) err_msg
399            write(*,*) 'q,deltaD=',q,deltaD(xt/q)
400            write(*,*) 'deltaDmax=',deltaDmax
401            iso_verif_aberrant_choix_nostop=1
402            if (abs(xt-q)/q.lt.errmax) then
403                write(*,*) 'attention, n''a-t-on pas confondu'// &
404     &           ' iso_HDO et iso_eau dans l''appel à la verif?'
405            endif
406        endif
407        endif
408
409#ifdef ISOVERIF
410        if (.not.(iso_HDO.gt.0)) then   
411         write(*,*) 'iso205: err_msg=',err_msg,': pourquoi verif?'
412          stop
413        endif
414#else
415        write(*,*) 'iso 209: err_msg=',err_msg,': pourquoi verif?'
416        stop
417#endif                   
418       
419        return
420        end function iso_verif_aberrant_choix_nostop
421
422        function iso_verif_aberrant_enc_choix_nostop(xt,q,   &
423     &            qmin,deltaDmax,err_msg)
424        use isotopes_mod, ONLY: iso_HDO
425        implicit none
426        ! si le rapprot iso R est plus grand que deltaDlim, on affiche message
427        ! d'erreur et stoppe
428
429        ! input:
430        real xt,q,qmin,deltaDmax
431        character*(*) err_msg ! message d''erreur à afficher
432        ! output
433        integer iso_verif_aberrant_enc_choix_nostop
434        ! locals
435        !real deltaD
436
437        iso_verif_aberrant_enc_choix_nostop=0
438        if (q.gt.qmin) then
439        if ((deltaD(xt/q).gt.deltaDmax).or. &
440     &          (deltaD(xt/q).lt.deltaDmin)) then                 
441            write(*,*) 'erreur detectee par '// &
442     &             'iso_verif_aberrant_choix_nostop:'
443            write(*,*) err_msg
444            write(*,*) 'q,deltaD=',q,deltaD(xt/q)
445            iso_verif_aberrant_enc_choix_nostop=1
446            if (abs(xt-q)/q.lt.errmax) then
447                write(*,*) 'attention, n''a-t-on pas confondu'// &
448     &           ' iso_HDO et iso_eau dans l''appel à la verif?'
449            endif
450        endif
451        endif
452
453#ifdef ISOVERIF
454        if (.not.(iso_HDO.gt.0)) then   
455         write(*,*) 'iso205: err_msg=',err_msg,': pourquoi verif?'
456          stop
457        endif
458#else
459        write(*,*) 'iso 209: err_msg=',err_msg,': pourquoi verif?'
460        stop
461#endif                   
462       
463        return
464        end function iso_verif_aberrant_enc_choix_nostop
465
466        !*******************
467
468        subroutine iso_verif_aberrant_o17(R17,R18,err_msg)
469        implicit none
470        ! si l'O17-excess est aberrant, on affiche un message
471        ! d'erreur et stoppe
472
473        ! input:
474        real R17,R18
475        character*(*) err_msg ! message d''erreur à afficher
476        !real o17excess
477
478        ! locals
479        !integer iso_verif_aberrant_o17_nostop
480
481!        write(*,*) 'O17_verif=',O17_verif
482        if (O17_verif) then
483            if (iso_verif_aberrant_o17_nostop(R17,R18,err_msg) &
484     &           .eq.1) then
485                stop
486            endif
487        endif !if (O17_verif) then
488
489#ifdef ISOVERIF
490#else
491        write(*,*) 'err_msg=',err_msg,': pourquoi verif?'
492        stop
493#endif                   
494       
495        return
496        end subroutine iso_verif_aberrant_o17
497
498        !*******************
499
500        function iso_verif_aberrant_o17_nostop(R17,R18,err_msg)
501        USE isotopes_mod, ONLY: tnat,iso_O17,iso_O18
502        implicit none
503        ! si l'O17-excess est aberrant, on affiche un message
504        ! d'erreur et renvoit 1
505
506        ! input:
507        real R17,R18
508        character*(*) err_msg ! message d''erreur à afficher
509        !local
510        !real o17excess
511        ! output
512        integer iso_verif_aberrant_o17_nostop
513
514        if (O17_verif) then
515
516        iso_verif_aberrant_o17_nostop=0
517        if ((o17excess(R17,R18).gt.o17excess_haut).or. &
518     &            (o17excess(R17,R18).lt.o17excess_bas)) then 
519            write(*,*) 'erreur detectee par iso_verif_aberrant_O17:'
520            write(*,*) err_msg
521            write(*,*) 'o17excess=',o17excess(R17,R18)
522            write(*,*) 'deltaO17=',(R17/tnat(iso_o17)-1.0)*1000.0
523            write(*,*) 'deltaO18=',(R18/tnat(iso_O18)-1.0)*1000.0
524            ! attention, vérifier que la ligne suivante est bien activée
525            iso_verif_aberrant_o17_nostop=1
526        endif
527
528        endif  !if (O17_verif) then
529
530#ifdef ISOVERIF
531#else
532        write(*,*) 'err_msg=',err_msg,': pourquoi verif?'
533        stop
534#endif                   
535       
536        return
537        end function iso_verif_aberrant_o17_nostop
538
539
540        !************************************
541
542        subroutine iso_verif_noNaN(x,err_msg)
543        implicit none
544        ! si x est NaN, on affiche message
545        ! d'erreur et stoppe
546
547        ! input:
548        real x
549        character*(*) err_msg ! message d''erreur à afficher
550
551        ! locals
552        !integer iso_verif_noNAN_nostop
553
554        if (iso_verif_noNAN_nostop(x,err_msg).eq.1) then
555            stop
556        endif
557#ifdef ISOVERIF
558#else
559        write(*,*) 'err_msg iso443=',err_msg,': pourquoi verif?'
560        stop
561#endif           
562       
563        end subroutine iso_verif_noNaN
564
565                !************************************
566
567        function iso_verif_noNaN_nostop(x,err_msg)
568        implicit none
569        ! si x est NaN, on affiche message
570        ! d'erreur et return 1 si erreur
571
572        ! input:
573        real x
574        character*(*) err_msg ! message d''erreur à afficher
575
576        ! output
577        integer iso_verif_noNAN_nostop
578
579        if ((x.gt.-borne).and.(x.lt.borne)) then
580                iso_verif_noNAN_nostop=0
581        else
582            write(*,*) 'erreur detectee par iso_verif_nonNAN:'
583            write(*,*) err_msg
584            write(*,*) 'x=',x
585            iso_verif_noNAN_nostop=1
586        endif     
587
588#ifdef ISOVERIF
589#else
590        write(*,*) 'err_msg iso 482=',err_msg,': pourquoi verif?'
591        stop
592#endif           
593
594        return
595        end function iso_verif_noNaN_nostop
596
597        subroutine iso_verif_noNaN_vect2D(x,err_msg,ni,n,m)
598        implicit none
599        ! si x est NaN, on affiche message
600        ! d'erreur et return 1 si erreur
601
602        ! input:
603          integer n,m,ni
604        real x(ni,n,m)
605        character*(*) err_msg ! message d''erreur à afficher
606
607        ! output
608
609        ! locals       
610        integer i,j,ixt
611
612      do i=1,n
613       do j=1,m
614        do ixt=1,ni
615         if ((x(ixt,i,j).gt.-borne).and. &
616     &            (x(ixt,i,j).lt.borne)) then
617         else !if ((x(ixt,i,j).gt.-borne).and.
618            write(*,*) 'erreur detectee par iso_verif_nonNAN:'
619            write(*,*) err_msg
620            write(*,*) 'x,ixt,i,j=',x(ixt,i,j),ixt,i,j
621            stop
622         endif  !if ((x(ixt,i,j).gt.-borne).and.
623        enddo !do ixt=1,ni
624       enddo !do j=1,m
625      enddo !do i=1,n     
626
627#ifdef ISOVERIF
628#else
629        write(*,*) 'err_msg iso525=',err_msg,': pourquoi verif?'
630        stop
631#endif           
632
633        end subroutine iso_verif_noNaN_vect2D
634
635
636        subroutine iso_verif_noNaN_vect1D(x,err_msg,ni,n)
637        implicit none
638        ! si x est NaN, on affiche message
639        ! d'erreur et return 1 si erreur
640
641        ! input:
642          integer n,ni
643        real x(ni,n)
644        character*(*) err_msg ! message d''erreur à afficher
645
646        ! output
647
648        ! locals       
649        integer i,ixt
650
651      do i=1,n
652        do ixt=1,ni
653         if ((x(ixt,i).gt.-borne).and. &
654     &            (x(ixt,i).lt.borne)) then
655         else !if ((x(ixt,i,j).gt.-borne).and.
656            write(*,*) 'erreur detectee par iso_verif_nonNAN:'
657            write(*,*) err_msg
658            write(*,*) 'x,ixt,i=',x(ixt,i),ixt,i
659            stop
660         endif  !if ((x(ixt,i,j).gt.-borne).and.
661        enddo !do ixt=1,ni
662      enddo !do i=1,n     
663
664#ifdef ISOVERIF
665#else
666        write(*,*) 'err_msg iso525=',err_msg,': pourquoi verif?'
667        stop
668#endif           
669
670        end subroutine iso_verif_noNaN_vect1D
671
672
673
674        !************************
675        subroutine iso_verif_egalite_choix(a,b,err_msg,erabs,errel)
676        implicit none
677        ! compare a et b. Si pas egal, on affiche message
678        ! d'erreur et stoppe
679        ! pour egalite, on verifie erreur absolue et arreur relative
680
681        ! input:
682        real a, b
683        real erabs,errel !erreur absolue et relative
684        character*(*) err_msg ! message d''erreur à afficher
685
686        ! locals
687        !integer iso_verif_egalite_choix_nostop
688
689        if (iso_verif_egalite_choix_nostop( &
690     &            a,b,err_msg,erabs,errel).eq.1) then
691           stop
692        endif
693
694#ifdef ISOVERIF
695#else
696        write(*,*) 'err_msg=',err_msg,': pourquoi verif?'
697        stop
698#endif           
699
700        end subroutine iso_verif_egalite_choix
701
702!************************
703        function iso_verif_egalite_choix_nostop &
704     &           (a,b,err_msg,erabs,errel)
705        implicit none
706        ! compare a et b. Si pas egal, on affiche message
707        ! d'erreur et stoppe
708        ! pour egalite, on verifie erreur absolue et arreur relative
709
710        ! input:
711        real a, b
712        real erabs,errel !erreur absolue et relative
713        character*(*) err_msg ! message d''erreur à afficher
714
715        ! output
716        integer iso_verif_egalite_choix_nostop
717        ! locals
718        !integer iso_verif_noNaN_nostop
719
720        iso_verif_egalite_choix_nostop=0
721
722#ifdef ISOVERIF
723        if (iso_verif_noNaN_nostop(a,err_msg).eq.1) then
724            write(*,*) 'a=',a
725            iso_verif_egalite_choix_nostop=1
726        endif     
727        if (iso_verif_noNaN_nostop(b,err_msg).eq.1) then
728            write(*,*) 'b=',b
729            iso_verif_egalite_choix_nostop=1
730        endif     
731#endif
732
733        if (abs(a-b).gt.erabs) then
734          if (abs((a-b)/max(max(abs(b),abs(a)),1e-18)) &
735     &            .gt.errel) then
736            write(*,*) 'erreur detectee par iso_verif_egalite:'
737            write(*,*) err_msg
738            write(*,*) 'a=',a
739            write(*,*) 'b=',b
740            write(*,*) 'erabs,errel=',erabs,errel
741            iso_verif_egalite_choix_nostop=1
742!            stop
743          endif
744        endif
745
746#ifdef ISOVERIF
747#else
748        write(*,*) 'err_msg=',err_msg,': pourquoi verif?'
749        stop
750#endif           
751       
752        return
753        end function iso_verif_egalite_choix_nostop       
754
755
756
757        !******************
758        subroutine iso_verif_positif(x,err_msg)
759        use isotopes_mod, ONLY: ridicule
760        implicit none
761        ! si x<0, on plante.
762        ! si très limite, on le met à 0.
763
764        ! input:
765        real x
766        character*(*) err_msg ! message d''erreur à afficher
767
768        ! locals
769        !integer iso_verif_positif_choix_nostop
770
771        if (iso_verif_positif_choix_nostop(x,ridicule,err_msg) &
772     &           .eq.1) then
773           stop
774        endif
775     
776#ifdef ISOVERIF
777#else
778        write(*,*) 'iso_verif 637: err_msg=',err_msg, &
779     &           ': pourquoi verif?'
780        stop
781#endif
782       
783        end subroutine iso_verif_positif
784
785        !******************
786        subroutine iso_verif_positif_vect(x,n,err_msg)
787        use isotopes_mod, ONLY: ridicule
788        implicit none
789        ! si x<0, on plante.
790
791        ! input:
792        integer n
793        real x(n)
794        character*(*) err_msg ! message d''erreur à afficher
795
796        ! locals
797        integer i
798        integer iso_verif_positif_nostop
799        integer ifaux
800
801        iso_verif_positif_nostop=0
802        do i=1,n
803          if (x(i).lt.-ridicule) then
804            iso_verif_positif_nostop=1
805            ifaux=i
806          endif
807        enddo
808        if (iso_verif_positif_nostop.eq.1) then
809          write(*,*) 'erreur detectee par iso_verif_positif_vect:'
810          write(*,*) err_msg
811          write(*,*) 'i,x=',ifaux,x(ifaux)
812          stop
813        endif   
814       
815        end subroutine iso_verif_positif_vect
816
817        subroutine iso_verif_positif_choix_vect(x,n,ridic,err_msg)
818        implicit none
819        ! si x<0, on plante.
820
821        ! input:
822        integer n
823        real x(n)
824        real ridic
825        character*(*) err_msg ! message d''erreur à afficher
826        ! locals
827        integer i
828        integer iso_verif_positif_nostop
829        integer ifaux
830
831        iso_verif_positif_nostop=0
832        do i=1,n
833          if (x(i).lt.-ridic) then
834                iso_verif_positif_nostop=1
835                ifaux=i
836          endif
837        enddo
838        if (iso_verif_positif_nostop.eq.1) then                       
839         write(*,*) 'erreur detectee par iso_verif_positif_choix_vect:'
840         write(*,*) err_msg
841         write(*,*) 'i,x=',ifaux,x(ifaux)
842         stop
843        endif   
844       
845        end subroutine iso_verif_positif_choix_vect
846
847
848!******************
849        subroutine iso_verif_positif_strict(x,err_msg)
850        implicit none
851        ! si x<0, on plante.
852        ! si très limite, on le met à 0.
853
854        ! input:
855        real x
856        character*(*) err_msg ! message d''erreur à afficher
857
858        ! locals
859        !integer iso_verif_positif_strict_nostop
860
861        if (iso_verif_positif_strict_nostop(x,err_msg) &
862     &           .eq.1) then
863           stop
864        endif           
865       
866        end subroutine iso_verif_positif_strict
867
868!******************
869
870        function iso_verif_positif_strict_nostop(x,err_msg)
871        implicit none
872        ! si x<0, on plante.
873        ! si très limite, on le met à 0.
874
875        ! input:
876        real x
877        character*(*) err_msg ! message d''erreur à afficher*
878
879        ! output
880        integer iso_verif_positif_strict_nostop
881
882        if (x.gt.0.0) then
883            iso_verif_positif_strict_nostop=0
884        else     
885            write(*,*) 'erreur detectee par iso_verif_positif_strict:'
886            write(*,*) err_msg
887            write(*,*) 'x=',x 
888            iso_verif_positif_strict_nostop=1   
889!            stop 
890        endif   
891       
892        return
893        end function iso_verif_positif_strict_nostop
894
895!******************
896
897        subroutine iso_verif_positif_choix(x,ridic,err_msg)
898        implicit none
899        ! si x<0, on plante.
900        ! si très limite, on le met à 0.
901
902        ! input:
903        real x
904        real ridic
905        character*(*) err_msg ! message d''erreur à afficher
906
907        ! locals
908        !integer iso_verif_positif_choix_nostop
909
910        if (iso_verif_positif_choix_nostop(x,ridic,err_msg) &
911     &           .eq.1) then
912           stop
913        endif
914     
915#ifdef ISOVERIF
916#else
917        write(*,*) 'iso_verif 801: err_msg=',err_msg, &
918     &           ': pourquoi verif?'
919        stop
920#endif           
921       
922        end subroutine iso_verif_positif_choix       
923
924        !******************
925        function iso_verif_positif_nostop(x,err_msg)
926        use isotopes_mod, ONLY: ridicule
927        implicit none
928        ! si x<0, on plante.
929        ! si très limite, on le met à 0.
930
931        ! input:
932        real x
933        character*(*) err_msg ! message d''erreur à afficher
934
935        ! output
936        integer iso_verif_positif_nostop
937
938        ! locals
939        !integer iso_verif_positif_choix_nostop
940
941        iso_verif_positif_nostop=iso_verif_positif_choix_nostop &
942     &          (x,ridicule,err_msg)
943
944#ifdef ISOVERIF
945#else
946        write(*,*) 'iso_verif 837: err_msg=',err_msg, &
947     &           ': pourquoi verif?'
948        stop
949
950#endif         
951       
952        return
953        end function iso_verif_positif_nostop
954
955        !******************
956        function iso_verif_positif_choix_nostop(x,ridic,err_msg)
957        implicit none
958        ! si x<0, on plante.
959        ! si très limite, on le met à 0.
960
961        ! input:
962        real x
963        real ridic
964        character*(*) err_msg ! message d''erreur à afficher
965
966        ! output
967        integer iso_verif_positif_choix_nostop
968
969        if (x.lt.-ridic) then
970            write(*,*) 'erreur detectee par iso_verif_positif_nostop:'
971            write(*,*) err_msg
972            write(*,*) 'x=',x
973            iso_verif_positif_choix_nostop=1
974        else
975!            x=max(x,0.0)
976            iso_verif_positif_choix_nostop=0
977        endif
978
979#ifdef ISOVERIF
980#else
981        write(*,*) 'iso_verif 877: err_msg=',err_msg, &
982     &           ': pourquoi verif?'
983        stop
984#endif
985       
986        return
987        end function iso_verif_positif_choix_nostop
988
989
990        !**************
991
992       
993        subroutine iso_verif_O18_aberrant(Rd,Ro,err_msg)
994        implicit none
995
996        ! vérifie que:
997        ! 1) deltaD et deltaO18 dans bonne gamme
998        ! 2) dexcess dans bonne gamme
999        ! input:
1000        real Rd,Ro
1001        character*(*) err_msg ! message d''erreur à afficher
1002
1003        ! local
1004        !integer iso_verif_O18_aberrant_nostop
1005
1006        if (iso_verif_O18_aberrant_nostop(Rd,Ro,err_msg).eq.1) then
1007            stop
1008        endif
1009
1010        end subroutine iso_verif_O18_aberrant
1011       
1012        function iso_verif_O18_aberrant_nostop(Rd,Ro,err_msg)
1013        USE isotopes_mod, ONLY: tnat, iso_HDO, iso_O18
1014        implicit none
1015
1016        ! vérifie que:
1017        ! 1) deltaD et deltaO18 dans bonne gamme
1018        ! 2) dexcess dans bonne gamme
1019
1020        ! input:
1021        real Rd,Ro
1022        character*(*) err_msg ! message d''erreur à afficher
1023
1024        ! outputs
1025        integer iso_verif_O18_aberrant_nostop
1026
1027        !locals
1028        real deltaD,deltao,dexcess
1029
1030        deltaD=(Rd/tnat(iso_hdo)-1)*1000
1031        deltao=(Ro/tnat(iso_O18)-1)*1000
1032        dexcess=deltaD-8*deltao
1033
1034        iso_verif_O18_aberrant_nostop=0
1035        if ((deltaD.lt.deltaDmin).or.(deltao.lt.deltaDmin/2.0).or. &
1036     &        (deltaD.gt.deltalim).or.(deltao.gt.deltalim/8.0).or. &
1037     &        ((deltaD.gt.-500.0).and.((dexcess.lt.dexcess_min) &
1038     &        .or.(dexcess.gt.dexcess_max)))) then
1039            write(*,*) 'erreur detectee par iso_verif_O18_aberrant:'
1040            write(*,*) err_msg
1041            write(*,*) 'delta180=',deltao
1042            write(*,*) 'deltaD=',deltaD
1043            write(*,*) 'Dexcess=',dexcess
1044            write(*,*) 'tnat=',tnat
1045!            stop
1046            iso_verif_O18_aberrant_nostop=1
1047          endif
1048
1049#ifdef ISOVERIF
1050#else
1051        write(*,*) 'err_msg=',err_msg,': pourquoi verif?'
1052        stop
1053#endif                   
1054
1055        return
1056        end function iso_verif_O18_aberrant_nostop
1057
1058
1059        ! **********
1060        function deltaD(R)
1061        USE isotopes_mod, ONLY: tnat,iso_HDO
1062        implicit none
1063        real R,deltaD
1064
1065       
1066        if (iso_HDO.gt.0) then
1067           deltaD=(R/tnat(iso_HDO)-1)*1000.0
1068        else
1069            write(*,*) 'iso_verif_egalite>deltaD 260: iso_HDO.gt.0=', &
1070     &           iso_HDO.gt.0
1071        endif
1072        return
1073        end function deltaD
1074
1075        ! **********
1076        function deltaO(R)
1077        USE isotopes_mod, ONLY: tnat,iso_O18
1078        implicit none
1079        real R,deltaO
1080       
1081        if (iso_O18.gt.0) then
1082           deltaO=(R/tnat(iso_O18)-1)*1000.0
1083        else
1084            write(*,*) 'iso_verif_egalite>deltaO18 260: iso_O18.gt.0=', &
1085     &           iso_O18.gt.0
1086        endif
1087        return
1088        end function deltaO
1089
1090        ! **********
1091        function dexcess(RD,RO)
1092        USE isotopes_mod, ONLY: tnat,iso_O18,iso_HDO
1093        implicit none
1094        real RD,RO,deltaD,deltaO,dexcess
1095       
1096        if ((iso_O18.gt.0).and.(iso_HDO.gt.0)) then
1097           deltaD=(RD/tnat(iso_HDO)-1)*1000.0
1098           deltaO=(RO/tnat(iso_O18)-1)*1000.0
1099           dexcess=deltaD-8*deltaO
1100        else
1101            write(*,*) 'iso_verif_egalite 1109: iso_O18,iso_HDO=',iso_O18,iso_HDO
1102        endif
1103        return
1104        end function dexcess
1105
1106
1107        ! **********
1108        function delta_all(R,ixt)
1109        USE isotopes_mod, ONLY: tnat
1110        implicit none
1111        real R,delta_all
1112        integer ixt
1113       
1114        delta_all=(R/tnat(ixt)-1)*1000.0
1115        return
1116        end function delta_all
1117
1118        ! **********
1119        function delta_to_R(delta,ixt)
1120        USE isotopes_mod, ONLY: tnat
1121        implicit none
1122        real delta,delta_to_R
1123        integer ixt
1124       
1125        delta_to_R=(delta/1000.0+1.0)*tnat(ixt)
1126        return
1127        end function delta_to_R
1128
1129         ! **********
1130        function o17excess(R17,R18)
1131        USE isotopes_mod, ONLY: tnat,iso_O18,iso_O17
1132        implicit none
1133        real R17,R18,o17excess
1134       
1135        if ((iso_O17.gt.0).and.(iso_O18.gt.0)) then
1136           
1137           o17excess=1e6*(log(R17/tnat(iso_o17)) &
1138     &           -0.528*log(R18/tnat(iso_O18)))
1139!           write(*,*) 'o17excess=',o17excess
1140        else
1141            write(*,*) 'iso_verif_egalite>deltaD 260: iso_O17.gt.0,18=', &
1142     &           iso_O17.gt.0,iso_O18.gt.0
1143        endif
1144        return
1145        end function o17excess
1146
1147        !       ****************
1148
1149          subroutine iso_verif_egalite_vect2D( &
1150     &           xt,q,err_msg,ni,n,m)
1151       
1152        USE isotopes_mod, ONLY: iso_eau
1153          implicit none
1154
1155          ! inputs
1156          integer n,m,ni
1157          real q(n,m)
1158          real xt(ni,n,m)
1159          character*(*) err_msg
1160
1161        ! locals
1162        integer iso_verif_egalite_nostop_loc
1163        integer i,j,ixt
1164        integer ifaux,jfaux
1165
1166        !write(*,*) 'iso_verif_egalite_vect2D 1099 tmp: q(2,1),xt(iso_eau,2,1)=',q(2,1),xt(iso_eau,2,1)
1167        !write(*,*) 'ni,n,m=',ni,n,m,errmax,errmaxrel
1168        if (iso_eau.gt.0) then
1169        iso_verif_egalite_nostop_loc=0
1170        do i=1,n
1171         do j=1,m
1172          if (abs(q(i,j)-xt(iso_eau,i,j)).gt.errmax) then
1173           if (abs((q(i,j)-xt(iso_eau,i,j))/ &
1174     &           max(max(abs(q(i,j)),abs(xt(iso_eau,i,j))),1e-18)) &
1175     &           .gt.errmaxrel) then
1176              iso_verif_egalite_nostop_loc=1
1177              ifaux=i
1178              jfaux=j
1179           endif
1180          endif
1181         enddo !do j=1,m
1182        enddo !do i=1,n
1183
1184        if (iso_verif_egalite_nostop_loc.eq.1) then
1185          write(*,*) 'erreur detectee par iso_verif_egalite_vect2D:'
1186          write(*,*) err_msg
1187          write(*,*) 'i,j=',ifaux,jfaux
1188          write(*,*) 'xt,q=',xt(iso_eau,ifaux,jfaux),q(ifaux,jfaux)
1189          stop
1190        endif
1191        endif
1192       
1193#ifdef ISOVERIF
1194        call iso_verif_noNaN_vect2D(xt,err_msg,ni,n,m)
1195#endif         
1196
1197        return
1198        end subroutine iso_verif_egalite_vect2D
1199
1200        subroutine iso_verif_egalite_vect1D( &
1201     &           xt,q,err_msg,ni,n)
1202
1203        USE isotopes_mod, ONLY: iso_eau
1204        implicit none
1205
1206        ! inputs
1207        integer n,ni
1208        real q(n)
1209        real xt(ni,n)
1210        character*(*) err_msg
1211
1212        ! locals
1213        integer iso_verif_egalite_nostop_loc
1214        integer i
1215        integer ifaux
1216
1217        if (iso_eau.gt.0) then
1218        iso_verif_egalite_nostop_loc=0
1219        do i=1,n
1220          if (abs(q(i)-xt(iso_eau,i)).gt.errmax) then
1221           if (abs((q(i)-xt(iso_eau,i))/ &
1222     &           max(max(abs(q(i)),abs(xt(iso_eau,i))),1e-18)) &
1223     &           .gt.errmaxrel) then
1224              iso_verif_egalite_nostop_loc=1
1225              ifaux=i
1226           endif !if (abs((q(i)-xt(iso_eau,i))/
1227          endif !if (abs(q(i)-xt(iso_eau,i)).gt.errmax) then
1228        enddo !do i=1,n
1229
1230        if (iso_verif_egalite_nostop_loc.eq.1) then
1231          write(*,*) 'erreur detectee par iso_verif_egalite_vect2D:'
1232          write(*,*) err_msg
1233          write(*,*) 'i=',ifaux
1234          write(*,*) 'xt,q=',xt(iso_eau,ifaux),q(ifaux)
1235          stop
1236        endif  !if (iso_verif_egalite_nostop.eq.1) then
1237        endif !if (iso_eau.gt.0) then
1238
1239        end subroutine iso_verif_egalite_vect1D       
1240
1241        subroutine iso_verif_egalite_std_vect( &
1242     &           a,b,err_msg,n,m,errmax,errmaxrel)
1243
1244          implicit none
1245
1246          ! inputs
1247          integer n,m,ni
1248          real a(n,m)
1249          real b(n,m)
1250          character*(*) err_msg
1251          real errmax,errmaxrel
1252
1253        ! locals
1254        integer iso_verif_egalite_nostop_loc
1255        integer i,j
1256        integer ifaux,jfaux
1257
1258        iso_verif_egalite_nostop_loc=0
1259        do i=1,n
1260         do j=1,m
1261          if (abs(a(i,j)-b(i,j)).gt.errmax) then
1262           if (abs((a(i,j)-b(i,j))/ &
1263     &           max(max(abs(a(i,j)),abs(b(i,j))),1e-18)) &
1264     &           .gt.errmaxrel) then
1265              iso_verif_egalite_nostop_loc=1
1266              ifaux=i
1267              jfaux=j
1268           endif
1269          endif
1270         enddo !do j=1,m
1271        enddo !do i=1,n
1272
1273        if (iso_verif_egalite_nostop_loc.eq.1) then
1274          write(*,*) 'erreur detectee par iso_verif_egalite_vect2D:'
1275          write(*,*) err_msg
1276          write(*,*) 'i,j=',ifaux,jfaux
1277          write(*,*) 'a,b=',a(ifaux,jfaux),b(ifaux,jfaux)
1278          stop
1279        endif
1280
1281        return
1282        end subroutine iso_verif_egalite_std_vect
1283
1284        subroutine iso_verif_aberrant_vect2D( &
1285     &           xt,q,err_msg,ni,n,m)
1286        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO
1287          implicit none
1288
1289          ! inputs   
1290          integer n,m,ni
1291          real q(n,m)
1292          real xt(ni,n,m)
1293          character*(*) err_msg
1294
1295        ! locals
1296        integer iso_verif_aberrant_nostop_loc
1297        integer i,j
1298        integer ifaux,jfaux
1299        !real deltaD
1300
1301        if (iso_HDO.gt.0) then
1302        iso_verif_aberrant_nostop_loc=0
1303        do i=1,n
1304         do j=1,m
1305          if (q(i,j).gt.ridicule) then
1306            if (((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1307     &                   .gt.deltalim).or. &
1308     &          ((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1309     &                   .lt.-borne)) then   
1310              iso_verif_aberrant_nostop_loc=1
1311              ifaux=i
1312              jfaux=j
1313           endif
1314          endif
1315         enddo !do j=1,m
1316        enddo !do i=1,n
1317
1318        if (iso_verif_aberrant_nostop_loc.eq.1) then
1319          write(*,*) 'erreur detectee par iso_verif_aberrant_vect2D:'
1320          write(*,*) err_msg
1321          write(*,*) 'i,j=',ifaux,jfaux
1322          write(*,*) 'deltaD=',deltaD(xt(iso_HDO,ifaux,jfaux) &
1323     &           /q(ifaux,jfaux))
1324          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1325          stop
1326        endif 
1327        endif !if (iso_HDO.gt.0) then
1328
1329        end subroutine iso_verif_aberrant_vect2D       
1330
1331        subroutine iso_verif_aberrant_enc_vect2D( &
1332     &           xt,q,err_msg,ni,n,m)
1333
1334        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO
1335          implicit none
1336
1337          ! inputs   
1338          integer n,m,ni
1339          real q(n,m)
1340          real xt(ni,n,m)
1341          character*(*) err_msg
1342
1343        ! locals
1344        integer iso_verif_aberrant_nostop_loc
1345        integer i,j
1346        integer ifaux,jfaux
1347        !real deltaD
1348
1349        if (iso_HDO.gt.0) then
1350        iso_verif_aberrant_nostop_loc=0
1351        do i=1,n
1352         do j=1,m
1353          if (q(i,j).gt.ridicule) then
1354            if (((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1355     &                   .gt.deltalim).or. &
1356     &          ((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1357     &                   .lt.deltaDmin).or. &
1358     &           (xt(iso_HDO,i,j).lt.-borne).or. &
1359     &           (xt(iso_HDO,i,j).gt.borne)) then     
1360              iso_verif_aberrant_nostop_loc=1
1361              ifaux=i
1362              jfaux=j
1363           endif
1364          endif
1365         enddo !do j=1,m
1366        enddo !do i=1,n
1367
1368        if (iso_verif_aberrant_nostop_loc.eq.1) then
1369          write(*,*) 'erreur detectee par ', &
1370     &           'iso_verif_aberrant_enc_vect2D:'
1371          write(*,*) err_msg
1372          write(*,*) 'i,j=',ifaux,jfaux
1373          write(*,*) 'deltaD=',deltaD(xt(iso_HDO,ifaux,jfaux) &
1374     &           /q(ifaux,jfaux))
1375          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1376          write(*,*) 'q(ifaux,jfaux)=',q(ifaux,jfaux)
1377          call abort_physic('isotopes_verif_mod','iso_verif_aberrant_enc_vect2D',1)
1378        endif 
1379        endif !if (iso_HDO.gt.0) then
1380
1381        end subroutine iso_verif_aberrant_enc_vect2D       
1382
1383       
1384        subroutine iso_verif_aberrant_enc_vect2D_ns( &
1385     &           xt,q,err_msg,ni,n,m)
1386
1387        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO
1388          implicit none
1389
1390          ! inputs   
1391          integer n,m,ni
1392          real q(n,m)
1393          real xt(ni,n,m)
1394          character*(*) err_msg
1395
1396        ! locals
1397        integer iso_verif_aberrant_nostop_loc
1398        integer i,j
1399        integer ifaux,jfaux
1400        !real deltaD
1401
1402        if (iso_HDO.gt.0) then
1403        iso_verif_aberrant_nostop_loc=0
1404        do i=1,n
1405         do j=1,m
1406          if (q(i,j).gt.ridicule) then
1407            if (((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1408     &                   .gt.deltalim).or. &
1409     &          ((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1410     &                   .lt.deltaDmin)) then   
1411              iso_verif_aberrant_nostop_loc=1
1412              ifaux=i
1413              jfaux=j
1414           endif
1415          endif
1416         enddo !do j=1,m
1417        enddo !do i=1,n
1418
1419        if (iso_verif_aberrant_nostop_loc.eq.1) then
1420          write(*,*) 'erreur detectee par ', &
1421     &           'iso_verif_aberrant_vect2D_ns:'
1422          write(*,*) err_msg
1423          write(*,*) 'i,j=',ifaux,jfaux
1424          write(*,*) 'deltaD=',deltaD(xt(iso_HDO,ifaux,jfaux) &
1425     &           /q(ifaux,jfaux))
1426          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1427!          stop
1428        endif 
1429        endif !if (iso_HDO.gt.0) then
1430
1431        end subroutine iso_verif_aberrant_enc_vect2D_ns       
1432
1433
1434         subroutine iso_verif_aberrant_vect2Dch( &
1435     &           xt,q,err_msg,ni,n,m,deltaDmax)
1436
1437        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO
1438          implicit none
1439
1440
1441          ! inputs   
1442          integer n,m,ni
1443          real q(n,m)
1444          real xt(ni,n,m)
1445          character*(*) err_msg
1446          real deltaDmax
1447
1448        ! locals
1449        integer iso_verif_aberrant_nostop_loc
1450        integer i,j
1451        integer ifaux,jfaux
1452        !real deltaD
1453
1454        if (iso_HDO.gt.0) then
1455        iso_verif_aberrant_nostop_loc=0
1456        do i=1,n
1457         do j=1,m
1458          if (q(i,j).gt.ridicule) then
1459            if (((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1460     &                   .gt.deltaDmax).or. &
1461     &          ((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1462     &                   .lt.-borne)) then   
1463              iso_verif_aberrant_nostop_loc=1
1464              ifaux=i
1465              jfaux=j
1466           endif
1467          endif
1468         enddo !do j=1,m
1469        enddo !do i=1,n
1470
1471        if (iso_verif_aberrant_nostop_loc.eq.1) then
1472          write(*,*) 'erreur detectee par iso_verif_aberrant_vect2D:'
1473          write(*,*) err_msg
1474          write(*,*) 'i,j=',ifaux,jfaux
1475          write(*,*) 'deltaD=',deltaD(xt(iso_HDO,ifaux,jfaux) &
1476     &           /q(ifaux,jfaux))
1477          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1478          stop
1479        endif 
1480        endif !if (iso_HDO.gt.0) then
1481
1482        end subroutine iso_verif_aberrant_vect2Dch     
1483
1484        subroutine iso_verif_O18_aberrant_enc_vect2D( &
1485     &           xt,q,err_msg,ni,n,m)
1486
1487        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO,iso_O18
1488          implicit none
1489
1490          ! inputs   
1491          integer n,m,ni
1492          real q(n,m)
1493          real xt(ni,n,m)
1494          character*(*) err_msg
1495
1496        ! locals
1497        integer iso_verif_aberrant_nostop_loc
1498        integer i,j
1499        integer ifaux,jfaux
1500        real deltaDloc,deltaoloc,dexcessloc
1501
1502        if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then
1503        iso_verif_aberrant_nostop_loc=0
1504        do i=1,n
1505         do j=1,m
1506          if (q(i,j).gt.ridicule) then
1507
1508            deltaDloc=(xt(iso_HDO,i,j)/q(i,j)/tnat(iso_hdo)-1)*1000
1509            deltaoloc=(xt(iso_O18,i,j)/q(i,j)/tnat(iso_O18)-1)*1000
1510            dexcessloc=deltaDloc-8*deltaoloc
1511            if ((deltaDloc.lt.deltaDmin).or.(deltaoloc.lt.deltaDmin/2.0).or. &
1512     &        (deltaDloc.gt.deltalim).or.(deltaoloc.gt.deltalim/8.0).or. &
1513     &        ((deltaDloc.gt.-500.0).and.((dexcessloc.lt.dexcess_min) &
1514     &        .or.(dexcessloc.gt.dexcess_max)))) then
1515              iso_verif_aberrant_nostop_loc=1
1516              ifaux=i
1517              jfaux=j
1518              write(*,*) 'deltaD,deltao,dexcess=',deltaDloc,deltaoloc,dexcessloc
1519           endif
1520          endif
1521         enddo !do j=1,m
1522        enddo !do i=1,n
1523
1524        if (iso_verif_aberrant_nostop_loc.eq.1) then
1525          write(*,*) 'erreur detectee par ', &
1526     &           'iso_verif_aberrant_enc_vect2D:'
1527          write(*,*) err_msg
1528          write(*,*) 'i,j=',ifaux,jfaux
1529          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1530          write(*,*) 'q(ifaux,jfaux)=',q(ifaux,jfaux)
1531          call abort_physic('isotopes_verif_mod','iso_verif_aberrant_enc_vect2D',1)
1532        endif 
1533        endif !if (iso_HDO.gt.0) then
1534
1535        end subroutine iso_verif_O18_aberrant_enc_vect2D   
1536
1537
1538        subroutine select_dim23_from4D(n1,n2,n3,n4, &
1539     &          var,vec,i1,i4)
1540        implicit none
1541
1542        ! inputs
1543        integer n1,n2,n3,n4
1544        real var(n1,n2,n3,n4)
1545        integer i1,i4
1546        ! outputs
1547        real vec(n2,n3)
1548        ! locals
1549        integer i2,i3
1550
1551        do i2=1,n2
1552         do i3=1,n3
1553          vec(i2,i3)=var(i1,i2,i3,i4)
1554         enddo
1555        enddo
1556
1557        return
1558        end subroutine select_dim23_from4D
1559
1560       
1561        subroutine select_dim4_from4D(ntime,nlev,nlat,nlon, &
1562     &          var,vec,itime,ilev,ilat)
1563        implicit none
1564
1565        ! inputs
1566        integer ntime,nlev,nlat,nlon
1567        real var(ntime,nlev,nlat,nlon)
1568        integer itime,ilev,ilat
1569        ! outputs
1570        real vec(nlon)
1571        ! locals
1572        integer ilon
1573
1574        do ilon=1,nlon
1575          vec(ilon)=var(itime,ilev,ilat,ilon)
1576        enddo
1577
1578        return
1579        end subroutine select_dim4_from4D
1580
1581        subroutine select_dim5_from5D(n1,n2,n3,n4,n5, &
1582     &          var,vec,i1,i2,i3,i4)
1583        implicit none
1584
1585        ! inputs
1586        integer n1,n2,n3,n4,n5
1587        real var(n1,n2,n3,n4,n5)
1588        integer i1,i2,i3,i4
1589        ! outputs
1590        real vec(n5)
1591        ! locals
1592        integer i5
1593
1594        do i5=1,n5
1595          vec(i5)=var(i1,i2,i3,i4,i5)
1596        enddo
1597
1598        end subroutine select_dim5_from5D
1599
1600       
1601        subroutine select_dim3_from3D(ntime,nlat,nlon, &
1602     &          var,vec,itime,ilat)
1603        implicit none
1604
1605        ! inputs
1606        integer ntime,nlat,nlon
1607        real var(ntime,nlat,nlon)
1608        integer itime,ilat
1609        ! outputs
1610        real vec(nlon)
1611        ! locals
1612        integer ilon
1613
1614        do ilon=1,nlon
1615          vec(ilon)=var(itime,ilat,ilon)
1616        enddo
1617
1618        end subroutine select_dim3_from3D
1619
1620       
1621        subroutine select_dim23_from3D(n1,n2,n3, &
1622     &          var,vec,i1)
1623        implicit none
1624
1625        ! inputs
1626        integer n1,n2,n3
1627        real var(n1,n2,n3)
1628        integer i1
1629        ! outputs
1630        real vec(n2,n3)
1631        ! locals
1632        integer i2,i3
1633
1634        do i2=1,n2
1635         do i3=1,n3
1636          vec(i2,i3)=var(i1,i2,i3)
1637         enddo
1638        enddo
1639
1640        end subroutine select_dim23_from3D
1641
1642        subroutine putinto_dim23_from4D(n1,n2,n3,n4, &
1643     &          var,vec,i1,i4)
1644        implicit none
1645
1646        ! inputs
1647        integer n1,n2,n3,n4
1648        real vec(n2,n3)
1649        integer i1,i4
1650        ! inout
1651        real var(n1,n2,n3,n4)
1652        ! locals
1653        integer i2,i3
1654
1655       do i2=1,n2
1656        do i3=1,n3
1657          var(i1,i2,i3,i4)=vec(i2,i3)
1658         enddo
1659        enddo
1660
1661        end subroutine putinto_dim23_from4D
1662
1663       
1664        subroutine putinto_dim12_from4D(n1,n2,n3,n4, &
1665     &          var,vec,i3,i4)
1666        implicit none
1667
1668        ! inputs
1669        integer n1,n2,n3,n4
1670        real vec(n1,n2)
1671        integer i3,i4
1672        ! inout
1673        real var(n1,n2,n3,n4)
1674        ! locals
1675        integer i1,i2
1676
1677       do i1=1,n1
1678        do i2=1,n2
1679          var(i1,i2,i3,i4)=vec(i1,i2)
1680         enddo
1681        enddo
1682
1683        end subroutine putinto_dim12_from4D
1684       
1685        subroutine putinto_dim23_from3D(n1,n2,n3, &
1686     &          var,vec,i1)
1687        implicit none
1688
1689        ! inputs
1690        integer n1,n2,n3
1691        real vec(n2,n3)
1692        integer i1
1693        ! inout
1694        real var(n1,n2,n3)
1695        ! locals
1696        integer i2,i3
1697
1698       do i2=1,n2
1699        do i3=1,n3
1700          var(i1,i2,i3)=vec(i2,i3)
1701         enddo
1702        enddo
1703
1704        end subroutine putinto_dim23_from3D
1705
1706       
1707
1708        subroutine iso_verif_noNaN_par2D(x,err_msg,ni,n,m,ib,ie)
1709        implicit none
1710        ! si x est NaN, on affiche message
1711        ! d'erreur et return 1 si erreur
1712
1713        ! input:
1714          integer n,m,ni,ib,ie
1715        real x(ni,n,m)
1716        character*(*) err_msg ! message d''erreur à afficher
1717
1718        ! output
1719
1720        ! locals       
1721        integer i,j,ixt
1722
1723      do i=ib,ie
1724       do j=1,m
1725        do ixt=1,ni
1726         if ((x(ixt,i,j).gt.-borne).and. &
1727     &            (x(ixt,i,j).lt.borne)) then
1728         else !if ((x(ixt,i,j).gt.-borne).and.
1729            write(*,*) 'erreur detectee par iso_verif_nonNAN:'
1730            write(*,*) err_msg
1731            write(*,*) 'x,ixt,i,j=',x(ixt,i,j),ixt,i,j
1732            stop
1733         endif  !if ((x(ixt,i,j).gt.-borne).and.
1734        enddo !do ixt=1,ni
1735       enddo !do j=1,m
1736      enddo !do i=1,n     
1737
1738#ifdef ISOVERIF
1739#else
1740        write(*,*) 'err_msg iso1772=',err_msg,': pourquoi verif?'
1741        stop
1742#endif           
1743
1744        return
1745        end subroutine iso_verif_noNaN_par2D
1746
1747       
1748        subroutine iso_verif_aberrant_enc_par2D( &
1749     &           xt,q,err_msg,ni,n,m,ib,ie)
1750
1751        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO
1752          implicit none
1753
1754          ! inputs   
1755          integer n,m,ni,ib,ie
1756          real q(n,m)
1757          real xt(ni,n,m)
1758          character*(*) err_msg
1759
1760        ! locals
1761        integer iso_verif_aberrant_nostop_loc
1762        integer i,j
1763        integer ifaux,jfaux
1764        !real deltaD
1765
1766        if (iso_HDO.gt.0) then
1767        iso_verif_aberrant_nostop_loc=0
1768        do i=ib,ie
1769         do j=1,m
1770          if (q(i,j).gt.ridicule) then
1771            if (((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1772     &                   .gt.deltalim).or. &
1773     &          ((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1774     &                   .lt.deltaDmin)) then   
1775              iso_verif_aberrant_nostop_loc=1
1776              ifaux=i
1777              jfaux=j
1778           endif
1779          endif
1780         enddo !do j=1,m
1781        enddo !do i=1,n
1782
1783        if (iso_verif_aberrant_nostop_loc.eq.1) then
1784          write(*,*) 'erreur detectee par iso_verif_aberrant_par2D:'
1785          write(*,*) err_msg
1786          write(*,*) 'i,j=',ifaux,jfaux
1787          write(*,*) 'deltaD=',deltaD(xt(iso_HDO,ifaux,jfaux) &
1788     &           /q(ifaux,jfaux))
1789          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1790          write(*,*) 'q(ifaux,jfaux)=',q(ifaux,jfaux)
1791          stop
1792        endif 
1793        endif !if (iso_HDO.gt.0) then
1794
1795        end subroutine iso_verif_aberrant_enc_par2D       
1796
1797       
1798          subroutine iso_verif_egalite_par2D( &
1799     &           xt,q,err_msg,ni,n,m,ib,ie)
1800       
1801        USE isotopes_mod, ONLY: iso_eau
1802          implicit none
1803
1804          ! inputs
1805          integer n,m,ni,ib,ie
1806          real q(n,m)
1807          real xt(ni,n,m)
1808          character*(*) err_msg
1809
1810        ! locals
1811        integer iso_verif_egalite_nostop_loc
1812        integer i,j
1813        integer ifaux,jfaux
1814
1815        if (iso_eau.gt.0) then
1816        iso_verif_egalite_nostop_loc=0
1817        do i=ib,ie
1818         do j=1,m
1819          if (abs(q(i,j)-xt(iso_eau,i,j)).gt.errmax) then
1820           if (abs((q(i,j)-xt(iso_eau,i,j))/ &
1821     &           max(max(abs(q(i,j)),abs(xt(iso_eau,i,j))),1e-18)) &
1822     &           .gt.errmaxrel) then
1823              iso_verif_egalite_nostop_loc=1
1824              ifaux=i
1825              jfaux=j
1826           endif
1827          endif
1828         enddo !do j=1,m
1829        enddo !do i=1,n
1830
1831        if (iso_verif_egalite_nostop_loc.eq.1) then
1832          write(*,*) 'erreur detectee par iso_verif_egalite_vect2D:'
1833          write(*,*) err_msg
1834          write(*,*) 'i,j=',ifaux,jfaux
1835          write(*,*) 'xt,q=',xt(iso_eau,ifaux,jfaux),q(ifaux,jfaux)
1836          stop
1837        endif
1838        endif
1839
1840        end subroutine iso_verif_egalite_par2D
1841
1842#ifdef ISOTRAC
1843
1844      function iso_verif_traceur_choix_nostop(x,err_msg, &
1845     &       errmax,errmaxrel,ridicule_trac,deltalimtrac)
1846        use isotopes_mod, ONLY: iso_HDO
1847        implicit none
1848        ! vérifier des choses sur les traceurs
1849        ! * toutes les zones donne t l'istope total
1850        ! * pas de deltaD aberrant
1851       
1852       ! inputs
1853       real x(ntraciso)
1854       character*(*) err_msg ! message d''erreur à afficher
1855       real errmax,errmaxrel,ridicule_trac,deltalimtrac
1856
1857       ! output
1858       integer iso_verif_traceur_choix_nostop
1859
1860       ! locals
1861       !integer iso_verif_traceur_noNaN_nostop
1862       !integer iso_verif_tracm_choix_nostop
1863       !integer iso_verif_tracdD_choix_nostop
1864       !integer iso_verif_tracpos_choix_nostop
1865
1866        iso_verif_traceur_choix_nostop=0 
1867
1868        ! verif noNaN
1869        if (iso_verif_traceur_noNaN_nostop(x,err_msg).eq.1) then
1870             iso_verif_traceur_choix_nostop=1
1871        endif
1872       
1873        ! verif masse
1874        if (iso_verif_tracm_choix_nostop(x,err_msg, &
1875     &           errmax,errmaxrel).eq.1) then
1876             iso_verif_traceur_choix_nostop=1
1877        endif             
1878
1879        ! verif deltaD
1880        if (iso_HDO.gt.0) then
1881        if (iso_verif_tracdD_choix_nostop(x,err_msg, &
1882     &           ridicule_trac,deltalimtrac).eq.1) then
1883             iso_verif_traceur_choix_nostop=1
1884        endif 
1885        endif !if (iso_HDO.gt.0) then 
1886
1887        ! verif pas aberramment negatif
1888        if (iso_verif_tracpos_choix_nostop(x,err_msg, &
1889     &           1e-5).eq.1) then
1890             iso_verif_traceur_choix_nostop=1
1891        endif
1892
1893        end function iso_verif_traceur_choix_nostop
1894
1895        function iso_verif_tracnps_choix_nostop(x,err_msg, &
1896     &       errmax,errmaxrel,ridicule_trac,deltalimtrac)
1897        USE isotopes_mod, ONLY: iso_HDO
1898        implicit none
1899        ! vérifier des choses sur les traceurs
1900        ! * toutes les zones donne t l'istope total
1901        ! * pas de deltaD aberrant
1902        ! on ne vérfie pas la positivité
1903       
1904       ! inputs
1905       real x(ntraciso)
1906       character*(*) err_msg ! message d''erreur à afficher
1907       real errmax,errmaxrel,ridicule_trac,deltalimtrac
1908
1909       ! output
1910       integer iso_verif_tracnps_choix_nostop
1911
1912       ! locals
1913       !integer iso_verif_traceur_noNaN_nostop
1914       !integer iso_verif_tracm_choix_nostop
1915       !integer iso_verif_tracdD_choix_nostop
1916
1917        iso_verif_tracnps_choix_nostop=0 
1918
1919        ! verif noNaN
1920        if (iso_verif_traceur_noNaN_nostop(x,err_msg).eq.1) then
1921             iso_verif_tracnps_choix_nostop=1
1922        endif
1923       
1924        ! verif masse
1925        if (iso_verif_tracm_choix_nostop(x,err_msg, &
1926     &           errmax,errmaxrel).eq.1) then
1927             iso_verif_tracnps_choix_nostop=1
1928        endif             
1929
1930        ! verif deltaD
1931        if (iso_HDO.gt.0) then
1932        if (iso_verif_tracdD_choix_nostop(x,err_msg, &
1933     &           ridicule_trac,deltalimtrac).eq.1) then
1934             iso_verif_tracnps_choix_nostop=1
1935        endif   
1936        endif ! if (iso_HDO.gt.0) then 
1937
1938        return
1939        end function iso_verif_tracnps_choix_nostop
1940
1941        function iso_verif_tracpos_choix_nostop(x,err_msg,seuil)
1942        use isotopes_mod, only: isoName
1943        implicit none
1944
1945        ! inputs
1946       real x(ntraciso)
1947       character*(*) err_msg ! message d''erreur à afficher
1948       real seuil
1949
1950       ! output
1951       integer iso_verif_tracpos_choix_nostop
1952
1953       ! locals
1954       integer lnblnk
1955       integer iiso,ixt
1956       !integer iso_verif_positif_choix_nostop
1957
1958       iso_verif_tracpos_choix_nostop=0
1959
1960       do ixt=niso+1,ntraciso
1961          if (iso_verif_positif_choix_nostop(x(ixt),seuil,err_msg// &
1962     &           ', verif positif, iso'//TRIM(isoName(ixt))).eq.1) then
1963            iso_verif_tracpos_choix_nostop=1
1964          endif
1965        enddo
1966
1967        end function iso_verif_tracpos_choix_nostop
1968
1969
1970        function iso_verif_traceur_noNaN_nostop(x,err_msg)
1971        use isotopes_mod, only: isoName
1972        implicit none
1973
1974        ! on vérifie juste que pas NaN
1975        ! inputs
1976       real x(ntraciso)
1977       character*(*) err_msg ! message d''erreur à afficher
1978
1979       ! output
1980       integer iso_verif_traceur_noNaN_nostop
1981
1982       ! locals
1983       integer lnblnk
1984       integer iiso,ixt
1985       !integer iso_verif_nonaN_nostop
1986
1987       iso_verif_traceur_noNaN_nostop=0
1988
1989        do ixt=niso+1,ntraciso
1990!          write(*,*) 'iso_verif_traceurs 154: iiso,ixt=',iiso,ixt
1991          if (iso_verif_noNaN_nostop(x(ixt),err_msg// &
1992     &           ', verif trac no NaN, iso'//TRIM(isoName(ixt))) &
1993     &           .eq.1) then
1994            iso_verif_traceur_noNaN_nostop=1
1995          endif
1996        enddo
1997
1998        end function iso_verif_traceur_noNaN_nostop
1999
2000        function iso_verif_tracm_choix_nostop(x,err_msg, &
2001     &           errmaxin,errmaxrelin)
2002
2003        use isotopes_mod, ONLY: ridicule,isoName
2004        ! on vérifie juste bilan de masse
2005        implicit none
2006       
2007        ! inputs
2008        real x(ntraciso)
2009        character*(*) err_msg ! message d''erreur à afficher
2010        real errmaxin,errmaxrelin
2011
2012        ! output
2013        integer iso_verif_tracm_choix_nostop
2014
2015       ! locals
2016       !integer iso_verif_egalite_choix_nostop
2017       integer iiso,izone,ixt
2018       real xtractot
2019
2020       iso_verif_tracm_choix_nostop=0
2021
2022        do iiso=1,niso
2023
2024          xtractot=0.0
2025          do izone=1,nzone 
2026            ixt=itZonIso(izone,iiso)
2027            xtractot=xtractot+x(ixt)
2028          enddo
2029
2030          if (iso_verif_egalite_choix_nostop(xtractot,x(iiso), &
2031     &        err_msg//', verif trac egalite1, iso '// &
2032     &        TRIM(isoName(iiso)), &
2033     &        errmaxin,errmaxrelin).eq.1) then
2034            write(*,*) 'iso_verif_traceur 202: x=',x
2035!            write(*,*) 'xtractot=',xtractot
2036            do izone=1,nzone 
2037              ixt=itZonIso(izone,iiso)
2038              write(*,*) 'izone,iiso,ixt=',izone,iiso,ixt
2039            enddo
2040            iso_verif_tracm_choix_nostop=1
2041          endif
2042
2043          ! verif ajoutée le 19 fev 2011
2044          if ((abs(xtractot).lt.ridicule**2).and. &
2045     &           (abs(x(iiso)).gt.ridicule)) then
2046            write(*,*) err_msg,', verif masse traceurs, iso ', &
2047     &          TRIM(isoName(iiso))
2048            write(*,*) 'iso_verif_traceur 209: x=',x
2049!            iso_verif_tracm_choix_nostop=1
2050          endif
2051
2052        enddo !do iiso=1,ntraceurs_iso 
2053
2054        end function iso_verif_tracm_choix_nostop
2055
2056        function iso_verif_tracdD_choix_nostop(x,err_msg, &
2057     &           ridicule_trac,deltalimtrac)
2058
2059        USE isotopes_mod, ONLY: iso_eau, iso_HDO
2060        use isotrac_mod, only: strtrac
2061        ! on vérifie juste deltaD
2062        implicit none
2063               
2064        ! inputs
2065        real x(ntraciso)
2066        character*(*) err_msg ! message d''erreur à afficher
2067        real ridicule_trac,deltalimtrac
2068
2069        ! output
2070        integer iso_verif_tracdD_choix_nostop       
2071
2072       ! locals
2073       integer izone,ieau,ixt
2074       !integer iso_verif_aberrant_choix_nostop
2075
2076        iso_verif_tracdD_choix_nostop=0
2077
2078        if ((iso_eau.gt.0).and.(iso_HDO.gt.0)) then
2079        do izone=1,nzone
2080             ieau=itZonIso(izone,iso_eau)
2081             ixt=itZonIso(izone,iso_HDO)
2082
2083             if (iso_verif_aberrant_choix_nostop(x(ixt),x(ieau), &
2084     &           ridicule_trac,deltalimtrac,err_msg// &
2085     &           ', verif trac no aberrant zone '//strtrac(izone)) &
2086     &           .eq.1) then
2087               iso_verif_tracdD_choix_nostop=1
2088             endif
2089!             if (x(ieau).gt.ridicule) then
2090!              call iso_verif_aberrant(x(ixt)/x(ieau),
2091!     :           err_msg//', verif trac no aberrant zone '
2092!     :           //strtrac(izone))
2093!             endif
2094        enddo !do izone=1,nzone
2095       endif ! if ((iso_eau.gt.0).and.(iso_HDO.gt.0)) then
2096
2097       end function iso_verif_tracdD_choix_nostop
2098
2099INTEGER FUNCTION iso_verif_tag17_q_deltaD_chns(x,err_msg) RESULT(res)
2100  USE isotopes_mod, ONLY: iso_HDO, iso_eau, ridicule
2101  USE isotrac_mod,  ONLY: nzone_temp, option_traceurs
2102  IMPLICIT NONE
2103  REAL,             INTENT(IN) :: x(ntraciso)
2104  CHARACTER(LEN=*), INTENT(IN) :: err_msg
2105  INTEGER :: ieau, ixt, ieau1
2106  res = 0
2107  IF(ALL([17,18]/=option_traceurs)) RETURN
2108  !--- Check whether * deltaD(highest tagging layer) < 200 permil
2109  !                  * q <
2110  ieau=itZonIso(nzone_temp,iso_eau)
2111  ixt=itZonIso(nzone_temp,iso_HDO)
2112  IF(x(ieau)>ridicule) THEN
2113    IF(iso_verif_positif_nostop(-200.0-deltaD(x(ixt)/x(ieau)), err_msg//': deltaDt05 trop fort')==1) THEN
2114      res=1; write(*,*) 'x=',x
2115    END IF
2116  END IF
2117  IF(iso_verif_positif_nostop(2.0e-3-x(ieau),err_msg//': qt05 trop fort')==1) THEN
2118    res=1; write(*,*) 'x=',x
2119  END IF
2120  !--- Check whether q is small ; then, qt01 < 10%
2121  IF(x(iso_eau)<2.0e-3) THEN
2122    ieau1= itZonIso(1,iso_eau)
2123    IF(iso_verif_positif_nostop(0.1-(x(ieau1)/x(iso_eau)),err_msg//': qt01 trop abondant')==1) THEN
2124      res=1; write(*,*) 'x=',x
2125    END IF
2126  END IF
2127END FUNCTION iso_verif_tag17_q_deltaD_chns
2128
2129SUBROUTINE iso_verif_trac17_q_deltaD(x,err_msg)
2130  USE isotrac_mod,  ONLY: nzone_temp, option_traceurs
2131  IMPLICIT NONE
2132  REAL,             INTENT(IN) :: x(ntraciso)
2133  CHARACTER(LEN=*), INTENT(IN) :: err_msg
2134  IF(ALL([17,18]/=option_traceurs)) RETURN
2135  IF(nzone_temp>=5) THEN
2136    IF(iso_verif_tag17_q_deltaD_chns(x,err_msg)==1) STOP
2137  END IF
2138END SUBROUTINE iso_verif_trac17_q_deltaD
2139
2140      subroutine iso_verif_traceur(x,err_msg)
2141        use isotrac_mod, only: ridicule_trac
2142        implicit none
2143        ! vérifier des choses sur les traceurs
2144        ! * toutes les zones donne t l'istope total
2145        ! * pas de deltaD aberrant
2146
2147        ! on prend les valeurs pas défaut pour
2148        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2149       
2150       ! inputs
2151       real x(ntraciso)
2152       character*(*) err_msg ! message d''erreur à afficher
2153
2154       ! locals
2155       !integer iso_verif_traceur_choix_nostop 
2156
2157        if (iso_verif_traceur_choix_nostop(x,err_msg, &
2158     &       errmax,errmaxrel,ridicule_trac,deltalimtrac) &
2159     &       .eq.1) then
2160                stop
2161        endif
2162
2163        end subroutine iso_verif_traceur
2164
2165       
2166      subroutine iso_verif_traceur_retourne3D(x,n1,n2,n3, &
2167     &           i1,i2,i3,err_msg)
2168        use isotrac_mod, only: ridicule_trac
2169
2170        implicit none
2171        ! vérifier des choses sur les traceurs
2172        ! * toutes les zones donne t l'istope total
2173        ! * pas de deltaD aberrant
2174
2175        ! on prend les valeurs pas défaut pour
2176        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2177       
2178       ! inputs
2179       integer n1,n2,n3
2180       real x(n1,n2,n3,ntraciso)
2181       character*(*) err_msg ! message d''erreur à afficher
2182       integer i1,i2,i3
2183
2184       ! locals
2185       !integer iso_verif_traceur_choix_nostop 
2186       real xiso(ntraciso)
2187
2188        call select_dim4_from4D(n1,n2,n3,ntraciso, &
2189     &      x,xiso,i1,i2,i3)
2190        if (iso_verif_traceur_choix_nostop(xiso,err_msg, &
2191     &       errmax,errmaxrel,ridicule_trac,deltalimtrac) &
2192     &       .eq.1) then
2193                stop
2194        endif
2195
2196        end subroutine iso_verif_traceur_retourne3D
2197
2198        subroutine iso_verif_traceur_retourne4D(x,n1,n2,n3,n4, &
2199     &           i1,i2,i3,i4,err_msg)
2200        use isotrac_mod, only: ridicule_trac
2201
2202        implicit none
2203        ! vérifier des choses sur les traceurs
2204        ! * toutes les zones donne t l'istope total
2205        ! * pas de deltaD aberrant
2206
2207        ! on prend les valeurs pas défaut pour
2208        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2209       
2210       ! inputs
2211       integer n1,n2,n3,n4
2212       real x(n1,n2,n3,n4,ntraciso)
2213       character*(*) err_msg ! message d''erreur à afficher
2214       integer i1,i2,i3,i4
2215
2216       ! locals
2217       !integer iso_verif_traceur_choix_nostop 
2218       real xiso(ntraciso)
2219
2220        call select_dim5_from5D(n1,n2,n3,n4,ntraciso, &
2221     &      x,xiso,i1,i2,i3,i4)
2222        if (iso_verif_traceur_choix_nostop(xiso,err_msg, &
2223     &       errmax,errmaxrel,ridicule_trac,deltalimtrac) &
2224     &       .eq.1) then
2225                stop
2226        endif
2227
2228        end subroutine iso_verif_traceur_retourne4D
2229
2230       
2231      subroutine iso_verif_traceur_retourne2D(x,n1,n2, &
2232     &           i1,i2,err_msg)
2233        use isotrac_mod, only: ridicule_trac
2234        implicit none
2235        ! vérifier des choses sur les traceurs
2236        ! * toutes les zones donne t l'istope total
2237        ! * pas de deltaD aberrant
2238
2239        ! on prend les valeurs pas défaut pour
2240        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2241       
2242       ! inputs
2243       integer n1,n2
2244       real x(n1,n2,ntraciso)
2245       character*(*) err_msg ! message d''erreur à afficher
2246       integer i1,i2
2247
2248       ! locals
2249       !integer iso_verif_traceur_choix_nostop 
2250       real xiso(ntraciso)
2251
2252        call select_dim3_from3D(n1,n2,ntraciso, &
2253     &      x,xiso,i1,i2)
2254        if (iso_verif_traceur_choix_nostop(xiso,err_msg, &
2255     &       errmax,errmaxrel,ridicule_trac,deltalimtrac) &
2256     &       .eq.1) then
2257                stop
2258        endif
2259
2260        end subroutine iso_verif_traceur_retourne2D
2261
2262        subroutine iso_verif_traceur_vect(x,n,m,err_msg)
2263        USE isotopes_mod, ONLY: iso_HDO
2264        implicit none
2265        ! vérifier des choses sur les traceurs
2266        ! * toutes les zones donne t l'istope total
2267        ! * pas de deltaD aberrant
2268
2269        ! on prend les valeurs pas défaut pour
2270        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2271       
2272       ! inputs
2273       integer n,m
2274       real x(ntraciso,n,m)
2275       character*(*) err_msg ! message d''erreur à afficher
2276
2277       ! locals
2278       logical iso_verif_traceur_nostop
2279       integer i,j,ixt,iiso,izone,ieau
2280       integer ifaux,jfaux,ixtfaux
2281       
2282       call iso_verif_traceur_noNaN_vect(x,n,m,err_msg)       
2283
2284        ! verif masse: iso_verif_tracm_choix_nostop
2285        call iso_verif_trac_masse_vect(x,n,m,err_msg,errmax,errmaxrel)
2286
2287        ! verif deltaD: iso_verif_tracdD_choix_nostop   
2288        if (iso_HDO.gt.0) then 
2289        call iso_verif_tracdd_vect(x,n,m,err_msg)     
2290        endif !if (iso_HDO.gt.0) then       
2291
2292        ! verif pas aberramment negatif: iso_verif_tracpos_choix_nostop
2293        call iso_verif_tracpos_vect(x,n,m,err_msg,1e-5) 
2294       
2295        end subroutine iso_verif_traceur_vect
2296
2297        subroutine iso_verif_tracnps_vect(x,n,m,err_msg)
2298        USE isotopes_mod, ONLY: iso_HDO
2299        implicit none
2300        ! vérifier des choses sur les traceurs
2301        ! * toutes les zones donne t l'istope total
2302        ! * pas de deltaD aberrant
2303
2304        ! on prend les valeurs pas défaut pour
2305        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2306       
2307       ! inputs
2308       integer n,m
2309       real x(ntraciso,n,m)
2310       character*(*) err_msg ! message d''erreur à afficher
2311
2312       ! locals
2313       logical iso_verif_traceur_nostop
2314       integer i,j,ixt,iiso,izone,ieau
2315       integer ifaux,jfaux,ixtfaux
2316       
2317       call iso_verif_traceur_noNaN_vect(x,n,m,err_msg)       
2318
2319        ! verif masse: iso_verif_tracm_choix_nostop
2320        call iso_verif_trac_masse_vect(x,n,m,err_msg,errmax,errmaxrel)
2321
2322        ! verif deltaD: iso_verif_tracdD_choix_nostop   
2323        if (iso_HDO.gt.0) then 
2324        call iso_verif_tracdd_vect(x,n,m,err_msg)     
2325        endif !if (iso_HDO.gt.0) then       
2326       
2327        end subroutine iso_verif_tracnps_vect
2328
2329
2330        subroutine iso_verif_traceur_noNaN_vect(x,n,m,err_msg)
2331        implicit none
2332       
2333       ! inputs
2334       integer n,m
2335       real x(ntraciso,n,m)
2336       character*(*) err_msg ! message d''erreur à afficher
2337
2338       ! locals
2339       logical iso_verif_traceur_nostop
2340       integer i,j,ixt,iiso
2341       integer ifaux,jfaux,ixtfaux
2342
2343
2344       iso_verif_traceur_nostop=.false.
2345        ! verif noNaN: iso_verif_traceur_noNaN_nostop       
2346        do j=1,m
2347        do i=1,n
2348        do ixt=niso+1,ntraciso
2349          if ((x(ixt,i,j).gt.-borne).and.(x(ixt,i,j).lt.borne)) then
2350          else !if ((x.gt.-borne).and.(x.lt.borne)) then
2351              iso_verif_traceur_nostop=.true.
2352              ifaux=i
2353              jfaux=i
2354          endif !if ((x.gt.-borne).and.(x.lt.borne)) then
2355        enddo !do ixt=niso+1,ntraciso
2356        enddo ! do i=1,n
2357        enddo ! do j=1,m
2358       
2359
2360        if (iso_verif_traceur_nostop) then
2361            write(*,*) 'erreur detectée par iso_verif_nonNAN ', &
2362     &           'dans iso_verif_traceur_vect'
2363            write(*,*) ''
2364            write(*,*) err_msg
2365            write(*,*) 'x=',x(:,ifaux,jfaux)
2366            stop
2367        endif
2368
2369        end subroutine iso_verif_traceur_noNaN_vect
2370
2371       
2372        subroutine iso_verif_trac_masse_vect(x,n,m,err_msg, &
2373     &            errmax,errmaxrel)
2374        use isotopes_mod, only: isoName
2375        implicit none
2376       
2377        ! inputs
2378       integer n,m
2379       real x(ntraciso,n,m)
2380       character*(*) err_msg ! message d''erreur à afficher
2381       real errmax,errmaxrel
2382
2383       ! locals
2384       logical iso_verif_traceur_nostop
2385       integer i,j,ixt,iiso,izone
2386       integer ifaux,jfaux,ixtfaux
2387       real xtractot(n,m)
2388       real xiiso(n,m)
2389
2390        do iiso=1,niso       
2391        do j=1,m
2392         do i=1,n       
2393          xtractot(i,j)=0.0
2394          xiiso(i,j)=x(iiso,i,j)
2395          do izone=1,nzone
2396            ixt=itZonIso(izone,iiso)
2397            xtractot(i,j)=xtractot(i,j)+x(ixt,i,j)           
2398          enddo !do izone=1,nzone
2399         enddo !do i=1,n
2400        enddo !do j=1,m
2401       
2402
2403        call iso_verif_egalite_std_vect( &
2404     &           xtractot,xiiso, &
2405     &           err_msg//', verif trac egalite2, iso ' &
2406     &           //TRIM(isoName(iiso)), &
2407     &           n,m,errmax,errmaxrel)
2408        enddo !do iiso=1,niso
2409
2410        end subroutine iso_verif_trac_masse_vect
2411
2412        subroutine iso_verif_tracdd_vect(x,n,m,err_msg)
2413        use isotopes_mod, only: iso_HDO,iso_eau
2414        use isotrac_mod, only: strtrac
2415        implicit none
2416       
2417        ! inputs
2418       integer n,m
2419       real x(ntraciso,n,m)
2420       character*(*) err_msg ! message d''erreur à afficher
2421
2422       ! locals
2423       integer i,j,iiso,izone,ieau,ixt
2424       real xiiso(niso,n,m)
2425       real xeau(n,m)
2426       integer lnblnk
2427
2428       if (iso_HDO.gt.0) then
2429        do izone=1,nzone
2430          ieau=itZonIso(izone,iso_eau)
2431          do iiso=1,niso
2432           ixt=itZonIso(izone,iiso)
2433           do j=1,m
2434            do i=1,n
2435             xiiso(iiso,i,j)=x(ixt,i,j)
2436            enddo !do i=1,n
2437           enddo !do j=1,m
2438          enddo !do iiso=1,niso
2439           
2440          do j=1,m
2441           do i=1,n
2442            xeau(i,j)=x(ieau,i,j)
2443           enddo !do i=1,n
2444          enddo !do j=1,m
2445         
2446          call iso_verif_aberrant_vect2Dch( &
2447     &           xiiso,xeau,err_msg//strtrac(izone),niso,n,m, &
2448     &           deltalimtrac)
2449         enddo !do izone=1,nzone
2450        endif !if (iso_HDO.gt.0) then
2451
2452        end subroutine iso_verif_tracdd_vect
2453
2454        subroutine iso_verif_tracpos_vect(x,n,m,err_msg,seuil)
2455        implicit none
2456
2457       ! inputs
2458       integer n,m
2459       real x(ntraciso,n,m)
2460       character*(*) err_msg ! message d''erreur à afficher
2461       real seuil
2462
2463       ! locals
2464       integer i,j,ixt
2465       logical iso_verif_traceur_nostop
2466       integer ifaux,jfaux,ixtfaux
2467
2468        iso_verif_traceur_nostop=.false.       
2469        do j=1,m
2470        do i=1,n
2471        do ixt=niso+1,ntraciso
2472          if (x(ixt,i,j).lt.-seuil) then
2473              ifaux=i
2474              jfaux=j
2475              ixtfaux=ixt
2476              iso_verif_traceur_nostop=.true.
2477          endif
2478        enddo !do ixt=niso+1,ntraciso
2479        enddo !do i=1,n
2480        enddo !do j=1,m       
2481
2482        if (iso_verif_traceur_nostop) then
2483            write(*,*) 'erreur detectée par verif positif ', &
2484     &           'dans iso_verif_traceur_vect'
2485            write(*,*) ''
2486            write(*,*) err_msg
2487            write(*,*) 'x=',x(:,ifaux,jfaux)
2488            stop
2489        endif
2490
2491        end subroutine iso_verif_tracpos_vect
2492
2493
2494
2495        subroutine iso_verif_tracnps(x,err_msg)
2496        use isotrac_mod, only: ridicule_trac
2497
2498        implicit none
2499        ! vérifier des choses sur les traceurs
2500        ! * toutes les zones donne t l'istope total
2501        ! * pas de deltaD aberrant
2502
2503        ! on prend les valeurs pas défaut pour
2504        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2505       
2506       ! inputs
2507       real x(ntraciso)
2508       character*(*) err_msg ! message d''erreur à afficher
2509
2510       ! locals
2511       !integer iso_verif_tracnps_choix_nostop 
2512
2513        if (iso_verif_tracnps_choix_nostop(x,err_msg, &
2514     &       errmax,errmaxrel,ridicule_trac,deltalimtrac) &
2515     &       .eq.1) then
2516                stop
2517        endif
2518
2519        end subroutine iso_verif_tracnps
2520
2521        subroutine iso_verif_tracpos_choix(x,err_msg,seuil)
2522        implicit none
2523        ! vérifier des choses sur les traceurs
2524        ! * toutes les zones donne t l'istope total
2525        ! * pas de deltaD aberrant
2526
2527        ! on prend les valeurs pas défaut pour
2528        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2529       
2530       ! inputs
2531       real x(ntraciso)
2532       character*(*) err_msg ! message d''erreur à afficher
2533       real seuil
2534
2535       ! locals
2536       !integer iso_verif_tracpos_choix_nostop 
2537
2538        if (iso_verif_tracpos_choix_nostop(x,err_msg,seuil) &
2539     &       .eq.1) then
2540                stop
2541        endif
2542
2543        end subroutine iso_verif_tracpos_choix
2544
2545        subroutine iso_verif_traceur_choix(x,err_msg, &
2546     &       errmax,errmaxrel,ridicule_trac_loc,deltalimtrac)
2547        implicit none
2548        ! vérifier des choses sur les traceurs
2549        ! * toutes les zones donne t l'istope total
2550        ! * pas de deltaD aberrant
2551       
2552       ! inputs
2553       real x(ntraciso)
2554       character*(*) err_msg ! message d''erreur à afficher
2555       real errmax,errmaxrel,ridicule_trac_loc,deltalimtrac
2556
2557       ! locals
2558       !integer iso_verif_traceur_choix_nostop 
2559
2560        if (iso_verif_traceur_choix_nostop(x,err_msg, &
2561     &       errmax,errmaxrel,ridicule_trac_loc,deltalimtrac) &
2562     &       .eq.1) then
2563                stop
2564        endif
2565
2566        end subroutine iso_verif_traceur_choix
2567
2568        function iso_verif_traceur_nostop(x,err_msg)
2569        use isotrac_mod, only: ridicule_trac
2570        !use isotopes_verif, only: errmax,errmaxrel,deltalimtrac
2571        implicit none
2572        ! vérifier des choses sur les traceurs
2573        ! * toutes les zones donne t l'istope total
2574        ! * pas de deltaD aberrant
2575
2576        ! on prend les valeurs pas défaut pour
2577        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2578       
2579       ! inputs
2580       real x(ntraciso)
2581       character*(*) err_msg ! message d''erreur à afficher
2582
2583       ! output
2584       integer iso_verif_traceur_nostop
2585
2586       ! locals
2587       !integer iso_verif_traceur_choix_nostop 
2588
2589        iso_verif_traceur_nostop= &
2590     &       iso_verif_traceur_choix_nostop(x,err_msg, &
2591     &       errmax,errmaxrel,ridicule_trac,deltalimtrac)
2592
2593        end function iso_verif_traceur_nostop
2594
2595
2596      subroutine iso_verif_traceur_justmass(x,err_msg)
2597        implicit none
2598        ! on vérifie que noNaN et masse
2599
2600        ! on prend les valeurs pas défaut pour
2601        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2602       
2603       ! inputs
2604       real x(ntraciso)
2605       character*(*) err_msg ! message d''erreur à afficher
2606
2607       ! locals
2608       !integer iso_verif_traceur_noNaN_nostop
2609       !integer iso_verif_tracm_choix_nostop
2610
2611        ! verif noNaN
2612        if (iso_verif_traceur_noNaN_nostop(x,err_msg).eq.1) then
2613             stop
2614        endif
2615       
2616        ! verif masse
2617        if (iso_verif_tracm_choix_nostop(x,err_msg, &
2618     &           errmax,errmaxrel).eq.1) then
2619             stop
2620        endif   
2621       
2622        end subroutine iso_verif_traceur_justmass
2623
2624        function iso_verif_traceur_jm_nostop(x,err_msg)
2625        implicit none
2626        ! on vérifie que noNaN et masse
2627
2628        ! on prend les valeurs pas défaut pour
2629        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2630       
2631       ! inputs
2632       real x(ntraciso)
2633       character*(*) err_msg ! message d''erreur à afficher
2634
2635       ! output
2636       integer iso_verif_traceur_jm_nostop
2637
2638       ! locals
2639!       integer iso_verif_traceur_noNaN_nostop
2640       !integer iso_verif_tracm_choix_nostop
2641
2642        iso_verif_traceur_jm_nostop=0
2643!        ! verif noNaN
2644!        if (iso_verif_traceur_noNaN_nostop(x,err_msg).eq.1) then
2645!             iso_verif_traceur_jm_nostop=1
2646!        endif
2647       
2648        ! verif masse
2649        if (iso_verif_tracm_choix_nostop(x,err_msg, &
2650     &           errmax,errmaxrel).eq.1) then
2651             iso_verif_traceur_jm_nostop=1
2652        endif   
2653       
2654        end function iso_verif_traceur_jm_nostop
2655
2656        subroutine iso_verif_tag17_q_deltaD_vect(x,n,m,err_msg)
2657        USE isotopes_mod, ONLY: tnat,iso_eau, ridicule,iso_HDO
2658        use isotrac_mod, only: option_traceurs,nzone_temp
2659        implicit none
2660
2661        ! inputs
2662        integer n,m
2663        real x(ntraciso,n,m)
2664        character*(*) err_msg
2665
2666        ! locals
2667        !integer iso_verif_positif_nostop
2668        !real deltaD
2669        integer ieau,ixt,ieau1
2670        integer i,k
2671
2672        if ((option_traceurs.eq.17).or. &
2673     &           (option_traceurs.eq.18)) then
2674        ! verifier que deltaD du tag de la couche la plus haute <
2675        ! 200 permil, et vérifier que son q est inférieur à
2676        ieau=itZonIso(nzone_temp,iso_eau)
2677        ixt=itZonIso(nzone_temp,iso_HDO)
2678        ieau1=itZonIso(1,iso_eau)
2679        do i=1,n
2680         do k=1,m
2681           if (x(ieau,i,k).gt.ridicule) then
2682             if ((x(ixt,i,k)/x(ieau,i,k)/tnat(iso_HDO)-1)*1000 &
2683     &            .gt.-200.0) then
2684                write(*,*) err_msg,', vect:deltaDt05 trop fort'
2685                write(*,*) 'i,k=',i,k
2686                write(*,*) 'x(:,i,k)=',x(:,i,k)
2687                stop
2688             endif !if (iso_verif_positif_nostop(-200.0-deltaD(x(ixt),x(ieau)),
2689           endif !if (x(ieau).gt.ridicule) then
2690           if (x(ieau,i,k).gt.2.0e-3) then
2691                write(*,*) err_msg,', vect:qt05 trop fort'
2692                write(*,*) 'i,k=',i,k
2693                write(*,*) 'x(:,i,k)=',x(:,i,k)
2694                stop
2695           endif !if (iso_verif_positif_nostop(1.0e-3-x(ieau),
2696           if (x(iso_eau,i,k).lt.2.0e-3) then
2697                if (x(ieau1,i,k)/x(iso_eau,i,k).gt.0.1) then
2698                   write(*,*) err_msg,', vect: qt01 trop abondant'
2699                   write(*,*) 'i,k=',i,k
2700                   write(*,*) 'ieau1,iso_eau,x(ieau1,i,k),', &
2701     &                 'x(iso_eau,i,k)=',ieau1,iso_eau, &
2702     &                  x(ieau1,i,k),x(iso_eau,i,k) 
2703                   write(*,*) 'x(:,i,k)=',x(:,i,k)
2704                   stop
2705                endif !if (x(ieau1,i,k)/x(iso_eau,i,k).gt.0.1) then
2706            endif
2707          enddo !do k=1,m
2708        enddo !do i=1,n
2709
2710        endif !if (option_traceurs.eq.17) then
2711
2712        end subroutine iso_verif_tag17_q_deltaD_vect
2713
2714
2715        subroutine iso_verif_tag17_q_deltaD_vect_ret3D(x,n,m,nq,err_msg)
2716        USE isotopes_mod, ONLY: tnat,iso_eau,iso_HDO,ridicule
2717        use isotrac_mod, only: option_traceurs,nzone_temp
2718        implicit none
2719
2720        ! inputs
2721        integer n,m,nq
2722        real x(n,m,nq,ntraciso)
2723        character*(*) err_msg
2724
2725        ! locals
2726        !integer iso_verif_positif_nostop
2727        !real deltaD
2728        integer ieau,ixt,ieau1
2729        integer i,k,iq
2730
2731        if ((option_traceurs.eq.17).or. &
2732     &           (option_traceurs.eq.18)) then
2733        ! verifier que deltaD du tag de la couche la plus haute <
2734        ! 200 permil, et vérifier que son q est inférieur à
2735        ieau=itZonIso(nzone_temp,iso_eau)
2736        ixt=itZonIso(nzone_temp,iso_HDO)
2737        ieau1=itZonIso(1,iso_eau)
2738        do iq=1,nq
2739        do i=1,n
2740         do k=1,m
2741           if (x(i,k,iq,ieau).gt.ridicule) then
2742             if ((x(i,k,iq,ixt)/x(i,k,iq,ieau)/tnat(iso_HDO)-1)*1000 &
2743     &            .gt.-200.0) then
2744                write(*,*) err_msg,', vect:deltaDt05 trop fort'
2745                write(*,*) 'i,k=',i,k
2746                write(*,*) 'x(i,k,iq,:)=',x(i,k,iq,:)
2747                stop
2748             endif !if (iso_verif_positif_nostop(-200.0-deltaD(x(ixt),x(ieau)),
2749           endif !if (x(ieau).gt.ridicule) then
2750           if (x(i,k,iq,ieau).gt.2.0e-3) then
2751                write(*,*) err_msg,', vect:qt05 trop fort'
2752                write(*,*) 'i,k=',i,k
2753                write(*,*) 'x(i,k,iq,:)=',x(i,k,iq,:)
2754                stop
2755           endif !if (iso_verif_positif_nostop(1.0e-3-x(ieau),
2756           if (x(i,k,iq,iso_eau).lt.2.0e-3) then
2757                if (x(i,k,iq,ieau1)/x(i,k,iq,iso_eau).gt.0.1) then
2758                   write(*,*) err_msg,', vect: qt01 trop abondant'
2759                   write(*,*) 'i,k=',i,k
2760                   write(*,*) 'ieau1,iso_eau,x(i,k,iq,ieau1),', &
2761     &                 'x(i,k,iq,ieau)=',ieau1,iso_eau, &
2762     &                  x(i,k,iq,ieau1),x(i,k,iq,iso_eau) 
2763                   write(*,*) 'x(i,k,iq,:)=',x(i,k,iq,:)
2764                   stop
2765                endif !if (x(ieau1,i,k)/x(iso_eau,i,k).gt.0.1) then
2766            endif
2767          enddo !do k=1,m
2768        enddo !do i=1,n
2769        enddo ! do iq=1,nq
2770
2771        endif !if (option_traceurs.eq.17) then
2772
2773        end subroutine iso_verif_tag17_q_deltaD_vect_ret3D
2774
2775
2776#endif
2777! endif ISOTRAC
2778
2779END MODULE isotopes_verif_mod
2780
2781#endif         
2782! endif ISOVERIF
2783
Note: See TracBrowser for help on using the repository browser.