The TELEMAC-MASCARET system  trunk
mt06cc.f
Go to the documentation of this file.
1 ! *****************
2  SUBROUTINE mt06cc
3 ! *****************
4 !
5  &( a11 , a12 , a13 , a14 , a15 , a16 ,
6  & a22 , a23 , a24 , a25 , a26 ,
7  & a33 , a34 , a35 , a36 ,
8  & a44 , a45 , a46 ,
9  & a55 , a56 ,
10  & a66 ,
11  & xmul,sf,f,surfac,
12  & ikle1,ikle2,ikle3,ikle4,ikle5,ikle6,nelem,nelmax)
13 !
14 !***********************************************************************
15 ! BIEF V6P1 21/08/2010
16 !***********************************************************************
17 !
18 !brief BUILDS THE FOLLOWING MATRIX:
19 !code
20 !+ SUM(F*PSII*PSIJ)
21 !+
22 !+ WITH: QUADRATIC P1
23 !+ QUADRATIC P2
24 !+ F P1 OR QUADRATIC
25 !
26 !warning THE JACOBIAN MUST BE POSITIVE
27 !
28 !history ALGIANE FROEHLY (MATMECA)
29 !+ 19/06/08
30 !+ V5P9
31 !+
32 !
33 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
34 !+ 13/07/2010
35 !+ V6P0
36 !+ Translation of French comments within the FORTRAN sources into
37 !+ English comments
38 !
39 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
40 !+ 21/08/2010
41 !+ V6P0
42 !+ Creation of DOXYGEN tags for automated documentation and
43 !+ cross-referencing of the FORTRAN sources
44 !
45 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
46 !| A11 |<--| ELEMENTS OF MATRIX
47 !| ... |<--| ELEMENTS OF MATRIX
48 !| A66 |<--| ELEMENTS OF MATRIX
49 !| F |-->| FUNCTION F USED IN THE FORMULA
50 !| IKLE1 |-->| FIRST POINTS OF TRIANGLES
51 !| IKLE2 |-->| SECOND POINTS OF TRIANGLES
52 !| IKLE3 |-->| THIRD POINTS OF TRIANGLES
53 !| IKLE4 |-->| FOURTH POINTS OF TRIANGLES (QUADRATIC)
54 !| IKLE5 |-->| FIFTH POINTS OF TRIANGLES (QUADRATIC)
55 !| IKLE6 |-->| SIXTH POINTS OF TRIANGLES (QUADRATIC)
56 !| NELEM |-->| NUMBER OF ELEMENTS
57 !| NELMAX |-->| MAXIMUM NUMBER OF ELEMENTS
58 !| SF |-->| BIEF_OBJ STRUCTURE OF F
59 !| SURFAC |-->| AREA OF TRIANGLES
60 !| XMUL |-->| MULTIPLICATION FACTOR
61 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
62 !
63  USE bief!, EX_MT06CC => MT06CC
64 !
65  USE declarations_special
66  IMPLICIT NONE
67 !
68 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
69 !
70  INTEGER, INTENT(IN) :: NELEM,NELMAX
71  INTEGER, INTENT(IN) :: IKLE1(nelmax),IKLE2(nelmax)
72  INTEGER, INTENT(IN) :: IKLE3(nelmax),IKLE4(nelmax)
73  INTEGER, INTENT(IN) :: IKLE5(nelmax),IKLE6(nelmax)
74 !
75  DOUBLE PRECISION, INTENT(INOUT) :: A11(*),A12(*),A13(*)
76  DOUBLE PRECISION, INTENT(INOUT) :: A14(*),A15(*),A16(*)
77  DOUBLE PRECISION, INTENT(INOUT) :: A22(*),A23(*)
78  DOUBLE PRECISION, INTENT(INOUT) :: A24(*),A25(*),A26(*)
79  DOUBLE PRECISION, INTENT(INOUT) :: A33(*)
80  DOUBLE PRECISION, INTENT(INOUT) :: A34(*),A35(*),A36(*)
81  DOUBLE PRECISION, INTENT(INOUT) :: A44(*),A45(*),A46(*)
82  DOUBLE PRECISION, INTENT(INOUT) :: A55(*),A56(*)
83  DOUBLE PRECISION, INTENT(INOUT) :: A66(*)
84 !
85  DOUBLE PRECISION, INTENT(IN) :: XMUL
86  DOUBLE PRECISION, INTENT(IN) :: F(*)
87 !
88 ! STRUCTURE OF F
89  TYPE(bief_obj), INTENT(IN) :: SF
90 !
91  DOUBLE PRECISION, INTENT(IN) :: SURFAC(nelmax)
92 !
93 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
94 !
95 ! DECLARATIONS SPECIFIC TO THIS SUBROUTINE
96 !
97  DOUBLE PRECISION F1,F2,F3,F4,F5,F6,XSUR030,XSUR045,XSUR180
98  DOUBLE PRECISION XSUR315,XSUR210,XSUR630,XSU1260
99  DOUBLE PRECISION AUX315,AUX210,AUX630,AUX1260
100  INTEGER IELMF,IELEM
101 !
102 !=======================================================================
103 !
104 ! EXTRACTS THE TYPE OF ELEMENT FOR F
105 !
106  ielmf = sf%ELM
107 !
108 ! CASE WHERE F IS P0
109 !
110  IF(ielmf.EQ.10) THEN
111 !
112  xsur030 = xmul / 30.d0
113  xsur045 = xmul / 45.d0
114  xsur180 = xmul /180.d0
115 !
116  DO ielem = 1 , nelem
117 !
118 ! INITIALISES THE GEOMETRICAL VARIABLES
119 !
120  f1 = f(ielem) * surfac(ielem)
121 !
122 ! DIAGONAL TERMS
123 !
124  a11(ielem) = xsur030 * f1
125  a22(ielem) = a11(ielem)
126  a33(ielem) = a11(ielem)
127  a44(ielem) = 8.d0 * xsur045 * f1
128  a55(ielem) = a44(ielem)
129  a66(ielem) = a44(ielem)
130 !
131 ! EXTRADIAGONAL TERMS
132 !
133  a12(ielem) = - xsur180 * f1
134  a13(ielem) = a12(ielem)
135  a14(ielem) = 0.d0
136  a15(ielem) = - xsur045*f1
137  a16(ielem) = 0.d0
138 !
139  a23(ielem) = a12(ielem)
140  a24(ielem) = 0.d0
141  a25(ielem) = 0.d0
142  a26(ielem) = a15(ielem)
143 !
144  a34(ielem) = a15(ielem)
145  a35(ielem) = 0.d0
146  a36(ielem) = 0.d0
147 !
148  a45(ielem) = 4.d0*xsur045*f1
149  a46(ielem) = a45(ielem)
150 !
151  a56(ielem) = a45(ielem)
152 !
153  ENDDO ! IELEM
154 !
155 !-----------------------------------------------------------------------
156 !
157 ! CASE WHERE F IS LINEAR
158 !
159  ELSEIF(ielmf.EQ.11) THEN
160 !
161  xsur210 = xmul / 210.d0
162  xsur315 = xmul / 315.d0
163  xsu1260 = xmul / 1260.d0
164 !
165  DO ielem = 1 , nelem
166  aux210 = surfac(ielem) * xsur210
167  aux315 = surfac(ielem) * xsur315
168  aux1260= surfac(ielem) * xsu1260
169 !
170 ! INITIALISES THE GEOMETRICAL VARIABLES
171 !
172  f1 = f(ikle1(ielem))
173  f2 = f(ikle2(ielem))
174  f3 = f(ikle3(ielem))
175 !
176 ! INITIALISES THE INTERMEDIATE VARIABLES
177 !
178 !
179 ! DIAGONAL TERMS
180 !
181  a11(ielem) = (5.d0*f1+ f2+ f3)*aux210
182  a22(ielem) = ( f1+5.d0*f2+ f3)*aux210
183  a33(ielem) = ( f1+ f2+5.d0*f3)*aux210
184  a44(ielem) = 8.d0*(3.d0*f1+3.d0*f2+ f3)*aux315
185  a55(ielem) = 8.d0*( f1+3.d0*f2+3.d0*f3)*aux315
186  a66(ielem) = 8.d0*(3.d0*f1+ f2+3.d0*f3)*aux315
187 !
188 ! EXTRADIAGONAL TERMS
189 !
190  a12(ielem) = -(4.d0*f1+4.d0*f2- f3)*aux1260
191  a13(ielem) = -(4.d0*f1- f2+4.d0*f3)*aux1260
192  a14(ielem) = (3.d0*f1-2.d0*f2- f3)*aux315
193  a15(ielem) = -( f1+3.d0*f2+3.d0*f3)*aux315
194  a16(ielem) = (3.d0*f1- f2-2.d0*f3)*aux315
195 !
196  a23(ielem) = ( f1-4.d0*f2-4.d0*f3)*aux1260
197  a24(ielem) = -(2.d0*f1-3.d0*f2+ f3)*aux315
198  a25(ielem) = -( f1-3.d0*f2+2.d0*f3)*aux315
199  a26(ielem) = -(3.d0*f1+ f2+3.d0*f3)*aux315
200 !
201  a34(ielem) = -(3.d0*f1+3.d0*f2+ f3)*aux315
202  a35(ielem) = -( f1+2.d0*f2-3.d0*f3)*aux315
203  a36(ielem) = -(2.d0*f1+ f2-3.d0*f3)*aux315
204 !
205  a45(ielem) = (2.d0*f1+3.d0*f2+2.d0*f3)*aux315*4.d0
206  a46(ielem) = (3.d0*f1+2.d0*f2+2.d0*f3)*aux315*4.d0
207 !
208  a56(ielem) = (2.d0*f1+3.d0*f3+2.d0*f2)*aux315*4.d0
209 !
210  ENDDO ! IELEM
211 !
212 !-----------------------------------------------------------------------
213 !
214  ELSEIF(ielmf.EQ.13) THEN
215 !
216 !-----------------------------------------------------------------------
217 !
218 ! QUADRATIC DISCRETISATION OF F:
219 !
220  xsur315 = xmul / 315.d0
221  xsur630 = xmul / 630.d0
222  xsu1260 = xmul / 1260.d0
223 !
224  DO ielem = 1 , nelem
225  aux315 = surfac(ielem) * xsur315
226  aux630 = surfac(ielem) * xsur630
227  aux1260= surfac(ielem) * xsu1260
228 !
229 ! INITIALISES THE GEOMETRICAL VARIABLES
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 ! DIAGONAL TERMS
239 !
240  a11(ielem) = (6.d0*(f4+f6)+9.d0*f1+2.d0*f5-f2-f3) * aux630
241  a22(ielem) = (6.d0*(f4+f5)+2.d0*f6+9.d0*f2-f1-f3) * aux630
242  a33(ielem) = (6.d0*(f6+f5)+9.d0*f3+2.d0*f4-f1-f2) * aux630
243  a44(ielem) = 4.d0*(3.d0*(f6+f5)+9.d0*f4-f3) * aux315
244  a55(ielem) = 4.d0*(3.d0*(f4+f6)+9.d0*f5-f1) * aux315
245  a66(ielem) = 4.d0*(3.d0*(f4+f5)+9.d0*f6-f2) * aux315
246 !
247 ! EXTRADIAGONAL TERMS
248 !
249  a12(ielem) = -(2.d0*(f1+f2)+4.d0*f4-f3) * aux1260
250  a13(ielem) = -(2.d0*(f1+f3)+4.d0*f6-f2) * aux1260
251  a14(ielem) = (3.d0* f1 -2.d0*f5-f2) * aux315
252  a15(ielem) = -(2.d0*(f4+f6)+4.d0*f5-f1) * aux315
253  a16(ielem) = (3.d0*f1 -2.d0*f5-f3) * aux315
254 !
255  a23(ielem) = (-2.d0*(f2+f3)-4.d0*f5+f1) * aux1260
256  a24(ielem) = (-2.d0*f6 +3.d0*f2-f1) * aux315
257  a25(ielem) = (-2.d0*f6 +3.d0*f2-f3) * aux315
258  a26(ielem) = (-2.d0*(f4+f5)-4.d0*f6+f2) * aux315
259 !
260  a34(ielem) = (-2.d0*(f6+f5)-4.d0*f4+f3) * aux315
261  a35(ielem) = (-2.d0*f4 +3.d0*f3-f2) * aux315
262  a36(ielem) = (-2.d0*f4 +3.d0*f3-f1) * aux315
263 !
264  a45(ielem) = 2.d0*(6.d0*(f4+f5)+4.d0*f6-f1-f3) * aux315
265  a46(ielem) = 2.d0*(6.d0*(f4+f6)+4.d0*f5-f2-f3) * aux315
266 !
267  a56(ielem) = 2.d0*(6.d0*(f6+f5)+4.d0*f4-f2-f1) * aux315
268 !
269  ENDDO ! IELEM
270 !
271  ELSE
272 !
273  WRITE(lu,101) ielmf,sf%NAME
274 101 FORMAT(1x,'MT06CC (BIEF) :',/,
275  & 1x,'DISCRETIZATION OF F NOT AVAILABLE:',1i6,
276  & 1x,'REAL NAME: ',a6)
277  CALL plante(1)
278  stop
279 !
280  ENDIF
281 !
282 !-----------------------------------------------------------------------
283 !
284  RETURN
285  END
declarations_special
Definition: declarations_special.F:3
mt06cc
subroutine mt06cc(A11, A12, A13, A14, A15, A16, A22, A23, A24, A25, A26, A33, A34, A35, A36, A44, A45, A46, A55, A56, A66, XMUL, SF, F, SURFAC, IKLE1, IKLE2, IKLE3, IKLE4, IKLE5, IKLE6, NELEM, NELMAX)
Definition: mt06cc.f:14
declarations_special::lu
integer lu
Definition: declarations_special.F:45
bief
Definition: bief.f:3
Generated on Tue Jan 5 2021 14:37:08 for The TELEMAC-MASCARET system by doxygen 1.8.13