The TELEMAC-MASCARET system  trunk
qnlin2.f
Go to the documentation of this file.
1 ! *****************
2  SUBROUTINE qnlin2
3 ! *****************
4 !
5  &( tstot , tsder , iangnl, coefnl, nf , ndire ,
6  & npoin2, f , xkmoy , taux1 , dfini , xcoef )
7 !
8 !***********************************************************************
9 ! TOMAWAC V6P1 24/06/2011
10 !***********************************************************************
11 !
12 !brief COMPUTES THE CONTRIBUTION OF THE NON-LINEAR INTERACTIONS
13 !+ SOURCE TERM (FREQUENCY QUADRUPLETS) USING THE MDIA METHOD
14 !+ ("MULTIPLE DISCRETE INTERACTION APPROXIMATION")
15 !+ PROPOSED BY TOLMAN (2004)
16 !
17 !+ PROCEDURE SPECIFIC TO THE CASE WHERE THE FREQUENCIES
18 !+ FOLLOW A GEOMETRICAL PROGRESSION AND THE DIRECTIONS
19 !+ ARE EVENLY DISTRIBUTED OVER [0;2.PI].
20 !
21 !note THIS SUBROUTINE USES THE OUTPUT FROM 'PRENL2' TO OPTIMISE
22 !+ THE COMPUTATIONS FOR MDIA.
23 !
24 !reference TOLMAN H.L. (2004):
25 !+ "INVERSE MODELING OF DISCRETE INTERACTION APPROXIMATIONS
26 !+ FOR NONLINEAR INTERACTIONS IN WIND WAVES". OCEAN
27 !+ MODELLING, 6, 405-422
28 !
29 !history E. GAGNAIRE-RENOU
30 !+ 04/2011
31 !+ V6P1
32 !+ CREATED
33 !
34 !history G.MATTAROLO (EDF - LNHE)
35 !+ 24/06/2011
36 !+ V6P1
37 !+ Translation of French names of the variables in argument
38 !
39 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
40 !| COEFNL |-->| COEFFICIENTS USED FOR DIA METHOD
41 !| DEPTH |-->| WATER DEPTH
42 !| DFINI |<->| WORK TABLE
43 !| F |-->| DIRECTIONAL SPECTRUM
44 !| F1 |-->| FIRST DISCRETIZED FREQUENCY
45 !| IANGNL |-->| ANGULAR INDICES TABLE
46 !| NF |-->| NUMBER OF FREQUENCIES
47 !| NDIRE |-->| NUMBER OF DIRECTIONS
48 !| NPOIN2 |-->| NUMBER OF POINTS IN 2D MESH
49 !| PROINF |-->| LOGICAL INDICATING INFINITE DEPTH ASSUMPTION
50 !| RAISF |-->| FREQUENTIAL RATIO
51 !| TAUX1 |<->| WORK TABLE
52 !| TSDER |<->| DERIVED PART OF THE SOURCE TERM CONTRIBUTION
53 !| TSTOT |<->| TOTAL PART OF THE SOURCE TERM CONTRIBUTION
54 !| XCOEF |-->| COEFFICIENT FOR MDIA METHOD
55 !| XKMOY |-->| AVERAGE WAVE NUMBER
56 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
57 ! !
58 ! APPELS : - PROGRAMME(S) APPELANT : SEMIMP !
59 ! ******** - PROGRAMME(S) APPELE(S) : CQUEUE !
60 !
62  & f1 , raisf
63 ! !
64  USE interface_tomawac, ex_qnlin2 => qnlin2
65  IMPLICIT NONE
66 !
67 !.....VARIABLES IN ARGUMENT
68 ! """"""""""""""""""""
69  INTEGER, INTENT(IN) :: NPOIN2, NDIRE , NF
70  INTEGER, INTENT(IN) :: IANGNL(ndire,16)
71  DOUBLE PRECISION, INTENT(IN) :: XCOEF
72  DOUBLE PRECISION, INTENT(IN) :: F(npoin2,ndire,nf), COEFNL(32)
73  DOUBLE PRECISION, INTENT(IN) :: XKMOY(npoin2)
74  DOUBLE PRECISION, INTENT(INOUT) :: TAUX1(npoin2), DFINI(npoin2)
75  DOUBLE PRECISION, INTENT(INOUT) :: TSTOT(npoin2,ndire,nf)
76  DOUBLE PRECISION, INTENT(INOUT) :: TSDER(npoin2,ndire,nf)
77 !
78 !.....LOCAL VARIABLES
79 ! """""""""""""""""
80  INTEGER IP , JP , JF , JFMIN , JFMAX
81  INTEGER JFP1 , JFM2 , JFP3 , JFM4 ,
82  & jf1_0 , jf1_1 , jf2_0 , jf2_1 , jf3_0 , jf3_1 , jf4_0 , jf4_1 ,
83  & jb1_0 , jb1_1 , jb2_0 , jb2_1 , jb3_0 , jb3_1 , jb4_0 , jb4_1 ,
84  & jp1p_0, jp1p_1, jp1m_0, jp1m_1, jp2p_0, jp2p_1, jp2m_0, jp2m_1,
85  & jp3p_0, jp3p_1, jp3m_0, jp3m_1, jp4p_0, jp4p_1, jp4m_0, jp4m_1
86  DOUBLE PRECISION US1PM4, US1MM4, US1PL4, US1ML4
87  DOUBLE PRECISION FREQ , XXFAC , TERM1 , W ,
88  & c01 , c02 , c03 , c04 , c05 , c06 , c07 , c08 ,
89  & c09 , c10 , c11 , c12 , c13 , c14 , c15 , c16 ,
90  & d01 , d02 , d03 , d04 , d05 , d06 , d07 , d08 ,
91  & d09 , d10 , d11 , d12 , d13 , d14 , d15 , d16 ,
92  & c01sq , c02sq , c03sq , c04sq , c05sq , c06sq , c07sq , c08sq ,
93  & c09sq , c10sq , c11sq , c12sq , c13sq , c14sq , c15sq , c16sq ,
94  & cf1_0 , cf1_1 , cf2_0 , cf2_1 , cf3_0 , cf3_1 , cf4_0 , cf4_1 ,
95  & qnl_a , qnl_b , qnl_c , qnl_d
96 
97  DOUBLE PRECISION :: TMP
98  DOUBLE PRECISION, DIMENSION(:), POINTER ::
99  & f1plus, f1moin, f2plus, f2moin, f3plus, f3moin, f4plus, f4moin,
100  & s_1p2m, s_1m2p, s_3p4m, s_3m4p, p_1p2m, p_1m2p, p_3p4m, p_3m4p
101  DOUBLE PRECISION, DIMENSION(:), POINTER ::
102  & q_apb , q_apc , q_bpd , q_cpd
103 
104  INTEGER :: IERR
105 
106 !TODO PUT IN POINT_TOMAWAC AND DECLARATIONS_TOMAWAC
107  LOGICAL, SAVE :: DEJA = .false.
108  DOUBLE PRECISION, ALLOCATABLE, TARGET, DIMENSION(:,:),SAVE ::
109  & trav_mdia
110 
111  IF (.NOT.deja) THEN
112  deja = .true.
113  ALLOCATE(trav_mdia(npoin2,22),stat=ierr)
114  CALL check_allocate(ierr,'TRAV_MDIA')
115  ENDIF
116 
117  f1plus => trav_mdia(:,1)
118  f1moin => trav_mdia(:,2)
119  f2plus => trav_mdia(:,3)
120  f2moin => trav_mdia(:,4)
121  f3plus => trav_mdia(:,5)
122  f3moin => trav_mdia(:,6)
123  f4plus => trav_mdia(:,7)
124  f4moin => trav_mdia(:,8)
125  s_1p2m => trav_mdia(:,9)
126  s_1m2p => trav_mdia(:,10)
127  s_3p4m => trav_mdia(:,11)
128  s_3m4p => trav_mdia(:,12)
129  p_1p2m => trav_mdia(:,13)
130  p_1m2p => trav_mdia(:,14)
131  p_3p4m => trav_mdia(:,15)
132  p_3m4p => trav_mdia(:,16)
133  q_apb => trav_mdia(:,17)
134  q_apc => trav_mdia(:,18)
135  q_bpd => trav_mdia(:,19)
136  q_cpd => trav_mdia(:,20)
137 
138 
139 !
140 !
141 !.....RECOVERS THE COEFFICIENTS COMPUTED IN 'PRENL2'
142 ! """"""""""""""""""""""""""""""""""""""""""""""""""
143  c01 = coefnl( 1)
144  c02 = coefnl( 2)
145  c03 = coefnl( 3)
146  c04 = coefnl( 4)
147  c05 = coefnl( 5)
148  c06 = coefnl( 6)
149  c07 = coefnl( 7)
150  c08 = coefnl( 8)
151  c09 = coefnl(17)
152  c10 = coefnl(18)
153  c11 = coefnl(19)
154  c12 = coefnl(20)
155  c13 = coefnl(21)
156  c14 = coefnl(22)
157  c15 = coefnl(23)
158  c16 = coefnl(24)
159 !
160  jfp1 = int(coefnl( 9)+1.d-7)
161  jfm2 = int(coefnl(10)-1.d-7)
162  jfp3 = int(coefnl(25)+1.d-7)
163  jfm4 = int(coefnl(26)-1.d-7)
164 !
165  us1pm4= coefnl(11)
166  us1mm4= coefnl(12)
167  us1pl4= coefnl(27)
168  us1ml4= coefnl(28)
169 !
170  jfmin = min(nint(coefnl(13)),nint(coefnl(29)))
171  jfmax = max(nint(coefnl(14)),nint(coefnl(30)))
172 !
173  c01sq = c01*c01
174  c02sq = c02*c02
175  c03sq = c03*c03
176  c04sq = c04*c04
177  c05sq = c05*c05
178  c06sq = c06*c06
179  c07sq = c07*c07
180  c08sq = c08*c08
181  c09sq = c09*c09
182  c10sq = c10*c10
183  c11sq = c11*c11
184  c12sq = c12*c12
185  c13sq = c13*c13
186  c14sq = c14*c14
187  c15sq = c15*c15
188  c16sq = c16*c16
189 !
190 !.....CORRECTION FACTOR FOR FINITE WATER DEPTH
191 ! """"""""""""""""""""""""""""""""""""""""""""
192  IF (.NOT.proinf) THEN
193 !
194  DO ip=1,npoin2
195  term1 = max(0.75d0*depth(ip)*xkmoy(ip),0.5d0)
196  dfini(ip) = 1.d0+(5.5d0/term1)*(1.d0-0.833d0*term1)
197  & *exp(-1.25d0*term1)
198  ENDDO
199  ENDIF
200 !
201 !.....FIRST LOOP ON THE FREQUENCIES
202 ! """""""""""""""""""""""""""""""""""
203  DO jf=jfmin,jfmax
204 !
205 !.......COMPUTES THE CONSIDERED FREQUENCY
206 ! """"""""""""""""""""""""""""""""""
207  freq = f1*raisf**(jf-1)
208 !
209 !.......GETS THE INDICES OF THE FREQUENCIES EITHER SIDE OF THE
210 !......."HIGHER" 1 FREQUENCY:
211 ! FREQ(JF1_0) < (1+MU).FREQ(JF) < FREQ(JF1_1)
212 ! """""""""""""""""""""""""""""""""""""""""""""""""""""""""""
213  jf1_0=jf+jfp1
214  jf1_1=jf1_0+1
215 !
216 !.......GETS THE INDICES OF THE FREQUENCIES EITHER SIDE OF THE
217 !......."LOWER" 2 FREQUENCY:
218 ! FREQ(JF2_0) < (1-MU).FREQ(JF) < FREQ(JF2_1)
219 ! """""""""""""""""""""""""""""""""""""""""""""""""""""""""""
220  jf2_0=jf+jfm2-1
221  jf2_1=jf2_0+1
222 !
223 !.......GETS THE INDICES OF THE FREQUENCIES EITHER SIDE OF THE
224 !......."HIGHER" 3 FREQUENCY:
225 ! FREQ(JF3_0) < (1+LAMBDA).FREQ(JF) < FREQ(JF3_1)
226 ! """""""""""""""""""""""""""""""""""""""""""""""""""""""""""
227  jf3_0=jf+jfp3
228  jf3_1=jf3_0+1
229 !
230 !.......GETS THE INDICES OF THE FREQUENCIES EITHER SIDE OF THE
231 !......."LOWER" 4 FREQUENCY:
232 ! FREQ(JF4_0) < (1-LAMBDA).FREQ(JF) < FREQ(JF4_1)
233 ! """""""""""""""""""""""""""""""""""""""""""""""""""""""""""
234  jf4_0=jf+jfm4-1
235  jf4_1=jf4_0+1
236 !
237 !.......LIMITS THE INDICES TO NF AND TAKES INTO ACCOUNT ANALYTICALLY
238 ! THE SPECTRUM TAIL (DECREASE IN -TAILF).
239 ! """"""""""""""""""""""""""""""""""""""""""""""""""""""""""
240  CALL cqueue( jf1_1 , jb1_1 , cf1_1 )
241  CALL cqueue( jf1_0 , jb1_0 , cf1_0 )
242  CALL cqueue( jf2_1 , jb2_1 , cf2_1 )
243  CALL cqueue( jf2_0 , jb2_0 , cf2_0 )
244  CALL cqueue( jf3_1 , jb3_1 , cf3_1 )
245  CALL cqueue( jf3_0 , jb3_0 , cf3_0 )
246  CALL cqueue( jf4_1 , jb4_1 , cf4_1 )
247  CALL cqueue( jf4_0 , jb4_0 , cf4_0 )
248 !
249 !.......INTERPOLATION COEFFICIENTS FOR THE MODIFIED SPECTRUM
250 ! """""""""""""""""""""""""""""""""""""""""""""""""
251  d01=c01*cf1_0*us1pm4
252  d02=c02*cf1_0*us1pm4
253  d03=c03*cf1_1*us1pm4
254  d04=c04*cf1_1*us1pm4
255  d05=c05*cf2_0*us1mm4
256  d06=c06*cf2_0*us1mm4
257  d07=c07*cf2_1*us1mm4
258  d08=c08*cf2_1*us1mm4
259  d09=c09*cf3_0*us1pl4
260  d10=c10*cf3_0*us1pl4
261  d11=c11*cf3_1*us1pl4
262  d12=c12*cf3_1*us1pl4
263  d13=c13*cf4_0*us1ml4
264  d14=c14*cf4_0*us1ml4
265  d15=c15*cf4_1*us1ml4
266  d16=c16*cf4_1*us1ml4
267 !
268 !.......COMPUTES THE MULTIPLICATIVE COEFFICIENT (IN F**11) AND TAKES
269 ! INTO ACCOUNT THE CORRECTION TERM IN FINITE DEPTH
270 ! """""""""""""""""""""""""""""""""""""""""""""""""""""""
271  xxfac= 0.5d0 * xcoef / gravit**4 * freq**11
272  IF (proinf) THEN
273  DO ip=1,npoin2
274  taux1(ip) = xxfac
275  ENDDO
276  ELSE
277  DO ip=1,npoin2
278  taux1(ip) = dfini(ip)*xxfac
279  ENDDO
280  ENDIF
281 !
282 !.......SECOND LOOP ON THE DIRECTIONS.
283 ! """"""""""""""""""""""""""""""""""""
284  DO jp=1,ndire
285  jp1p_0 = iangnl(jp, 1)
286  jp1p_1 = iangnl(jp, 2)
287  jp2m_0 = iangnl(jp, 3)
288  jp2m_1 = iangnl(jp, 4)
289  jp1m_0 = iangnl(jp, 5)
290  jp1m_1 = iangnl(jp, 6)
291  jp2p_0 = iangnl(jp, 7)
292  jp2p_1 = iangnl(jp, 8)
293  jp3p_0 = iangnl(jp, 9)
294  jp3p_1 = iangnl(jp,10)
295  jp4m_0 = iangnl(jp,11)
296  jp4m_1 = iangnl(jp,12)
297  jp3m_0 = iangnl(jp,13)
298  jp3m_1 = iangnl(jp,14)
299  jp4p_0 = iangnl(jp,15)
300  jp4p_1 = iangnl(jp,16)
301 !
302 !
303 !........./---------------------------------------------------------------/
304 !........./ FREQUENCIES F1, F2, F3 AND F4 MAY HAVE ENERGY
305 !........./---------------------------------------------------------------/
306 !
307  DO ip=1,npoin2
308  f1plus(ip) = f(ip,jp1p_0,jb1_0)*d01 + f(ip,jp1p_1,jb1_0)*d02
309  & + f(ip,jp1p_0,jb1_1)*d03 + f(ip,jp1p_1,jb1_1)*d04
310  ENDDO
311  DO ip=1,npoin2
312  f1moin(ip) = f(ip,jp1m_0,jb1_0)*d01 + f(ip,jp1m_1,jb1_0)*d02
313  & + f(ip,jp1m_0,jb1_1)*d03 + f(ip,jp1m_1,jb1_1)*d04
314  ENDDO
315  DO ip=1,npoin2
316  f2plus(ip) = f(ip,jp2p_0,jb2_0)*d05 + f(ip,jp2p_1,jb2_0)*d06
317  & + f(ip,jp2p_0,jb2_1)*d07 + f(ip,jp2p_1,jb2_1)*d08
318  ENDDO
319  DO ip=1,npoin2
320  f2moin(ip) = f(ip,jp2m_0,jb2_0)*d05 + f(ip,jp2m_1,jb2_0)*d06
321  & + f(ip,jp2m_0,jb2_1)*d07 + f(ip,jp2m_1,jb2_1)*d08
322  ENDDO
323  DO ip=1,npoin2
324  f3plus(ip) = f(ip,jp3p_0,jb3_0)*d09 + f(ip,jp3p_1,jb3_0)*d10
325  & + f(ip,jp3p_0,jb3_1)*d11 + f(ip,jp3p_1,jb3_1)*d12
326  ENDDO
327  DO ip=1,npoin2
328  f3moin(ip) = f(ip,jp3m_0,jb3_0)*d09 + f(ip,jp3m_1,jb3_0)*d10
329  & + f(ip,jp3m_0,jb3_1)*d11 + f(ip,jp3m_1,jb3_1)*d12
330  ENDDO
331  DO ip=1,npoin2
332  f4plus(ip) = f(ip,jp4p_0,jb4_0)*d13 + f(ip,jp4p_1,jb4_0)*d14
333  & + f(ip,jp4p_0,jb4_1)*d15 + f(ip,jp4p_1,jb4_1)*d16
334  ENDDO
335  DO ip=1,npoin2
336  f4moin(ip) = f(ip,jp4m_0,jb4_0)*d13 + f(ip,jp4m_1,jb4_0)*d14
337  & + f(ip,jp4m_0,jb4_1)*d15 + f(ip,jp4m_1,jb4_1)*d16
338  ENDDO
339 
340  DO ip=1,npoin2
341  s_1p2m(ip) = f1plus(ip) + f2moin(ip)
342  s_1m2p(ip) = f1moin(ip) + f2plus(ip)
343  s_3p4m(ip) = f3plus(ip) + f4moin(ip)
344  s_3m4p(ip) = f3moin(ip) + f4plus(ip)
345  p_1p2m(ip) = f1plus(ip) * f2moin(ip)
346  p_1m2p(ip) = f1moin(ip) * f2plus(ip)
347  p_3p4m(ip) = f3plus(ip) * f4moin(ip)
348  p_3m4p(ip) = f3moin(ip) * f4plus(ip)
349  ENDDO
350 
351  DO ip=1,npoin2
352  w=taux1(ip)
353  qnl_a=w*( p_1p2m(ip)*s_3p4m(ip) - p_3p4m(ip)*s_1p2m(ip) )
354  qnl_b=w*( p_1p2m(ip)*s_3m4p(ip) - p_3m4p(ip)*s_1p2m(ip) )
355  qnl_c=w*( p_1m2p(ip)*s_3p4m(ip) - p_3p4m(ip)*s_1m2p(ip) )
356  qnl_d=w*( p_1m2p(ip)*s_3m4p(ip) - p_3m4p(ip)*s_1m2p(ip) )
357 !
358  q_apb(ip) = qnl_a + qnl_b
359  q_apc(ip) = qnl_a + qnl_c
360  q_bpd(ip) = qnl_b + qnl_d
361  q_cpd(ip) = qnl_c + qnl_d
362  ENDDO
363 
364  !overwrite the pervious terms, in order to diminish the number of summations
365  DO ip=1,npoin2
366  s_1p2m(ip) = s_1p2m(ip) + s_1m2p(ip)
367  s_3p4m(ip) = s_3p4m(ip) + s_3m4p(ip)
368  p_1p2m(ip) = p_1p2m(ip) + p_1m2p(ip)
369  p_3p4m(ip) = p_3p4m(ip) + p_3m4p(ip)
370  ENDDO
371 
372 !
373 
374  IF (jb4_0.EQ.jf4_0) THEN
375  DO ip=1,npoin2
376 !...........For the frequency index before F4(-) : configurations A and C
377  tstot(ip,jp4m_0,jf4_0)=tstot(ip,jp4m_0,jf4_0)
378  & + q_apc(ip)*c13
379  tstot(ip,jp4m_1,jf4_0)=tstot(ip,jp4m_1,jf4_0)
380  & + q_apc(ip)*c14
381  w=taux1(ip)
382  tmp = (-f3plus(ip)*s_1p2m(ip) + p_1p2m(ip))*us1ml4*w
383  tsder(ip,jp4m_0,jf4_0)=tsder(ip,jp4m_0,jf4_0)
384  & + tmp* c13sq
385  tsder(ip,jp4m_1,jf4_0)=tsder(ip,jp4m_1,jf4_0)
386  & + tmp* c14sq
387  ENDDO
388  DO ip=1,npoin2
389 !...........For the frequency index before F4(+) : configurations B and D
390  tstot(ip,jp4p_0,jf4_0)=tstot(ip,jp4p_0,jf4_0)
391  & + q_bpd(ip)*c13
392  tstot(ip,jp4p_1,jf4_0)=tstot(ip,jp4p_1,jf4_0)
393  & + q_bpd(ip)*c14
394  w=taux1(ip)
395  tmp = (-f3moin(ip)*s_1p2m(ip)+ p_1p2m(ip))*us1ml4*w
396 
397  tsder(ip,jp4p_0,jf4_0)=tsder(ip,jp4p_0,jf4_0)
398  & + tmp* c13sq
399  tsder(ip,jp4p_1,jf4_0)=tsder(ip,jp4p_1,jf4_0)
400  & + tmp* c14sq
401  ENDDO
402  ENDIF
403 !
404  IF (jb4_1.EQ.jf4_1) THEN
405 !...........For the frequency index after F4(-) : configurations A and C
406  DO ip=1,npoin2
407  tstot(ip,jp4m_0,jf4_1)=tstot(ip,jp4m_0,jf4_1)
408  & + q_apc(ip)*c15
409  tstot(ip,jp4m_1,jf4_1)=tstot(ip,jp4m_1,jf4_1)
410  & + q_apc(ip)*c16
411  w=taux1(ip)
412  tmp = (-f3plus(ip)*s_1p2m(ip) + p_1p2m(ip))*us1ml4*w
413  tsder(ip,jp4m_0,jf4_1)=tsder(ip,jp4m_0,jf4_1)
414  & + tmp*c15sq
415  tsder(ip,jp4m_1,jf4_1)=tsder(ip,jp4m_1,jf4_1)
416  & + tmp*c16sq
417  ENDDO
418 !...........For the frequency index after F4(+) : configurations B and D
419  DO ip=1,npoin2
420  tstot(ip,jp4p_0,jf4_1)=tstot(ip,jp4p_0,jf4_1)
421  & + q_bpd(ip)*c15
422  tstot(ip,jp4p_1,jf4_1)=tstot(ip,jp4p_1,jf4_1)
423  & + q_bpd(ip)*c16
424  w=taux1(ip)
425  tmp = (-f3moin(ip)*s_1p2m(ip) + p_1p2m(ip))*us1ml4*w
426  tsder(ip,jp4p_0,jf4_1)=tsder(ip,jp4p_0,jf4_1)
427  & + tmp*c15sq
428  tsder(ip,jp4p_1,jf4_1)=tsder(ip,jp4p_1,jf4_1)
429  & + tmp*c16sq
430  ENDDO
431  ENDIF
432 !
433  IF (jb2_0.EQ.jf2_0) THEN
434 !...........For the frequency index before F2(-) : configurations A and B
435  DO ip=1,npoin2
436  tstot(ip,jp2m_0,jf2_0)=tstot(ip,jp2m_0,jf2_0)
437  & - q_apb(ip)*c05
438  tstot(ip,jp2m_1,jf2_0)=tstot(ip,jp2m_1,jf2_0)
439  & - q_apb(ip)*c06
440  w=taux1(ip)
441  tmp = ( f1plus(ip)*s_3p4m(ip) - p_3p4m(ip))*us1mm4*w
442  tsder(ip,jp2m_0,jf2_0)=tsder(ip,jp2m_0,jf2_0)
443  & - tmp*c05sq
444  tsder(ip,jp2m_1,jf2_0)=tsder(ip,jp2m_1,jf2_0)
445  & - tmp*c06sq
446  ENDDO
447 !...........For the frequency index before F2(+) : configurations C and D
448  DO ip=1,npoin2
449  tstot(ip,jp2p_0,jf2_0)=tstot(ip,jp2p_0,jf2_0)
450  & - q_cpd(ip)*c05
451  tstot(ip,jp2p_1,jf2_0)=tstot(ip,jp2p_1,jf2_0)
452  & - q_cpd(ip)*c06
453  w=taux1(ip)
454  tmp = ( f1moin(ip)*s_3p4m(ip) - p_3p4m(ip))*us1mm4*w
455  tsder(ip,jp2p_0,jf2_0)=tsder(ip,jp2p_0,jf2_0)
456  & - tmp*c05sq
457  tsder(ip,jp2p_1,jf2_0)=tsder(ip,jp2p_1,jf2_0)
458  & - tmp*c06sq
459  ENDDO
460  ENDIF
461 !
462  IF (jb2_1.EQ.jf2_1) THEN
463 !...........For the frequency index after F2(-) : configurations A and B
464  DO ip=1,npoin2
465  tstot(ip,jp2m_0,jf2_1)=tstot(ip,jp2m_0,jf2_1)
466  & - q_apb(ip)*c07
467  tstot(ip,jp2m_1,jf2_1)=tstot(ip,jp2m_1,jf2_1)
468  & - q_apb(ip)*c08
469  w=taux1(ip)
470  tmp = ( f1plus(ip)*s_3p4m(ip) - p_3p4m(ip))*us1mm4*w
471  tsder(ip,jp2m_0,jf2_1)=tsder(ip,jp2m_0,jf2_1)
472  & - tmp*c07sq
473  tsder(ip,jp2m_1,jf2_1)=tsder(ip,jp2m_1,jf2_1)
474  & - tmp*c08sq
475  ENDDO
476 !...........For the frequency index after F2(+) : configurations C and D
477  DO ip=1,npoin2
478  tstot(ip,jp2p_0,jf2_1)=tstot(ip,jp2p_0,jf2_1)
479  & - q_cpd(ip)*c07
480  tstot(ip,jp2p_1,jf2_1)=tstot(ip,jp2p_1,jf2_1)
481  & - q_cpd(ip)*c08
482  w=taux1(ip)
483  tmp = ( f1moin(ip)*s_3p4m(ip) - p_3p4m(ip))*us1mm4*w
484  tsder(ip,jp2p_0,jf2_1)=tsder(ip,jp2p_0,jf2_1)
485  & - tmp*c07sq
486  tsder(ip,jp2p_1,jf2_1)=tsder(ip,jp2p_1,jf2_1)
487  & - tmp*c08sq
488  ENDDO
489  ENDIF
490 !
491  IF (jb1_0.EQ.jf1_0) THEN
492 !...........For the frequency index before F1(+) : configurations A and B
493  DO ip=1,npoin2
494  tstot(ip,jp1p_0,jf1_0)=tstot(ip,jp1p_0,jf1_0)
495  & - q_apb(ip)*c01
496  tstot(ip,jp1p_1,jf1_0)=tstot(ip,jp1p_1,jf1_0)
497  & - q_apb(ip)*c02
498  w=taux1(ip)
499  tmp = ( f2moin(ip)*s_3p4m(ip) - p_3p4m(ip))*us1pm4*w
500  tsder(ip,jp1p_0,jf1_0)=tsder(ip,jp1p_0,jf1_0)
501  & -tmp*c01sq
502  tsder(ip,jp1p_1,jf1_0)=tsder(ip,jp1p_1,jf1_0)
503  & -tmp*c02sq
504  ENDDO
505 !...........For the frequency index before F1(-) : configurations C and D
506  DO ip=1,npoin2
507  tstot(ip,jp1m_0,jf1_0)=tstot(ip,jp1m_0,jf1_0)
508  & - q_cpd(ip)*c01
509  tstot(ip,jp1m_1,jf1_0)=tstot(ip,jp1m_1,jf1_0)
510  & - q_cpd(ip)*c02
511  w=taux1(ip)
512  tmp = ( f2plus(ip)*s_3p4m(ip) - p_3p4m(ip))*us1pm4*w
513  tsder(ip,jp1m_0,jf1_0)=tsder(ip,jp1m_0,jf1_0)
514  & -tmp*c01sq
515  tsder(ip,jp1m_1,jf1_0)=tsder(ip,jp1m_1,jf1_0)
516  & -tmp*c02sq
517  ENDDO
518  ENDIF
519 !
520  IF (jb1_1.EQ.jf1_1) THEN
521 !...........For the frequency index after F1(+) : configurations A and B
522  DO ip=1,npoin2
523  tstot(ip,jp1p_0,jf1_1)=tstot(ip,jp1p_0,jf1_1)
524  & - q_apb(ip)*c03
525  tstot(ip,jp1p_1,jf1_1)=tstot(ip,jp1p_1,jf1_1)
526  & - q_apb(ip)*c04
527  w=taux1(ip)
528  tmp = ( f2moin(ip)*s_3p4m(ip) - p_3p4m(ip))*us1pm4*w
529  tsder(ip,jp1p_0,jf1_1)=tsder(ip,jp1p_0,jf1_1)
530  & -tmp*c03sq
531  tsder(ip,jp1p_1,jf1_1)=tsder(ip,jp1p_1,jf1_1)
532  & -tmp*c04sq
533  ENDDO
534 !...........For the frequency index after F1(-) : configurations C and D
535  DO ip=1,npoin2
536  tstot(ip,jp1m_0,jf1_1)=tstot(ip,jp1m_0,jf1_1)
537  & - q_cpd(ip)*c03
538  tstot(ip,jp1m_1,jf1_1)=tstot(ip,jp1m_1,jf1_1)
539  & - q_cpd(ip)*c04
540  w=taux1(ip)
541  tmp = ( f2plus(ip)*s_3p4m(ip) - p_3p4m(ip))*us1pm4*w
542  tsder(ip,jp1m_0,jf1_1)=tsder(ip,jp1m_0,jf1_1)
543  & -tmp*c03sq
544  tsder(ip,jp1m_1,jf1_1)=tsder(ip,jp1m_1,jf1_1)
545  & -tmp*c04sq
546  ENDDO
547  ENDIF
548 !
549  IF (jb3_0.EQ.jf3_0) THEN
550 !...........For the frequency index before F3(+) : configurations A and C
551  DO ip=1,npoin2
552  tstot(ip,jp3p_0,jf3_0)=tstot(ip,jp3p_0,jf3_0)
553  & + q_apc(ip)*c09
554  tstot(ip,jp3p_1,jf3_0)=tstot(ip,jp3p_1,jf3_0)
555  & + q_apc(ip)*c10
556  w=taux1(ip)
557  tmp = (-f4moin(ip)*s_1p2m(ip) + p_1p2m(ip))*us1pl4*w
558  tsder(ip,jp3p_0,jf3_0)=tsder(ip,jp3p_0,jf3_0)
559  & + tmp * c09sq
560  tsder(ip,jp3p_1,jf3_0)=tsder(ip,jp3p_1,jf3_0)
561  & + tmp * c10sq
562  ENDDO
563 !...........For the frequency index before F3(-) : configurations B and D
564  DO ip=1,npoin2
565  tstot(ip,jp3m_0,jf3_0)=tstot(ip,jp3m_0,jf3_0)
566  & + q_bpd(ip)*c09
567  tstot(ip,jp3m_1,jf3_0)=tstot(ip,jp3m_1,jf3_0)
568  & + q_bpd(ip)*c10
569  w=taux1(ip)
570  tmp = (-f4plus(ip)*s_1p2m(ip) + p_1p2m(ip) )*us1pl4*w
571  tsder(ip,jp3m_0,jf3_0)=tsder(ip,jp3m_0,jf3_0)
572  & + tmp * c09sq
573  tsder(ip,jp3m_1,jf3_0)=tsder(ip,jp3m_1,jf3_0)
574  & + tmp * c10sq
575  ENDDO
576  ENDIF
577 !
578  IF (jb3_1.EQ.jf3_1) THEN
579 !...........For the frequency index after F3(+) : configurations A and C
580  DO ip=1,npoin2
581  tstot(ip,jp3p_0,jf3_1)=tstot(ip,jp3p_0,jf3_1)
582  & + q_apc(ip)*c11
583  tstot(ip,jp3p_1,jf3_1)=tstot(ip,jp3p_1,jf3_1)
584  & + q_apc(ip)*c12
585  w=taux1(ip)
586  tmp = (-f4moin(ip)*s_1p2m(ip) + p_1p2m(ip))*us1pl4*w
587  tsder(ip,jp3p_0,jf3_1)=tsder(ip,jp3p_0,jf3_1)
588  & + tmp*c11sq
589  tsder(ip,jp3p_1,jf3_1)=tsder(ip,jp3p_1,jf3_1)
590  & + tmp*c12sq
591  ENDDO
592 !...........For the frequency index after F3(-) : configurations B and D
593  DO ip=1,npoin2
594  tstot(ip,jp3m_0,jf3_1)=tstot(ip,jp3m_0,jf3_1)
595  & + q_bpd(ip)*c11
596  tstot(ip,jp3m_1,jf3_1)=tstot(ip,jp3m_1,jf3_1)
597  & + q_bpd(ip)*c12
598  w=taux1(ip)
599  tmp = (-f4plus(ip)*s_1p2m(ip) + p_1p2m(ip))*us1pl4*w
600  tsder(ip,jp3m_0,jf3_1)=tsder(ip,jp3m_0,jf3_1)
601  & + tmp*c11sq
602  tsder(ip,jp3m_1,jf3_1)=tsder(ip,jp3m_1,jf3_1)
603  & + tmp*c12sq
604  ENDDO
605  ENDIF
606 !
607  ENDDO
608 !
609  ENDDO
610 !
611  RETURN
612  END
613 !=======================================================================
subroutine qnlin2(TSTOT, TSDER, IANGNL, COEFNL, NF, NDIRE, NPOIN2, F, XKMOY, TAUX1, DFINI, XCOEF)
Definition: qnlin2.f:8
double precision, dimension(:), pointer depth
subroutine cqueue(JFRE, JBIS, COEF1)
Definition: cqueue.f:7