The TELEMAC-MASCARET system  trunk
flux_ef_vf.f
Go to the documentation of this file.
1 ! *********************
2  SUBROUTINE flux_ef_vf
3 ! *********************
4 !
5  &(flow,phiel,nseg,nelem,nelmax,eltseg,oriseg,ikle,iniflo,iopt,fn,
6  & yaflulim,flulim,yaflulimebe,flulimebe)
7 !
8 !***********************************************************************
9 ! BIEF V7P3
10 !***********************************************************************
11 !
12 !brief Computing fluxes between points within a triangle, with
13 !+ options.
14 !
15 !warning If YAFLULIM and YAFLULIMEBE are both true, only FLULIMEBE
16 !+ is taken into account !!!!!!!!!!!
17 !
18 !
19 !history JMH and Leo Postma
20 !+ 06/05/2009
21 !+ V5P9
22 !+ First version
23 !
24 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
25 !+ 13/07/2010
26 !+ V6P0
27 !+ Translation of French comments within the FORTRAN sources into
28 !+ English comments
29 !
30 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
31 !+ 21/08/2010
32 !+ V6P0
33 !+ Creation of DOXYGEN tags for automated documentation and
34 !+ cross-referencing of the FORTRAN sources
35 !
36 !history J-M HERVOUET & SARA PAVAN (EDF R&D, LNHE)
37 !+ 22/10/2013
38 !+ V7P0
39 !+ Correction of PSI scheme. Reduction of fluxes based on same
40 !+ contribution than classical distributive schemes (i.e. without
41 !+ integration by part).
42 !
43 !history J-M HERVOUET (EDF R&D, LNHE)
44 !+ 05/11/2016
45 !+ V7P3
46 !+ Adding flux limitation with two options, EBE or per segment
47 !
48 !history J-M HERVOUET (EDF R&D, LNHE)
49 !+ 09/09/2017
50 !+ V7P3
51 !+ Argument NELMAX added for dimensioning arrays, instead of NELEM.
52 !
53 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
54 !| ELTSEG |-->| SEGMENTS OF EVERY TRIANGLE
55 !| FLOW |<--| FLUXES PER SEGMENTS
56 !| FN |-->| OPTIONAL ARGUMENT FOR PSI SCHEME
57 !| IKLE |-->| CONNECTIVITY TABLE
58 !| INIFLO |-->| IF(YES) FLOW WILL BE INITIALISED AT 0
59 !| IOPT |-->| OPTION FOR THE CONSTANT PER ELEMENT
60 !| NELEM |-->| NUMBER OF ELEMENTS
61 !| NELMAX |-->| MAXIMUM NUMBER OF ELEMENTS
62 !| NSEG |-->| NUMBER OF SEGMENTS
63 !| ORISEG |-->| ORIENTATION OF SEGMENTS OF EVERY TRIANGLE
64 !| PHIEL |-->| PER ELEMENT, FLUXES LEAVING POINTS
65 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
66 !
67  USE bief, ex_flux_ef_vf => flux_ef_vf
68 !
70  IMPLICIT NONE
71 !
72 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
73 !
74  INTEGER, INTENT(IN) :: NSEG,IOPT,NELEM,NELMAX
75  INTEGER, INTENT(IN) :: ELTSEG(nelmax,3)
76  INTEGER, INTENT(IN) :: ORISEG(nelmax,3)
77  INTEGER, INTENT(IN) :: IKLE(nelmax,3)
78  DOUBLE PRECISION, INTENT(INOUT) :: FLOW(nseg)
79  DOUBLE PRECISION, INTENT(IN) :: PHIEL(nelmax,3)
80  LOGICAL, INTENT(IN) :: INIFLO
81  TYPE(bief_obj), INTENT(IN), OPTIONAL :: FN
82  DOUBLE PRECISION, INTENT(IN), OPTIONAL :: FLULIMEBE(*)
83  DOUBLE PRECISION, INTENT(IN), OPTIONAL :: FLULIM(*)
84  LOGICAL, INTENT(IN), OPTIONAL :: YAFLULIM,YAFLULIMEBE
85 !
86 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
87 !
88  INTEGER IELEM,ISEG,I1,I2,I3
89  DOUBLE PRECISION A1,A2,A3,K1,K2,K3,THIRD,F12,F23,F31,F1,F2,F3
90  DOUBLE PRECISION FN1,FN2,FN3,F21,F32,F13
91  DOUBLE PRECISION BETA1FI,BETA2FI,BETA3FI,FI
92  LOGICAL TFLULIM,TFLULIMEBE
93 !
94  INTRINSIC abs,min,max
95 !
96 !-----------------------------------------------------------------------
97 !
98 ! CHECKING OPTIONAL PARAMETERS
99 !
100  IF(PRESENT(yaflulim)) THEN
101  tflulim=yaflulim
102  IF(tflulim.AND..NOT.PRESENT(flulim)) THEN
103  WRITE(lu,*) 'FLUX_EF_VF:'
104  WRITE(lu,*) 'VARIABLE FLULIM MISSING'
105  CALL plante(1)
106  stop
107  ENDIF
108  ELSE
109  tflulim=.false.
110  ENDIF
111  IF(PRESENT(yaflulimebe)) THEN
112  tflulimebe=yaflulimebe
113  IF(tflulimebe.AND..NOT.PRESENT(flulimebe)) THEN
114  WRITE(lu,*) 'FLUX_EF_VF:'
115  WRITE(lu,*) 'VARIABLE FLULIMEBE MISSING'
116  CALL plante(1)
117  stop
118  ENDIF
119  ELSE
120  tflulimebe=.false.
121  ENDIF
122 !
123  third=1.d0/3.d0
124 !
125 !-----------------------------------------------------------------------
126 !
127 ! INITIALISES FLOW TO 0.D0
128 !
129  IF(iniflo) THEN
130  DO iseg = 1,nseg
131  flow(iseg) = 0.d0
132  ENDDO
133  ENDIF
134 !
135 !-----------------------------------------------------------------------
136 !
137  IF(iopt.EQ.-1) THEN
138 !
139 !-----------------------------------------------------------------------
140 !
141 ! FLUXES ALREADY COMPUTED BEFORE CALLING THIS SUBROUTINE
142 ! THEY ARE JUST ASSEMBLED HERE
143 !
144  IF(tflulimebe) THEN
145 ! WITH EBE FLUX LIMITATION
146  DO ielem = 1,nelem
147 ! SEGMENT 1
148  iseg = eltseg(ielem,1)
149  IF(oriseg(ielem,1).EQ.1) THEN
150  flow(iseg) = flow(iseg) +
151  & phiel(ielem,1)*flulimebe((1-1)*nelmax+ielem)
152  ELSE
153  flow(iseg) = flow(iseg) -
154  & phiel(ielem,1)*flulimebe((1-1)*nelmax+ielem)
155  ENDIF
156 ! SEGMENT 2
157  iseg = eltseg(ielem,2)
158  IF(oriseg(ielem,2).EQ.1) THEN
159  flow(iseg) = flow(iseg) +
160  & phiel(ielem,2)*flulimebe((2-1)*nelmax+ielem)
161  ELSE
162  flow(iseg) = flow(iseg) -
163  & phiel(ielem,2)*flulimebe((2-1)*nelmax+ielem)
164  ENDIF
165 ! SEGMENT 3
166  iseg = eltseg(ielem,3)
167  IF(oriseg(ielem,3).EQ.1) THEN
168  flow(iseg) = flow(iseg) +
169  & phiel(ielem,3)*flulimebe((3-1)*nelmax+ielem)
170  ELSE
171  flow(iseg) = flow(iseg) -
172  & phiel(ielem,3)*flulimebe((3-1)*nelmax+ielem)
173  ENDIF
174  ENDDO
175  ELSE
176 ! WITHOUT FLUX LIMITATION (OR LATER)
177  DO ielem = 1,nelem
178 ! SEGMENT 1
179  iseg = eltseg(ielem,1)
180  IF(oriseg(ielem,1).EQ.1) THEN
181  flow(iseg) = flow(iseg) + phiel(ielem,1)
182  ELSE
183  flow(iseg) = flow(iseg) - phiel(ielem,1)
184  ENDIF
185 ! SEGMENT 2
186  iseg = eltseg(ielem,2)
187  IF(oriseg(ielem,2).EQ.1) THEN
188  flow(iseg) = flow(iseg) + phiel(ielem,2)
189  ELSE
190  flow(iseg) = flow(iseg) - phiel(ielem,2)
191  ENDIF
192 ! SEGMENT 3
193  iseg = eltseg(ielem,3)
194  IF(oriseg(ielem,3).EQ.1) THEN
195  flow(iseg) = flow(iseg) + phiel(ielem,3)
196  ELSE
197  flow(iseg) = flow(iseg) - phiel(ielem,3)
198  ENDIF
199  ENDDO
200  ENDIF
201 !
202 !-----------------------------------------------------------------------
203 !
204  ELSEIF(iopt.EQ.0) THEN
205 !
206 !-----------------------------------------------------------------------
207 !
208 ! WITH NO CONSTANT
209 !
210  IF(tflulimebe) THEN
211 ! WITH EBE FLUX LIMITATION
212  DO ielem = 1,nelem
213  f1 = phiel(ielem,1)
214  f2 = phiel(ielem,2)
215  f3 = phiel(ielem,3)
216 ! SEGMENT 1
217  iseg = eltseg(ielem,1)
218  IF(oriseg(ielem,1).EQ.1) THEN
219  flow(iseg) = flow(iseg) +
220  & third*(f1-f2)*flulimebe((1-1)*nelmax + ielem)
221  ELSE
222  flow(iseg) = flow(iseg) -
223  & third*(f1-f2)*flulimebe((1-1)*nelmax + ielem)
224  ENDIF
225 ! SEGMENT 2
226  iseg = eltseg(ielem,2)
227  IF(oriseg(ielem,2).EQ.1) THEN
228  flow(iseg) = flow(iseg) +
229  & third*(f2-f3)*flulimebe((2-1)*nelmax + ielem)
230  ELSE
231  flow(iseg) = flow(iseg) -
232  & third*(f2-f3)*flulimebe((2-1)*nelmax + ielem)
233  ENDIF
234 ! SEGMENT 3
235  iseg = eltseg(ielem,3)
236  IF(oriseg(ielem,3).EQ.1) THEN
237  flow(iseg) = flow(iseg) +
238  & third*(f3-f1)*flulimebe((3-1)*nelmax + ielem)
239  ELSE
240  flow(iseg) = flow(iseg) -
241  & third*(f3-f1)*flulimebe((3-1)*nelmax + ielem)
242  ENDIF
243  ENDDO
244  ELSE
245 ! WITHOUT FLUX LIMITATION
246  DO ielem = 1,nelem
247  f1 = phiel(ielem,1)
248  f2 = phiel(ielem,2)
249  f3 = phiel(ielem,3)
250 ! SEGMENT 1
251  iseg = eltseg(ielem,1)
252  IF(oriseg(ielem,1).EQ.1) THEN
253  flow(iseg) = flow(iseg) + third*(f1-f2)
254  ELSE
255  flow(iseg) = flow(iseg) - third*(f1-f2)
256  ENDIF
257 ! SEGMENT 2
258  iseg = eltseg(ielem,2)
259  IF(oriseg(ielem,2).EQ.1) THEN
260  flow(iseg) = flow(iseg) + third*(f2-f3)
261  ELSE
262  flow(iseg) = flow(iseg) - third*(f2-f3)
263  ENDIF
264 ! SEGMENT 3
265  iseg = eltseg(ielem,3)
266  IF(oriseg(ielem,3).EQ.1) THEN
267  flow(iseg) = flow(iseg) + third*(f3-f1)
268  ELSE
269  flow(iseg) = flow(iseg) - third*(f3-f1)
270  ENDIF
271  ENDDO
272  ENDIF
273 !
274 !-----------------------------------------------------------------------
275 !
276  ELSEIF(iopt.EQ.2) THEN
277 !
278 !-----------------------------------------------------------------------
279 !
280 ! LEO POSTMA'S METHOD (EQUIVALENT TO FLUXES GIVEN BY N-SCHEME)
281 !
282  IF(tflulimebe) THEN
283  DO ielem = 1,nelem
284  f1 = phiel(ielem,1)
285  f2 = phiel(ielem,2)
286  f3 = phiel(ielem,3)
287  a1 = abs(f1)
288  a2 = abs(f2)
289  a3 = abs(f3)
290  IF(a1.GE.a2.AND.a1.GE.a3) THEN
291 ! ALL FLOW TO AND FROM NODE 1
292  iseg = eltseg(ielem,1)
293  IF(oriseg(ielem,1).EQ.1) THEN
294  flow(iseg) = flow(iseg) -
295  & f2*flulimebe((1-1)*nelmax + ielem)
296  ELSE
297  flow(iseg) = flow(iseg) +
298  & f2*flulimebe((1-1)*nelmax + ielem)
299  ENDIF
300  iseg = eltseg(ielem,3)
301  IF(oriseg(ielem,3).EQ.2) THEN
302  flow(iseg) = flow(iseg) -
303  & f3*flulimebe((3-1)*nelmax + ielem)
304  ELSE
305  flow(iseg) = flow(iseg) +
306  & f3*flulimebe((3-1)*nelmax + ielem)
307  ENDIF
308  ELSEIF(a2.GE.a1.AND.a2.GE.a3) THEN
309 ! ALL FLOW TO AND FROM NODE 2
310  iseg = eltseg(ielem,1)
311  IF(oriseg(ielem,1).EQ.2) THEN
312  flow(iseg) = flow(iseg) -
313  & f1*flulimebe((1-1)*nelmax + ielem)
314  ELSE
315  flow(iseg) = flow(iseg) +
316  & f1*flulimebe((1-1)*nelmax + ielem)
317  ENDIF
318  iseg = eltseg(ielem,2)
319  IF(oriseg(ielem,2).EQ.1) THEN
320  flow(iseg) = flow(iseg) -
321  & f3*flulimebe((2-1)*nelmax + ielem)
322  ELSE
323  flow(iseg) = flow(iseg) +
324  & f3*flulimebe((2-1)*nelmax + ielem)
325  ENDIF
326  ELSE
327 ! ALL FLOW TO AND FROM NODE 3
328  iseg = eltseg(ielem,2)
329  IF(oriseg(ielem,2).EQ.2) THEN
330  flow(iseg) = flow(iseg) -
331  & f2*flulimebe((2-1)*nelmax + ielem)
332  ELSE
333  flow(iseg) = flow(iseg) +
334  & f2*flulimebe((2-1)*nelmax + ielem)
335  ENDIF
336  iseg = eltseg(ielem,3)
337  IF(oriseg(ielem,3).EQ.1) THEN
338  flow(iseg) = flow(iseg) -
339  & f1*flulimebe((3-1)*nelmax + ielem)
340  ELSE
341  flow(iseg) = flow(iseg) +
342  & f1*flulimebe((3-1)*nelmax + ielem)
343  ENDIF
344  ENDIF
345  ENDDO
346  ELSE
347  DO ielem = 1,nelem
348  f1 = phiel(ielem,1)
349  f2 = phiel(ielem,2)
350  f3 = phiel(ielem,3)
351  a1 = abs(f1)
352  a2 = abs(f2)
353  a3 = abs(f3)
354  IF(a1.GE.a2.AND.a1.GE.a3) THEN
355 ! ALL FLOW TO AND FROM NODE 1
356  iseg = eltseg(ielem,1)
357  IF(oriseg(ielem,1).EQ.1) THEN
358  flow(iseg) = flow(iseg) - f2
359  ELSE
360  flow(iseg) = flow(iseg) + f2
361  ENDIF
362  iseg = eltseg(ielem,3)
363  IF(oriseg(ielem,3).EQ.2) THEN
364  flow(iseg) = flow(iseg) - f3
365  ELSE
366  flow(iseg) = flow(iseg) + f3
367  ENDIF
368  ELSEIF(a2.GE.a1.AND.a2.GE.a3) THEN
369 ! ALL FLOW TO AND FROM NODE 2
370  iseg = eltseg(ielem,1)
371  IF(oriseg(ielem,1).EQ.2) THEN
372  flow(iseg) = flow(iseg) - f1
373  ELSE
374  flow(iseg) = flow(iseg) + f1
375  ENDIF
376  iseg = eltseg(ielem,2)
377  IF(oriseg(ielem,2).EQ.1) THEN
378  flow(iseg) = flow(iseg) - f3
379  ELSE
380  flow(iseg) = flow(iseg) + f3
381  ENDIF
382  ELSE
383 ! ALL FLOW TO AND FROM NODE 3
384  iseg = eltseg(ielem,2)
385  IF(oriseg(ielem,2).EQ.2) THEN
386  flow(iseg) = flow(iseg) - f2
387  ELSE
388  flow(iseg) = flow(iseg) + f2
389  ENDIF
390  iseg = eltseg(ielem,3)
391  IF(oriseg(ielem,3).EQ.1) THEN
392  flow(iseg) = flow(iseg) - f1
393  ELSE
394  flow(iseg) = flow(iseg) + f1
395  ENDIF
396  ENDIF
397  ENDDO
398  ENDIF
399 !
400 !-----------------------------------------------------------------------
401 !
402  ELSEIF(iopt.EQ.3.AND.PRESENT(fn)) THEN
403 !
404 !-----------------------------------------------------------------------
405 !
406 ! PSI-SCHEME
407 !
408  IF(tflulimebe) THEN
409  DO ielem = 1,nelem
410 !
411  k1 = -phiel(ielem,1)
412  k2 = -phiel(ielem,2)
413  k3 = -phiel(ielem,3)
414 !
415 ! STARTS WITH N-SCHEME (EQUIVALENT TO LEO POSTMA'S IMPLEMENTATION)
416 ! FIJ HERE LIKE LAMBDA(I,J) IN BOOK, SO FLUXES FROM J TO I.
417 !
418  f12=max(min(k1,-k2),0.d0)
419  f23=max(min(k2,-k3),0.d0)
420  f31=max(min(k3,-k1),0.d0)
421  f21=max(min(k2,-k1),0.d0)
422  f32=max(min(k3,-k2),0.d0)
423  f13=max(min(k1,-k3),0.d0)
424 !
425  i1=ikle(ielem,1)
426  i2=ikle(ielem,2)
427  i3=ikle(ielem,3)
428  fn1=fn%R(i1)
429  fn2=fn%R(i2)
430  fn3=fn%R(i3)
431 !
432  beta1fi=f12*(fn1-fn2)+f13*(fn1-fn3)
433  beta2fi=f21*(fn2-fn1)+f23*(fn2-fn3)
434  beta3fi=f31*(fn3-fn1)+f32*(fn3-fn2)
435 !
436  fi=beta1fi+beta2fi+beta3fi
437 !
438 ! NOW PSI-SCHEME
439 !
440 ! WHAT FOLLOWS IS INSPIRED FROM SUBROUTINE VC08AA
441 ! WHERE FIJ IS LIJ
442 !
443  IF(fi.GT.0.d0) THEN
444  IF(beta1fi.GT.fi) THEN
445  f12=f12*fi/beta1fi
446  f13=f13*fi/beta1fi
447  ELSEIF(beta1fi.LT.0.d0) THEN
448  f12=0.d0
449  f13=0.d0
450  ENDIF
451  IF(beta2fi.GT.fi) THEN
452  f21=f21*fi/beta2fi
453  f23=f23*fi/beta2fi
454  ELSEIF(beta2fi.LT.0.d0) THEN
455  f21=0.d0
456  f23=0.d0
457  ENDIF
458  IF(beta3fi.GT.fi) THEN
459  f31=f31*fi/beta3fi
460  f32=f32*fi/beta3fi
461  ELSEIF(beta3fi.LT.0.d0) THEN
462  f31=0.d0
463  f32=0.d0
464  ENDIF
465  ELSEIF(fi.LT.0.d0) THEN
466  IF(beta1fi.LT.fi) THEN
467  f12=f12*fi/beta1fi
468  f13=f13*fi/beta1fi
469  ELSEIF(beta1fi.GT.0.d0) THEN
470  f12=0.d0
471  f13=0.d0
472  ENDIF
473  IF(beta2fi.LT.fi) THEN
474  f21=f21*fi/beta2fi
475  f23=f23*fi/beta2fi
476  ELSEIF(beta2fi.GT.0.d0) THEN
477  f21=0.d0
478  f23=0.d0
479  ENDIF
480  IF(beta3fi.LT.fi) THEN
481  f31=f31*fi/beta3fi
482  f32=f32*fi/beta3fi
483  ELSEIF(beta3fi.GT.0.d0) THEN
484  f31=0.d0
485  f32=0.d0
486  ENDIF
487  ELSE
488  f12=0.d0
489  f23=0.d0
490  f31=0.d0
491  f21=0.d0
492  f32=0.d0
493  f13=0.d0
494  ENDIF
495 !
496 ! ASSEMBLES FLUXES
497 ! A DIFFERENCE WITH STANDARD DISTRIBUTIVE SCHEMES
498 !
499 ! HERE BEWARE CONVENTION ON FLOW
500 !
501 ! SEGMENT 1
502  iseg = eltseg(ielem,1)
503  IF(oriseg(ielem,1).EQ.1) THEN
504  flow(iseg) = flow(iseg) +
505  & (+ f21 - f12)*flulimebe((1-1)*nelmax + ielem)
506  ELSE
507  flow(iseg) = flow(iseg) +
508  & (- f21 + f12)*flulimebe((1-1)*nelmax + ielem)
509  ENDIF
510 ! SEGMENT 2
511  iseg = eltseg(ielem,2)
512  IF(oriseg(ielem,2).EQ.1) THEN
513  flow(iseg) = flow(iseg) +
514  & (+ f32 - f23)*flulimebe((2-1)*nelmax + ielem)
515  ELSE
516  flow(iseg) = flow(iseg) +
517  & (- f32 + f23)*flulimebe((2-1)*nelmax + ielem)
518  ENDIF
519 ! SEGMENT 3
520  iseg = eltseg(ielem,3)
521  IF(oriseg(ielem,3).EQ.1) THEN
522  flow(iseg) = flow(iseg) +
523  & (+ f13 - f31)*flulimebe((3-1)*nelmax + ielem)
524  ELSE
525  flow(iseg) = flow(iseg) +
526  & (- f13 + f31)*flulimebe((3-1)*nelmax + ielem)
527  ENDIF
528  ENDDO
529  ELSE
530  DO ielem = 1,nelem
531 !
532  k1 = -phiel(ielem,1)
533  k2 = -phiel(ielem,2)
534  k3 = -phiel(ielem,3)
535 !
536 ! STARTS WITH N-SCHEME (EQUIVALENT TO LEO POSTMA'S IMPLEMENTATION)
537 ! FIJ HERE LIKE LAMBDA(I,J) IN BOOK, SO FLUXES FROM J TO I.
538 !
539  f12=max(min(k1,-k2),0.d0)
540  f23=max(min(k2,-k3),0.d0)
541  f31=max(min(k3,-k1),0.d0)
542  f21=max(min(k2,-k1),0.d0)
543  f32=max(min(k3,-k2),0.d0)
544  f13=max(min(k1,-k3),0.d0)
545 !
546  i1=ikle(ielem,1)
547  i2=ikle(ielem,2)
548  i3=ikle(ielem,3)
549  fn1=fn%R(i1)
550  fn2=fn%R(i2)
551  fn3=fn%R(i3)
552 !
553  beta1fi=f12*(fn1-fn2)+f13*(fn1-fn3)
554  beta2fi=f21*(fn2-fn1)+f23*(fn2-fn3)
555  beta3fi=f31*(fn3-fn1)+f32*(fn3-fn2)
556 !
557  fi=beta1fi+beta2fi+beta3fi
558 !
559 ! NOW PSI-SCHEME
560 !
561 ! WHAT FOLLOWS IS INSPIRED FROM SUBROUTINE VC08AA
562 ! WHERE FIJ IS LIJ
563 !
564  IF(fi.GT.0.d0) THEN
565  IF(beta1fi.GT.fi) THEN
566  f12=f12*fi/beta1fi
567  f13=f13*fi/beta1fi
568  ELSEIF(beta1fi.LT.0.d0) THEN
569  f12=0.d0
570  f13=0.d0
571  ENDIF
572  IF(beta2fi.GT.fi) THEN
573  f21=f21*fi/beta2fi
574  f23=f23*fi/beta2fi
575  ELSEIF(beta2fi.LT.0.d0) THEN
576  f21=0.d0
577  f23=0.d0
578  ENDIF
579  IF(beta3fi.GT.fi) THEN
580  f31=f31*fi/beta3fi
581  f32=f32*fi/beta3fi
582  ELSEIF(beta3fi.LT.0.d0) THEN
583  f31=0.d0
584  f32=0.d0
585  ENDIF
586  ELSEIF(fi.LT.0.d0) THEN
587  IF(beta1fi.LT.fi) THEN
588  f12=f12*fi/beta1fi
589  f13=f13*fi/beta1fi
590  ELSEIF(beta1fi.GT.0.d0) THEN
591  f12=0.d0
592  f13=0.d0
593  ENDIF
594  IF(beta2fi.LT.fi) THEN
595  f21=f21*fi/beta2fi
596  f23=f23*fi/beta2fi
597  ELSEIF(beta2fi.GT.0.d0) THEN
598  f21=0.d0
599  f23=0.d0
600  ENDIF
601  IF(beta3fi.LT.fi) THEN
602  f31=f31*fi/beta3fi
603  f32=f32*fi/beta3fi
604  ELSEIF(beta3fi.GT.0.d0) THEN
605  f31=0.d0
606  f32=0.d0
607  ENDIF
608  ELSE
609  f12=0.d0
610  f23=0.d0
611  f31=0.d0
612  f21=0.d0
613  f32=0.d0
614  f13=0.d0
615  ENDIF
616 !
617 ! ASSEMBLES FLUXES
618 ! A DIFFERENCE WITH STANDARD DISTRIBUTIVE SCHEMES
619 !
620 ! HERE BEWARE CONVENTION ON FLOW
621 !
622 ! SEGMENT 1
623  iseg = eltseg(ielem,1)
624  IF(oriseg(ielem,1).EQ.1) THEN
625  flow(iseg) = flow(iseg) + f21 - f12
626  ELSE
627  flow(iseg) = flow(iseg) - f21 + f12
628  ENDIF
629 ! SEGMENT 2
630  iseg = eltseg(ielem,2)
631  IF(oriseg(ielem,2).EQ.1) THEN
632  flow(iseg) = flow(iseg) + f32 - f23
633  ELSE
634  flow(iseg) = flow(iseg) - f32 + f23
635  ENDIF
636 ! SEGMENT 3
637  iseg = eltseg(ielem,3)
638  IF(oriseg(ielem,3).EQ.1) THEN
639  flow(iseg) = flow(iseg) + f13 - f31
640  ELSE
641  flow(iseg) = flow(iseg) - f13 + f31
642  ENDIF
643  ENDDO
644  ENDIF
645 !
646 !-----------------------------------------------------------------------
647 !
648  ELSE
649 !
650  WRITE(lu,*) 'FLUX_EF_VF:'
651  WRITE(lu,*) 'UNKNOWN OPTION: ',iopt
652 ! IF(IOPT.EQ.3.AND..NOT.PRESENT(FN)) THEN
653  IF(iopt.EQ.3) THEN
654  WRITE(lu,*) 'OPTION 3: ADVECTED FUNCTION REQUIRED'
655  ENDIF
656  CALL plante(1)
657  stop
658 !
659  ENDIF
660 !
661 !-----------------------------------------------------------------------
662 !
663 ! OPTIONAL FLUX LIMITATION BY SEGMENT
664 !
665  IF(tflulim.AND..NOT.tflulimebe) THEN
666  DO iseg=1,nseg
667  flow(iseg) = flow(iseg)*flulim(iseg)
668  ENDDO
669  ENDIF
670 !
671 !-----------------------------------------------------------------------
672 !
673  RETURN
674  END
675 
subroutine flux_ef_vf(FLOW, PHIEL, NSEG, NELEM, NELMAX, ELTSEG, ORISEG, IKLE, INIFLO, IOPT, FN, YAFLULIM, FLULIM, YAFLULIMEBE, FLULIMEBE)
Definition: flux_ef_vf.f:8
Definition: bief.f:3