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

Last change on this file since 5456 was 4491, checked in by crisi, 21 months ago

Bug corrections in LMDZiso, especially for water tagging

  • 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!            stop
1045            iso_verif_O18_aberrant_nostop=1
1046          endif
1047
1048#ifdef ISOVERIF
1049#else
1050        write(*,*) 'err_msg=',err_msg,': pourquoi verif?'
1051        stop
1052#endif                   
1053
1054        return
1055        end function iso_verif_O18_aberrant_nostop
1056
1057
1058        ! **********
1059        function deltaD(R)
1060        USE isotopes_mod, ONLY: tnat,iso_HDO
1061        implicit none
1062        real R,deltaD
1063
1064       
1065        if (iso_HDO.gt.0) then
1066           deltaD=(R/tnat(iso_HDO)-1)*1000.0
1067        else
1068            write(*,*) 'iso_verif_egalite>deltaD 260: iso_HDO.gt.0=', &
1069     &           iso_HDO.gt.0
1070        endif
1071        return
1072        end function deltaD
1073
1074        ! **********
1075        function deltaO(R)
1076        USE isotopes_mod, ONLY: tnat,iso_O18
1077        implicit none
1078        real R,deltaO
1079       
1080        if (iso_O18.gt.0) then
1081           deltaO=(R/tnat(iso_O18)-1)*1000.0
1082        else
1083            write(*,*) 'iso_verif_egalite>deltaO18 260: iso_O18.gt.0=', &
1084     &           iso_O18.gt.0
1085        endif
1086        return
1087        end function deltaO
1088
1089        ! **********
1090        function dexcess(RD,RO)
1091        USE isotopes_mod, ONLY: tnat,iso_O18,iso_HDO
1092        implicit none
1093        real RD,RO,deltaD,deltaO,dexcess
1094       
1095        if ((iso_O18.gt.0).and.(iso_HDO.gt.0)) then
1096           deltaD=(RD/tnat(iso_HDO)-1)*1000.0
1097           deltaO=(RO/tnat(iso_O18)-1)*1000.0
1098           dexcess=deltaD-8*deltaO
1099        else
1100            write(*,*) 'iso_verif_egalite 1109: iso_O18,iso_HDO=',iso_O18,iso_HDO
1101        endif
1102        return
1103        end function dexcess
1104
1105
1106        ! **********
1107        function delta_all(R,ixt)
1108        USE isotopes_mod, ONLY: tnat
1109        implicit none
1110        real R,delta_all
1111        integer ixt
1112       
1113        delta_all=(R/tnat(ixt)-1)*1000.0
1114        return
1115        end function delta_all
1116
1117        ! **********
1118        function delta_to_R(delta,ixt)
1119        USE isotopes_mod, ONLY: tnat
1120        implicit none
1121        real delta,delta_to_R
1122        integer ixt
1123       
1124        delta_to_R=(delta/1000.0+1.0)*tnat(ixt)
1125        return
1126        end function delta_to_R
1127
1128         ! **********
1129        function o17excess(R17,R18)
1130        USE isotopes_mod, ONLY: tnat,iso_O18,iso_O17
1131        implicit none
1132        real R17,R18,o17excess
1133       
1134        if ((iso_O17.gt.0).and.(iso_O18.gt.0)) then
1135           
1136           o17excess=1e6*(log(R17/tnat(iso_o17)) &
1137     &           -0.528*log(R18/tnat(iso_O18)))
1138!           write(*,*) 'o17excess=',o17excess
1139        else
1140            write(*,*) 'iso_verif_egalite>deltaD 260: iso_O17.gt.0,18=', &
1141     &           iso_O17.gt.0,iso_O18.gt.0
1142        endif
1143        return
1144        end function o17excess
1145
1146        !       ****************
1147
1148          subroutine iso_verif_egalite_vect2D( &
1149     &           xt,q,err_msg,ni,n,m)
1150       
1151        USE isotopes_mod, ONLY: iso_eau
1152          implicit none
1153
1154          ! inputs
1155          integer n,m,ni
1156          real q(n,m)
1157          real xt(ni,n,m)
1158          character*(*) err_msg
1159
1160        ! locals
1161        integer iso_verif_egalite_nostop_loc
1162        integer i,j,ixt
1163        integer ifaux,jfaux
1164
1165        !write(*,*) 'iso_verif_egalite_vect2D 1099 tmp: q(2,1),xt(iso_eau,2,1)=',q(2,1),xt(iso_eau,2,1)
1166        !write(*,*) 'ni,n,m=',ni,n,m,errmax,errmaxrel
1167        if (iso_eau.gt.0) then
1168        iso_verif_egalite_nostop_loc=0
1169        do i=1,n
1170         do j=1,m
1171          if (abs(q(i,j)-xt(iso_eau,i,j)).gt.errmax) then
1172           if (abs((q(i,j)-xt(iso_eau,i,j))/ &
1173     &           max(max(abs(q(i,j)),abs(xt(iso_eau,i,j))),1e-18)) &
1174     &           .gt.errmaxrel) then
1175              iso_verif_egalite_nostop_loc=1
1176              ifaux=i
1177              jfaux=j
1178           endif
1179          endif
1180         enddo !do j=1,m
1181        enddo !do i=1,n
1182
1183        if (iso_verif_egalite_nostop_loc.eq.1) then
1184          write(*,*) 'erreur detectee par iso_verif_egalite_vect2D:'
1185          write(*,*) err_msg
1186          write(*,*) 'i,j=',ifaux,jfaux
1187          write(*,*) 'xt,q=',xt(iso_eau,ifaux,jfaux),q(ifaux,jfaux)
1188          stop
1189        endif
1190        endif
1191       
1192#ifdef ISOVERIF
1193        call iso_verif_noNaN_vect2D(xt,err_msg,ni,n,m)
1194#endif         
1195
1196        return
1197        end subroutine iso_verif_egalite_vect2D
1198
1199        subroutine iso_verif_egalite_vect1D( &
1200     &           xt,q,err_msg,ni,n)
1201
1202        USE isotopes_mod, ONLY: iso_eau
1203        implicit none
1204
1205        ! inputs
1206        integer n,ni
1207        real q(n)
1208        real xt(ni,n)
1209        character*(*) err_msg
1210
1211        ! locals
1212        integer iso_verif_egalite_nostop_loc
1213        integer i
1214        integer ifaux
1215
1216        if (iso_eau.gt.0) then
1217        iso_verif_egalite_nostop_loc=0
1218        do i=1,n
1219          if (abs(q(i)-xt(iso_eau,i)).gt.errmax) then
1220           if (abs((q(i)-xt(iso_eau,i))/ &
1221     &           max(max(abs(q(i)),abs(xt(iso_eau,i))),1e-18)) &
1222     &           .gt.errmaxrel) then
1223              iso_verif_egalite_nostop_loc=1
1224              ifaux=i
1225           endif !if (abs((q(i)-xt(iso_eau,i))/
1226          endif !if (abs(q(i)-xt(iso_eau,i)).gt.errmax) then
1227        enddo !do i=1,n
1228
1229        if (iso_verif_egalite_nostop_loc.eq.1) then
1230          write(*,*) 'erreur detectee par iso_verif_egalite_vect2D:'
1231          write(*,*) err_msg
1232          write(*,*) 'i=',ifaux
1233          write(*,*) 'xt,q=',xt(iso_eau,ifaux),q(ifaux)
1234          stop
1235        endif  !if (iso_verif_egalite_nostop.eq.1) then
1236        endif !if (iso_eau.gt.0) then
1237
1238        end subroutine iso_verif_egalite_vect1D       
1239
1240        subroutine iso_verif_egalite_std_vect( &
1241     &           a,b,err_msg,n,m,errmax,errmaxrel)
1242
1243          implicit none
1244
1245          ! inputs
1246          integer n,m,ni
1247          real a(n,m)
1248          real b(n,m)
1249          character*(*) err_msg
1250          real errmax,errmaxrel
1251
1252        ! locals
1253        integer iso_verif_egalite_nostop_loc
1254        integer i,j
1255        integer ifaux,jfaux
1256
1257        iso_verif_egalite_nostop_loc=0
1258        do i=1,n
1259         do j=1,m
1260          if (abs(a(i,j)-b(i,j)).gt.errmax) then
1261           if (abs((a(i,j)-b(i,j))/ &
1262     &           max(max(abs(a(i,j)),abs(b(i,j))),1e-18)) &
1263     &           .gt.errmaxrel) then
1264              iso_verif_egalite_nostop_loc=1
1265              ifaux=i
1266              jfaux=j
1267           endif
1268          endif
1269         enddo !do j=1,m
1270        enddo !do i=1,n
1271
1272        if (iso_verif_egalite_nostop_loc.eq.1) then
1273          write(*,*) 'erreur detectee par iso_verif_egalite_vect2D:'
1274          write(*,*) err_msg
1275          write(*,*) 'i,j=',ifaux,jfaux
1276          write(*,*) 'a,b=',a(ifaux,jfaux),b(ifaux,jfaux)
1277          stop
1278        endif
1279
1280        return
1281        end subroutine iso_verif_egalite_std_vect
1282
1283        subroutine iso_verif_aberrant_vect2D( &
1284     &           xt,q,err_msg,ni,n,m)
1285        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO
1286          implicit none
1287
1288          ! inputs   
1289          integer n,m,ni
1290          real q(n,m)
1291          real xt(ni,n,m)
1292          character*(*) err_msg
1293
1294        ! locals
1295        integer iso_verif_aberrant_nostop_loc
1296        integer i,j
1297        integer ifaux,jfaux
1298        !real deltaD
1299
1300        if (iso_HDO.gt.0) then
1301        iso_verif_aberrant_nostop_loc=0
1302        do i=1,n
1303         do j=1,m
1304          if (q(i,j).gt.ridicule) then
1305            if (((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1306     &                   .gt.deltalim).or. &
1307     &          ((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1308     &                   .lt.-borne)) then   
1309              iso_verif_aberrant_nostop_loc=1
1310              ifaux=i
1311              jfaux=j
1312           endif
1313          endif
1314         enddo !do j=1,m
1315        enddo !do i=1,n
1316
1317        if (iso_verif_aberrant_nostop_loc.eq.1) then
1318          write(*,*) 'erreur detectee par iso_verif_aberrant_vect2D:'
1319          write(*,*) err_msg
1320          write(*,*) 'i,j=',ifaux,jfaux
1321          write(*,*) 'deltaD=',deltaD(xt(iso_HDO,ifaux,jfaux) &
1322     &           /q(ifaux,jfaux))
1323          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1324          stop
1325        endif 
1326        endif !if (iso_HDO.gt.0) then
1327
1328        end subroutine iso_verif_aberrant_vect2D       
1329
1330        subroutine iso_verif_aberrant_enc_vect2D( &
1331     &           xt,q,err_msg,ni,n,m)
1332
1333        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO
1334          implicit none
1335
1336          ! inputs   
1337          integer n,m,ni
1338          real q(n,m)
1339          real xt(ni,n,m)
1340          character*(*) err_msg
1341
1342        ! locals
1343        integer iso_verif_aberrant_nostop_loc
1344        integer i,j
1345        integer ifaux,jfaux
1346        !real deltaD
1347
1348        if (iso_HDO.gt.0) then
1349        iso_verif_aberrant_nostop_loc=0
1350        do i=1,n
1351         do j=1,m
1352          if (q(i,j).gt.ridicule) then
1353            if (((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1354     &                   .gt.deltalim).or. &
1355     &          ((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1356     &                   .lt.deltaDmin).or. &
1357     &           (xt(iso_HDO,i,j).lt.-borne).or. &
1358     &           (xt(iso_HDO,i,j).gt.borne)) then     
1359              iso_verif_aberrant_nostop_loc=1
1360              ifaux=i
1361              jfaux=j
1362           endif
1363          endif
1364         enddo !do j=1,m
1365        enddo !do i=1,n
1366
1367        if (iso_verif_aberrant_nostop_loc.eq.1) then
1368          write(*,*) 'erreur detectee par ', &
1369     &           'iso_verif_aberrant_enc_vect2D:'
1370          write(*,*) err_msg
1371          write(*,*) 'i,j=',ifaux,jfaux
1372          write(*,*) 'deltaD=',deltaD(xt(iso_HDO,ifaux,jfaux) &
1373     &           /q(ifaux,jfaux))
1374          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1375          write(*,*) 'q(ifaux,jfaux)=',q(ifaux,jfaux)
1376          call abort_physic('isotopes_verif_mod','iso_verif_aberrant_enc_vect2D',1)
1377        endif 
1378        endif !if (iso_HDO.gt.0) then
1379
1380        end subroutine iso_verif_aberrant_enc_vect2D       
1381
1382       
1383        subroutine iso_verif_aberrant_enc_vect2D_ns( &
1384     &           xt,q,err_msg,ni,n,m)
1385
1386        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO
1387          implicit none
1388
1389          ! inputs   
1390          integer n,m,ni
1391          real q(n,m)
1392          real xt(ni,n,m)
1393          character*(*) err_msg
1394
1395        ! locals
1396        integer iso_verif_aberrant_nostop_loc
1397        integer i,j
1398        integer ifaux,jfaux
1399        !real deltaD
1400
1401        if (iso_HDO.gt.0) then
1402        iso_verif_aberrant_nostop_loc=0
1403        do i=1,n
1404         do j=1,m
1405          if (q(i,j).gt.ridicule) then
1406            if (((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1407     &                   .gt.deltalim).or. &
1408     &          ((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1409     &                   .lt.deltaDmin)) then   
1410              iso_verif_aberrant_nostop_loc=1
1411              ifaux=i
1412              jfaux=j
1413           endif
1414          endif
1415         enddo !do j=1,m
1416        enddo !do i=1,n
1417
1418        if (iso_verif_aberrant_nostop_loc.eq.1) then
1419          write(*,*) 'erreur detectee par ', &
1420     &           'iso_verif_aberrant_vect2D_ns:'
1421          write(*,*) err_msg
1422          write(*,*) 'i,j=',ifaux,jfaux
1423          write(*,*) 'deltaD=',deltaD(xt(iso_HDO,ifaux,jfaux) &
1424     &           /q(ifaux,jfaux))
1425          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1426!          stop
1427        endif 
1428        endif !if (iso_HDO.gt.0) then
1429
1430        end subroutine iso_verif_aberrant_enc_vect2D_ns       
1431
1432
1433         subroutine iso_verif_aberrant_vect2Dch( &
1434     &           xt,q,err_msg,ni,n,m,deltaDmax)
1435
1436        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO
1437          implicit none
1438
1439
1440          ! inputs   
1441          integer n,m,ni
1442          real q(n,m)
1443          real xt(ni,n,m)
1444          character*(*) err_msg
1445          real deltaDmax
1446
1447        ! locals
1448        integer iso_verif_aberrant_nostop_loc
1449        integer i,j
1450        integer ifaux,jfaux
1451        !real deltaD
1452
1453        if (iso_HDO.gt.0) then
1454        iso_verif_aberrant_nostop_loc=0
1455        do i=1,n
1456         do j=1,m
1457          if (q(i,j).gt.ridicule) then
1458            if (((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1459     &                   .gt.deltaDmax).or. &
1460     &          ((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1461     &                   .lt.-borne)) then   
1462              iso_verif_aberrant_nostop_loc=1
1463              ifaux=i
1464              jfaux=j
1465           endif
1466          endif
1467         enddo !do j=1,m
1468        enddo !do i=1,n
1469
1470        if (iso_verif_aberrant_nostop_loc.eq.1) then
1471          write(*,*) 'erreur detectee par iso_verif_aberrant_vect2D:'
1472          write(*,*) err_msg
1473          write(*,*) 'i,j=',ifaux,jfaux
1474          write(*,*) 'deltaD=',deltaD(xt(iso_HDO,ifaux,jfaux) &
1475     &           /q(ifaux,jfaux))
1476          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1477          stop
1478        endif 
1479        endif !if (iso_HDO.gt.0) then
1480
1481        end subroutine iso_verif_aberrant_vect2Dch     
1482
1483        subroutine iso_verif_O18_aberrant_enc_vect2D( &
1484     &           xt,q,err_msg,ni,n,m)
1485
1486        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO,iso_O18
1487          implicit none
1488
1489          ! inputs   
1490          integer n,m,ni
1491          real q(n,m)
1492          real xt(ni,n,m)
1493          character*(*) err_msg
1494
1495        ! locals
1496        integer iso_verif_aberrant_nostop_loc
1497        integer i,j
1498        integer ifaux,jfaux
1499        real deltaDloc,deltaoloc,dexcessloc
1500
1501        if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then
1502        iso_verif_aberrant_nostop_loc=0
1503        do i=1,n
1504         do j=1,m
1505          if (q(i,j).gt.ridicule) then
1506
1507            deltaDloc=(xt(iso_HDO,i,j)/q(i,j)/tnat(iso_hdo)-1)*1000
1508            deltaoloc=(xt(iso_O18,i,j)/q(i,j)/tnat(iso_O18)-1)*1000
1509            dexcessloc=deltaDloc-8*deltaoloc
1510            if ((deltaDloc.lt.deltaDmin).or.(deltaoloc.lt.deltaDmin/2.0).or. &
1511     &        (deltaDloc.gt.deltalim).or.(deltaoloc.gt.deltalim/8.0).or. &
1512     &        ((deltaDloc.gt.-500.0).and.((dexcessloc.lt.dexcess_min) &
1513     &        .or.(dexcessloc.gt.dexcess_max)))) then
1514              iso_verif_aberrant_nostop_loc=1
1515              ifaux=i
1516              jfaux=j
1517              write(*,*) 'deltaD,deltao,dexcess=',deltaDloc,deltaoloc,dexcessloc
1518           endif
1519          endif
1520         enddo !do j=1,m
1521        enddo !do i=1,n
1522
1523        if (iso_verif_aberrant_nostop_loc.eq.1) then
1524          write(*,*) 'erreur detectee par ', &
1525     &           'iso_verif_aberrant_enc_vect2D:'
1526          write(*,*) err_msg
1527          write(*,*) 'i,j=',ifaux,jfaux
1528          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1529          write(*,*) 'q(ifaux,jfaux)=',q(ifaux,jfaux)
1530          call abort_physic('isotopes_verif_mod','iso_verif_aberrant_enc_vect2D',1)
1531        endif 
1532        endif !if (iso_HDO.gt.0) then
1533
1534        end subroutine iso_verif_O18_aberrant_enc_vect2D   
1535
1536
1537        subroutine select_dim23_from4D(n1,n2,n3,n4, &
1538     &          var,vec,i1,i4)
1539        implicit none
1540
1541        ! inputs
1542        integer n1,n2,n3,n4
1543        real var(n1,n2,n3,n4)
1544        integer i1,i4
1545        ! outputs
1546        real vec(n2,n3)
1547        ! locals
1548        integer i2,i3
1549
1550        do i2=1,n2
1551         do i3=1,n3
1552          vec(i2,i3)=var(i1,i2,i3,i4)
1553         enddo
1554        enddo
1555
1556        return
1557        end subroutine select_dim23_from4D
1558
1559       
1560        subroutine select_dim4_from4D(ntime,nlev,nlat,nlon, &
1561     &          var,vec,itime,ilev,ilat)
1562        implicit none
1563
1564        ! inputs
1565        integer ntime,nlev,nlat,nlon
1566        real var(ntime,nlev,nlat,nlon)
1567        integer itime,ilev,ilat
1568        ! outputs
1569        real vec(nlon)
1570        ! locals
1571        integer ilon
1572
1573        do ilon=1,nlon
1574          vec(ilon)=var(itime,ilev,ilat,ilon)
1575        enddo
1576
1577        return
1578        end subroutine select_dim4_from4D
1579
1580        subroutine select_dim5_from5D(n1,n2,n3,n4,n5, &
1581     &          var,vec,i1,i2,i3,i4)
1582        implicit none
1583
1584        ! inputs
1585        integer n1,n2,n3,n4,n5
1586        real var(n1,n2,n3,n4,n5)
1587        integer i1,i2,i3,i4
1588        ! outputs
1589        real vec(n5)
1590        ! locals
1591        integer i5
1592
1593        do i5=1,n5
1594          vec(i5)=var(i1,i2,i3,i4,i5)
1595        enddo
1596
1597        end subroutine select_dim5_from5D
1598
1599       
1600        subroutine select_dim3_from3D(ntime,nlat,nlon, &
1601     &          var,vec,itime,ilat)
1602        implicit none
1603
1604        ! inputs
1605        integer ntime,nlat,nlon
1606        real var(ntime,nlat,nlon)
1607        integer itime,ilat
1608        ! outputs
1609        real vec(nlon)
1610        ! locals
1611        integer ilon
1612
1613        do ilon=1,nlon
1614          vec(ilon)=var(itime,ilat,ilon)
1615        enddo
1616
1617        end subroutine select_dim3_from3D
1618
1619       
1620        subroutine select_dim23_from3D(n1,n2,n3, &
1621     &          var,vec,i1)
1622        implicit none
1623
1624        ! inputs
1625        integer n1,n2,n3
1626        real var(n1,n2,n3)
1627        integer i1
1628        ! outputs
1629        real vec(n2,n3)
1630        ! locals
1631        integer i2,i3
1632
1633        do i2=1,n2
1634         do i3=1,n3
1635          vec(i2,i3)=var(i1,i2,i3)
1636         enddo
1637        enddo
1638
1639        end subroutine select_dim23_from3D
1640
1641        subroutine putinto_dim23_from4D(n1,n2,n3,n4, &
1642     &          var,vec,i1,i4)
1643        implicit none
1644
1645        ! inputs
1646        integer n1,n2,n3,n4
1647        real vec(n2,n3)
1648        integer i1,i4
1649        ! inout
1650        real var(n1,n2,n3,n4)
1651        ! locals
1652        integer i2,i3
1653
1654       do i2=1,n2
1655        do i3=1,n3
1656          var(i1,i2,i3,i4)=vec(i2,i3)
1657         enddo
1658        enddo
1659
1660        end subroutine putinto_dim23_from4D
1661
1662       
1663        subroutine putinto_dim12_from4D(n1,n2,n3,n4, &
1664     &          var,vec,i3,i4)
1665        implicit none
1666
1667        ! inputs
1668        integer n1,n2,n3,n4
1669        real vec(n1,n2)
1670        integer i3,i4
1671        ! inout
1672        real var(n1,n2,n3,n4)
1673        ! locals
1674        integer i1,i2
1675
1676       do i1=1,n1
1677        do i2=1,n2
1678          var(i1,i2,i3,i4)=vec(i1,i2)
1679         enddo
1680        enddo
1681
1682        end subroutine putinto_dim12_from4D
1683       
1684        subroutine putinto_dim23_from3D(n1,n2,n3, &
1685     &          var,vec,i1)
1686        implicit none
1687
1688        ! inputs
1689        integer n1,n2,n3
1690        real vec(n2,n3)
1691        integer i1
1692        ! inout
1693        real var(n1,n2,n3)
1694        ! locals
1695        integer i2,i3
1696
1697       do i2=1,n2
1698        do i3=1,n3
1699          var(i1,i2,i3)=vec(i2,i3)
1700         enddo
1701        enddo
1702
1703        end subroutine putinto_dim23_from3D
1704
1705       
1706
1707        subroutine iso_verif_noNaN_par2D(x,err_msg,ni,n,m,ib,ie)
1708        implicit none
1709        ! si x est NaN, on affiche message
1710        ! d'erreur et return 1 si erreur
1711
1712        ! input:
1713          integer n,m,ni,ib,ie
1714        real x(ni,n,m)
1715        character*(*) err_msg ! message d''erreur à afficher
1716
1717        ! output
1718
1719        ! locals       
1720        integer i,j,ixt
1721
1722      do i=ib,ie
1723       do j=1,m
1724        do ixt=1,ni
1725         if ((x(ixt,i,j).gt.-borne).and. &
1726     &            (x(ixt,i,j).lt.borne)) then
1727         else !if ((x(ixt,i,j).gt.-borne).and.
1728            write(*,*) 'erreur detectee par iso_verif_nonNAN:'
1729            write(*,*) err_msg
1730            write(*,*) 'x,ixt,i,j=',x(ixt,i,j),ixt,i,j
1731            stop
1732         endif  !if ((x(ixt,i,j).gt.-borne).and.
1733        enddo !do ixt=1,ni
1734       enddo !do j=1,m
1735      enddo !do i=1,n     
1736
1737#ifdef ISOVERIF
1738#else
1739        write(*,*) 'err_msg iso1772=',err_msg,': pourquoi verif?'
1740        stop
1741#endif           
1742
1743        return
1744        end subroutine iso_verif_noNaN_par2D
1745
1746       
1747        subroutine iso_verif_aberrant_enc_par2D( &
1748     &           xt,q,err_msg,ni,n,m,ib,ie)
1749
1750        use isotopes_mod, ONLY: ridicule,tnat,iso_HDO
1751          implicit none
1752
1753          ! inputs   
1754          integer n,m,ni,ib,ie
1755          real q(n,m)
1756          real xt(ni,n,m)
1757          character*(*) err_msg
1758
1759        ! locals
1760        integer iso_verif_aberrant_nostop_loc
1761        integer i,j
1762        integer ifaux,jfaux
1763        !real deltaD
1764
1765        if (iso_HDO.gt.0) then
1766        iso_verif_aberrant_nostop_loc=0
1767        do i=ib,ie
1768         do j=1,m
1769          if (q(i,j).gt.ridicule) then
1770            if (((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1771     &                   .gt.deltalim).or. &
1772     &          ((xt(iso_HDO,i,j)/q(i,j)/tnat(iso_HDO)-1)*1000.0 &
1773     &                   .lt.deltaDmin)) then   
1774              iso_verif_aberrant_nostop_loc=1
1775              ifaux=i
1776              jfaux=j
1777           endif
1778          endif
1779         enddo !do j=1,m
1780        enddo !do i=1,n
1781
1782        if (iso_verif_aberrant_nostop_loc.eq.1) then
1783          write(*,*) 'erreur detectee par iso_verif_aberrant_par2D:'
1784          write(*,*) err_msg
1785          write(*,*) 'i,j=',ifaux,jfaux
1786          write(*,*) 'deltaD=',deltaD(xt(iso_HDO,ifaux,jfaux) &
1787     &           /q(ifaux,jfaux))
1788          write(*,*) 'xt(:,ifaux,jfaux)=',xt(:,ifaux,jfaux)
1789          write(*,*) 'q(ifaux,jfaux)=',q(ifaux,jfaux)
1790          stop
1791        endif 
1792        endif !if (iso_HDO.gt.0) then
1793
1794        end subroutine iso_verif_aberrant_enc_par2D       
1795
1796       
1797          subroutine iso_verif_egalite_par2D( &
1798     &           xt,q,err_msg,ni,n,m,ib,ie)
1799       
1800        USE isotopes_mod, ONLY: iso_eau
1801          implicit none
1802
1803          ! inputs
1804          integer n,m,ni,ib,ie
1805          real q(n,m)
1806          real xt(ni,n,m)
1807          character*(*) err_msg
1808
1809        ! locals
1810        integer iso_verif_egalite_nostop_loc
1811        integer i,j
1812        integer ifaux,jfaux
1813
1814        if (iso_eau.gt.0) then
1815        iso_verif_egalite_nostop_loc=0
1816        do i=ib,ie
1817         do j=1,m
1818          if (abs(q(i,j)-xt(iso_eau,i,j)).gt.errmax) then
1819           if (abs((q(i,j)-xt(iso_eau,i,j))/ &
1820     &           max(max(abs(q(i,j)),abs(xt(iso_eau,i,j))),1e-18)) &
1821     &           .gt.errmaxrel) then
1822              iso_verif_egalite_nostop_loc=1
1823              ifaux=i
1824              jfaux=j
1825           endif
1826          endif
1827         enddo !do j=1,m
1828        enddo !do i=1,n
1829
1830        if (iso_verif_egalite_nostop_loc.eq.1) then
1831          write(*,*) 'erreur detectee par iso_verif_egalite_vect2D:'
1832          write(*,*) err_msg
1833          write(*,*) 'i,j=',ifaux,jfaux
1834          write(*,*) 'xt,q=',xt(iso_eau,ifaux,jfaux),q(ifaux,jfaux)
1835          stop
1836        endif
1837        endif
1838
1839        end subroutine iso_verif_egalite_par2D
1840
1841#ifdef ISOTRAC
1842
1843      function iso_verif_traceur_choix_nostop(x,err_msg, &
1844     &       errmax,errmaxrel,ridicule_trac,deltalimtrac)
1845        use isotopes_mod, ONLY: iso_HDO
1846        implicit none
1847        ! vérifier des choses sur les traceurs
1848        ! * toutes les zones donne t l'istope total
1849        ! * pas de deltaD aberrant
1850       
1851       ! inputs
1852       real x(ntraciso)
1853       character*(*) err_msg ! message d''erreur à afficher
1854       real errmax,errmaxrel,ridicule_trac,deltalimtrac
1855
1856       ! output
1857       integer iso_verif_traceur_choix_nostop
1858
1859       ! locals
1860       !integer iso_verif_traceur_noNaN_nostop
1861       !integer iso_verif_tracm_choix_nostop
1862       !integer iso_verif_tracdD_choix_nostop
1863       !integer iso_verif_tracpos_choix_nostop
1864
1865        iso_verif_traceur_choix_nostop=0 
1866
1867        ! verif noNaN
1868        if (iso_verif_traceur_noNaN_nostop(x,err_msg).eq.1) then
1869             iso_verif_traceur_choix_nostop=1
1870        endif
1871       
1872        ! verif masse
1873        if (iso_verif_tracm_choix_nostop(x,err_msg, &
1874     &           errmax,errmaxrel).eq.1) then
1875             iso_verif_traceur_choix_nostop=1
1876        endif             
1877
1878        ! verif deltaD
1879        if (iso_HDO.gt.0) then
1880        if (iso_verif_tracdD_choix_nostop(x,err_msg, &
1881     &           ridicule_trac,deltalimtrac).eq.1) then
1882             iso_verif_traceur_choix_nostop=1
1883        endif 
1884        endif !if (iso_HDO.gt.0) then 
1885
1886        ! verif pas aberramment negatif
1887        if (iso_verif_tracpos_choix_nostop(x,err_msg, &
1888     &           1e-5).eq.1) then
1889             iso_verif_traceur_choix_nostop=1
1890        endif
1891
1892        end function iso_verif_traceur_choix_nostop
1893
1894        function iso_verif_tracnps_choix_nostop(x,err_msg, &
1895     &       errmax,errmaxrel,ridicule_trac,deltalimtrac)
1896        USE isotopes_mod, ONLY: iso_HDO
1897        implicit none
1898        ! vérifier des choses sur les traceurs
1899        ! * toutes les zones donne t l'istope total
1900        ! * pas de deltaD aberrant
1901        ! on ne vérfie pas la positivité
1902       
1903       ! inputs
1904       real x(ntraciso)
1905       character*(*) err_msg ! message d''erreur à afficher
1906       real errmax,errmaxrel,ridicule_trac,deltalimtrac
1907
1908       ! output
1909       integer iso_verif_tracnps_choix_nostop
1910
1911       ! locals
1912       !integer iso_verif_traceur_noNaN_nostop
1913       !integer iso_verif_tracm_choix_nostop
1914       !integer iso_verif_tracdD_choix_nostop
1915
1916        iso_verif_tracnps_choix_nostop=0 
1917
1918        ! verif noNaN
1919        if (iso_verif_traceur_noNaN_nostop(x,err_msg).eq.1) then
1920             iso_verif_tracnps_choix_nostop=1
1921        endif
1922       
1923        ! verif masse
1924        if (iso_verif_tracm_choix_nostop(x,err_msg, &
1925     &           errmax,errmaxrel).eq.1) then
1926             iso_verif_tracnps_choix_nostop=1
1927        endif             
1928
1929        ! verif deltaD
1930        if (iso_HDO.gt.0) then
1931        if (iso_verif_tracdD_choix_nostop(x,err_msg, &
1932     &           ridicule_trac,deltalimtrac).eq.1) then
1933             iso_verif_tracnps_choix_nostop=1
1934        endif   
1935        endif ! if (iso_HDO.gt.0) then 
1936
1937        return
1938        end function iso_verif_tracnps_choix_nostop
1939
1940        function iso_verif_tracpos_choix_nostop(x,err_msg,seuil)
1941        use isotopes_mod, only: isoName
1942        implicit none
1943
1944        ! inputs
1945       real x(ntraciso)
1946       character*(*) err_msg ! message d''erreur à afficher
1947       real seuil
1948
1949       ! output
1950       integer iso_verif_tracpos_choix_nostop
1951
1952       ! locals
1953       integer lnblnk
1954       integer iiso,ixt
1955       !integer iso_verif_positif_choix_nostop
1956
1957       iso_verif_tracpos_choix_nostop=0
1958
1959       do ixt=niso+1,ntraciso
1960          if (iso_verif_positif_choix_nostop(x(ixt),seuil,err_msg// &
1961     &           ', verif positif, iso'//TRIM(isoName(ixt))).eq.1) then
1962            iso_verif_tracpos_choix_nostop=1
1963          endif
1964        enddo
1965
1966        end function iso_verif_tracpos_choix_nostop
1967
1968
1969        function iso_verif_traceur_noNaN_nostop(x,err_msg)
1970        use isotopes_mod, only: isoName
1971        implicit none
1972
1973        ! on vérifie juste que pas NaN
1974        ! inputs
1975       real x(ntraciso)
1976       character*(*) err_msg ! message d''erreur à afficher
1977
1978       ! output
1979       integer iso_verif_traceur_noNaN_nostop
1980
1981       ! locals
1982       integer lnblnk
1983       integer iiso,ixt
1984       !integer iso_verif_nonaN_nostop
1985
1986       iso_verif_traceur_noNaN_nostop=0
1987
1988        do ixt=niso+1,ntraciso
1989!          write(*,*) 'iso_verif_traceurs 154: iiso,ixt=',iiso,ixt
1990          if (iso_verif_noNaN_nostop(x(ixt),err_msg// &
1991     &           ', verif trac no NaN, iso'//TRIM(isoName(ixt))) &
1992     &           .eq.1) then
1993            iso_verif_traceur_noNaN_nostop=1
1994          endif
1995        enddo
1996
1997        end function iso_verif_traceur_noNaN_nostop
1998
1999        function iso_verif_tracm_choix_nostop(x,err_msg, &
2000     &           errmaxin,errmaxrelin)
2001
2002        use isotopes_mod, ONLY: ridicule,isoName
2003        ! on vérifie juste bilan de masse
2004        implicit none
2005       
2006        ! inputs
2007        real x(ntraciso)
2008        character*(*) err_msg ! message d''erreur à afficher
2009        real errmaxin,errmaxrelin
2010
2011        ! output
2012        integer iso_verif_tracm_choix_nostop
2013
2014       ! locals
2015       !integer iso_verif_egalite_choix_nostop
2016       integer iiso,izone,ixt
2017       real xtractot
2018
2019       iso_verif_tracm_choix_nostop=0
2020
2021        do iiso=1,niso
2022
2023          xtractot=0.0
2024          do izone=1,nzone 
2025            ixt=itZonIso(izone,iiso)
2026            xtractot=xtractot+x(ixt)
2027          enddo
2028
2029          if (iso_verif_egalite_choix_nostop(xtractot,x(iiso), &
2030     &        err_msg//', verif trac egalite1, iso '// &
2031     &        TRIM(isoName(iiso)), &
2032     &        errmaxin,errmaxrelin).eq.1) then
2033            write(*,*) 'iso_verif_traceur 202: x=',x
2034!            write(*,*) 'xtractot=',xtractot
2035            do izone=1,nzone 
2036              ixt=itZonIso(izone,iiso)
2037              write(*,*) 'izone,iiso,ixt=',izone,iiso,ixt
2038            enddo
2039            iso_verif_tracm_choix_nostop=1
2040          endif
2041
2042          ! verif ajoutée le 19 fev 2011
2043          if ((abs(xtractot).lt.ridicule**2).and. &
2044     &           (abs(x(iiso)).gt.ridicule)) then
2045            write(*,*) err_msg,', verif masse traceurs, iso ', &
2046     &          TRIM(isoName(iiso))
2047            write(*,*) 'iso_verif_traceur 209: x=',x
2048!            iso_verif_tracm_choix_nostop=1
2049          endif
2050
2051        enddo !do iiso=1,ntraceurs_iso 
2052
2053        end function iso_verif_tracm_choix_nostop
2054
2055        function iso_verif_tracdD_choix_nostop(x,err_msg, &
2056     &           ridicule_trac,deltalimtrac)
2057
2058        USE isotopes_mod, ONLY: iso_eau, iso_HDO
2059        use isotrac_mod, only: strtrac
2060        ! on vérifie juste deltaD
2061        implicit none
2062               
2063        ! inputs
2064        real x(ntraciso)
2065        character*(*) err_msg ! message d''erreur à afficher
2066        real ridicule_trac,deltalimtrac
2067
2068        ! output
2069        integer iso_verif_tracdD_choix_nostop       
2070
2071       ! locals
2072       integer izone,ieau,ixt
2073       !integer iso_verif_aberrant_choix_nostop
2074
2075        iso_verif_tracdD_choix_nostop=0
2076
2077        if ((iso_eau.gt.0).and.(iso_HDO.gt.0)) then
2078        do izone=1,nzone
2079             ieau=itZonIso(izone,iso_eau)
2080             ixt=itZonIso(izone,iso_HDO)
2081
2082             if (iso_verif_aberrant_choix_nostop(x(ixt),x(ieau), &
2083     &           ridicule_trac,deltalimtrac,err_msg// &
2084     &           ', verif trac no aberrant zone '//strtrac(izone)) &
2085     &           .eq.1) then
2086               iso_verif_tracdD_choix_nostop=1
2087             endif
2088!             if (x(ieau).gt.ridicule) then
2089!              call iso_verif_aberrant(x(ixt)/x(ieau),
2090!     :           err_msg//', verif trac no aberrant zone '
2091!     :           //strtrac(izone))
2092!             endif
2093        enddo !do izone=1,nzone
2094       endif ! if ((iso_eau.gt.0).and.(iso_HDO.gt.0)) then
2095
2096       end function iso_verif_tracdD_choix_nostop
2097
2098INTEGER FUNCTION iso_verif_tag17_q_deltaD_chns(x,err_msg) RESULT(res)
2099  USE isotopes_mod, ONLY: iso_HDO, iso_eau, ridicule
2100  USE isotrac_mod,  ONLY: nzone_temp, option_traceurs
2101  IMPLICIT NONE
2102  REAL,             INTENT(IN) :: x(ntraciso)
2103  CHARACTER(LEN=*), INTENT(IN) :: err_msg
2104  INTEGER :: ieau, ixt, ieau1
2105  res = 0
2106  IF(ALL([17,18]/=option_traceurs)) RETURN
2107  !--- Check whether * deltaD(highest tagging layer) < 200 permil
2108  !                  * q <
2109  ieau=itZonIso(nzone_temp,iso_eau)
2110  ixt=itZonIso(nzone_temp,iso_HDO)
2111  IF(x(ieau)>ridicule) THEN
2112    IF(iso_verif_positif_nostop(-200.0-deltaD(x(ixt)/x(ieau)), err_msg//': deltaDt05 trop fort')==1) THEN
2113      res=1; write(*,*) 'x=',x
2114    END IF
2115  END IF
2116  IF(iso_verif_positif_nostop(2.0e-3-x(ieau),err_msg//': qt05 trop fort')==1) THEN
2117    res=1; write(*,*) 'x=',x
2118  END IF
2119  !--- Check whether q is small ; then, qt01 < 10%
2120  IF(x(iso_eau)<2.0e-3) THEN
2121    ieau1= itZonIso(1,iso_eau)
2122    IF(iso_verif_positif_nostop(0.1-(x(ieau1)/x(iso_eau)),err_msg//': qt01 trop abondant')==1) THEN
2123      res=1; write(*,*) 'x=',x
2124    END IF
2125  END IF
2126END FUNCTION iso_verif_tag17_q_deltaD_chns
2127
2128SUBROUTINE iso_verif_trac17_q_deltaD(x,err_msg)
2129  USE isotrac_mod,  ONLY: nzone_temp, option_traceurs
2130  IMPLICIT NONE
2131  REAL,             INTENT(IN) :: x(ntraciso)
2132  CHARACTER(LEN=*), INTENT(IN) :: err_msg
2133  IF(ALL([17,18]/=option_traceurs)) RETURN
2134  IF(nzone_temp>=5) THEN
2135    IF(iso_verif_tag17_q_deltaD_chns(x,err_msg)==1) STOP
2136  END IF
2137END SUBROUTINE iso_verif_trac17_q_deltaD
2138
2139      subroutine iso_verif_traceur(x,err_msg)
2140        use isotrac_mod, only: ridicule_trac
2141        implicit none
2142        ! vérifier des choses sur les traceurs
2143        ! * toutes les zones donne t l'istope total
2144        ! * pas de deltaD aberrant
2145
2146        ! on prend les valeurs pas défaut pour
2147        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2148       
2149       ! inputs
2150       real x(ntraciso)
2151       character*(*) err_msg ! message d''erreur à afficher
2152
2153       ! locals
2154       !integer iso_verif_traceur_choix_nostop 
2155
2156        if (iso_verif_traceur_choix_nostop(x,err_msg, &
2157     &       errmax,errmaxrel,ridicule_trac,deltalimtrac) &
2158     &       .eq.1) then
2159                stop
2160        endif
2161
2162        end subroutine iso_verif_traceur
2163
2164       
2165      subroutine iso_verif_traceur_retourne3D(x,n1,n2,n3, &
2166     &           i1,i2,i3,err_msg)
2167        use isotrac_mod, only: ridicule_trac
2168
2169        implicit none
2170        ! vérifier des choses sur les traceurs
2171        ! * toutes les zones donne t l'istope total
2172        ! * pas de deltaD aberrant
2173
2174        ! on prend les valeurs pas défaut pour
2175        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2176       
2177       ! inputs
2178       integer n1,n2,n3
2179       real x(n1,n2,n3,ntraciso)
2180       character*(*) err_msg ! message d''erreur à afficher
2181       integer i1,i2,i3
2182
2183       ! locals
2184       !integer iso_verif_traceur_choix_nostop 
2185       real xiso(ntraciso)
2186
2187        call select_dim4_from4D(n1,n2,n3,ntraciso, &
2188     &      x,xiso,i1,i2,i3)
2189        if (iso_verif_traceur_choix_nostop(xiso,err_msg, &
2190     &       errmax,errmaxrel,ridicule_trac,deltalimtrac) &
2191     &       .eq.1) then
2192                stop
2193        endif
2194
2195        end subroutine iso_verif_traceur_retourne3D
2196
2197        subroutine iso_verif_traceur_retourne4D(x,n1,n2,n3,n4, &
2198     &           i1,i2,i3,i4,err_msg)
2199        use isotrac_mod, only: ridicule_trac
2200
2201        implicit none
2202        ! vérifier des choses sur les traceurs
2203        ! * toutes les zones donne t l'istope total
2204        ! * pas de deltaD aberrant
2205
2206        ! on prend les valeurs pas défaut pour
2207        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2208       
2209       ! inputs
2210       integer n1,n2,n3,n4
2211       real x(n1,n2,n3,n4,ntraciso)
2212       character*(*) err_msg ! message d''erreur à afficher
2213       integer i1,i2,i3,i4
2214
2215       ! locals
2216       !integer iso_verif_traceur_choix_nostop 
2217       real xiso(ntraciso)
2218
2219        call select_dim5_from5D(n1,n2,n3,n4,ntraciso, &
2220     &      x,xiso,i1,i2,i3,i4)
2221        if (iso_verif_traceur_choix_nostop(xiso,err_msg, &
2222     &       errmax,errmaxrel,ridicule_trac,deltalimtrac) &
2223     &       .eq.1) then
2224                stop
2225        endif
2226
2227        end subroutine iso_verif_traceur_retourne4D
2228
2229       
2230      subroutine iso_verif_traceur_retourne2D(x,n1,n2, &
2231     &           i1,i2,err_msg)
2232        use isotrac_mod, only: ridicule_trac
2233        implicit none
2234        ! vérifier des choses sur les traceurs
2235        ! * toutes les zones donne t l'istope total
2236        ! * pas de deltaD aberrant
2237
2238        ! on prend les valeurs pas défaut pour
2239        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2240       
2241       ! inputs
2242       integer n1,n2
2243       real x(n1,n2,ntraciso)
2244       character*(*) err_msg ! message d''erreur à afficher
2245       integer i1,i2
2246
2247       ! locals
2248       !integer iso_verif_traceur_choix_nostop 
2249       real xiso(ntraciso)
2250
2251        call select_dim3_from3D(n1,n2,ntraciso, &
2252     &      x,xiso,i1,i2)
2253        if (iso_verif_traceur_choix_nostop(xiso,err_msg, &
2254     &       errmax,errmaxrel,ridicule_trac,deltalimtrac) &
2255     &       .eq.1) then
2256                stop
2257        endif
2258
2259        end subroutine iso_verif_traceur_retourne2D
2260
2261        subroutine iso_verif_traceur_vect(x,n,m,err_msg)
2262        USE isotopes_mod, ONLY: iso_HDO
2263        implicit none
2264        ! vérifier des choses sur les traceurs
2265        ! * toutes les zones donne t l'istope total
2266        ! * pas de deltaD aberrant
2267
2268        ! on prend les valeurs pas défaut pour
2269        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2270       
2271       ! inputs
2272       integer n,m
2273       real x(ntraciso,n,m)
2274       character*(*) err_msg ! message d''erreur à afficher
2275
2276       ! locals
2277       logical iso_verif_traceur_nostop
2278       integer i,j,ixt,iiso,izone,ieau
2279       integer ifaux,jfaux,ixtfaux
2280       
2281       call iso_verif_traceur_noNaN_vect(x,n,m,err_msg)       
2282
2283        ! verif masse: iso_verif_tracm_choix_nostop
2284        call iso_verif_trac_masse_vect(x,n,m,err_msg,errmax,errmaxrel)
2285
2286        ! verif deltaD: iso_verif_tracdD_choix_nostop   
2287        if (iso_HDO.gt.0) then 
2288        call iso_verif_tracdd_vect(x,n,m,err_msg)     
2289        endif !if (iso_HDO.gt.0) then       
2290
2291        ! verif pas aberramment negatif: iso_verif_tracpos_choix_nostop
2292        call iso_verif_tracpos_vect(x,n,m,err_msg,1e-5) 
2293       
2294        end subroutine iso_verif_traceur_vect
2295
2296        subroutine iso_verif_tracnps_vect(x,n,m,err_msg)
2297        USE isotopes_mod, ONLY: iso_HDO
2298        implicit none
2299        ! vérifier des choses sur les traceurs
2300        ! * toutes les zones donne t l'istope total
2301        ! * pas de deltaD aberrant
2302
2303        ! on prend les valeurs pas défaut pour
2304        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2305       
2306       ! inputs
2307       integer n,m
2308       real x(ntraciso,n,m)
2309       character*(*) err_msg ! message d''erreur à afficher
2310
2311       ! locals
2312       logical iso_verif_traceur_nostop
2313       integer i,j,ixt,iiso,izone,ieau
2314       integer ifaux,jfaux,ixtfaux
2315       
2316       call iso_verif_traceur_noNaN_vect(x,n,m,err_msg)       
2317
2318        ! verif masse: iso_verif_tracm_choix_nostop
2319        call iso_verif_trac_masse_vect(x,n,m,err_msg,errmax,errmaxrel)
2320
2321        ! verif deltaD: iso_verif_tracdD_choix_nostop   
2322        if (iso_HDO.gt.0) then 
2323        call iso_verif_tracdd_vect(x,n,m,err_msg)     
2324        endif !if (iso_HDO.gt.0) then       
2325       
2326        end subroutine iso_verif_tracnps_vect
2327
2328
2329        subroutine iso_verif_traceur_noNaN_vect(x,n,m,err_msg)
2330        implicit none
2331       
2332       ! inputs
2333       integer n,m
2334       real x(ntraciso,n,m)
2335       character*(*) err_msg ! message d''erreur à afficher
2336
2337       ! locals
2338       logical iso_verif_traceur_nostop
2339       integer i,j,ixt,iiso
2340       integer ifaux,jfaux,ixtfaux
2341
2342
2343       iso_verif_traceur_nostop=.false.
2344        ! verif noNaN: iso_verif_traceur_noNaN_nostop       
2345        do j=1,m
2346        do i=1,n
2347        do ixt=niso+1,ntraciso
2348          if ((x(ixt,i,j).gt.-borne).and.(x(ixt,i,j).lt.borne)) then
2349          else !if ((x.gt.-borne).and.(x.lt.borne)) then
2350              iso_verif_traceur_nostop=.true.
2351              ifaux=i
2352              jfaux=i
2353          endif !if ((x.gt.-borne).and.(x.lt.borne)) then
2354        enddo !do ixt=niso+1,ntraciso
2355        enddo ! do i=1,n
2356        enddo ! do j=1,m
2357       
2358
2359        if (iso_verif_traceur_nostop) then
2360            write(*,*) 'erreur detectée par iso_verif_nonNAN ', &
2361     &           'dans iso_verif_traceur_vect'
2362            write(*,*) ''
2363            write(*,*) err_msg
2364            write(*,*) 'x=',x(:,ifaux,jfaux)
2365            stop
2366        endif
2367
2368        end subroutine iso_verif_traceur_noNaN_vect
2369
2370       
2371        subroutine iso_verif_trac_masse_vect(x,n,m,err_msg, &
2372     &            errmax,errmaxrel)
2373        use isotopes_mod, only: isoName
2374        implicit none
2375       
2376        ! inputs
2377       integer n,m
2378       real x(ntraciso,n,m)
2379       character*(*) err_msg ! message d''erreur à afficher
2380       real errmax,errmaxrel
2381
2382       ! locals
2383       logical iso_verif_traceur_nostop
2384       integer i,j,ixt,iiso,izone
2385       integer ifaux,jfaux,ixtfaux
2386       real xtractot(n,m)
2387       real xiiso(n,m)
2388
2389        do iiso=1,niso       
2390        do j=1,m
2391         do i=1,n       
2392          xtractot(i,j)=0.0
2393          xiiso(i,j)=x(iiso,i,j)
2394          do izone=1,nzone
2395            ixt=itZonIso(izone,iiso)
2396            xtractot(i,j)=xtractot(i,j)+x(ixt,i,j)           
2397          enddo !do izone=1,nzone
2398         enddo !do i=1,n
2399        enddo !do j=1,m
2400       
2401
2402        call iso_verif_egalite_std_vect( &
2403     &           xtractot,xiiso, &
2404     &           err_msg//', verif trac egalite2, iso ' &
2405     &           //TRIM(isoName(iiso)), &
2406     &           n,m,errmax,errmaxrel)
2407        enddo !do iiso=1,niso
2408
2409        end subroutine iso_verif_trac_masse_vect
2410
2411        subroutine iso_verif_tracdd_vect(x,n,m,err_msg)
2412        use isotopes_mod, only: iso_HDO,iso_eau
2413        use isotrac_mod, only: strtrac
2414        implicit none
2415       
2416        ! inputs
2417       integer n,m
2418       real x(ntraciso,n,m)
2419       character*(*) err_msg ! message d''erreur à afficher
2420
2421       ! locals
2422       integer i,j,iiso,izone,ieau,ixt
2423       real xiiso(niso,n,m)
2424       real xeau(n,m)
2425       integer lnblnk
2426
2427       if (iso_HDO.gt.0) then
2428        do izone=1,nzone
2429          ieau=itZonIso(izone,iso_eau)
2430          do iiso=1,niso
2431           ixt=itZonIso(izone,iiso)
2432           do j=1,m
2433            do i=1,n
2434             xiiso(iiso,i,j)=x(ixt,i,j)
2435            enddo !do i=1,n
2436           enddo !do j=1,m
2437          enddo !do iiso=1,niso
2438           
2439          do j=1,m
2440           do i=1,n
2441            xeau(i,j)=x(ieau,i,j)
2442           enddo !do i=1,n
2443          enddo !do j=1,m
2444         
2445          call iso_verif_aberrant_vect2Dch( &
2446     &           xiiso,xeau,err_msg//strtrac(izone),niso,n,m, &
2447     &           deltalimtrac)
2448         enddo !do izone=1,nzone
2449        endif !if (iso_HDO.gt.0) then
2450
2451        end subroutine iso_verif_tracdd_vect
2452
2453        subroutine iso_verif_tracpos_vect(x,n,m,err_msg,seuil)
2454        implicit none
2455
2456       ! inputs
2457       integer n,m
2458       real x(ntraciso,n,m)
2459       character*(*) err_msg ! message d''erreur à afficher
2460       real seuil
2461
2462       ! locals
2463       integer i,j,ixt
2464       logical iso_verif_traceur_nostop
2465       integer ifaux,jfaux,ixtfaux
2466
2467        iso_verif_traceur_nostop=.false.       
2468        do j=1,m
2469        do i=1,n
2470        do ixt=niso+1,ntraciso
2471          if (x(ixt,i,j).lt.-seuil) then
2472              ifaux=i
2473              jfaux=j
2474              ixtfaux=ixt
2475              iso_verif_traceur_nostop=.true.
2476          endif
2477        enddo !do ixt=niso+1,ntraciso
2478        enddo !do i=1,n
2479        enddo !do j=1,m       
2480
2481        if (iso_verif_traceur_nostop) then
2482            write(*,*) 'erreur detectée par verif positif ', &
2483     &           'dans iso_verif_traceur_vect'
2484            write(*,*) ''
2485            write(*,*) err_msg
2486            write(*,*) 'x=',x(:,ifaux,jfaux)
2487            stop
2488        endif
2489
2490        end subroutine iso_verif_tracpos_vect
2491
2492
2493
2494        subroutine iso_verif_tracnps(x,err_msg)
2495        use isotrac_mod, only: ridicule_trac
2496
2497        implicit none
2498        ! vérifier des choses sur les traceurs
2499        ! * toutes les zones donne t l'istope total
2500        ! * pas de deltaD aberrant
2501
2502        ! on prend les valeurs pas défaut pour
2503        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2504       
2505       ! inputs
2506       real x(ntraciso)
2507       character*(*) err_msg ! message d''erreur à afficher
2508
2509       ! locals
2510       !integer iso_verif_tracnps_choix_nostop 
2511
2512        if (iso_verif_tracnps_choix_nostop(x,err_msg, &
2513     &       errmax,errmaxrel,ridicule_trac,deltalimtrac) &
2514     &       .eq.1) then
2515                stop
2516        endif
2517
2518        end subroutine iso_verif_tracnps
2519
2520        subroutine iso_verif_tracpos_choix(x,err_msg,seuil)
2521        implicit none
2522        ! vérifier des choses sur les traceurs
2523        ! * toutes les zones donne t l'istope total
2524        ! * pas de deltaD aberrant
2525
2526        ! on prend les valeurs pas défaut pour
2527        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2528       
2529       ! inputs
2530       real x(ntraciso)
2531       character*(*) err_msg ! message d''erreur à afficher
2532       real seuil
2533
2534       ! locals
2535       !integer iso_verif_tracpos_choix_nostop 
2536
2537        if (iso_verif_tracpos_choix_nostop(x,err_msg,seuil) &
2538     &       .eq.1) then
2539                stop
2540        endif
2541
2542        end subroutine iso_verif_tracpos_choix
2543
2544        subroutine iso_verif_traceur_choix(x,err_msg, &
2545     &       errmax,errmaxrel,ridicule_trac_loc,deltalimtrac)
2546        implicit none
2547        ! vérifier des choses sur les traceurs
2548        ! * toutes les zones donne t l'istope total
2549        ! * pas de deltaD aberrant
2550       
2551       ! inputs
2552       real x(ntraciso)
2553       character*(*) err_msg ! message d''erreur à afficher
2554       real errmax,errmaxrel,ridicule_trac_loc,deltalimtrac
2555
2556       ! locals
2557       !integer iso_verif_traceur_choix_nostop 
2558
2559        if (iso_verif_traceur_choix_nostop(x,err_msg, &
2560     &       errmax,errmaxrel,ridicule_trac_loc,deltalimtrac) &
2561     &       .eq.1) then
2562                stop
2563        endif
2564
2565        end subroutine iso_verif_traceur_choix
2566
2567        function iso_verif_traceur_nostop(x,err_msg)
2568        use isotrac_mod, only: ridicule_trac
2569        !use isotopes_verif, only: errmax,errmaxrel,deltalimtrac
2570        implicit none
2571        ! vérifier des choses sur les traceurs
2572        ! * toutes les zones donne t l'istope total
2573        ! * pas de deltaD aberrant
2574
2575        ! on prend les valeurs pas défaut pour
2576        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2577       
2578       ! inputs
2579       real x(ntraciso)
2580       character*(*) err_msg ! message d''erreur à afficher
2581
2582       ! output
2583       integer iso_verif_traceur_nostop
2584
2585       ! locals
2586       !integer iso_verif_traceur_choix_nostop 
2587
2588        iso_verif_traceur_nostop= &
2589     &       iso_verif_traceur_choix_nostop(x,err_msg, &
2590     &       errmax,errmaxrel,ridicule_trac,deltalimtrac)
2591
2592        end function iso_verif_traceur_nostop
2593
2594
2595      subroutine iso_verif_traceur_justmass(x,err_msg)
2596        implicit none
2597        ! on vérifie que noNaN et masse
2598
2599        ! on prend les valeurs pas défaut pour
2600        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2601       
2602       ! inputs
2603       real x(ntraciso)
2604       character*(*) err_msg ! message d''erreur à afficher
2605
2606       ! locals
2607       !integer iso_verif_traceur_noNaN_nostop
2608       !integer iso_verif_tracm_choix_nostop
2609
2610        ! verif noNaN
2611        if (iso_verif_traceur_noNaN_nostop(x,err_msg).eq.1) then
2612             stop
2613        endif
2614       
2615        ! verif masse
2616        if (iso_verif_tracm_choix_nostop(x,err_msg, &
2617     &           errmax,errmaxrel).eq.1) then
2618             stop
2619        endif   
2620       
2621        end subroutine iso_verif_traceur_justmass
2622
2623        function iso_verif_traceur_jm_nostop(x,err_msg)
2624        implicit none
2625        ! on vérifie que noNaN et masse
2626
2627        ! on prend les valeurs pas défaut pour
2628        ! errmax,errmaxrel,ridicule_trac,deltalimtrac
2629       
2630       ! inputs
2631       real x(ntraciso)
2632       character*(*) err_msg ! message d''erreur à afficher
2633
2634       ! output
2635       integer iso_verif_traceur_jm_nostop
2636
2637       ! locals
2638!       integer iso_verif_traceur_noNaN_nostop
2639       !integer iso_verif_tracm_choix_nostop
2640
2641        iso_verif_traceur_jm_nostop=0
2642!        ! verif noNaN
2643!        if (iso_verif_traceur_noNaN_nostop(x,err_msg).eq.1) then
2644!             iso_verif_traceur_jm_nostop=1
2645!        endif
2646       
2647        ! verif masse
2648        if (iso_verif_tracm_choix_nostop(x,err_msg, &
2649     &           errmax,errmaxrel).eq.1) then
2650             iso_verif_traceur_jm_nostop=1
2651        endif   
2652       
2653        end function iso_verif_traceur_jm_nostop
2654
2655        subroutine iso_verif_tag17_q_deltaD_vect(x,n,m,err_msg)
2656        USE isotopes_mod, ONLY: tnat,iso_eau, ridicule,iso_HDO
2657        use isotrac_mod, only: option_traceurs,nzone_temp
2658        implicit none
2659
2660        ! inputs
2661        integer n,m
2662        real x(ntraciso,n,m)
2663        character*(*) err_msg
2664
2665        ! locals
2666        !integer iso_verif_positif_nostop
2667        !real deltaD
2668        integer ieau,ixt,ieau1
2669        integer i,k
2670
2671        if ((option_traceurs.eq.17).or. &
2672     &           (option_traceurs.eq.18)) then
2673        ! verifier que deltaD du tag de la couche la plus haute <
2674        ! 200 permil, et vérifier que son q est inférieur à
2675        ieau=itZonIso(nzone_temp,iso_eau)
2676        ixt=itZonIso(nzone_temp,iso_HDO)
2677        ieau1=itZonIso(1,iso_eau)
2678        do i=1,n
2679         do k=1,m
2680           if (x(ieau,i,k).gt.ridicule) then
2681             if ((x(ixt,i,k)/x(ieau,i,k)/tnat(iso_HDO)-1)*1000 &
2682     &            .gt.-200.0) then
2683                write(*,*) err_msg,', vect:deltaDt05 trop fort'
2684                write(*,*) 'i,k=',i,k
2685                write(*,*) 'x(:,i,k)=',x(:,i,k)
2686                stop
2687             endif !if (iso_verif_positif_nostop(-200.0-deltaD(x(ixt),x(ieau)),
2688           endif !if (x(ieau).gt.ridicule) then
2689           if (x(ieau,i,k).gt.2.0e-3) then
2690                write(*,*) err_msg,', vect:qt05 trop fort'
2691                write(*,*) 'i,k=',i,k
2692                write(*,*) 'x(:,i,k)=',x(:,i,k)
2693                stop
2694           endif !if (iso_verif_positif_nostop(1.0e-3-x(ieau),
2695           if (x(iso_eau,i,k).lt.2.0e-3) then
2696                if (x(ieau1,i,k)/x(iso_eau,i,k).gt.0.1) then
2697                   write(*,*) err_msg,', vect: qt01 trop abondant'
2698                   write(*,*) 'i,k=',i,k
2699                   write(*,*) 'ieau1,iso_eau,x(ieau1,i,k),', &
2700     &                 'x(iso_eau,i,k)=',ieau1,iso_eau, &
2701     &                  x(ieau1,i,k),x(iso_eau,i,k) 
2702                   write(*,*) 'x(:,i,k)=',x(:,i,k)
2703                   stop
2704                endif !if (x(ieau1,i,k)/x(iso_eau,i,k).gt.0.1) then
2705            endif
2706          enddo !do k=1,m
2707        enddo !do i=1,n
2708
2709        endif !if (option_traceurs.eq.17) then
2710
2711        end subroutine iso_verif_tag17_q_deltaD_vect
2712
2713
2714        subroutine iso_verif_tag17_q_deltaD_vect_ret3D(x,n,m,nq,err_msg)
2715        USE isotopes_mod, ONLY: tnat,iso_eau,iso_HDO,ridicule
2716        use isotrac_mod, only: option_traceurs,nzone_temp
2717        implicit none
2718
2719        ! inputs
2720        integer n,m,nq
2721        real x(n,m,nq,ntraciso)
2722        character*(*) err_msg
2723
2724        ! locals
2725        !integer iso_verif_positif_nostop
2726        !real deltaD
2727        integer ieau,ixt,ieau1
2728        integer i,k,iq
2729
2730        if ((option_traceurs.eq.17).or. &
2731     &           (option_traceurs.eq.18)) then
2732        ! verifier que deltaD du tag de la couche la plus haute <
2733        ! 200 permil, et vérifier que son q est inférieur à
2734        ieau=itZonIso(nzone_temp,iso_eau)
2735        ixt=itZonIso(nzone_temp,iso_HDO)
2736        ieau1=itZonIso(1,iso_eau)
2737        do iq=1,nq
2738        do i=1,n
2739         do k=1,m
2740           if (x(i,k,iq,ieau).gt.ridicule) then
2741             if ((x(i,k,iq,ixt)/x(i,k,iq,ieau)/tnat(iso_HDO)-1)*1000 &
2742     &            .gt.-200.0) then
2743                write(*,*) err_msg,', vect:deltaDt05 trop fort'
2744                write(*,*) 'i,k=',i,k
2745                write(*,*) 'x(i,k,iq,:)=',x(i,k,iq,:)
2746                stop
2747             endif !if (iso_verif_positif_nostop(-200.0-deltaD(x(ixt),x(ieau)),
2748           endif !if (x(ieau).gt.ridicule) then
2749           if (x(i,k,iq,ieau).gt.2.0e-3) then
2750                write(*,*) err_msg,', vect:qt05 trop fort'
2751                write(*,*) 'i,k=',i,k
2752                write(*,*) 'x(i,k,iq,:)=',x(i,k,iq,:)
2753                stop
2754           endif !if (iso_verif_positif_nostop(1.0e-3-x(ieau),
2755           if (x(i,k,iq,iso_eau).lt.2.0e-3) then
2756                if (x(i,k,iq,ieau1)/x(i,k,iq,iso_eau).gt.0.1) then
2757                   write(*,*) err_msg,', vect: qt01 trop abondant'
2758                   write(*,*) 'i,k=',i,k
2759                   write(*,*) 'ieau1,iso_eau,x(i,k,iq,ieau1),', &
2760     &                 'x(i,k,iq,ieau)=',ieau1,iso_eau, &
2761     &                  x(i,k,iq,ieau1),x(i,k,iq,iso_eau) 
2762                   write(*,*) 'x(i,k,iq,:)=',x(i,k,iq,:)
2763                   stop
2764                endif !if (x(ieau1,i,k)/x(iso_eau,i,k).gt.0.1) then
2765            endif
2766          enddo !do k=1,m
2767        enddo !do i=1,n
2768        enddo ! do iq=1,nq
2769
2770        endif !if (option_traceurs.eq.17) then
2771
2772        end subroutine iso_verif_tag17_q_deltaD_vect_ret3D
2773
2774
2775#endif
2776! endif ISOTRAC
2777
2778END MODULE isotopes_verif_mod
2779
2780#endif         
2781! endif ISOVERIF
2782
Note: See TracBrowser for help on using the repository browser.