The TELEMAC-MASCARET system  trunk
om1302.f
Go to the documentation of this file.
1 ! *****************
2  SUBROUTINE om1302
3 ! *****************
4 !
5  &(op , dm,typdim,xm,typexm, dn,typdin,xn,typexn, c,
6  & nulone,nelbor,nbor,nelmax,ndiag,nptfr,nelebx,neleb)
7 !
8 !***********************************************************************
9 ! BIEF V7P0 21/08/2010
10 !***********************************************************************
11 !
12 !brief OPERATIONS ON MATRICES
13 !code
14 !+ M: P2 TRIANGLE
15 !+ N: BOUNDARY MATRIX
16 !+ D: DIAGONAL MATRIX
17 !+ C: CONSTANT
18 !+
19 !+ OP IS A STRING OF 8 CHARACTERS, WHICH INDICATES THE OPERATION TO BE
20 !+ PERFORMED ON MATRICES M AND N, D AND C.
21 !+
22 !+ THE RESULT IS MATRIX M.
23 !+
24 !+ OP = 'M=M+N ' : ADDS N TO M
25 !
26 !code
27 !+ CONVENTION FOR THE STORAGE OF EXTRA-DIAGONAL TERMS:
28 !+
29 !+ XM( , 1) ----> M(1,2) XN(1) ----> M(1,2) BOUNDARY SEG 1
30 !+ XM( , 2) ----> M(1,3) XN(2) ----> M(1,2) BOUNDARY SEG 2
31 !+ XM( , 3) ----> M(1,4) XN(3) ----> M(1,2) BOUNDARY SEG 3
32 !+ XM( , 4) ----> M(1,5) XN(4) ----> M(1,2) BOUNDARY SEG 4
33 !+ XM( , 5) ----> M(1,6) XN(5) ----> M(1,2) BOUNDARY SEG 5
34 !+ XM( , 6) ----> M(2,3) | | |
35 !+ XM( , 7) ----> M(2,4) | | |
36 !+ XM( , 8) ----> M(2,5) | | |
37 !+ XM( , 9) ----> M(2,6)
38 !+ XM( ,10) ----> M(3,4)
39 !+ XM( ,11) ----> M(3,5) XN(1+NPTFRX) ----> M(1,3) BOUNDARY SEG 1
40 !+ XM( ,12) ----> M(3,6) XN(2+NPTFRX) ----> M(1,3) BOUNDARY SEG 2
41 !+ XM( ,13) ----> M(4,5) XN(3+NPTFRX) ----> M(1,3) BOUNDARY SEG 3
42 !+ XM( ,14) ----> M(4,6) | | |
43 !+ XM( ,15) ----> M(5,6) | | |
44 !+ XM( ,16) ----> M(2,1) | | |
45 !+ XM( ,17) ----> M(3,1)
46 !+ XM( ,18) ----> M(4,1)
47 !+ XM( ,19) ----> M(5,1)
48 !+ XM( ,20) ----> M(6,1)
49 !+ XM( ,21) ----> M(3,2)
50 !+ XM( ,22) ----> M(4,2)
51 !+ XM( ,23) ----> M(5,2)
52 !+ XM( ,24) ----> M(6,2)
53 !+ XM( ,25) ----> M(4,3) XN(6*NPTFRX-2) ----> M(3,2) NPTFRX-2
54 !+ XM( ,26) ----> M(5,3) XN(6*NPTFRX-1) ----> M(3,2) NPTFRX-1
55 !+ XM( ,27) ----> M(6,3) XN(6*NPTFRX ) ----> M(3,2) NPTFRX
56 !+ XM( ,28) ----> M(5,4)
57 !+ XM( ,29) ----> M(6,4) IF N IS SYMMETRICAL, NO STORAGE OF
58 !+ XM( ,30) ----> M(6,5) SUB-DIAGONAL TERMS
59 !
60 !history J-M HERVOUET (LNHE)
61 !+ 09/07/2008
62 !+ V5P9
63 !+
64 !
65 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
66 !+ 13/07/2010
67 !+ V6P0
68 !+ Translation of French comments within the FORTRAN sources into
69 !+ English comments
70 !
71 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
72 !+ 21/08/2010
73 !+ V6P0
74 !+ Creation of DOXYGEN tags for automated documentation and
75 !+ cross-referencing of the FORTRAN sources
76 !
77 !history J-M HERVOUET (EDF LAB, LNHE)
78 !+ 13/03/2014
79 !+ V7P0
80 !+ Now written to enable different numbering of boundary points and
81 !+ boundary segments.
82 !
83 !history S.E.BOURBAN (HRW)
84 !+ 21/03/2017
85 !+ V7P3
86 !+ Replacement of the DATA declarations by the PARAMETER associates
87 !
88 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
89 !| C |-->| A GIVEN CONSTANT USED IN OPERATION OP
90 !| DM |<->| DIAGONAL OF M
91 !| DN |-->| DIAGONAL OF N
92 !| NBOR |-->| GLOBAL NUMBER OF BOUNDARY POINTS
93 !| NDIAG |-->| NUMBER OF TERMS IN THE DIAGONAL
94 !| NELBOR |-->| FOR THE KTH BOUNDARY EDGE, GIVES THE CORRESPONDING
95 !| | | ELEMENT.
96 !| NELEB |-->| NUMBER OF BOUNDARY ELEMENTS
97 !| NELEBX |-->| MAXIMUM NUMBER OF BOUNDARY ELEMENTS
98 !| NELMAX |-->| MAXIMUM NUMBER OF ELEMENTS
99 !| NPTFR |-->| NUMBER OF BOUNDARY POINTS
100 !| NULONE |-->| GOES WITH ARRAY NELBOR. NELBOR GIVES THE
101 !| | | ADJACENT ELEMENT, NULONE GIVES THE LOCAL
102 !| | | NUMBER OF THE FIRST NODE OF THE BOUNDARY EDGE
103 !| | | I.E. 1, 2 OR 3 FOR TRIANGLES.
104 !| OP |-->| OPERATION TO BE DONE (SEE ABOVE)
105 !| TYPDIM |<->| TYPE OF DIAGONAL OF M:
106 !| | | TYPDIM = 'Q' : ANY VALUE
107 !| | | TYPDIM = 'I' : IDENTITY
108 !| | | TYPDIM = '0' : ZERO
109 !| TYPDIN |<->| TYPE OF DIAGONAL OF N:
110 !| | | TYPDIN = 'Q' : ANY VALUE
111 !| | | TYPDIN = 'I' : IDENTITY
112 !| | | TYPDIN = '0' : ZERO
113 !| TYPEXM |-->| TYPE OF OFF-DIAGONAL TERMS OF M:
114 !| | | TYPEXM = 'Q' : ANY VALUE
115 !| | | TYPEXM = 'S' : SYMMETRIC
116 !| | | TYPEXM = '0' : ZERO
117 !| TYPEXN |-->| TYPE OF OFF-DIAGONAL TERMS OF N:
118 !| | | TYPEXN = 'Q' : ANY VALUE
119 !| | | TYPEXN = 'S' : SYMMETRIC
120 !| | | TYPEXN = '0' : ZERO
121 !| XM |-->| OFF-DIAGONAL TERMS OF M
122 !| XN |-->| OFF-DIAGONAL TERMS OF N
123 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
124 !
125  USE bief, ex_om1302 => om1302
126 !
127  USE declarations_special
128  IMPLICIT NONE
129 !
130 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
131 !
132  INTEGER, INTENT(IN) :: NELMAX,NDIAG,NPTFR,NELEBX,NELEB
133  CHARACTER(LEN=8), INTENT(IN) :: OP
134  INTEGER, INTENT(IN) :: NULONE(nelebx),NELBOR(nelebx)
135  INTEGER, INTENT(IN) :: NBOR(*)
136  DOUBLE PRECISION, INTENT(IN) :: DN(*),XN(*)
137  DOUBLE PRECISION, INTENT(INOUT) :: DM(*),XM(nelmax,*)
138  CHARACTER(LEN=1), INTENT(INOUT) :: TYPDIM,TYPEXM,TYPDIN,TYPEXN
139  DOUBLE PRECISION, INTENT(IN) :: C
140 !
141 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
142 !
143  INTEGER K,IEL
144 !
145  DOUBLE PRECISION Z(1)
146 !
147 !-----------------------------------------------------------------------
148 ! BEWARE: ONLY WORKS FOR P2
149  INTEGER :: CORNSY(6,6)
150  parameter( cornsy = reshape( (/
151  & 1 , 6 , 17, 0 , 0 , 0 ,
152  & 3 , 8 , 12, 0 , 0 , 0 ,
153  & 7 , 11, 5 , 0 , 0 , 0 ,
154  & 16, 21, 2 , 0 , 0 , 0 ,
155  & 18, 23, 27, 0 , 0 , 0 ,
156  & 22, 26, 20, 0 , 0 , 0 /), shape=(/ 6,6 /) ) )
157  INTEGER :: CORSYM(6,3)
158  parameter( corsym = reshape( (/
159  & 1 , 6 , 2 , 0 , 0 , 0 ,
160  & 3 , 8 , 12, 0 , 0 , 0 ,
161  & 7 , 11, 5 , 0 , 0 , 0 /), shape=(/ 6,3 /) ) )
162 !
163 !-----------------------------------------------------------------------
164 !
165  IF(op(1:8).EQ.'M=M+N ') THEN
166 !
167  IF(typdim.EQ.'Q'.AND.typdim.EQ.'Q'.AND.ndiag.GE.nptfr) THEN
168  CALL ovdb( 'X=X+Y ' , dm , dn , z , c , nbor , nptfr )
169  ELSE
170  WRITE(lu,199) typdim(1:1),op(1:8),typdin(1:1)
171 199 FORMAT(1x,'OM1302 (BIEF) : TYPDIM = ',a1,' NOT IMPLEMENTED',
172  & /,1x,'FOR THE OPERATION : ',a8,' WITH TYPDIN = ',a1)
173  CALL plante(0)
174  stop
175  ENDIF
176 !
177  IF(typexm(1:1).EQ.'Q'.AND.typexn(1:1).EQ.'Q') THEN
178 !
179 ! CASE WHERE BOTH MATRICES ARE NONSYMMETRICAL
180 !
181  DO k = 1 , neleb
182  iel = nelbor(k)
183  xm( iel , cornsy(nulone(k),1) ) =
184  & xm( iel , cornsy(nulone(k),1) ) + xn(k)
185  xm( iel , cornsy(nulone(k),2) ) =
186  & xm( iel , cornsy(nulone(k),2) ) + xn(k+nelebx)
187  xm( iel , cornsy(nulone(k),3) ) =
188  & xm( iel , cornsy(nulone(k),3) ) + xn(k+2*nelebx)
189  xm( iel , cornsy(nulone(k),4) ) =
190  & xm( iel , cornsy(nulone(k),4) ) + xn(k+3*nelebx)
191  xm( iel , cornsy(nulone(k),5) ) =
192  & xm( iel , cornsy(nulone(k),5) ) + xn(k+4*nelebx)
193  xm( iel , cornsy(nulone(k),6) ) =
194  & xm( iel , cornsy(nulone(k),6) ) + xn(k+5*nelebx)
195 !
196  ENDDO ! K
197 !
198  ELSEIF(typexm(1:1).EQ.'Q'.AND.typexn(1:1).EQ.'S') THEN
199 !
200 ! CASE WHERE M CAN BE ANYTHING AND N IS SYMMETRICAL
201 !
202  DO k = 1 , neleb
203  iel = nelbor(k)
204  xm( iel , cornsy(nulone(k),1) ) =
205  & xm( iel , cornsy(nulone(k),1) ) + xn(k)
206  xm( iel , cornsy(nulone(k),2) ) =
207  & xm( iel , cornsy(nulone(k),2) ) + xn(k+nelebx)
208  xm( iel , cornsy(nulone(k),3) ) =
209  & xm( iel , cornsy(nulone(k),3) ) + xn(k+2*nelebx)
210  xm( iel , cornsy(nulone(k),4) ) =
211  & xm( iel , cornsy(nulone(k),4) ) + xn(k)
212  xm( iel , cornsy(nulone(k),5) ) =
213  & xm( iel , cornsy(nulone(k),5) ) + xn(k+nelebx)
214  xm( iel , cornsy(nulone(k),6) ) =
215  & xm( iel , cornsy(nulone(k),6) ) + xn(k+2*nelebx)
216  ENDDO ! K
217 !
218  ELSEIF(typexm(1:1).EQ.'S'.AND.typexn(1:1).EQ.'S') THEN
219 !
220 ! CASE WHERE BOTH MATRICES ARE SYMMETRICAL
221 !
222  DO k = 1 , neleb
223  iel = nelbor(k)
224  xm( iel , corsym(nulone(k),1) ) =
225  & xm( iel , corsym(nulone(k),1) ) + xn(k)
226  xm( iel , corsym(nulone(k),2) ) =
227  & xm( iel , corsym(nulone(k),2) ) + xn(k+nelebx)
228  xm( iel , corsym(nulone(k),3) ) =
229  & xm( iel , corsym(nulone(k),3) ) + xn(k+2*nelebx)
230  ENDDO ! K
231  ELSE
232  WRITE(lu,99) typexm(1:1),op(1:8),typexn(1:1)
233 99 FORMAT(1x,'OM1302 (BIEF) : TYPEXM = ',a1,' DOES NOT GO',
234  & /,1x,'FOR THE OPERATION : ',a8,' WITH TYPEXN = ',a1)
235  CALL plante(1)
236  stop
237  ENDIF
238 !
239 !-----------------------------------------------------------------------
240 !
241  ELSEIF(op(1:8).EQ.'M=M+TN ') THEN
242 !
243  CALL ovdb( 'X=X+Y ' , dm , dn , z , c , nbor , nptfr )
244 !
245  IF(typexm(1:1).EQ.'Q'.AND.typexn(1:1).EQ.'Q') THEN
246 !
247 ! CASE WHERE BOTH MATRICES ARE NONSYMMETRICAL
248 !
249  DO k = 1 , neleb
250  iel = nelbor(k)
251  xm( iel , cornsy(nulone(k),1) ) =
252  & xm( iel , cornsy(nulone(k),1) ) + xn(k+3*nelebx)
253  xm( iel , cornsy(nulone(k),2) ) =
254  & xm( iel , cornsy(nulone(k),2) ) + xn(k+4*nelebx)
255  xm( iel , cornsy(nulone(k),3) ) =
256  & xm( iel , cornsy(nulone(k),3) ) + xn(k+5*nelebx)
257  xm( iel , cornsy(nulone(k),4) ) =
258  & xm( iel , cornsy(nulone(k),4) ) + xn(k)
259  xm( iel , cornsy(nulone(k),5) ) =
260  & xm( iel , cornsy(nulone(k),5) ) + xn(k+nelebx)
261  xm( iel , cornsy(nulone(k),6) ) =
262  & xm( iel , cornsy(nulone(k),6) ) + xn(k+2*nelebx)
263  ENDDO ! K
264 !
265  ELSEIF(typexm(1:1).EQ.'Q'.AND.typexn(1:1).EQ.'S') THEN
266 !
267 ! CASE WHERE M CAN BE ANYTHING AND N IS SYMMETRICAL
268 !
269  DO k = 1 , neleb
270  iel = nelbor(k)
271  xm( iel , cornsy(nulone(k),1) ) =
272  & xm( iel , cornsy(nulone(k),1) ) + xn(k)
273  xm( iel , cornsy(nulone(k),2) ) =
274  & xm( iel , cornsy(nulone(k),2) ) + xn(k+nelebx)
275  xm( iel , cornsy(nulone(k),3) ) =
276  & xm( iel , cornsy(nulone(k),3) ) + xn(k+2*nelebx)
277  xm( iel , cornsy(nulone(k),4) ) =
278  & xm( iel , cornsy(nulone(k),4) ) + xn(k)
279  xm( iel , cornsy(nulone(k),5) ) =
280  & xm( iel , cornsy(nulone(k),5) ) + xn(k+nelebx)
281  xm( iel , cornsy(nulone(k),6) ) =
282  & xm( iel , cornsy(nulone(k),6) ) + xn(k+2*nelebx)
283  ENDDO ! K
284 !
285  ELSEIF(typexm(1:1).EQ.'S'.AND.typexn(1:1).EQ.'S') THEN
286 !
287 ! CASE WHERE BOTH MATRICES ARE SYMMETRICAL
288 !
289  DO k = 1 , neleb
290  iel = nelbor(k)
291  xm( iel , corsym(nulone(k),1) ) =
292  & xm( iel , corsym(nulone(k),1) ) + xn(k)
293  xm( iel , corsym(nulone(k),2) ) =
294  & xm( iel , corsym(nulone(k),2) ) + xn(k+nelebx)
295  xm( iel , corsym(nulone(k),3) ) =
296  & xm( iel , corsym(nulone(k),3) ) + xn(k+2*nelebx)
297  ENDDO ! K
298 !
299  ELSE
300  WRITE(lu,99) typexm(1:1),op(1:8),typexn(1:1)
301  CALL plante(1)
302  stop
303  ENDIF
304 !
305 !-----------------------------------------------------------------------
306 !
307  ELSE
308 !
309  WRITE(lu,141) op
310 141 FORMAT(1x,'OM1302 (BIEF) : UNKNOWN OPERATION : ',a8)
311  CALL plante(1)
312  stop
313 !
314  ENDIF
315 !
316 !-----------------------------------------------------------------------
317 !
318  RETURN
319  END
declarations_special
Definition: declarations_special.F:3
declarations_special::lu
integer lu
Definition: declarations_special.F:45
ovdb
subroutine ovdb(OP, X, Y, Z, C, NBOR, NPTFR)
Definition: ovdb.f:7
om1302
subroutine om1302(OP, DM, TYPDIM, XM, TYPEXM, DN, TYPDIN, XN, TYPEXN, C, NULONE, NELBOR, NBOR, NELMAX, NDIAG, NPTFR, NELEBX, NELEB)
Definition: om1302.f:8
bief
Definition: bief.f:3
Generated on Tue Jan 5 2021 14:37:20 for The TELEMAC-MASCARET system by doxygen 1.8.13