The TELEMAC-MASCARET system  trunk
mvseg.f
Go to the documentation of this file.
1 ! ****************
2  SUBROUTINE mvseg
3 ! ****************
4 !
5  &(op, x , da,typdia,xa1,xa2,
6  & typext, y,c,npoin,nelem,nseg1,nseg2,gloseg1,gloseg2,ielm1,ielm2)
7 !
8 !***********************************************************************
9 ! BIEF V6P1 21/08/2010
10 !***********************************************************************
11 !
12 !brief MATRIX VECTOR PRODUCT FOR EDGE-BASED STORAGE.
13 !code
14 !+ OP IS A STRING OF 8 CHARACTERS, WHICH INDICATES THE OPERATION TO BE
15 !+ PERFORMED ON VECTORS X,Y AND MATRIX M.
16 !+
17 !+ THE RESULT IS VECTOR X.
18 !+
19 !+ THESE OPERATIONS ARE DIFFERENT DEPENDING ON THE DIAGONAL TYPE
20 !+ AND THE TYPE OF EXTRADIAGONAL TERMS.
21 !+
22 !+ IMPLEMENTED OPERATIONS:
23 !+
24 !+ OP = 'X=AY ' : X = AY
25 !+ OP = 'X=CAY ' : X = CAY
26 !+ OP = 'X=-AY ' : X = -AY
27 !+ OP = 'X=X+AY ' : X = X + AY
28 !+ OP = 'X=X-AY ' : X = X - AY
29 !+ OP = 'X=X+CAY ' : X = X + C AY
30 !+ OP = 'X=TAY ' : X = TA Y (TRANSPOSE OF A)
31 !+ OP = 'X=-TAY ' : X = - TA Y (- TRANSPOSE OF A)
32 !+ OP = 'X=X+TAY ' : X = X + TA Y
33 !+ OP = 'X=X-TAY ' : X = X - TA Y
34 !+ OP = 'X=X+CTAY' : X = X + C TA Y
35 !
36 !history J-M HERVOUET (LNH) ; F LEPEINTRE (LNH)
37 !+ 05/02/91
38 !+ V5P1
39 !+
40 !
41 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
42 !+ 13/07/2010
43 !+ V6P0
44 !+ Translation of French comments within the FORTRAN sources into
45 !+ English comments
46 !
47 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
48 !+ 21/08/2010
49 !+ V6P0
50 !+ Creation of DOXYGEN tags for automated documentation and
51 !+ cross-referencing of the FORTRAN sources
52 !
53 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
54 !| C |-->| A GIVEN CONSTANT
55 !| DA |-->| MATRIX DIAGONAL
56 !| GLOSEG1 |-->| FIRST POINT OF SEGMENTS
57 !| GLOSEG2 |-->| SECOND POINT OF SEGMENTS
58 !| IELM1 |-->| TYPE OF LINE ELEMENT.
59 !| IELM2 |-->| TYPE OF COLUMN ELEMENT.
60 !| NELEM |-->| NUMBER OF ELEMENTS
61 !| NPOIN |-->| NUMBER OF LINEAR POINTS
62 !| NSEG1 |-->| NUMBER OF SEGMENTS OF THE LINE ELEMENT
63 !| NSEG2 |-->| NUMBER OF SEGMENTS OF THE COLUMN ELEMENT
64 !| OP |-->| OPERATION TO BE DONE (SEE ABOVE)
65 !| TYPDIA |-->| TYPE OF DIAGONAL:
66 !| | | TYPDIA = 'Q' : ANY VALUE
67 !| | | TYPDIA = 'I' : IDENTITY
68 !| | | TYPDIA = '0' : ZERO
69 !| TYPEXT |-->| TYPE OF OFF-DIAGONAL TERMS
70 !| | | TYPEXT = 'Q' : ANY VALUE
71 !| | | TYPEXT = 'S' : SYMMETRIC
72 !| | | TYPEXT = '0' : ZERO
73 !| X |<->| RESULT IN ASSEMBLED FORM
74 !| XA1 |-->| OFF-DIAGONAL TERM OF MATRIX
75 !| XA2 |-->| OFF-DIAGONAL TERM OF MATRIX
76 !| Y |-->| VECTOR USED IN THE OPERATION
77 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
78 !
79  USE bief,ex_mvseg => mvseg
80 !
82  IMPLICIT NONE
83 !
84 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
85 !
86  INTEGER, INTENT(IN) :: NPOIN
87  INTEGER, INTENT(IN) :: GLOSEG1(*),GLOSEG2(*)
88  INTEGER, INTENT(IN) :: NSEG1,NSEG2,NELEM,IELM1,IELM2
89 !
90  DOUBLE PRECISION, INTENT(IN) :: Y(*),DA(*)
91  DOUBLE PRECISION, INTENT(INOUT) :: X(*)
92  DOUBLE PRECISION, INTENT(IN) :: XA1(*),XA2(*)
93  DOUBLE PRECISION, INTENT(IN) :: C
94 !
95  CHARACTER(LEN=8),INTENT(IN) :: OP
96  CHARACTER(LEN=1),INTENT(IN) :: TYPDIA,TYPEXT
97 !
98 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
99 !
100  INTEGER ISEG,I,MINSEG,MAXSEG
101  DOUBLE PRECISION Z(1)
102 !
103  INTRINSIC min,max
104 !
105 !-----------------------------------------------------------------------
106 !
107  minseg = min(nseg1,nseg2)
108  maxseg = max(nseg1,nseg2)
109 !
110  IF(op(1:8).EQ.'X=AY ') THEN
111 !
112 ! CONTRIBUTION OF THE DIAGONAL:
113 !
114  IF(typdia(1:1).EQ.'Q') THEN
115  CALL ov ('X=YZ ', x , y , da , c , npoin )
116  ELSEIF(typdia(1:1).EQ.'I') THEN
117  CALL ov ('X=Y ', x , y , z , c , npoin )
118  ELSEIF(typdia(1:1).EQ.'0') THEN
119  CALL ov ('X=C ', x , y , da , 0.d0 , npoin )
120  ELSE
121  WRITE(lu,2001) typdia
122  CALL plante(1)
123  stop
124  ENDIF
125 !
126 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
127 !
128  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
129 !
130 ! SQUARE PART
131 !
132  DO iseg = 1 , minseg
133  x(gloseg1(iseg))=
134  & x(gloseg1(iseg))+xa1(iseg)*y(gloseg2(iseg))
135  x(gloseg2(iseg))=
136  & x(gloseg2(iseg))+xa2(iseg)*y(gloseg1(iseg))
137  ENDDO
138 !
139 ! THE REST OF THE RECTANGULAR MATRIX
140 !
141  IF(nseg1.GT.nseg2) THEN
142 ! PART OF X HAS NOT BEEN INITIALISED
143 ! BY THE CONTRIBUTION OF THE DIAGONAL
144  IF(ielm1.EQ.12.OR.ielm2.EQ.12) THEN
145 ! OPTIMISATION FOR QUASI-BUBBLE ELEMENT
146  DO i = npoin+1,npoin+nelem
147  x(i)=0.d0
148  ENDDO
149  ELSE
150  DO iseg = minseg+1,maxseg
151  x(gloseg2(iseg))=0.d0
152  ENDDO
153  ENDIF
154  DO iseg = minseg+1,maxseg
155  x(gloseg2(iseg))=
156  & x(gloseg2(iseg))+xa2(iseg)*y(gloseg1(iseg))
157  ENDDO
158  ELSEIF(nseg2.GT.nseg1) THEN
159  DO iseg = minseg+1,maxseg
160  x(gloseg1(iseg))=
161  & x(gloseg1(iseg))+xa2(iseg)*y(gloseg2(iseg))
162  ENDDO
163  ENDIF
164 !
165  ELSEIF(typext(1:1).NE.'0') THEN
166 !
167  WRITE(lu,1001) typext
168  CALL plante(1)
169  stop
170 !
171  ENDIF
172 !
173 !-----------------------------------------------------------------------
174 !
175  ELSEIF(op(1:8).EQ.'X=CAY ') THEN
176 !
177 ! CONTRIBUTION OF THE DIAGONAL:
178 !
179  IF(typdia(1:1).EQ.'Q') THEN
180  CALL ov ('X=CYZ ', x , y , da , c , npoin )
181  ELSEIF(typdia(1:1).EQ.'I') THEN
182  CALL ov ('X=CY ', x , y , z , c , npoin )
183  ELSEIF(typdia(1:1).EQ.'0') THEN
184  CALL ov ('X=C ', x , y , da , 0.d0 , npoin )
185  ELSE
186  WRITE(lu,2001) typdia
187  CALL plante(1)
188  stop
189  ENDIF
190 !
191 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
192 !
193  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
194 !
195  DO iseg = 1 , minseg
196  x(gloseg1(iseg))=
197  & x(gloseg1(iseg))+c*xa1(iseg)*y(gloseg2(iseg))
198  x(gloseg2(iseg))=
199  & x(gloseg2(iseg))+c*xa2(iseg)*y(gloseg1(iseg))
200  ENDDO
201  IF(nseg1.GT.nseg2) THEN
202  IF(ielm1.EQ.12.OR.ielm2.EQ.12) THEN
203 ! OPTIMISATION FOR QUASI-BUBBLE ELEMENT
204  DO i = npoin+1,npoin+nelem
205  x(i)=0.d0
206  ENDDO
207  ELSE
208  DO iseg = minseg+1,maxseg
209  x(gloseg2(iseg))=0.d0
210  ENDDO
211  ENDIF
212  DO iseg = minseg+1,maxseg
213  x(gloseg2(iseg))=
214  & x(gloseg2(iseg))+c*xa2(iseg)*y(gloseg1(iseg))
215  ENDDO
216  ELSEIF(nseg2.GT.nseg1) THEN
217  DO iseg = minseg+1,maxseg
218  x(gloseg1(iseg))=
219  & x(gloseg1(iseg))+c*xa2(iseg)*y(gloseg2(iseg))
220  ENDDO
221  ENDIF
222 !
223  ELSEIF(typext(1:1).NE.'0') THEN
224 !
225  WRITE(lu,1001) typext
226  CALL plante(1)
227  stop
228 !
229  ENDIF
230 !
231 !-----------------------------------------------------------------------
232 !
233  ELSEIF(op(1:8).EQ.'X=-AY ') THEN
234 !
235 ! CONTRIBUTION OF THE DIAGONAL:
236 !
237  IF(typdia(1:1).EQ.'Q') THEN
238  CALL ov ('X=-YZ ', x , y , da , c , npoin )
239  ELSEIF(typdia(1:1).EQ.'I') THEN
240  CALL ov ('X=-Y ', x , y , z , c , npoin )
241  ELSEIF(typdia(1:1).EQ.'0') THEN
242  CALL ov ('X=C ', x , y , da , 0.d0 , npoin )
243  ELSE
244  WRITE(lu,2001) typdia
245  CALL plante(1)
246  stop
247  ENDIF
248 !
249 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
250 !
251  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
252 !
253 !
254  DO iseg = 1 , minseg
255  x(gloseg1(iseg))=
256  & x(gloseg1(iseg))-xa1(iseg)*y(gloseg2(iseg))
257  x(gloseg2(iseg))=
258  & x(gloseg2(iseg))-xa2(iseg)*y(gloseg1(iseg))
259  ENDDO
260  IF(nseg1.GT.nseg2) THEN
261  IF(ielm1.EQ.12.OR.ielm2.EQ.12) THEN
262 ! OPTIMISATION FOR QUASI-BUBBLE ELEMENT
263  DO i = npoin+1,npoin+nelem
264  x(i)=0.d0
265  ENDDO
266  ELSE
267  DO iseg = minseg+1,maxseg
268  x(gloseg2(iseg))=0.d0
269  ENDDO
270  ENDIF
271  DO iseg = minseg+1,maxseg
272  x(gloseg2(iseg))=
273  & x(gloseg2(iseg))-xa2(iseg)*y(gloseg1(iseg))
274  ENDDO
275  ELSEIF(nseg2.GT.nseg1) THEN
276  DO iseg = minseg+1,maxseg
277  x(gloseg1(iseg))=
278  & x(gloseg1(iseg))-xa2(iseg)*y(gloseg2(iseg))
279  ENDDO
280  ENDIF
281 !
282  ELSEIF(typext(1:1).NE.'0') THEN
283 !
284  WRITE(lu,1001) typext
285  CALL plante(1)
286  stop
287 !
288  ENDIF
289 !
290 !-----------------------------------------------------------------------
291 !
292  ELSEIF(op(1:8).EQ.'X=X+AY ') THEN
293 !
294 ! CONTRIBUTION OF THE DIAGONAL:
295 !
296  IF(typdia(1:1).EQ.'Q') THEN
297  CALL ov ('X=X+YZ ', x , y , da , c , npoin )
298  ELSEIF(typdia(1:1).EQ.'I') THEN
299  CALL ov ('X=X+Y ', x , y , z , c , npoin )
300  ELSEIF(typdia(1:1).NE.'0') THEN
301  WRITE(lu,2001) typdia
302  CALL plante(1)
303  stop
304  ENDIF
305 !
306 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
307 !
308  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
309 !
310  DO iseg = 1 , minseg
311  x(gloseg1(iseg))=
312  & x(gloseg1(iseg))+xa1(iseg)*y(gloseg2(iseg))
313  x(gloseg2(iseg))=
314  & x(gloseg2(iseg))+xa2(iseg)*y(gloseg1(iseg))
315  ENDDO
316  IF(nseg1.GT.nseg2) THEN
317  DO iseg = minseg+1,maxseg
318  x(gloseg2(iseg))=
319  & x(gloseg2(iseg))+xa2(iseg)*y(gloseg1(iseg))
320  ENDDO
321  ELSEIF(nseg2.GT.nseg1) THEN
322  DO iseg = minseg+1,maxseg
323  x(gloseg1(iseg))=
324  & x(gloseg1(iseg))+xa2(iseg)*y(gloseg2(iseg))
325  ENDDO
326  ENDIF
327 !
328  ELSEIF(typext(1:1).NE.'0') THEN
329 !
330  WRITE(lu,1001) typext
331  CALL plante(1)
332  stop
333 !
334  ENDIF
335 !
336 !-----------------------------------------------------------------------
337 !
338  ELSEIF(op(1:8).EQ.'X=X-AY ') THEN
339 !
340 ! CONTRIBUTION OF THE DIAGONAL:
341 !
342  IF(typdia(1:1).EQ.'Q') THEN
343  CALL ov ('X=X-YZ ', x , y , da , c , npoin )
344  ELSEIF(typdia(1:1).EQ.'I') THEN
345  CALL ov ('X=X-Y ', x , y , z , c , npoin )
346  ELSEIF(typdia(1:1).NE.'0') THEN
347  WRITE(lu,2001) typdia
348  CALL plante(1)
349  stop
350  ENDIF
351 !
352 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
353 !
354  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
355 !
356  DO iseg = 1 , minseg
357  x(gloseg1(iseg))=
358  & x(gloseg1(iseg))-xa1(iseg)*y(gloseg2(iseg))
359  x(gloseg2(iseg))=
360  & x(gloseg2(iseg))-xa2(iseg)*y(gloseg1(iseg))
361  ENDDO
362  IF(nseg1.GT.nseg2) THEN
363  DO iseg = minseg+1,maxseg
364  x(gloseg2(iseg))=
365  & x(gloseg2(iseg))-xa2(iseg)*y(gloseg1(iseg))
366  ENDDO
367  ELSEIF(nseg2.GT.nseg1) THEN
368  DO iseg = minseg+1,maxseg
369  x(gloseg1(iseg))=
370  & x(gloseg1(iseg))-xa2(iseg)*y(gloseg2(iseg))
371  ENDDO
372  ENDIF
373 !
374  ELSEIF(typext(1:1).NE.'0') THEN
375 !
376  WRITE(lu,1001) typext
377  CALL plante(1)
378  stop
379 !
380  ENDIF
381 !
382 !-----------------------------------------------------------------------
383 !
384  ELSEIF(op(1:8).EQ.'X=X+CAY ') THEN
385 !
386 ! CONTRIBUTION OF THE DIAGONAL:
387 !
388  IF(typdia(1:1).EQ.'Q') THEN
389  CALL ov ('X=X+CYZ ', x , y , da , c , npoin )
390  ELSEIF(typdia(1:1).EQ.'I') THEN
391  CALL ov ('X=X+CY ', x , y , z , c , npoin )
392  ELSEIF(typdia(1:1).NE.'0') THEN
393  WRITE(lu,2001) typdia
394  CALL plante(1)
395  stop
396  ENDIF
397 !
398 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
399 !
400  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
401 !
402  DO iseg = 1 , minseg
403  x(gloseg1(iseg))=
404  & x(gloseg1(iseg))+c*xa1(iseg)*y(gloseg2(iseg))
405  x(gloseg2(iseg))=
406  & x(gloseg2(iseg))+c*xa2(iseg)*y(gloseg1(iseg))
407  ENDDO
408  IF(nseg1.GT.nseg2) THEN
409  DO iseg = minseg+1,maxseg
410  x(gloseg2(iseg))=
411  & x(gloseg2(iseg))+c*xa2(iseg)*y(gloseg1(iseg))
412  ENDDO
413  ELSEIF(nseg2.GT.nseg1) THEN
414  DO iseg = minseg+1,maxseg
415  x(gloseg1(iseg))=
416  & x(gloseg1(iseg))+c*xa2(iseg)*y(gloseg2(iseg))
417  ENDDO
418  ENDIF
419 !
420  ELSEIF(typext(1:1).NE.'0') THEN
421 !
422  WRITE(lu,1001) typext
423  CALL plante(1)
424  stop
425 !
426  ENDIF
427 !
428 !-----------------------------------------------------------------------
429 !
430  ELSEIF(op(1:8).EQ.'X=TAY ') THEN
431 !
432 ! CONTRIBUTION OF THE DIAGONAL:
433 !
434  IF(typdia(1:1).EQ.'Q') THEN
435  CALL ov ('X=YZ ', x , y , da , c , npoin )
436  ELSEIF(typdia(1:1).EQ.'I') THEN
437  CALL ov ('X=Y ', x , y , z , c , npoin )
438  ELSEIF(typdia(1:1).EQ.'0') THEN
439  CALL ov ('X=C ', x , y , da , 0.d0 , npoin )
440  ELSE
441  WRITE(lu,2001) typdia
442  CALL plante(1)
443  stop
444  ENDIF
445 !
446 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
447 !
448  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
449 !
450  DO iseg = 1 , minseg
451  x(gloseg1(iseg))=
452  & x(gloseg1(iseg))+xa2(iseg)*y(gloseg2(iseg))
453  x(gloseg2(iseg))=
454  & x(gloseg2(iseg))+xa1(iseg)*y(gloseg1(iseg))
455  ENDDO
456  IF(nseg1.GT.nseg2) THEN
457  DO iseg = minseg+1,maxseg
458  x(gloseg1(iseg))=
459  & x(gloseg1(iseg))+xa2(iseg)*y(gloseg2(iseg))
460  ENDDO
461  ELSEIF(nseg2.GT.nseg1) THEN
462  IF(ielm1.EQ.12.OR.ielm2.EQ.12) THEN
463 ! OPTIMISATION FOR QUASI-BUBBLE ELEMENT
464  DO i = npoin+1,npoin+nelem
465  x(i)=0.d0
466  ENDDO
467  ELSE
468  DO iseg = minseg+1,maxseg
469  x(gloseg2(iseg))=0.d0
470  ENDDO
471  ENDIF
472  DO iseg = minseg+1,maxseg
473  x(gloseg2(iseg))=
474  & x(gloseg2(iseg))+xa2(iseg)*y(gloseg1(iseg))
475  ENDDO
476  ENDIF
477 !
478  ELSEIF(typext(1:1).NE.'0') THEN
479 !
480  WRITE(lu,1001) typext
481  CALL plante(1)
482  stop
483 !
484  ENDIF
485 !
486 !-----------------------------------------------------------------------
487 !
488  ELSEIF(op(1:8).EQ.'X=-TAY ') THEN
489 !
490 ! CONTRIBUTION OF THE DIAGONAL
491 !
492  IF(typdia(1:1).EQ.'Q') THEN
493  CALL ov ('X=-YZ ', x , y , da , c , npoin )
494  ELSEIF(typdia(1:1).EQ.'I') THEN
495  CALL ov ('X=-Y ', x , y , z , c , npoin )
496  ELSEIF(typdia(1:1).EQ.'0') THEN
497  CALL ov ('X=C ', x , y , da , 0.d0 , npoin )
498  ELSE
499  WRITE(lu,2001) typdia
500  CALL plante(1)
501  stop
502  ENDIF
503 !
504 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
505 !
506  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
507 !
508  DO iseg = 1 , minseg
509  x(gloseg1(iseg))=
510  & x(gloseg1(iseg))-xa2(iseg)*y(gloseg2(iseg))
511  x(gloseg2(iseg))=
512  & x(gloseg2(iseg))-xa1(iseg)*y(gloseg1(iseg))
513  ENDDO
514  IF(nseg1.GT.nseg2) THEN
515  DO iseg = minseg+1,maxseg
516  x(gloseg1(iseg))=
517  & x(gloseg1(iseg))-xa2(iseg)*y(gloseg2(iseg))
518  ENDDO
519  ELSEIF(nseg2.GT.nseg1) THEN
520  IF(ielm1.EQ.12.OR.ielm2.EQ.12) THEN
521 ! OPTIMISATION FOR QUASI-BUBBLE ELEMENT
522  DO i = npoin+1,npoin+nelem
523  x(i)=0.d0
524  ENDDO
525  ELSE
526  DO iseg = minseg+1,maxseg
527  x(gloseg2(iseg))=0.d0
528  ENDDO
529  ENDIF
530  DO iseg = minseg+1,maxseg
531  x(gloseg2(iseg))=
532  & x(gloseg2(iseg))-xa2(iseg)*y(gloseg1(iseg))
533  ENDDO
534  ENDIF
535 !
536  ELSEIF(typext(1:1).NE.'0') THEN
537 !
538  WRITE(lu,1001) typext
539  CALL plante(1)
540  stop
541 !
542  ENDIF
543 !
544 !-----------------------------------------------------------------------
545 !
546  ELSEIF(op(1:8).EQ.'X=X+TAY ') THEN
547 !
548 ! CONTRIBUTION OF THE DIAGONAL
549 !
550  IF(typdia(1:1).EQ.'Q') THEN
551  CALL ov ('X=X+YZ ', x , y , da , c , npoin )
552  ELSEIF(typdia(1:1).EQ.'I') THEN
553  CALL ov ('X=X+Y ', x , y , z , c , npoin )
554  ELSEIF(typdia(1:1).NE.'0') THEN
555  WRITE(lu,2001) typdia
556  CALL plante(1)
557  stop
558  ENDIF
559 !
560 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
561 !
562  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
563 !
564  DO iseg = 1 , minseg
565  x(gloseg1(iseg))=
566  & x(gloseg1(iseg))+xa2(iseg)*y(gloseg2(iseg))
567  x(gloseg2(iseg))=
568  & x(gloseg2(iseg))+xa1(iseg)*y(gloseg1(iseg))
569  ENDDO
570  IF(nseg1.GT.nseg2) THEN
571  DO iseg = minseg+1,maxseg
572  x(gloseg1(iseg))=
573  & x(gloseg1(iseg))+xa2(iseg)*y(gloseg2(iseg))
574  ENDDO
575  ELSEIF(nseg2.GT.nseg1) THEN
576  DO iseg = minseg+1,maxseg
577  x(gloseg2(iseg))=
578  & x(gloseg2(iseg))+xa2(iseg)*y(gloseg1(iseg))
579  ENDDO
580  ENDIF
581 !
582  ELSEIF(typext(1:1).NE.'0') THEN
583 !
584  WRITE(lu,1001) typext
585  CALL plante(1)
586  stop
587 !
588  ENDIF
589 !
590 !-----------------------------------------------------------------------
591 !
592  ELSEIF(op(1:8).EQ.'X=X-TAY ') THEN
593 !
594 ! CONTRIBUTION OF THE DIAGONAL
595 !
596  IF(typdia(1:1).EQ.'Q') THEN
597  CALL ov ('X=X-YZ ', x , y , da , c , npoin )
598  ELSEIF(typdia(1:1).EQ.'I') THEN
599  CALL ov ('X=X-Y ', x , y , z , c , npoin )
600  ELSEIF(typdia(1:1).NE.'0') THEN
601  WRITE(lu,2001) typdia
602  CALL plante(1)
603  stop
604  ENDIF
605 !
606 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
607 !
608  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
609 !
610  DO iseg = 1 , minseg
611  x(gloseg1(iseg))=
612  & x(gloseg1(iseg))-xa2(iseg)*y(gloseg2(iseg))
613  x(gloseg2(iseg))=
614  & x(gloseg2(iseg))-xa1(iseg)*y(gloseg1(iseg))
615  ENDDO
616  IF(nseg1.GT.nseg2) THEN
617  DO iseg = minseg+1,maxseg
618  x(gloseg1(iseg))=
619  & x(gloseg1(iseg))-xa2(iseg)*y(gloseg2(iseg))
620  ENDDO
621  ELSEIF(nseg2.GT.nseg1) THEN
622  DO iseg = minseg+1,maxseg
623  x(gloseg2(iseg))=
624  & x(gloseg2(iseg))-xa2(iseg)*y(gloseg1(iseg))
625  ENDDO
626  ENDIF
627 !
628  ELSEIF(typext(1:1).NE.'0') THEN
629 !
630  WRITE(lu,1001) typext
631  CALL plante(1)
632  stop
633 !
634  ENDIF
635 !
636 !-----------------------------------------------------------------------
637 !
638  ELSEIF(op(1:8).EQ.'X=X+CTAY') THEN
639 !
640 ! CONTRIBUTION OF THE DIAGONAL
641 !
642  IF(typdia(1:1).EQ.'Q') THEN
643  CALL ov ('X=X+CYZ ', x , y , da , c , npoin )
644  ELSEIF(typdia(1:1).EQ.'I') THEN
645  CALL ov ('X=X+CY ', x , y , z , c , npoin )
646  ELSEIF(typdia(1:1).NE.'0') THEN
647  WRITE(lu,2001) typdia
648  CALL plante(1)
649  stop
650  ENDIF
651 !
652 ! CONTRIBUTION OF EXTRADIAGONAL TERMS:
653 !
654  IF(typext(1:1).EQ.'Q'.OR.typext(1:1).EQ.'S') THEN
655 !
656  DO iseg = 1 , minseg
657  x(gloseg1(iseg))=
658  & x(gloseg1(iseg))+c*xa2(iseg)*y(gloseg2(iseg))
659  x(gloseg2(iseg))=
660  & x(gloseg2(iseg))+c*xa1(iseg)*y(gloseg1(iseg))
661  ENDDO
662  IF(nseg1.GT.nseg2) THEN
663  DO iseg = minseg+1,maxseg
664  x(gloseg1(iseg))=
665  & x(gloseg1(iseg))+c*xa2(iseg)*y(gloseg2(iseg))
666  ENDDO
667  ELSEIF(nseg2.GT.nseg1) THEN
668  DO iseg = minseg+1,maxseg
669  x(gloseg2(iseg))=
670  & x(gloseg2(iseg))+c*xa2(iseg)*y(gloseg1(iseg))
671  ENDDO
672  ENDIF
673 !
674  ELSEIF(typext(1:1).NE.'0') THEN
675 !
676  WRITE(lu,1001) typext
677  CALL plante(1)
678  stop
679 !
680  ENDIF
681 !
682 !-----------------------------------------------------------------------
683 !
684  ELSE
685 !
686  WRITE(lu,3001) op
687  CALL plante(1)
688  stop
689 !
690 !-----------------------------------------------------------------------
691 !
692  ENDIF
693 !
694 !-----------------------------------------------------------------------
695 !
696  RETURN
697 !
698 1001 FORMAT(1x,'MVSEG (BIEF) : EXTRADIAG. TERMS UNKNOWN TYPE : ',a1)
699 2001 FORMAT(1x,'MVSEG (BIEF) : DIAGONAL : UNKNOWN TYPE : ',a1)
700 3001 FORMAT(1x,'MVSEG (BIEF) : UNKNOWN OPERATION : ',a8)
701 !
702  END
subroutine ov(OP, X, Y, Z, C, DIM1)
Definition: ov.f:7
subroutine mvseg(OP, X, DA, TYPDIA, XA1, XA2, TYPEXT, Y, C, NPOIN, NELEM, NSEG1, NSEG2, GLOSEG1, GLOSEG2, IELM1, IELM2)
Definition: mvseg.f:8
Definition: bief.f:3