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