The TELEMAC-MASCARET system  trunk
mt08ba.f
Go to the documentation of this file.
1 ! *****************
2  SUBROUTINE mt08ba
3 ! *****************
4 !
5  &( a11 , a12 , a13 ,
6  & a21 , a22 , a23 ,
7  & a31 , a32 , a33 ,
8  & a41 , a42 , a43 ,
9  & xmul,sf,f,xel,yel,ikle1,ikle2,ikle3,ikle4,
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 QUASI-BUBBLE TRIANGLE
27 !+ PSI2: BASES OF TYPE P1 TRIANGLE
28 !+
29 !+ IT WOULD BE A DERIVATIVE WRT Y WITH ICOORD=2
30 !
31 !warning THE JACOBIAN MUST BE POSITIVE
32 !
33 !history J-M HERVOUET (LNH) ; C MOULIN (LNH)
34 !+ 09/12/94
35 !+ V5P1
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 !| A12 |<--| ELEMENTS OF MATRIX
53 !| A13 |<--| ELEMENTS OF MATRIX
54 !| A21 |<--| ELEMENTS OF MATRIX
55 !| A22 |<--| ELEMENTS OF MATRIX
56 !| A23 |<--| ELEMENTS OF MATRIX
57 !| A31 |<--| ELEMENTS OF MATRIX
58 !| A32 |<--| ELEMENTS OF MATRIX
59 !| A33 |<--| ELEMENTS OF MATRIX
60 !| A41 |<--| ELEMENTS OF MATRIX
61 !| A42 |<--| ELEMENTS OF MATRIX
62 !| A43 |<--| ELEMENTS OF MATRIX
63 !| F |-->| FUNCTION USED IN THE FORMULA
64 !| ICOORD |-->| 1: DERIVATIVE ALONG X, 2: ALONG Y
65 !| IKLE1 |-->| FIRST POINTS OF TRIANGLES
66 !| IKLE2 |-->| SECOND POINTS OF TRIANGLES
67 !| IKLE3 |-->| THIRD POINTS OF TRIANGLES
68 !| NELEM |-->| NUMBER OF ELEMENTS
69 !| NELMAX |-->| MAXIMUM NUMBER OF ELEMENTS
70 !| SF |-->| STRUCTURE OF FUNCTIONS F
71 !| XEL |-->| ABSCISSAE OF POINTS IN THE MESH, PER ELEMENT
72 !| YEL |-->| ORDINATES OF POINTS IN THE MESH, PER ELEMENT
73 !| XMUL |-->| MULTIPLICATION FACTOR
74 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
75 !
76  USE bief, ex_mt08ba => mt08ba
77 !
78  USE declarations_special
79  IMPLICIT NONE
80 !
81 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
82 !
83  INTEGER, INTENT(IN) :: NELEM,NELMAX,ICOORD
84  INTEGER, INTENT(IN) :: IKLE1(nelmax),IKLE2(nelmax)
85  INTEGER, INTENT(IN) :: IKLE3(nelmax),IKLE4(nelmax)
86 !
87  DOUBLE PRECISION, INTENT(INOUT) :: A11(*),A12(*),A13(*)
88  DOUBLE PRECISION, INTENT(INOUT) :: A21(*),A22(*),A23(*)
89  DOUBLE PRECISION, INTENT(INOUT) :: A31(*),A32(*),A33(*)
90  DOUBLE PRECISION, INTENT(INOUT) :: A41(*),A42(*),A43(*)
91 !
92  DOUBLE PRECISION, INTENT(IN) :: XMUL
93  DOUBLE PRECISION, INTENT(IN) :: F(*)
94 !
95 ! STRUCTURE OF F
96  TYPE(bief_obj), INTENT(IN) :: SF
97 !
98  DOUBLE PRECISION, INTENT(IN) :: XEL(nelmax,3),YEL(nelmax,3)
99 !
100 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
101 !
102  INTEGER IELEM,IELMF
103  DOUBLE PRECISION X2,X3,Y2,Y3,F1,F2,F3,F4,XSUR18,XSUR72,XSU216
104 !
105 !-----------------------------------------------------------------------
106 !
107  ielmf=sf%ELM
108 !
109  xsur18 = xmul / 18.d0
110  xsur72 = xmul / 72.d0
111  xsu216 = xmul / 216.d0
112 !
113 !-----------------------------------------------------------------------
114 ! CASE WHERE F IS OF P1 DISCRETISATION
115 !-----------------------------------------------------------------------
116 !
117  IF(ielmf.EQ.11) THEN
118 !
119 !================================
120 ! CASE OF DERIVATIVE WRT X =
121 !================================
122 !
123  IF(icoord.EQ.1) THEN
124 !
125 ! LOOP ON THE ELEMENTS
126 !
127  DO ielem = 1 , nelem
128 !
129 ! INITIALISES THE GEOMETRICAL VARIABLES
130 !
131  y2 = yel(ielem,2)
132  y3 = yel(ielem,3)
133 !
134  f1 = f(ikle1(ielem))
135  f2 = f(ikle2(ielem))
136  f3 = f(ikle3(ielem))
137 !
138 ! EXTRADIAGONAL TERMS
139 !
140  a12(ielem)=(5*y3*f3+10*y3*f2+9*y3*f1-5*y2*f3-18*y2*
141  & f2-13*y2*f1)*xsu216
142  a13(ielem)=(18*y3*f3+5*y3*f2+13*y3*f1-10*y2*f3-5*y2*
143  & f2-9*y2*f1)*xsu216
144  a21(ielem)=-(5*y3*f3+9*y3*f2+10*y3*f1+4*y2*f2+8*y2*f1)*xsu216
145  a23(ielem)=-(18*y3*f3+13*y3*f2+5*y3*f1-8*y2*f3-4*y2*f2)*xsu216
146  a31(ielem)=(4*y3*f3+8*y3*f1+9*y2*f3+5*y2*f2+10*y2*f1)*xsu216
147  a32(ielem)=-(4*y3*f3+8*y3*f2-13*y2*f3-18*y2*f2-5*y2*f1)*xsu216
148  a41(ielem)=-(y3*f3+2*y3*f1-y2*f2-2*y2*f1)*xsur18
149  a42(ielem)=(y3*f3+2*y3*f2-y2*f3+y2*f1)*xsur18
150  a43(ielem)=(y3*f2-y3*f1-2*y2*f3-y2*f2)*xsur18
151 !
152 ! DIAGONAL TERMS
153 ! (SUM OF EACH COLUMN IN THE MATRIX IS 0)
154 !
155  a11(ielem) = - a21(ielem) - a31(ielem) - a41(ielem)
156  a22(ielem) = - a12(ielem) - a32(ielem) - a42(ielem)
157  a33(ielem) = - a13(ielem) - a23(ielem) - a43(ielem)
158 !
159  ENDDO ! IELEM
160 !
161  ELSEIF(icoord.EQ.2) THEN
162 !
163 !================================
164 ! CASE OF DERIVATIVE WRT Y =
165 !================================
166 !
167  DO ielem = 1 , nelem
168 !
169 ! INITIALISES THE GEOMETRICAL VARIABLES
170 !
171  x2 = xel(ielem,2)
172  x3 = xel(ielem,3)
173 !
174  f1 = f(ikle1(ielem))
175  f2 = f(ikle2(ielem))
176  f3 = f(ikle3(ielem))
177 !
178 ! EXTRADIAGONAL TERMS
179 !
180  a12(ielem)=(5*x2*f3+18*x2*f2+13*x2*f1-5*x3*f3-10*x3*
181  & f2-9*x3*f1)*xsu216
182  a13(ielem)=(10*x2*f3+5*x2*f2+9*x2*f1-18*x3*f3-5*x3*
183  & f2-13*x3*f1)*xsu216
184  a21(ielem)=(4*x2*f2+8*x2*f1+5*x3*f3+9*x3*f2+10*x3*f1 )*xsu216
185  a23(ielem)=-(8*x2*f3+4*x2*f2-18*x3*f3-13*x3*f2-5*x3*f1)*xsu216
186  a31(ielem)=-(9*x2*f3+5*x2*f2+10*x2*f1+4*x3*f3+8*x3*f1)*xsu216
187  a32(ielem)=-(13*x2*f3+18*x2*f2+5*x2*f1-4*x3*f3-8*x3*f2)*xsu216
188  a41(ielem)=-(x2*f2+2*x2*f1-x3*f3-2*x3*f1)*xsur18
189  a42(ielem)=(x2*f3-x2*f1-x3*f3-2*x3*f2)*xsur18
190  a43(ielem)=(2*x2*f3+x2*f2-x3*f2+x3*f1)*xsur18
191 !
192 ! DIAGONAL TERMS
193 ! (SUM OF EACH COLUMN IN THE MATRIX IS 0)
194 !
195  a11(ielem) = - a21(ielem) - a31(ielem) - a41(ielem)
196  a22(ielem) = - a12(ielem) - a32(ielem) - a42(ielem)
197  a33(ielem) = - a13(ielem) - a23(ielem) - a43(ielem)
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)=((2*y3-y2)*(f3+2*f4+f1)+(y3-2*y2)*(5*f4+7
238  & *f2+4*f1))*xsu216
239  a13(ielem)=((2*y3-y2)*(7*f3+5*f4+4*f1)+(y3-2*y2)*(2
240  & *f4+f2+f1))*xsu216
241  a21(ielem)=-((2*y3-y2)*(f3+2*f4+f2)+(y3+y2)*(5*f4+4*
242  & f2+7*f1))*xsu216
243  a23(ielem)=-((2*y3-y2)*(7*f3+5*f4+4*f2)+(y3+y2)*(2*
244  & f4+f2+f1))*xsu216
245  a31(ielem)=((y3+y2)*(4*f3+5*f4+7*f1)-(y3-2*y2)*(f3+2
246  & *f4+f2))*xsu216
247  a32(ielem)=((y3+y2)*(f3+2*f4+f1)-(y3-2*y2)*(4*f3+5*f4
248  & +7*f2))*xsu216
249  a41(ielem)=-(3*y3*f3+3*y3*f4-y3*f2+7*y3*f1+y2*f3-3*y2
250  & *f4-3*y2*f2-7*y2*f1)*xsur72
251  a42(ielem)=(3*y3*f3+3*y3*f4+7*y3*f2-y3*f1-4*y2*f3+4*
252  & y2*f1)*xsur72
253  a43(ielem)=(4*y3*f2-4*y3*f1-7*y2*f3-3*y2*f4-3*y2*f2+
254  & y2*f1)*xsur72
255 !
256 ! DIAGONAL TERMS
257 ! (SUM OF EACH COLUMN IN THE MATRIX IS 0)
258 !
259  a11(ielem) = - a21(ielem) - a31(ielem) - a41(ielem)
260  a22(ielem) = - a12(ielem) - a32(ielem) - a42(ielem)
261  a33(ielem) = - a13(ielem) - a23(ielem) - a43(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 !
283 ! EXTRADIAGONAL TERMS
284 !
285  a12(ielem)=((f3+2*f4+f1)*x2-2*(f3+2*f4+f1)*x3+2*(5*
286  & f4+7*f2+4*f1)*x2-(5*f4+7*f2+4*f1)*x3)*xsu216
287  a13(ielem)=((7*f3+5*f4+4*f1)*x2-2*(7*f3+5*f4+4*f1)
288  & *x3+2*(2*f4+f2+f1)*x2-(2*f4+f2+f1)*x3)*xsu216
289  a21(ielem)=-((f3+2*f4+f2)*x2-2*(f3+2*f4+f2)*x3-(5*f4+
290  & 4*f2+7*f1)*x2-(5*f4+4*f2+7*f1)*x3)*xsu216
291  a23(ielem)=-((7*f3+5*f4+4*f2)*x2-2*(7*f3+5*f4+4*f2
292  & )*x3-(2*f4+f2+f1)*x2-(2*f4+f2+f1)*x3)*xsu216
293  a31(ielem)=-((4*f3+5*f4+7*f1)*x2+(4*f3+5*f4+7*f1)*
294  & x3+2*(f3+2*f4+f2)*x2-(f3+2*f4+f2)*x3)*xsu216
295  a32(ielem)=-(2*(4*f3+5*f4+7*f2)*x2-(4*f3+5*f4+7*f2
296  & )*x3+(f3+2*f4+f1)*x2+(f3+2*f4+f1)*x3)*xsu216
297  a41(ielem)=(x2*f3-3*x2*f4-3*x2*f2-7*x2*f1+3*x3*f3+3*
298  & x3*f4-x3*f2+7*x3*f1)*xsur72
299  a42(ielem)=(4*x2*f3-4*x2*f1-3*x3*f3-3*x3*f4-7*x3*f2+
300  & x3*f1)*xsur72
301  a43(ielem)=(7*x2*f3+3*x2*f4+3*x2*f2-x2*f1-4*x3*f2+4*
302  & x3*f1)*xsur72
303 !
304 ! DIAGONAL TERMS
305 ! (SUM OF EACH COLUMN IN THE MATRIX IS 0)
306 !
307  a11(ielem) = - a21(ielem) - a31(ielem) - a41(ielem)
308  a22(ielem) = - a12(ielem) - a32(ielem) - a42(ielem)
309  a33(ielem) = - a13(ielem) - a23(ielem) - a43(ielem)
310 !
311  ENDDO ! IELEM
312 !
313  ELSE
314 !
315  WRITE(lu,201) icoord
316  CALL plante(0)
317  stop
318  ENDIF
319 !
320 !-----------------------------------------------------------------------
321 !
322  ELSE
323  WRITE(lu,101) ielmf
324 101 FORMAT(1x,'MT08BA (BIEF) :',/,
325  & 1x,'DISCRETIZATION OF F : ',1i6,' NOT AVAILABLE')
326  CALL plante(0)
327  stop
328  ENDIF
329 !
330 201 FORMAT(1x,'MT08BA (BIEF) : IMPOSSIBLE COMPONENT ',
331  & 1i6,' CHECK ICOORD')
332 !
333 !-----------------------------------------------------------------------
334 !
335  RETURN
336  END
declarations_special
Definition: declarations_special.F:3
declarations_special::lu
integer lu
Definition: declarations_special.F:45
mt08ba
subroutine mt08ba(A11, A12, A13, A21, A22, A23, A31, A32, A33, A41, A42, A43, XMUL, SF, F, XEL, YEL, IKLE1, IKLE2, IKLE3, IKLE4, NELEM, NELMAX, ICOORD)
Definition: mt08ba.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