The TELEMAC-MASCARET system  trunk
mt08ac.f
Go to the documentation of this file.
1 ! *****************
2  SUBROUTINE mt08ac
3 ! *****************
4 !
5  &( a11 , a12 , a13 , a14 , a15, a16,
6  & a21 , a22 , a23 , a24 , a25, a26,
7  & a31 , a32 , a33 , a34 , a35, a36,
8  & xmul,sf,f,xel,yel,ikle1,ikle2,ikle3,
9  & ikle4,ikle5,ikle6,
10  & nelem,nelmax,icoord)
11 !
12 !***********************************************************************
13 ! BIEF V6P1 21/08/2010
14 !***********************************************************************
15 !
16 !brief COMPUTES THE COEFFICIENTS OF THE FOLLOWING MATRIX:
17 !code
18 !+ EXAMPLE WITH ICOORD = 1
19 !+
20 !+ / D
21 !+ A(I,J)= -XMUL / PSI2(J) * F * --( PSI1(I) ) D(OMEGA)
22 !+ /OMEGA DX
23 !+
24 !+ BEWARE THE MINUS SIGN !!
25 !+
26 !+ PSI1: BASES OF TYPE P1 TRIANGLE
27 !+ PSI2: BASES OF TYPE P2 TRIANGLE
28 !+
29 !+ IT WOULD BE A DERIVATIVE WRT Y WITH ICOORD=2
30 !
31 !warning THE JACOBIAN MUST BE POSITIVE
32 !
33 !history ALGIANE FROEHLY (MATMECA)
34 !+ 19/06/08
35 !+ V5P9
36 !+
37 !
38 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
39 !+ 13/07/2010
40 !+ V6P0
41 !+ Translation of French comments within the FORTRAN sources into
42 !+ English comments
43 !
44 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
45 !+ 21/08/2010
46 !+ V6P0
47 !+ Creation of DOXYGEN tags for automated documentation and
48 !+ cross-referencing of the FORTRAN sources
49 !
50 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
51 !| A11 |<--| ELEMENTS OF MATRIX
52 !| ... |<--| ELEMENTS OF MATRIX
53 !| A36 |<--| ELEMENTS OF MATRIX
54 !| F |-->| FUNCTION USED IN THE FORMULA
55 !| ICOORD |-->| 1: DERIVATIVE ALONG X, 2: ALONG Y
56 !| IKLE1 |-->| FIRST POINTS OF TRIANGLES
57 !| IKLE2 |-->| SECOND POINTS OF TRIANGLES
58 !| IKLE3 |-->| THIRD POINTS OF TRIANGLES
59 !| IKLE4 |-->| FOURTH POINTS OF TRIANGLES (QUADRATIC)
60 !| IKLE5 |-->| FIFTH POINTS OF TRIANGLES (QUADRATIC)
61 !| IKLE6 |-->| SIXTH POINTS OF TRIANGLES (QUADRATIC)
62 !| NELEM |-->| NUMBER OF ELEMENTS
63 !| NELMAX |-->| MAXIMUM NUMBER OF ELEMENTS
64 !| SF |-->| STRUCTURE OF FUNCTIONS F
65 !| XEL |-->| ABSCISSAE OF POINTS IN THE MESH, PER ELEMENT
66 !| YEL |-->| ORDINATES OF POINTS IN THE MESH, PER ELEMENT
67 !| XMUL |-->| MULTIPLICATION FACTOR
68 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
69 !
70  USE bief!, EX_MT08AC => MT08AC
71 !
72  USE declarations_special
73  IMPLICIT NONE
74 !
75 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
76 !
77  INTEGER, INTENT(IN) :: NELEM,NELMAX,ICOORD
78  INTEGER, INTENT(IN) :: IKLE1(nelmax),IKLE2(nelmax)
79  INTEGER, INTENT(IN) :: IKLE3(nelmax),IKLE4(nelmax)
80  INTEGER, INTENT(IN) :: IKLE5(nelmax),IKLE6(nelmax)
81 !
82  DOUBLE PRECISION, INTENT(INOUT) :: A11(*),A12(*),A13(*)
83  DOUBLE PRECISION, INTENT(INOUT) :: A14(*),A15(*),A16(*)
84  DOUBLE PRECISION, INTENT(INOUT) :: A21(*),A22(*),A23(*)
85  DOUBLE PRECISION, INTENT(INOUT) :: A24(*),A25(*),A26(*)
86  DOUBLE PRECISION, INTENT(INOUT) :: A31(*),A32(*),A33(*)
87  DOUBLE PRECISION, INTENT(INOUT) :: A34(*),A35(*),A36(*)
88 !
89  DOUBLE PRECISION, INTENT(IN) :: XMUL
90  DOUBLE PRECISION, INTENT(IN) :: F(*)
91 !
92 ! STRUCTURE OF F
93  TYPE(bief_obj), INTENT(IN) :: SF
94 !
95  DOUBLE PRECISION, INTENT(IN) :: XEL(nelmax,3),YEL(nelmax,3)
96 !
97 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
98 !
99  INTEGER IELEM,IELMF
100  DOUBLE PRECISION X2,X3,Y2,Y3,F1,F2,F3,F4,F5,F6
101 !
102 !-----------------------------------------------------------------------
103 !
104  ielmf=sf%ELM
105 !
106 !-----------------------------------------------------------------------
107 ! CASE WHERE F IS OF P1 DISCRETISATION
108 !-----------------------------------------------------------------------
109 !
110  IF(ielmf.EQ.11) THEN
111 !
112 !================================
113 ! CASE OF DERIVATIVE WRT X =
114 !================================
115 !
116  IF(icoord.EQ.1) THEN
117 !
118 ! LOOP ON THE ELEMENTS
119 !
120  DO ielem = 1 , nelem
121 !
122 ! INITIALISES THE GEOMETRICAL VARIABLES
123 !
124  y2 = yel(ielem,2)
125  y3 = yel(ielem,3)
126 !
127  f1 = f(ikle1(ielem))
128  f2 = f(ikle2(ielem))
129  f3 = f(ikle3(ielem))
130 !
131 ! EXTRADIAGONAL TERMS
132 !
133  a12(ielem) = (y3-y2) * (-f1+2.d0*f2-f3) * (xmul/120.d0)
134  a13(ielem) = -(y3-y2) * (f1+f2-2.d0*f3) * (xmul/120.d0)
135  a14(ielem) = (y3-y2) * (2.d0*f1+f3+2.d0*f2) * (xmul/30.d0)
136  a15(ielem) = (y3-y2) * (f1+2.d0*f3+2.d0*f2) * (xmul/30.d0)
137  a16(ielem) = (y3-y2) * (f2+2.d0*f1+2.d0*f3) * (xmul/30.d0)
138  a21(ielem) = y3 * (f2-2.d0*f1+f3) * (xmul/120.d0)
139  a23(ielem) = -y3 * (-f1-f2+2.d0*f3) * (xmul/120.d0)
140  a24(ielem) = -y3 * (2.d0*f1+f3+2.d0*f2) * (xmul/30.d0)
141  a25(ielem) = -y3 * (f1+2.d0*f3+2.d0*f2) * (xmul/30.d0)
142  a26(ielem) = -y3 * (f2+2.d0*f1+2.d0*f3) * (xmul/30.d0)
143  a31(ielem) = -y2 * (f2-2.d0*f1+f3) * (xmul/120.d0)
144  a32(ielem) = y2 * (-f1+2.d0*f2-f3) * (xmul/120.d0)
145  a34(ielem) = y2 * (2.d0*f1+f3+2.d0*f2) * (xmul/30.d0)
146  a35(ielem) = y2 * (f1+2.d0*f3+2.d0*f2) * (xmul/30.d0)
147  a36(ielem) = y2 * (f2+2.d0*f1+2.d0*f3) * (xmul/30.d0)
148 !
149 ! DIAGONAL TERMS
150 ! (SUM OF EACH LINE IN THE MATRIX IS 0)
151 !
152  a11(ielem) = - a21(ielem) - a31(ielem)
153  a22(ielem) = - a12(ielem) - a32(ielem)
154  a33(ielem) = - a13(ielem) - a23(ielem)
155 !
156  ENDDO ! IELEM
157 !
158  ELSEIF(icoord.EQ.2) THEN
159 !
160 !================================
161 ! CASE OF DERIVATIVE WRT Y =
162 !================================
163 !
164  DO ielem = 1 , nelem
165 !
166 ! INITIALISES THE GEOMETRICAL VARIABLES
167 !
168  x2 = xel(ielem,2)
169  x3 = xel(ielem,3)
170 !
171  f1 = f(ikle1(ielem))
172  f2 = f(ikle2(ielem))
173  f3 = f(ikle3(ielem))
174 !
175 ! EXTRADIAGONAL TERMS
176 !
177  a12(ielem) = (x3-x2 ) * (f1-2.d0*f2+f3 ) * (xmul/120.d0)
178  a13(ielem) = (-x3+x2) * (-f1-f2+2.d0*f3) * (xmul/120.d0)
179  a14(ielem) = (-x3+x2) * (2.d0*f1+f3+2.d0*f2) * (xmul/30.d0)
180  a15(ielem) = (-x3+x2) * (f1+2.d0*f3+2.d0*f2) * (xmul/30.d0)
181  a16(ielem) = (-x3+x2) * (2.d0*f3+2.d0*f1+f2) * (xmul/30.d0)
182  a21(ielem) = -x3 * (f3-2.d0*f1+f2 ) * (xmul/120.d0)
183  a23(ielem) = x3 * (-f1-f2+2.d0*f3) * (xmul/120.d0)
184  a24(ielem) = x3 * (2.d0*f1+f3+2.d0*f2) * (xmul/30.d0)
185  a25(ielem) = x3 * (f1+2.d0*f3+2.d0*f2) * (xmul/30.d0)
186  a26(ielem) = x3 * (2.d0*f3+2.d0*f1+f2) * (xmul/30.d0)
187  a31(ielem) = x2 * (f3-2.d0*f1+f2 ) * (xmul/120.d0)
188  a32(ielem) = x2 * (f1-2.d0*f2+f3 ) * (xmul/120.d0)
189  a34(ielem) = -x2 * (2.d0*f1+f3+2.d0*f2) * (xmul/30.d0)
190  a35(ielem) = -x2 * (f1+2.d0*f3+2.d0*f2) * (xmul/30.d0)
191  a36(ielem) = -x2 * (2.d0*f3+2.d0*f1+f2) * (xmul/30.d0)
192 !
193 ! DIAGONAL TERMS
194 ! (SUM OF EACH LINE IN THE MATRIX IS 0)
195 !
196  a11(ielem) = - a21(ielem) - a31(ielem)
197  a22(ielem) = - a12(ielem) - a32(ielem)
198  a33(ielem) = - a13(ielem) - a23(ielem)
199 !
200  ENDDO ! IELEM
201 !
202  ELSE
203 !
204  WRITE(lu,201) icoord
205  CALL plante(0)
206  stop
207  ENDIF
208 !
209 !
210 !-----------------------------------------------------------------------
211 ! CASE WHERE F IS OF P2 DISCRETISATION
212 !-----------------------------------------------------------------------
213 !
214  ELSEIF(ielmf.EQ.13) THEN
215 !
216 !================================
217 ! CASE OF DERIVATIVE WRT X =
218 !================================
219 !
220  IF(icoord.EQ.1) THEN
221 !
222 ! LOOP ON THE ELEMENTS
223 !
224  DO ielem = 1 , nelem
225 !
226 ! INITIALISES THE GEOMETRICAL VARIABLES
227 !
228  y2 = yel(ielem,2)
229  y3 = yel(ielem,3)
230 !
231  f1 = f(ikle1(ielem))
232  f2 = f(ikle2(ielem))
233  f3 = f(ikle3(ielem))
234  f4 = f(ikle4(ielem))
235  f5 = f(ikle5(ielem))
236  f6 = f(ikle6(ielem))
237 !
238 ! EXTRADIAGONAL TERMS
239 !
240  a12(ielem) = (-y3+y2) *(f1-6.d0*f2+f3+4.d0*f6) * (xmul/360.d0)
241  a13(ielem) = (-y3+y2) *(f1+f2-6.d0*f3+4.d0*f4) * (xmul/360.d0)
242  a14(ielem) = (-y3+y2) *(f3-8.d0*f4-4.d0*f6-4.d0*f5) *(xmul/90.d0)
243  a15(ielem) = (-y3+y2) *(f1-4.d0*f4-4.d0*f6-8.d0*f5) *(xmul/90.d0)
244  a16(ielem) = (f2-4.d0*f4-8.d0*f6-4.d0*f5) *(-y3+y2) *(xmul/90.d0)
245  a21(ielem) =-y3 *(6.d0*f1-f2-f3-4.d0*f5) * (xmul/360.d0)
246  a23(ielem) = y3 *(f1+f2-6.d0*f3+4.d0*f4) * (xmul/360.d0)
247  a24(ielem) = y3 *(f3-8.d0*f4-4.d0*f6-4.d0*f5) *(xmul/90.d0)
248  a25(ielem) = y3 *(f1-4.d0*f4-4.d0*f6-8.d0*f5) *(xmul/90.d0)
249  a26(ielem) = y3 *(f2-4.d0*f4-8.d0*f6-4.d0*f5) *(xmul/90.d0)
250  a31(ielem) = y2 *(6.d0*f1-f2-f3-4.d0*f5) * (xmul/360.d0)
251  a32(ielem) =-y2 *(f1-6.d0*f2+f3+4.d0*f6) * (xmul/360.d0)
252  a34(ielem) =-y2 *(f3-8.d0*f4-4.d0*f6-4.d0*f5) *(xmul/90.d0)
253  a35(ielem) =-y2*(f1-4.d0*f4-4.d0*f6-8.d0*f5) *( xmul/90.d0)
254  a36(ielem) =-y2*(f2-4.d0*f4-8.d0*f6-4.d0*f5) * (xmul/90.d0)
255 !
256 ! DIAGONAL TERMS
257 ! (SUM OF EACH LINE IN THE MATRIX IS 0)
258 !
259  a11(ielem) = - a21(ielem) - a31(ielem)
260  a22(ielem) = - a12(ielem) - a32(ielem)
261  a33(ielem) = - a13(ielem) - a23(ielem)
262 !
263  ENDDO ! IELEM
264 !
265  ELSEIF(icoord.EQ.2) THEN
266 !
267 !================================
268 ! CASE OF DERIVATIVE WRT Y =
269 !================================
270 !
271  DO ielem = 1 , nelem
272 !
273 ! INITIALISES THE GEOMETRICAL VARIABLES
274 !
275  x2 = xel(ielem,2)
276  x3 = xel(ielem,3)
277 !
278  f1 = f(ikle1(ielem))
279  f2 = f(ikle2(ielem))
280  f3 = f(ikle3(ielem))
281  f4 = f(ikle4(ielem))
282  f5 = f(ikle5(ielem))
283  f6 = f(ikle6(ielem))
284 !
285 ! EXTRADIAGONAL TERMS
286 !
287  a12(ielem) = (x3-x2) *(f1-6.d0*f2+f3+4.d0*f6) * (xmul/360.d0)
288  a13(ielem) = (x3-x2) *(f1+f2-6.d0*f3+4.d0*f4) * (xmul/360.d0)
289  a14(ielem) = (x3-x2) *(f3-8.d0*f4-4.d0*f5-4.d0*f6) *(xmul/90.d0)
290  a15(ielem) = (x3-x2) *(f1-4.d0*f4-8.d0*f5-4.d0*f6) *(xmul/90.d0)
291  a16(ielem) = (x3-x2) *(f2-4.d0*f4-4.d0*f5-8.d0*f6) *(xmul/90.d0)
292  a21(ielem) = x3 *(6.d0*f1-f2-f3-4.d0*f5) * (xmul/360.d0)
293  a23(ielem) =-x3 *(f1+f2-6.d0*f3+4.d0*f4) * (xmul/360.d0)
294  a24(ielem) =-x3 *(f3-8.d0*f4-4.d0*f5-4.d0*f6) * (xmul/90.d0)
295  a25(ielem) =-x3 *(f1-4.d0*f4-8.d0*f5-4.d0*f6) * (xmul/90.d0)
296  a26(ielem) =-x3 *(f2-4.d0*f4-4.d0*f5-8.d0*f6) * (xmul/90.d0)
297  a31(ielem) =-x2 *(6.d0*f1-f2-f3-4.d0*f5) * (xmul/360.d0)
298  a32(ielem) = x2 *(f1-6.d0*f2+f3+4.d0*f6) * (xmul/360.d0)
299  a34(ielem) = x2 *(f3-8.d0*f4-4.d0*f5-4.d0*f6) * (xmul/90.d0)
300  a35(ielem) = x2 *(f1-4.d0*f4-8.d0*f5-4.d0*f6) * (xmul/90.d0)
301  a36(ielem) = x2 *(f2-4.d0*f4-4.d0*f5-8.d0*f6) * (xmul/90.d0)
302 !
303 ! DIAGONAL TERMS
304 ! (SUM OF EACH LINE IN THE MATRIX IS 0)
305 !
306  a11(ielem) = - a21(ielem) - a31(ielem)
307  a22(ielem) = - a12(ielem) - a32(ielem)
308  a33(ielem) = - a13(ielem) - a23(ielem)
309 !
310  ENDDO ! IELEM
311 !
312  ELSE
313 !
314  WRITE(lu,201) icoord
315  CALL plante(1)
316  stop
317  ENDIF
318 !
319 !-----------------------------------------------------------------------
320 !
321  ELSE
322  WRITE(lu,101) ielmf
323 101 FORMAT(1x,'MT08AC (BIEF) :',/,
324  & 1x,'DISCRETIZATION OF F : ',1i6,' NOT AVAILABLE')
325  CALL plante(1)
326  stop
327  ENDIF
328 !
329 201 FORMAT(1x,'MT08AC (BIEF) : IMPOSSIBLE COMPONENT ',
330  & 1i6,' CHECK ICOORD')
331 !
332 !-----------------------------------------------------------------------
333 !
334  RETURN
335  END
declarations_special
Definition: declarations_special.F:3
declarations_special::lu
integer lu
Definition: declarations_special.F:45
mt08ac
subroutine mt08ac(A11, A12, A13, A14, A15, A16, A21, A22, A23, A24, A25, A26, A31, A32, A33, A34, A35, A36, XMUL, SF, F, XEL, YEL, IKLE1, IKLE2, IKLE3, IKLE4, IKLE5, IKLE6, NELEM, NELMAX, ICOORD)
Definition: mt08ac.f:12
bief
Definition: bief.f:3
Generated on Tue Jan 5 2021 14:37:11 for The TELEMAC-MASCARET system by doxygen 1.8.13