The TELEMAC-MASCARET system  trunk
om2121.f
Go to the documentation of this file.
1 ! *****************
2  SUBROUTINE om2121
3 ! *****************
4 !
5  &(op , dm,typdim,xm,typexm, dn,typdin,xn,typexn, d,c,
6  & ikle,nelem,nelmax,ndiag)
7 !
8 !***********************************************************************
9 ! BIEF V6P1 21/08/2010
10 !***********************************************************************
11 !
12 !brief OPERATIONS ON MATRICES WITH Q1 QUADRILATERAL
13 !+ OR ANY OTHER ELEMENT WITH THE SAME NUMBER OF POINTS.
14 !code
15 !+ D: DIAGONAL MATRIX
16 !+ C: CONSTANT
17 !+
18 !+ OP IS A STRING OF 8 CHARACTERS, WHICH INDICATES THE OPERATION TO BE
19 !+ PERFORMED ON MATRICES M AND N, D AND C.
20 !+
21 !+ THE RESULT IS MATRIX M.
22 !+
23 !+ OP = 'M=N ' : COPIES N IN M
24 !+ OP = 'M=TN ' : COPIES TRANSPOSE OF N IN M
25 !+ OP = 'M=CN ' : MULTIPLIES N BY C
26 !+ OP = 'M=CM ' : MULTIPLIES M BY C
27 !+ OP = 'M=M+CN ' : ADDS CN TO M
28 !+ OP = 'M=M+CTN ' : ADDS C TRANSPOSE(N) TO M
29 !+ OP = 'M=M+TN ' : ADDS TRANSPOSE(N) TO M
30 !+ OP = 'M=M+N ' : ADDS N TO M
31 !+ OP = 'M=MD ' : M X D
32 !+ OP = 'M=DM ' : D X M
33 !+ OP = 'M=DMD ' : D X M X D
34 !+ OP = 'M=M+D ' : ADDS D TO M
35 !+ OP = 'M=X(M) ' : NOT SYMMETRICAL FORM OF M
36 !+ (OLD MATSNS)
37 !+ OP = 'M=MSK(M)' : MASKS M EXTRADIAGONAL TERMS
38 !+ (OLD MASKEX)
39 !+ THE MASK IS TAKEN FROM D
40 !
41 !code
42 !+ CONVENTION FOR THE STORAGE OF EXTRA-DIAGONAL TERMS:
43 !+
44 !+ XM(IELEM, 1) ----> M(1,2)
45 !+ XM(IELEM, 2) ----> M(1,3)
46 !+ XM(IELEM, 3) ----> M(1,4)
47 !+ XM(IELEM, 4) ----> M(2,3)
48 !+ XM(IELEM, 5) ----> M(2,4)
49 !+ XM(IELEM, 6) ----> M(3,4)
50 !+ XM(IELEM, 7) ----> M(2,1)
51 !+ XM(IELEM, 8) ----> M(3,1)
52 !+ XM(IELEM, 9) ----> M(4,1)
53 !+ XM(IELEM,10) ----> M(3,2)
54 !+ XM(IELEM,11) ----> M(4,2)
55 !+ XM(IELEM,12) ----> M(4,3)
56 !
57 !history J-M HERVOUET (LNHE) ; F LEPEINTRE (LNH)
58 !+ 21/09/93
59 !+ V5P6
60 !+
61 !
62 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
63 !+ 13/07/2010
64 !+ V6P0
65 !+ Translation of French comments within the FORTRAN sources into
66 !+ English comments
67 !
68 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
69 !+ 21/08/2010
70 !+ V6P0
71 !+ Creation of DOXYGEN tags for automated documentation and
72 !+ cross-referencing of the FORTRAN sources
73 !
74 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
75 !| C |-->| A GIVEN CONSTANT USED IN OPERATION OP
76 !| D |-->| A DIAGONAL MATRIX
77 !| DM |<->| DIAGONAL OF M
78 !| DN |-->| DIAGONAL OF N
79 !| IKLE |-->| CONNECTIVITY TABLE.
80 !| NDIAG |-->| NUMBER OF TERMS IN THE DIAGONAL
81 !| NELEM |-->| NUMBER OF ELEMENTS
82 !| NELMAX |-->| MAXIMUM NUMBER OF ELEMENTS
83 !| OP |-->| OPERATION TO BE DONE (SEE ABOVE)
84 !| TYPDIM |<->| TYPE OF DIAGONAL OF M:
85 !| | | TYPDIM = 'Q' : ANY VALUE
86 !| | | TYPDIM = 'I' : IDENTITY
87 !| | | TYPDIM = '0' : ZERO
88 !| TYPDIN |<->| TYPE OF DIAGONAL OF N:
89 !| | | TYPDIN = 'Q' : ANY VALUE
90 !| | | TYPDIN = 'I' : IDENTITY
91 !| | | TYPDIN = '0' : ZERO
92 !| TYPEXM |-->| TYPE OF OFF-DIAGONAL TERMS OF M:
93 !| | | TYPEXM = 'Q' : ANY VALUE
94 !| | | TYPEXM = 'S' : SYMMETRIC
95 !| | | TYPEXM = '0' : ZERO
96 !| TYPEXN |-->| TYPE OF OFF-DIAGONAL TERMS OF N:
97 !| | | TYPEXN = 'Q' : ANY VALUE
98 !| | | TYPEXN = 'S' : SYMMETRIC
99 !| | | TYPEXN = '0' : ZERO
100 !| XM |-->| OFF-DIAGONAL TERMS OF M
101 !| XN |-->| OFF-DIAGONAL TERMS OF N
102 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
103 !
104  USE bief, ex_om2121 => om2121
105 !
106  USE declarations_special
107  IMPLICIT NONE
108 !
109 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
110 !
111  INTEGER, INTENT(IN) :: NELEM,NELMAX,NDIAG
112  INTEGER, INTENT(IN) :: IKLE(nelmax,4)
113  CHARACTER(LEN=8), INTENT(IN) :: OP
114  DOUBLE PRECISION, INTENT(IN) :: DN(*),D(*),XN(nelmax,*)
115  DOUBLE PRECISION, INTENT(INOUT) :: DM(*),XM(nelmax,*)
116  CHARACTER(LEN=1), INTENT(INOUT) :: TYPDIM,TYPEXM,TYPDIN,TYPEXN
117  DOUBLE PRECISION, INTENT(IN) :: C
118 !
119 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
120 !
121  INTEGER I,J,IELEM
122 !
123 !-----------------------------------------------------------------------
124 !
125  IF(op(3:8).EQ.'N ') THEN
126 !
127  IF(typdin(1:1).EQ.'Q') THEN
128  CALL ov('X=Y ', x=dm, y=dn, dim1=ndiag)
129  ELSEIF(typdin(1:1).EQ.'I'.OR.typdin(1:1).EQ.'0') THEN
130 ! NOTHING TO DO, ONLY NEEDS TO COPY TYPDIN
131  ELSE
132  WRITE(lu,6) typdin(1:1)
133 6 FORMAT(1x,'OM2121 (BIEF) : TYPDIN UNKNOWN :',a1)
134  CALL plante(1)
135  stop
136  ENDIF
137  typdim(1:1)=typdin(1:1)
138 !
139  IF(typexn(1:1).EQ.'S') THEN
140  DO i=1,6
141  CALL ov('X=Y ', x=xm(1,i), y=xn(1,i), dim1=nelem)
142  ENDDO ! I
143  ELSEIF(typexn(1:1).EQ.'Q') THEN
144  DO i=1,12
145  CALL ov('X=Y ', x=xm(1,i), y=xn(1,i), dim1=nelem)
146  ENDDO ! I
147  ELSEIF(typexn(1:1).NE.'0') THEN
148  WRITE(lu,40) typexn(1:1)
149 40 FORMAT(1x,'OM2121 (BIEF) : TYPEXN UNKNOWN :',a1)
150  CALL plante(1)
151  stop
152  ENDIF
153 !
154  typexm(1:1)=typexn(1:1)
155 !
156 !-----------------------------------------------------------------------
157 !
158  ELSEIF(op(3:8).EQ.'TN ') THEN
159 !
160  CALL ov('X=Y ', x=dm, y=dn, dim1=ndiag)
161 !
162  IF(typexn(1:1).EQ.'S') THEN
163  DO i=1,6
164  CALL ov('X=Y ', x=xm(1,i), y=xn(1,i), dim1=nelem)
165  ENDDO ! I
166  ELSEIF(typexn(1:1).EQ.'Q') THEN
167  DO i=1,6
168  CALL ov('X=Y ', x=xm(1,i), y=xn(1,i+6), dim1=nelem)
169  CALL ov('X=Y ', x=xm(1,i+6), y=xn(1,i), dim1=nelem)
170  ENDDO ! I
171  ELSEIF(typexn(1:1).NE.'0') THEN
172  WRITE(lu,40) typexn(1:1)
173  CALL plante(1)
174  stop
175  ENDIF
176 !
177  typdim(1:1)=typdin(1:1)
178  typexm(1:1)=typexn(1:1)
179 !
180 !-----------------------------------------------------------------------
181 !
182  ELSEIF(op(3:8).EQ.'CN ') THEN
183 !
184  CALL ov('X=CY ', x=dm, y=dn, c=c, dim1=ndiag)
185 !
186  IF(typexn(1:1).EQ.'S') THEN
187  DO i=1,6
188  CALL ov('X=CY ', x=xm(1,i), y=xn(1,i), c=c, dim1=nelem)
189  ENDDO ! I
190  ELSEIF(typexn(1:1).EQ.'Q') THEN
191  DO i=1,12
192  CALL ov('X=CY ', x=xm(1,i), y=xn(1,i), c=c, dim1=nelem)
193  ENDDO ! I
194  ELSEIF(typexn(1:1).NE.'0') THEN
195  WRITE(lu,40) typexn(1:1)
196  CALL plante(1)
197  stop
198  ENDIF
199 !
200  typdim(1:1)=typdin(1:1)
201  typexm(1:1)=typexn(1:1)
202 !
203 !-----------------------------------------------------------------------
204 !
205  ELSEIF(op(3:8).EQ.'CM ') THEN
206 !
207  CALL ov('X=CX ', x=dm, c=c, dim1=ndiag)
208 !
209  IF(typexm(1:1).EQ.'S') THEN
210  DO i=1,6
211  CALL ov('X=CX ', x=xm(1,i), c=c, dim1=nelem)
212  ENDDO ! I
213  ELSEIF(typexm(1:1).EQ.'Q') THEN
214  DO i=1,12
215  CALL ov('X=CX ', x=xm(1,i), c=c, dim1=nelem)
216  ENDDO ! I
217  ELSEIF(typexm(1:1).NE.'0') THEN
218  WRITE(lu,40) typexm(1:1)
219  CALL plante(1)
220  stop
221  ENDIF
222 !
223 !-----------------------------------------------------------------------
224 !
225  ELSEIF(op(3:8).EQ.'M+CN ' .OR.
226  & (op(3:8).EQ.'M+CTN '.AND.typexn(1:1).NE.'Q') ) THEN
227 !
228  IF(typdin(1:1).EQ.'I') THEN
229  CALL ov('X=X+C ', x=dm, c=c, dim1=ndiag)
230  ELSEIF(typdin(1:1).NE.'0') THEN
231  CALL ov('X=X+CY ', x=dm, y=dn, c=c, dim1=ndiag)
232  ENDIF
233 !
234  IF(typexn(1:1).EQ.'S') THEN
235  DO i=1,6
236  CALL ov('X=X+CY ', x=xm(1,i), y=xn(1,i), c=c, dim1=nelem)
237  ENDDO ! I
238  IF(typexm(1:1).EQ.'Q') THEN
239  DO i=1,6
240  CALL ov('X=X+CY ', x=xm(1,i+6), y=xn(1,i), c=c, dim1=nelem)
241  ENDDO ! I
242  ENDIF
243  ELSEIF(typexn(1:1).EQ.'Q') THEN
244  IF(typexm(1:1).NE.'Q') THEN
245  WRITE(lu,98) typexm(1:1),op(1:8),typexn(1:1)
246 98 FORMAT(1x,'OM2121 (BIEF) : TYPEXM = ',a1,' DOES NOT GO ',
247  & /,1x,'FOR THE OPERATION : ',a8,' WITH TYPEXN = ',a1)
248  CALL plante(1)
249  stop
250  ENDIF
251  DO i=1,12
252  CALL ov('X=X+CY ', x=xm(1,i), y=xn(1,i), c=c, dim1=nelem)
253  ENDDO ! I
254  ELSEIF(typexn(1:1).NE.'0') THEN
255  WRITE(lu,40) typexn(1:1)
256  CALL plante(1)
257  stop
258  ENDIF
259 !
260 !-----------------------------------------------------------------------
261 !
262  ELSEIF(op(3:8).EQ.'M+CTN ') THEN
263 !
264 ! THE CASES WHERE N IS SYMMETRICAL ARE TREATED WITH M=M+CN
265 !
266  CALL ov('X=X+CY ', x=dm, y=dn, c=c, dim1=ndiag)
267 !
268  IF(typexn(1:1).EQ.'Q') THEN
269  IF(typexm(1:1).NE.'Q') THEN
270  WRITE(lu,98) typexm(1:1),op(1:8),typexn(1:1)
271  CALL plante(1)
272  stop
273  ENDIF
274  DO i=1,6
275  CALL ov('X=X+CY ', x=xm(1,i), y=xn(1,i+6), c=c, dim1=nelem)
276  CALL ov('X=X+CY ', x=xm(1,i+6), y=xn(1,i ),
277  & c=c, dim1=nelem )
278  ENDDO ! I
279  ELSE
280  WRITE(lu,40) typexn(1:1)
281  CALL plante(1)
282  stop
283  ENDIF
284 !
285 !-----------------------------------------------------------------------
286 !
287  ELSEIF(op(3:8).EQ.'M+N '.OR.
288  & (op(3:8).EQ.'M+TN ').AND.typexn(1:1).NE.'Q') THEN
289 !
290  CALL ov('X=X+Y ', x=dm, y=dn, dim1=ndiag)
291 !
292  IF(typexn(1:1).EQ.'S') THEN
293  DO i=1,6
294  CALL ov('X=X+Y ', x=xm(1,i), y=xn(1,i), dim1=nelem)
295  ENDDO ! I
296  IF(typexm(1:1).EQ.'Q') THEN
297  DO i=1,6
298  CALL ov('X=X+Y ', x=xm(1,i+6), y=xn(1,i), dim1=nelem)
299  ENDDO ! I
300  ENDIF
301  ELSEIF(typexn(1:1).EQ.'Q') THEN
302  IF(typexm(1:1).NE.'Q') THEN
303  WRITE(lu,98) typexm(1:1),op(1:8),typexn(1:1)
304  CALL plante(1)
305  stop
306  ENDIF
307  DO i=1,12
308  CALL ov('X=X+Y ', x=xm(1,i), y=xn(1,i), dim1=nelem)
309  ENDDO ! I
310  ELSEIF(typexn(1:1).NE.'0') THEN
311  WRITE(lu,40) typexn(1:1)
312  CALL plante(1)
313  stop
314  ENDIF
315 !
316 !-----------------------------------------------------------------------
317 !
318  ELSEIF(op(3:8).EQ.'M+TN ') THEN
319 !
320 ! THE CASES WHERE N IS SYMMETRICAL ARE TREATED WITH M=M+N
321 !
322  CALL ov('X=X+Y ', x=dm, y=dn, dim1=ndiag)
323 !
324  IF(typexn(1:1).EQ.'Q') THEN
325  IF(typexm(1:1).NE.'Q') THEN
326  WRITE(lu,98) typexm(1:1),op(1:8),typexn(1:1)
327  CALL plante(1)
328  stop
329  ENDIF
330  DO i=1,6
331  CALL ov('X=X+Y ', x=xm(1,i), y=xn(1,i+6), dim1=nelem)
332  CALL ov('X=X+Y ', x=xm(1,i+6), y=xn(1,i ), dim1=nelem)
333  ENDDO ! I
334  ELSE
335  WRITE(lu,40) typexn(1:1)
336  CALL plante(1)
337  stop
338  ENDIF
339 !
340 !-----------------------------------------------------------------------
341 !
342  ELSEIF(op(3:8).EQ.'MD ') THEN
343 !
344 ! DIAGONAL TERMS
345 !
346  IF(typdim(1:1).EQ.'Q') THEN
347  CALL ov('X=XY ', x=dm, y=d, dim1=ndiag)
348  ELSEIF(typdim(1:1).EQ.'I') THEN
349  CALL ov('X=Y ', x=dm, y=d, dim1=ndiag)
350  typdim(1:1)='Q'
351  ELSEIF(typdim(1:1).NE.'0') THEN
352  WRITE(lu,13) typdim(1:1)
353  CALL plante(1)
354  stop
355  ENDIF
356 !
357 ! EXTRADIAGONAL TERMS
358 !
359 !
360  IF(typexm(1:1).EQ.'Q') THEN
361 !
362  DO ielem = 1 , nelem
363 !
364  xm(ielem, 1) = xm(ielem, 1) * d(ikle(ielem,2))
365  xm(ielem, 2) = xm(ielem, 2) * d(ikle(ielem,3))
366  xm(ielem, 3) = xm(ielem, 3) * d(ikle(ielem,4))
367 !
368  xm(ielem, 4) = xm(ielem, 4) * d(ikle(ielem,3))
369  xm(ielem, 5) = xm(ielem, 5) * d(ikle(ielem,4))
370  xm(ielem, 6) = xm(ielem, 6) * d(ikle(ielem,4))
371 !
372  xm(ielem, 7) = xm(ielem, 7) * d(ikle(ielem,1))
373  xm(ielem, 8) = xm(ielem, 8) * d(ikle(ielem,1))
374  xm(ielem, 9) = xm(ielem, 9) * d(ikle(ielem,1))
375 !
376  xm(ielem,10) = xm(ielem,10) * d(ikle(ielem,2))
377  xm(ielem,11) = xm(ielem,11) * d(ikle(ielem,2))
378  xm(ielem,12) = xm(ielem,12) * d(ikle(ielem,3))
379 !
380  ENDDO ! IELEM
381 !
382  ELSEIF(typexm(1:1).EQ.'S') THEN
383  WRITE(lu,180)
384 180 FORMAT(1x,
385  & 'OM2121 (BIEF) : M=MD NOT AVAILABLE IF M SYMMETRIC')
386  CALL plante(1)
387  stop
388  ELSEIF(typexm(1:1).NE.'0') THEN
389  WRITE(lu,200)
390  CALL plante(1)
391  stop
392  ENDIF
393 !
394 !-----------------------------------------------------------------------
395 !
396  ELSEIF(op(3:8).EQ.'DM ') THEN
397 !
398 ! DIAGONAL TERMS
399 !
400  IF(typdim(1:1).EQ.'Q') THEN
401  CALL ov('X=XY ', x=dm, y=d, dim1=ndiag)
402  ELSEIF(typdim(1:1).EQ.'I') THEN
403  CALL ov('X=Y ', x=dm, y=d, dim1=ndiag)
404  typdim(1:1)='Q'
405  ELSEIF(typdim(1:1).NE.'0') THEN
406  WRITE(lu,13) typdim(1:1)
407  CALL plante(1)
408  stop
409  ENDIF
410 !
411 ! EXTRADIAGONAL TERMS
412 !
413  IF(typexm(1:1).EQ.'Q') THEN
414 !
415  DO ielem = 1 , nelem
416 !
417  xm(ielem, 7) = xm(ielem, 7) * d(ikle(ielem,2))
418  xm(ielem, 8) = xm(ielem, 8) * d(ikle(ielem,3))
419  xm(ielem, 9) = xm(ielem, 9) * d(ikle(ielem,4))
420 !
421  xm(ielem, 1) = xm(ielem, 1) * d(ikle(ielem,1))
422  xm(ielem,10) = xm(ielem,10) * d(ikle(ielem,3))
423  xm(ielem,11) = xm(ielem,11) * d(ikle(ielem,4))
424 !
425  xm(ielem, 2) = xm(ielem, 2) * d(ikle(ielem,1))
426  xm(ielem, 4) = xm(ielem, 4) * d(ikle(ielem,2))
427  xm(ielem,12) = xm(ielem,12) * d(ikle(ielem,4))
428 !
429  xm(ielem, 3) = xm(ielem, 3) * d(ikle(ielem,1))
430  xm(ielem, 5) = xm(ielem, 5) * d(ikle(ielem,2))
431  xm(ielem, 6) = xm(ielem, 6) * d(ikle(ielem,3))
432 !
433  ENDDO ! IELEM
434 !
435  ELSEIF(typexm(1:1).EQ.'S') THEN
436  WRITE(lu,230)
437 230 FORMAT(1x,
438  & 'OM2121 (BIEF) : M=MD NOT AVAILABLE IF M SYMMETRIC')
439  CALL plante(1)
440  stop
441  ELSEIF(typexm(1:1).NE.'0') THEN
442  WRITE(lu,200)
443 200 FORMAT(1x,'OM2121 (BIEF) : TYPEXM NOT AVAILABLE : ',a1)
444  CALL plante(1)
445  stop
446  ENDIF
447 !
448 !-----------------------------------------------------------------------
449 !
450  ELSEIF(op(3:8).EQ.'DMD ') THEN
451 !
452 ! DIAGONAL TERMS
453 !
454  IF(typdim(1:1).EQ.'Q') THEN
455  CALL ov('X=XY ', x=dm, y=d, dim1=ndiag)
456  CALL ov('X=XY ', x=dm, y=d, dim1=ndiag)
457  ELSEIF(typdim(1:1).EQ.'I') THEN
458  CALL ov('X=YZ ', x=dm, y=d, z=d, dim1=ndiag)
459  typdim(1:1)='Q'
460  ELSEIF(typdim(1:1).NE.'0') THEN
461  WRITE(lu,13) typdim(1:1)
462 13 FORMAT(1x,'OM2121 (BIEF) : TYPDIM UNKNOWN :',a1)
463  CALL plante(1)
464  stop
465  ENDIF
466 !
467 ! EXTRADIAGONAL TERMS
468 !
469  IF(typexm(1:1).EQ.'S') THEN
470 !
471  DO ielem = 1 , nelem
472  xm(ielem, 1)=xm(ielem, 1) * d(ikle(ielem,2))*d(ikle(ielem,1))
473  xm(ielem, 2)=xm(ielem, 2) * d(ikle(ielem,3))*d(ikle(ielem,1))
474  xm(ielem, 3)=xm(ielem, 3) * d(ikle(ielem,4))*d(ikle(ielem,1))
475  xm(ielem, 4)=xm(ielem, 4) * d(ikle(ielem,3))*d(ikle(ielem,2))
476  xm(ielem, 5)=xm(ielem, 5) * d(ikle(ielem,4))*d(ikle(ielem,2))
477  xm(ielem, 6)=xm(ielem, 6) * d(ikle(ielem,4))*d(ikle(ielem,3))
478  ENDDO ! IELEM
479 !
480  ELSEIF(typexm(1:1).EQ.'Q') THEN
481 !
482  DO ielem = 1 , nelem
483  xm(ielem, 1)=xm(ielem, 1) * d(ikle(ielem,2))*d(ikle(ielem,1))
484  xm(ielem, 2)=xm(ielem, 2) * d(ikle(ielem,3))*d(ikle(ielem,1))
485  xm(ielem, 3)=xm(ielem, 3) * d(ikle(ielem,4))*d(ikle(ielem,1))
486  xm(ielem, 4)=xm(ielem, 4) * d(ikle(ielem,3))*d(ikle(ielem,2))
487  xm(ielem, 5)=xm(ielem, 5) * d(ikle(ielem,4))*d(ikle(ielem,2))
488  xm(ielem, 6)=xm(ielem, 6) * d(ikle(ielem,4))*d(ikle(ielem,3))
489  xm(ielem, 7)=xm(ielem, 7) * d(ikle(ielem,2))*d(ikle(ielem,1))
490  xm(ielem, 8)=xm(ielem, 8) * d(ikle(ielem,3))*d(ikle(ielem,1))
491  xm(ielem, 9)=xm(ielem, 9) * d(ikle(ielem,4))*d(ikle(ielem,1))
492  xm(ielem,10)=xm(ielem,10) * d(ikle(ielem,3))*d(ikle(ielem,2))
493  xm(ielem,11)=xm(ielem,11) * d(ikle(ielem,4))*d(ikle(ielem,2))
494  xm(ielem,12)=xm(ielem,12) * d(ikle(ielem,4))*d(ikle(ielem,3))
495  ENDDO ! IELEM
496 !
497  ELSEIF(typexm(1:1).NE.'0') THEN
498  WRITE(lu,241) typexm(1:1)
499 241 FORMAT(1x,'OM2121 (BIEF) : TYPEXM UNKNOWN :',a1)
500  CALL plante(1)
501  stop
502  ENDIF
503 !
504 !-----------------------------------------------------------------------
505 !
506  ELSEIF(op(3:8).EQ.'M+D ') THEN
507 !
508  IF(typdim(1:1).EQ.'Q') THEN
509  CALL ov('X=X+Y ', x=dm, y=d, dim1=ndiag)
510  ELSE
511  WRITE(lu,13) typdim(1:1)
512  CALL plante(1)
513  stop
514  ENDIF
515 !
516 !-----------------------------------------------------------------------
517 !
518  ELSEIF(op(3:8).EQ.'0 ') THEN
519 !
520  CALL ov('X=C ', x=dm, c=0.d0, dim1=ndiag)
521 !
522  IF(typexm(1:1).EQ.'S') THEN
523  DO i=1,6
524  CALL ov('X=C ', x=xm(1,i), c=0.d0, dim1=nelem)
525  ENDDO
526  ELSEIF(typexm(1:1).EQ.'Q') THEN
527  DO i=1,12
528  CALL ov('X=C ', x=xm(1,i), c=0.d0, dim1=nelem)
529  ENDDO
530  ELSEIF(typexm(1:1).NE.'0') THEN
531  WRITE(lu,711) typexm(1:1)
532 711 FORMAT(1x,'OM2121 (BIEF) : TYPEXM UNKNOWN :',a1)
533  CALL plante(1)
534  stop
535  ENDIF
536 !
537  typdim(1:1)='0'
538 ! TYPEXM(1:1) IS NOT CHANGED
539 !
540 !-----------------------------------------------------------------------
541 !
542  ELSEIF(op(3:8).EQ.'X(M) ') THEN
543 !
544  IF(typexm(1:1).EQ.'S') THEN
545  CALL ov('X=Y ', x=xm(1, 7), y=xm(1,1), dim1=nelem)
546  CALL ov('X=Y ', x=xm(1, 8), y=xm(1,2), dim1=nelem)
547  CALL ov('X=Y ', x=xm(1, 9), y=xm(1,3), dim1=nelem)
548  CALL ov('X=Y ', x=xm(1,10), y=xm(1,4), dim1=nelem)
549  CALL ov('X=Y ', x=xm(1,11), y=xm(1,5), dim1=nelem)
550  CALL ov('X=Y ', x=xm(1,12), y=xm(1,6), dim1=nelem)
551  ELSEIF(typexm(1:1).NE.'0') THEN
552  WRITE(lu,811) typexm(1:1)
553 811 FORMAT(1x,'OM2121 (BIEF) : MATRIX ALREADY NON SYMMETRICAL: ',
554  & a1)
555  CALL plante(1)
556  stop
557  ENDIF
558  typexm(1:1)='Q'
559 !
560 !-----------------------------------------------------------------------
561 !
562  ELSEIF(op(3:8).EQ.'MSK(M)') THEN
563 !
564  IF(typexm(1:1).EQ.'S') THEN
565  j = 6
566  ELSEIF(typexm(1:1).EQ.'Q') THEN
567  j = 12
568  ELSEIF(typexm(1:1).EQ.'0') THEN
569  j = 0
570  ELSE
571  WRITE(lu,711) typexm
572  j = 0
573  CALL plante(1)
574  stop
575  ENDIF
576 !
577  IF(j.GT.0) THEN
578  DO i = 1,j
579  CALL ov('X=XY ', x=xm(1,i), y=d, dim1=nelem)
580  ENDDO ! I
581  ENDIF
582 !
583 !-----------------------------------------------------------------------
584 !
585  ELSE
586 !
587  WRITE(lu,42) op
588 42 FORMAT(1x,'OM2121 (BIEF) : UNKNOWN OPERATION : ',a8)
589  CALL plante(1)
590  stop
591 !
592  ENDIF
593 !
594 !-----------------------------------------------------------------------
595 !
596  RETURN
597  END
ov
subroutine ov(OP, X, Y, Z, C, DIM1)
Definition: ov.f:7
declarations_special
Definition: declarations_special.F:3
om2121
subroutine om2121(OP, DM, TYPDIM, XM, TYPEXM, DN, TYPDIN, XN, TYPEXN, D, C, IKLE, NELEM, NELMAX, NDIAG)
Definition: om2121.f:8
declarations_special::lu
integer lu
Definition: declarations_special.F:45
bief
Definition: bief.f:3
Generated on Tue Jan 5 2021 14:37:20 for The TELEMAC-MASCARET system by doxygen 1.8.13