source: LMDZ6/trunk/libf/dyn3d/check_isotopes.F @ 4098

Last change on this file since 4098 was 4050, checked in by dcugnet, 3 years ago

Second commit for new tracers.

  • include most of the keys in the tracers descriptor vector "tracers(:)".
  • fix in phylmdiso/cv3_routines: fq_* variables were used where their fxt_* counterparts were expected.
  • multiple IF(nqdesc(iq)>0) and IF(nqfils(iq)>0) tests suppressed, because they are not needed: "do ... enddo" loops with 0 upper bound are not executed.
  • remove French accents from comments (encoding problem) in phylmdiso/cv3_routines and phylmdiso/cv30_routines.
  • modifications in "isotopes_verif_mod", where the call to function "iso_verif_tag17_q_deltad_chn" in "iso_verif_tag17_q_deltad_chn" was not detected at linking stage, although defined in the same module (?).
File size: 8.0 KB
RevLine 
[2270]1        subroutine check_isotopes_seq(q,ip1jmp1,err_msg)
[4050]2        USE infotrac, ONLY: nqtot, nqo, niso, ntraciso, ntraceurs_zone,
3     &                     ok_isotopes, ok_isotrac, use_iso,
4     &                     iqiso, index_trac,indnum_fn_num, tnat
[2270]5        implicit none
6
7#include "dimensions.h"
8
9        ! inputs
10        integer ip1jmp1
11        real q(ip1jmp1,llm,nqtot)
12        character*(*) err_msg ! message d''erreur à afficher
13
14        ! locals
15        integer ixt,phase,k,i,iq,iiso,izone,ieau,iqeau
16        real xtractot,xiiso
17        real borne
18        real qmin
19        real errmax ! erreur maximale en absolu.
20        real errmaxrel ! erreur maximale en relatif autorisée
21        real deltaDmax,deltaDmin
22        real ridicule
23        parameter (borne=1e19)
24        parameter (errmax=1e-8)
25        parameter (errmaxrel=1e-3)
26        parameter (qmin=1e-11)
[4037]27        parameter (deltaDmax=1000.0,deltaDmin=-999.0)
[2270]28        parameter (ridicule=1e-12)
29        real deltaD
30
31        if (ok_isotopes) then
32
[2281]33        write(*,*) 'check_isotopes 31: err_msg=',err_msg
[2270]34        ! verifier que rien n'est NaN
35        do ixt=1,ntraciso
36          do phase=1,nqo
37            iq=iqiso(ixt,phase)
38            do k=1,llm
39              DO i = 1,ip1jmp1
40                if ((q(i,k,iq).gt.-borne).and.
41     :            (q(i,k,iq).lt.borne)) then
42                else !if ((x(ixt,i,j).gt.-borne).and.
43                  write(*,*) 'erreur detectee par iso_verif_noNaN:'
44                  write(*,*) err_msg
45                  write(*,*) 'q,i,k,iq=',q(i,k,iq),i,k,iq
46                  write(*,*) 'borne=',borne
[4037]47                  call abort_gcm('check_isotopes_loc','plantage iso',0)
[2270]48                endif  !if ((x(ixt,i,j).gt.-borne).and.
49              enddo !DO i = 1,ip1jmp1
50            enddo !do k=1,llm
51          enddo !do phase=1,nqo
52        enddo !do ixt=1,ntraciso
53
[2281]54        !write(*,*) 'check_isotopes 52'
[2270]55        ! verifier que l'eau normale est OK
56        if (use_iso(1)) then
57          ixt=indnum_fn_num(1)
58          do phase=1,nqo
59            iq=iqiso(ixt,phase)
60            do k=1,llm
61            DO i = 1,ip1jmp1 
62              if ((abs(q(i,k,phase)-q(i,k,iq)).gt.errmax).and.
63     :          (abs((q(i,k,phase)-q(i,k,iq))/
64     :           max(max(abs(q(i,k,phase)),abs(q(i,k,iq))),1e-18))
65     :           .gt.errmaxrel)) then
66                  write(*,*) 'erreur detectee par iso_verif_egalite:'
67                  write(*,*) err_msg
68                  write(*,*) 'ixt,phase=',ixt,phase
69                  write(*,*) 'q,iq,i,k=',q(i,k,iq),iq,i,k
70                  write(*,*) 'q(i,k,phase)=',q(i,k,phase)
[4037]71                  call abort_gcm('check_isotopes_loc','plantage iso',0)
[2270]72              endif !if ((abs(q(i,k,phase)-q(i,k,iq)).gt.errmax).and.
73              ! bidouille pour éviter divergence:
74              q(i,k,iq)= q(i,k,phase)
75            enddo ! DO i = 1,ip1jmp1
76            enddo !do k=1,llm
77          enddo ! do phase=1,nqo
78        endif !if (use_iso(1)) then
79       
[2281]80        !write(*,*) 'check_isotopes 78'
[2270]81        ! verifier que HDO est raisonable
82        if (use_iso(2)) then
83          ixt=indnum_fn_num(2)
84          do phase=1,nqo
85            iq=iqiso(ixt,phase)
86            do k=1,llm
87            DO i = 1,ip1jmp1
88            if (q(i,k,iq).gt.qmin) then
89             deltaD=(q(i,k,iq)/q(i,k,phase)/tnat(2)-1)*1000
90             if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then
91                  write(*,*) 'erreur detectee par iso_verif_aberrant:'
92                  write(*,*) err_msg
93                  write(*,*) 'ixt,phase=',ixt,phase
94                  write(*,*) 'q,iq,i,k,=',q(i,k,iq),iq,i,k
95                  write(*,*) 'q=',q(i,k,:)
96                  write(*,*) 'deltaD=',deltaD
[4037]97                  call abort_gcm('check_isotopes_loc','plantage iso',0)
[2270]98             endif !if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then
99            endif !if (q(i,k,iq).gt.qmin) then
100            enddo !DO i = 1,ip1jmp1
101            enddo !do k=1,llm
102          enddo ! do phase=1,nqo
103        endif !if (use_iso(2)) then
104
[2281]105        !write(*,*) 'check_isotopes 103'
[2270]106        ! verifier que O18 est raisonable
107        if (use_iso(3)) then
108          ixt=indnum_fn_num(3)
109          do phase=1,nqo
110            iq=iqiso(ixt,phase)
111            do k=1,llm
112            DO i = 1,ip1jmp1
113            if (q(i,k,iq).gt.qmin) then
114             deltaD=(q(i,k,iq)/q(i,k,phase)/tnat(3)-1)*1000
115             if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then
116                  write(*,*) 'erreur detectee iso_verif_aberrant O18:'
117                  write(*,*) err_msg
118                  write(*,*) 'ixt,phase=',ixt,phase
119                  write(*,*) 'q,iq,i,k,=',q(i,k,phase),iq,i,k
120                  write(*,*) 'xt=',q(i,k,:)
121                  write(*,*) 'deltaO18=',deltaD
[4037]122                  call abort_gcm('check_isotopes_loc','plantage iso',0)
[2270]123             endif !if ((deltaD.gt.deltaDmax).or.(deltaD.lt.deltaDmin)) then
124            endif !if (q(i,k,iq).gt.qmin) then
125            enddo !DO i = 1,ip1jmp1
126            enddo !do k=1,llm
127          enddo ! do phase=1,nqo
128        endif !if (use_iso(2)) then
129
130
[2281]131        !write(*,*) 'check_isotopes 129'
[2270]132        if (ok_isotrac) then
133
134          if (use_iso(2).and.use_iso(1)) then
135            do izone=1,ntraceurs_zone
136             ixt=index_trac(izone,indnum_fn_num(2))
137             ieau=index_trac(izone,indnum_fn_num(1))
138             do phase=1,nqo
139               iq=iqiso(ixt,phase)
140               iqeau=iqiso(ieau,phase)
141               do k=1,llm
142                DO i = 1,ip1jmp1
143                if (q(i,k,iq).gt.qmin) then
144                 deltaD=(q(i,k,iq)/q(i,k,iqeau)/tnat(2)-1)*1000
145                 if ((deltaD.gt.deltaDmax).or.
146     &                   (deltaD.lt.deltaDmin)) then
147                  write(*,*) 'erreur dans iso_verif_aberrant trac:'
148                  write(*,*) err_msg
149                  write(*,*) 'izone,phase=',izone,phase
150                  write(*,*) 'ixt,ieau=',ixt,ieau
151                  write(*,*) 'q,iq,i,k,=',q(i,k,iq),iq,i,k
152                  write(*,*) 'deltaD=',deltaD
[4037]153                  call abort_gcm('check_isotopes_loc','plantage iso',0)
[2270]154                 endif !if ((deltaD.gt.deltaDmax).or.
155                endif !if (q(i,k,iq).gt.qmin) then
156                enddo !DO i = 1,ip1jmp1
157                enddo  ! do k=1,llm
158              enddo ! do phase=1,nqo   
159            enddo !do izone=1,ntraceurs_zone
160          endif !if (use_iso(2).and.use_iso(1)) then
161
162          do iiso=1,niso
163           do phase=1,nqo
164              iq=iqiso(iiso,phase)
165              do k=1,llm
166                DO i = 1,ip1jmp1
167                   xtractot=0.0
168                   xiiso=q(i,k,iq)
169                   do izone=1,ntraceurs_zone
170                      iq=iqiso(index_trac(izone,iiso),phase)
171                      xtractot=xtractot+ q(i,k,iq)
172                   enddo !do izone=1,ntraceurs_zone
173                   if ((abs(xtractot-xiiso).gt.errmax).and.
174     :                  (abs(xtractot-xiiso)/
175     :                  max(max(abs(xtractot),abs(xiiso)),1e-18)
176     :                  .gt.errmaxrel)) then
177                  write(*,*) 'erreur detectee par iso_verif_traceurs:'
178                  write(*,*) err_msg
179                  write(*,*) 'iiso,phase=',iiso,phase
180                  write(*,*) 'i,k,=',i,k
181                  write(*,*) 'q(i,k,:)=',q(i,k,:)
[4037]182                  call abort_gcm('check_isotopes_loc','plantage iso',0)
[2270]183                 endif !if ((abs(q(i,k,phase)-q(i,k,iq)).gt.errmax).and.
184                 
185                 ! bidouille pour éviter divergence:
186                 if (abs(xtractot).gt.ridicule) then
187                   do izone=1,ntraceurs_zone
188                     ixt=index_trac(izone,iiso)
189                     q(i,k,iq)=q(i,k,iq)/xtractot*xiiso
190                   enddo !do izone=1,ntraceurs_zone               
191                  endif !if ((abs(xtractot).gt.ridicule) then
192                enddo !DO i = 1,ip1jmp1
193              enddo !do k=1,llm
194           enddo !do phase=1,nqo
195          enddo !do iiso=1,niso
196
197        endif !if (ok_isotrac) then
198
199        endif ! if (ok_isotopes)
[2281]200        !write(*,*) 'check_isotopes 198'
[2270]201       
202        end
203
Note: See TracBrowser for help on using the repository browser.