The TELEMAC-MASCARET system  trunk
vc03bb.f
Go to the documentation of this file.
1 ! *****************
2  SUBROUTINE vc03bb
3 ! *****************
4 !
5  &(xmul,sf,sg,sh,su,sv,f,g,h,u,v,xel,yel,
6  & ikle1,ikle2,ikle3,ikle4,nelem,nelmax,w1,w2,w3,w4 )
7 !
8 !***********************************************************************
9 ! BIEF V6P1 21/08/2010
10 !***********************************************************************
11 !
12 !brief COMPUTES THE FOLLOWING VECTOR IN FINITE ELEMENTS:
13 !code
14 !+ / DF DF
15 !+ V = XMUL / K GRAD(PSII) * ( U -- + V -- ) D(OMEGA)
16 !+ I /OMEGA DX DY
17 !+
18 !+ PSI(I) IS A BASE OF TYPE QUASI-BUBBLE TRIANGLE
19 !+
20 !+ F, U AND V ARE VECTORS
21 !+ K IS A VECTOR WITH COMPONENTS G AND H
22 !
23 !warning THE JACOBIAN MUST BE POSITIVE
24 !warning THE RESULT IS IN W IN NOT ASSEMBLED FORM
25 !
26 !history J-M HERVOUET (LNH) ; C MOULIN (LNH)
27 !+ 13/01/95
28 !+ V5P1
29 !+
30 !
31 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
32 !+ 13/07/2010
33 !+ V6P0
34 !+ Translation of French comments within the FORTRAN sources into
35 !+ English comments
36 !
37 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
38 !+ 21/08/2010
39 !+ V6P0
40 !+ Creation of DOXYGEN tags for automated documentation and
41 !+ cross-referencing of the FORTRAN sources
42 !
43 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
44 !| F |-->| FUNCTION USED IN THE VECTOR FORMULA
45 !| G |-->| FUNCTION USED IN THE VECTOR FORMULA
46 !| H |-->| FUNCTION USED IN THE VECTOR FORMULA
47 !| IKLE1 |-->| FIRST POINT OF TRIANGLES
48 !| IKLE2 |-->| SECOND POINT OF TRIANGLES
49 !| IKLE3 |-->| THIRD POINT OF TRIANGLES
50 !| IKLE4 |-->| QUASI-BUBBLE POINT OF TRIANGLES
51 !| NELEM |-->| NUMBER OF ELEMENTS
52 !| NELMAX |-->| MAXIMUM NUMBER OF ELEMENTS
53 !| SF |-->| BIEF_OBJ STRUCTURE OF F
54 !| SG |-->| BIEF_OBJ STRUCTURE OF G
55 !| SH |-->| BIEF_OBJ STRUCTURE OF H
56 !| SU |-->| BIEF_OBJ STRUCTURE OF U
57 !| SV |-->| BIEF_OBJ STRUCTURE OF V
58 !| U |-->| FUNCTION USED IN THE VECTOR FORMULA
59 !| V |-->| FUNCTION USED IN THE VECTOR FORMULA
60 !| W1 |<--| RESULT IN NON ASSEMBLED FORM
61 !| W2 |<--| RESULT IN NON ASSEMBLED FORM
62 !| W3 |<--| RESULT IN NON ASSEMBLED FORM
63 !| W4 |<--| RESULT IN NON ASSEMBLED FORM
64 !| XEL |-->| ABSCISSAE OF POINTS IN THE MESH, PER ELEMENT
65 !| XMUL |-->| MULTIPLICATION COEFFICIENT
66 !| YEL |-->| ORDINATES OF POINTS IN THE MESH, PER ELEMENT
67 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
68 !
69  USE bief, ex_vc03bb => vc03bb
70 !
71  USE declarations_special
72  IMPLICIT NONE
73 !
74 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
75 !
76  INTEGER, INTENT(IN) :: NELEM,NELMAX
77  INTEGER, INTENT(IN) :: IKLE1(nelmax),IKLE2(nelmax)
78  INTEGER, INTENT(IN) :: IKLE3(nelmax),IKLE4(nelmax)
79 !
80  DOUBLE PRECISION, INTENT(IN) :: XEL(nelmax,*),YEL(nelmax,*)
81  DOUBLE PRECISION, INTENT(INOUT) :: W1(nelmax),W2(nelmax)
82  DOUBLE PRECISION, INTENT(INOUT) :: W3(nelmax),W4(nelmax)
83  DOUBLE PRECISION, INTENT(IN) :: XMUL
84 !
85 ! STRUCTURES OF F, G, H, U, V AND REAL DATA
86 !
87  TYPE(bief_obj), INTENT(IN) :: SF,SG,SH,SU,SV
88  DOUBLE PRECISION, INTENT(IN) :: F(*),G(*),H(*),U(*),V(*)
89 !
90 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
91 !
92  INTEGER IELEM,IELMF,IELMG,IELMU,IELMV,IELMH
93  DOUBLE PRECISION X2,Y2,X3,Y3,F1,F2,F3,F4
94  DOUBLE PRECISION U1,U2,U3,U4,V1,V2,V3,V4,GEL,HEL
95 !
96 !-----------------------------------------------------------------------
97 !
98  ielmf=sf%ELM
99  ielmg=sg%ELM
100  ielmh=sh%ELM
101  ielmu=su%ELM
102  ielmv=sv%ELM
103 !
104 !-----------------------------------------------------------------------
105 !
106 ! F IS QUASI-BUBBLE; G AND H (NOT CHECKED) P0; AND U, V (NOT CHECKED) P1
107 !
108  IF( ielmf.EQ.12.AND.ielmg.EQ.10.AND.ielmh.EQ.10
109  & .AND.ielmu.EQ.11.AND.ielmv.EQ.11 ) THEN
110 !
111  DO ielem = 1 , nelem
112 !
113  x2 = xel(ielem,2)
114  x3 = xel(ielem,3)
115  y2 = yel(ielem,2)
116  y3 = yel(ielem,3)
117 !
118  f1 = f(ikle1(ielem))
119  f2 = f(ikle2(ielem)) - f1
120  f3 = f(ikle3(ielem)) - f1
121  f4 = f(ikle4(ielem)) - f1
122 !
123  u1 = u(ikle1(ielem))
124  u2 = u(ikle2(ielem))
125  u3 = u(ikle3(ielem))
126 !
127  v1 = v(ikle1(ielem))
128  v2 = v(ikle2(ielem))
129  v3 = v(ikle3(ielem))
130 !
131  gel = g(ielem)
132  hel = h(ielem)
133 !
134  w1(ielem) = (-(3*((2*x2*hel-x3*hel+gel*y3-2*gel*y2)*(x2*v3+4
135  & *x2*v2+4*x2*v1-u3*y2-4*u2*y2-4*u1*y2)-(x2*hel-2*x3*hel+
136  & 2*gel*y3-gel*y2)*(4*x3*v3+x3*v2+4*x3*v1-4*u3*y3-u2*y3-4
137  & *u1*y3))*f4-(2*x2*hel-x3*hel+gel*y3-2*gel*y2)*(x2*v3+4*x2*
138  & v2+4*x2*v1+x3*v3+4*x3*v2+4*x3*v1-u3*y3-u3*y2-4*u2*y3-
139  & 4*u2*y2-4*u1*y3-4*u1*y2)*f2+(x2*hel-2*x3*hel+2*gel*y3-gel
140  & *y2)*(4*x2*v3+x2*v2+4*x2*v1+4*x3*v3+x3*v2+4*x3*v1-4*
141  & u3*y3-4*u3*y2-u2*y3-u2*y2-4*u1*y3-4*u1*y2)*f3))
142  & *xmul/(54*(x2*y3-x3*y2))
143 !
144  w2(ielem) = (f2*((x2*hel+x3*hel-gel*y3-gel*y2)*(x2*v3+4*x2*v2+
145  & 4*x2*v1+x3*v3+4*x3*v2+4*x3*v1-u3*y3-u3*y2-4*u2*y3-4*
146  & u2*y2-4*u1*y3-4*u1*y2)+(x2*hel-2*x3*hel+2*gel*y3-gel*y2)*
147  & (4*x2*v3+4*x2*v2+x2*v1-8*x3*v3-8*x3*v2-2*x3*v1+8*u3
148  & *y3-4*u3*y2+8*u2*y3-4*u2*y2+2*u1*y3-u1*y2))+f3*(x2*hel
149  & -2*x3*hel+2*gel*y3-gel*y2)*(8*x2*v3+8*x2*v2+2*x2*v1-4*
150  & x3*v3-4*x3*v2-x3*v1+4*u3*y3-8*u3*y2+4*u2*y3-8*u2*y2+
151  & u1*y3-2*u1*y2)-3*f4*((x2*hel+x3*hel-gel*y3-gel*y2)*(x2*v3+
152  & 4*x2*v2+4*x2*v1-u3*y2-4*u2*y2-4*u1*y2)+(x2*hel-2*x3*hel
153  & +2*gel*y3-gel*y2)*(4*x2*v3+4*x2*v2+x2*v1-4*x3*v3-4*x3*
154  & v2-x3*v1+4*u3*y3-4*u3*y2+4*u2*y3-4*u2*y2+u1*y3-u1*y2)
155  & ))*xmul/(54*(x2*y3-x3*y2))
156 !
157  w3(ielem) = (f2*(2*x2*hel-x3*hel+gel*y3-2*gel*y2)*(4*x2*v3+4
158  & *x2*v2+x2*v1-8*x3*v3-8*x3*v2-2*x3*v1+8*u3*y3-4*u3*y2
159  & +8*u2*y3-4*u2*y2+2*u1*y3-u1*y2)+f3*((2*x2*hel-x3*hel+gel
160  & *y3-2*gel*y2)*(8*x2*v3+8*x2*v2+2*x2*v1-4*x3*v3-4*x3*
161  & v2-x3*v1+4*u3*y3-8*u3*y2+4*u2*y3-8*u2*y2+u1*y3-2*u1*
162  & y2)+(x2*hel+x3*hel-gel*y3-gel*y2)*(4*x2*v3+x2*v2+4*x2*v1+4
163  & *x3*v3+x3*v2+4*x3*v1-4*u3*y3-4*u3*y2-u2*y3-u2*y2-4*u1
164  & *y3-4*u1*y2))-3*f4*((2*x2*hel-x3*hel+gel*y3-2*gel*y2)*(4
165  & *x2*v3+4*x2*v2+x2*v1-4*x3*v3-4*x3*v2-x3*v1+4*u3*y3-4
166  & *u3*y2+4*u2*y3-4*u2*y2+u1*y3-u1*y2)+(x2*hel+x3*hel-gel*y3-
167  & gel*y2)*(4*x3*v3+x3*v2+4*x3*v1-4*u3*y3-u2*y3-4*u1*y3))
168  & )*xmul/(54*(x2*y3-x3*y2))
169 !
170  w4(ielem) = (-(((x2*hel-x3*hel+gel*y3-gel*y2)*(8*x2*v3+8*x2*v2
171  & +2*x2*v1-4*x3*v3-4*x3*v2-x3*v1+4*u3*y3-8*u3*y2+4*u2
172  & *y3-8*u2*y2+u1*y3-2*u1*y2)+(4*x2*v3+x2*v2+4*x2*v1+4*
173  & x3*v3+x3*v2+4*x3*v1-4*u3*y3-4*u3*y2-u2*y3-u2*y2-4*u1*
174  & y3-4*u1*y2)*(x3*hel-gel*y3))*f3-3*((x2*hel-x3*hel+gel*y3-gel*
175  & y2)*(4*x2*v3+4*x2*v2+x2*v1-4*x3*v3-4*x3*v2-x3*v1+4*
176  & u3*y3-4*u3*y2+4*u2*y3-4*u2*y2+u1*y3-u1*y2)+(x2*hel-gel*
177  & y2)*(x2*v3+4*x2*v2+4*x2*v1-u3*y2-4*u2*y2-4*u1*y2)+(x3
178  & *hel-gel*y3)*(4*x3*v3+x3*v2+4*x3*v1-4*u3*y3-u2*y3-4*u1*
179  & y3))*f4+((x2*hel-x3*hel+gel*y3-gel*y2)*(4*x2*v3+4*x2*v2+x2*
180  & v1-8*x3*v3-8*x3*v2-2*x3*v1+8*u3*y3-4*u3*y2+8*u2*y3-
181  & 4*u2*y2+2*u1*y3-u1*y2)+(x2*hel-gel*y2)*(x2*v3+4*x2*v2+4*
182  & x2*v1+x3*v3+4*x3*v2+4*x3*v1-u3*y3-u3*y2-4*u2*y3-4*u2*
183  & y2-4*u1*y3-4*u1*y2))*f2))*xmul/(18*(x2*y3-x3*y2))
184 !
185  ENDDO ! IELEM
186 !
187 !-----------------------------------------------------------------------
188 !
189 ! F IS QUASI-B; G AND H (NOT CHECKED) P0; AND U, V (NOT CHECKED) QUASI-B
190 !
191  ELSEIF(ielmf.EQ.12.AND.ielmg.EQ.10.AND.ielmh.EQ.10.
192  & and.ielmu.EQ.12.AND.ielmv.EQ.12 ) THEN
193 !
194  DO ielem = 1 , nelem
195 !
196  x2 = xel(ielem,2)
197  x3 = xel(ielem,3)
198  y2 = yel(ielem,2)
199  y3 = yel(ielem,3)
200 !
201  f1 = f(ikle1(ielem))
202  f2 = f(ikle2(ielem)) - f1
203  f3 = f(ikle3(ielem)) - f1
204  f4 = f(ikle4(ielem)) - f1
205 !
206  u1 = u(ikle1(ielem))
207  u2 = u(ikle2(ielem))
208  u3 = u(ikle3(ielem))
209  u4 = u(ikle4(ielem))
210 !
211  v1 = v(ikle1(ielem))
212  v2 = v(ikle2(ielem))
213  v3 = v(ikle3(ielem))
214  v4 = v(ikle4(ielem))
215 !
216  gel = g(ielem)
217  hel = h(ielem)
218 !
219  w1(ielem) = (-(3*((2*x2*hel-x3*hel+gel*y3-2*gel*y2)*(x2*v4+x2
220  & *v2+x2*v1-u4*y2-u2*y2-u1*y2)-(x2*hel-2*x3*hel+2*gel*y3-gel*
221  & y2)*(x3*v3+x3*v4+x3*v1-u3*y3-u4*y3-u1*y3))*f4-(2*x2*hel-
222  & x3*hel+gel*y3-2*gel*y2)*(x2*v4+x2*v2+x2*v1+x3*v4+x3*v2+x3*
223  & v1-u4*y3-u4*y2-u2*y3-u2*y2-u1*y3-u1*y2)*f2+(x2*hel-2*x3*
224  & hel+2*gel*y3-gel*y2)*(x2*v3+x2*v4+x2*v1+x3*v3+x3*v4+x3*v1-
225  & u3*y3-u3*y2-u4*y3-u4*y2-u1*y3-u1*y2)*f3))
226  & *xmul/(18*(x2*y3-x3*y2))
227 !
228  w2(ielem) = (f2*((x2*hel+x3*hel-gel*y3-gel*y2)*(x2*v4+x2*v2+x2*
229  & v1+x3*v4+x3*v2+x3*v1-u4*y3-u4*y2-u2*y3-u2*y2-u1*y3-u1*y2)
230  & +(x2*hel-2*x3*hel+2*gel*y3-gel*y2)*(x2*v3+x2*v4+x2*v2-2*x3
231  & *v3-2*x3*v4-2*x3*v2+2*u3*y3-u3*y2+2*u4*y3-u4*y2+2*u2
232  & *y3-u2*y2))+f3*(x2*hel-2*x3*hel+2*gel*y3-gel*y2)*(2*x2*v3+
233  & 2*x2*v4+2*x2*v2-x3*v3-x3*v4-x3*v2+u3*y3-2*u3*y2+u4*y3-
234  & 2*u4*y2+u2*y3-2*u2*y2)-3*f4*((x2*hel+x3*hel-gel*y3-gel*y2)*
235  & (x2*v4+x2*v2+x2*v1-u4*y2-u2*y2-u1*y2)+(x2*hel-2*x3*hel+2*
236  & gel*y3-gel*y2)*(x2*v3+x2*v4+x2*v2-x3*v3-x3*v4-x3*v2+u3*y3-
237  & u3*y2+u4*y3-u4*y2+u2*y3-u2*y2)))*xmul/(18*(x2*y3-x3*y2))
238 !
239  w3(ielem) = (f2*(2*x2*hel-x3*hel+gel*y3-2*gel*y2)*(x2*v3+x2*v4
240  & +x2*v2-2*x3*v3-2*x3*v4-2*x3*v2+2*u3*y3-u3*y2+2*u4*y3
241  & -u4*y2+2*u2*y3-u2*y2)+f3*((2*x2*hel-x3*hel+gel*y3-2*gel*y2
242  & )*(2*x2*v3+2*x2*v4+2*x2*v2-x3*v3-x3*v4-x3*v2+u3*y3-2*
243  & u3*y2+u4*y3-2*u4*y2+u2*y3-2*u2*y2)+(x2*hel+x3*hel-gel*y3-
244  & gel*y2)*(x2*v3+x2*v4+x2*v1+x3*v3+x3*v4+x3*v1-u3*y3-u3*y2-
245  & u4*y3-u4*y2-u1*y3-u1*y2))-3*f4*((2*x2*hel-x3*hel+gel*y3-2
246  & *gel*y2)*(x2*v3+x2*v4+x2*v2-x3*v3-x3*v4-x3*v2+u3*y3-u3*y2+
247  & u4*y3-u4*y2+u2*y3-u2*y2)+(x2*hel+x3*hel-gel*y3-gel*y2)*
248  & (x3*v3+x3*v4+x3*v1-u3*y3-u4*y3-u1*y3)))
249  & *xmul/(18*(x2*y3-x3*y2))
250 !
251  w4(ielem) = (-(((x2*hel-x3*hel+gel*y3-gel*y2)*(2*x2*v3+2*x2*v4
252  & +2*x2*v2-x3*v3-x3*v4-x3*v2+u3*y3-2*u3*y2+u4*y3-2*u4*y2
253  & +u2*y3-2*u2*y2)+(x2*v3+x2*v4+x2*v1+x3*v3+x3*v4+x3*v1-u3*
254  & y3-u3*y2-u4*y3-u4*y2-u1*y3-u1*y2)*(x3*hel-gel*y3))*f3-3*((
255  & x2*hel-x3*hel+gel*y3-gel*y2)*(x2*v3+x2*v4+x2*v2-x3*v3-x3*v4-
256  & x3*v2+u3*y3-u3*y2+u4*y3-u4*y2+u2*y3-u2*y2)+(x2*hel-gel*y2)*
257  & (x2*v4+x2*v2+x2*v1-u4*y2-u2*y2-u1*y2)+(x3*hel-gel*y3)*(x3*
258  & v3+x3*v4+x3*v1-u3*y3-u4*y3-u1*y3))*f4+((x2*hel-x3*hel+gel*y3
259  & -gel*y2)*(x2*v3+x2*v4+x2*v2-2*x3*v3-2*x3*v4-2*x3*v2+2*
260  & u3*y3-u3*y2+2*u4*y3-u4*y2+2*u2*y3-u2*y2)+(x2*hel-gel*y2)*
261  & (x2*v4+x2*v2+x2*v1+x3*v4+x3*v2+x3*v1-u4*y3-u4*y2-u2*y3-u2
262  & *y2-u1*y3-u1*y2))*f2))*xmul/(6*(x2*y3-x3*y2))
263 !
264  ENDDO ! IELEM
265 !
266 !-----------------------------------------------------------------------
267  ELSE
268 !-----------------------------------------------------------------------
269 !
270  WRITE(lu,101) ielmf,sf%NAME
271  WRITE(lu,111) ielmg,sg%NAME
272  WRITE(lu,201) ielmu,su%NAME
273  WRITE(lu,301)
274 101 FORMAT(1x,'VC03BB (BIEF) :',/,
275  & 1x,'DISCRETIZATION OF F:',1i6,
276  & 1x,'REAL NAME: ',a6)
277 111 FORMAT(1x,'DISCRETIZATION OF G:',1i6,
278  & 1x,'REAL NAME: ',a6)
279 201 FORMAT(1x,'DISCRETIZATION OF U:',1i6,
280  & 1x,'REAL NAME: ',a6)
281 301 FORMAT(1x,'CASE NOT IMPLEMENTED')
282  CALL plante(0)
283  stop
284 !
285  ENDIF
286 !
287 !-----------------------------------------------------------------------
288 !
289  RETURN
290  END
declarations_special
Definition: declarations_special.F:3
declarations_special::lu
integer lu
Definition: declarations_special.F:45
vc03bb
subroutine vc03bb(XMUL, SF, SG, SH, SU, SV, F, G, H, U, V, XEL, YEL, IKLE1, IKLE2, IKLE3, IKLE4, NELEM, NELMAX, W1, W2, W3, W4)
Definition: vc03bb.f:8
bief
Definition: bief.f:3
Generated on Tue Jan 5 2021 14:37:50 for The TELEMAC-MASCARET system by doxygen 1.8.13