The TELEMAC-MASCARET system  trunk
matvct.f
Go to the documentation of this file.
1 ! *****************
2  SUBROUTINE matvct
3 ! *****************
4 !
5  &(op, x , da,typdia,xa,typext, y ,
6  & c,ikle,npt,nelem,nelmax,w,lego,ielm1,ielm2,ielmx,lv,
7  & s,p,iklem1,dimikm,limvoi,mxptvs,npmax,npoin,
8  & gloseg,sizglo,sizxa,ndp,mesh,stox
9  & ,x_err,y_err,da_err)
10 !
11 !***********************************************************************
12 ! BIEF V7P2
13 !***********************************************************************
14 !
15 !brief MATRIX VECTOR OPERATIONS.
16 !code
17 !+ OP IS A STRING OF 8 CHARACTERS, WHICH INDICATES THE OPERATION TO BE
18 !+ PERFORMED ON VECTORS X,Y AND MATRIX M.
19 !+
20 !+ THE RESULT IS VECTOR X (A NON-ASSEMBLED PART OF WHICH CAN BE IN
21 !+ ARRAY W IF LEGO = .FALSE.)
22 !+
23 !+ THESE OPERATIONS ARE DIFFERENTS DEPENDING ON THE DIAGONAL TYPE
24 !+ AND THE OFF-DIAGONAL TERMS TYPE.
25 !+
26 !+ IMPLEMENTED OPERATIONS :
27 !+
28 !+ OP = 'X=AY ' : X = AY
29 !+ OP = 'X=X+AY ' : X = X + AY
30 !+ OP = 'X=X-AY ' : X = X - AY
31 !+ OP = 'X=X+CAY ' : X = X + C AY
32 !+ OP = 'X=TAY ' : X = TA Y (TA: TRANSPOSE OF A)
33 !+ OP = 'X=X+TAY ' : X = X + TA Y
34 !+ OP = 'X=X-TAY ' : X = X - TA Y
35 !+ OP = 'X=X+CTAY' : X = X + C TA Y
36 !
37 !history ALGIANE FROEHLY (MATMECA)
38 !+ 31/03/2008
39 !+
40 !+ QUADRATIC TRIANGLE
41 !
42 !history J-M HERVOUET (LNHE)
43 !+ 05/02/2010
44 !+ V6P0
45 !+ NSEG1 AND NSEG2 MODIFIED BEFORE CALLING MVSEG
46 !
47 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
48 !+ 13/07/2010
49 !+ V6P0
50 !+ Translation of French comments within the FORTRAN sources into
51 !+ English comments
52 !
53 !history N.DURAND (HRW), S.E.BOURBAN (HRW)
54 !+ 21/08/2010
55 !+ V6P0
56 !+ Creation of DOXYGEN tags for automated documentation and
57 !+ cross-referencing of the FORTRAN sources
58 !history R.NHEILI (Univerte de Perpignan, DALI)
59 !+ 24/02/2016
60 !+ V7P3
61 !+ ADD MODASS=3
62 !
63 !history S.E.BOURBAN (HRW)
64 !+ 21/03/2017
65 !+ V7P3
66 !+ Replacement of the DATA declarations by the PARAMETER associates
67 !
68 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
69 !| C |-->| A GIVEN CONSTANT
70 !| DA |-->| MATRIX DIAGONAL
71 !| DIMIKM |-->| FIRST DIMENSION OF IKLEM1.
72 !| GLOSEG |-->| FIRST AND SECOND POINT OF SEGMENTS
73 !| IELM1 |-->| TYPE OF ELEMENT FOR LINES
74 !| IELM2 |-->| TYPE OF ELEMENT FOR COLUMNS
75 !| IELMX |-->| TYPE OF ELEMENT OF RESULT
76 !| | | CAN BE IELM1 OR IELM2 DEPENDING ON OP
77 !| IKLE |-->| CONNECTIVITY TABLE.
78 !| IKLEM1 |-->| CONNECTIVITY TABLE USED FOR MATRIX-VECTOR 2
79 !| LEGO |-->| = .TRUE. W1,2,... ARE ASSEMBLED ON X
80 !| | | =.FALSE. W1,2,... ARE NOT ASSEMBLED
81 !| LIMVOI |-->| ARRAY USED FOR MATRIX-VECTOR 2
82 !| LV |-->| VECTOR LENGTH OF THE MACHINE
83 !| MXPTVS |-->| MAXIMUM NUMBER OF NEIGHBOURS OF A POINT
84 !| NELEM |-->| NUMBER OF ELEMENTS
85 !| NELMAX |-->| MAXIMUM NUMBER OF ELEMENTS
86 !| NPMAX |-->| MAXIMUM NUMBER OF POINTS IN THE MESH
87 !| NPOIN |-->| NUMBER OF POINTS
88 !| NPT |-->| DIMENSION OF DIAGONAL
89 !| OP |-->| OPERATION TO BE DONE
90 !| P |-->| TYPE OF MATRIX X VECTOR PRODUCT.
91 !| S |-->| TYPE OF STORAGE.
92 !| SIZGLO |-->| FIRST DIMENSION OF GLOSEG
93 !| SIZXA |-->| FIRST DIMENSION OF ARRAY XA
94 !| STOX |-->| STORAGE OPTION OF OFF-DIAGONAL TERMS
95 !| | | 1=(NELMAX,NDP) 2=(NDP,NELMAX)
96 !| TYPDIA |-->| TYPE OF DIAGONAL:
97 !| | | TYPDIA = 'Q' : ANY VALUE
98 !| | | TYPDIA = 'I' : IDENTITY
99 !| | | TYPDIA = '0' : ZERO
100 !| TYPEXT |-->| TYPE OF OFF-DIAGONAL TERMS
101 !| | | TYPEXT = 'Q' : ANY VALUE
102 !| | | TYPEXT = 'S' : SYMMETRIC
103 !| | | TYPEXT = '0' : ZERO
104 !| W |<--| WORK ARRAY WITH NON ASSEMBLED RESULT
105 !| X |<--| RESULTING VECTOR
106 !| XA |-->| OFF-DIAGONAL TERMS IN THE MATRIX A
107 !| Y |-->| A GIVEN VECTOR USED IN OPERATION OP
108 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
109 !
110  USE bief, ex_matvct => matvct
111  USE declarations_telemac, ONLY : modass
112 !
114  IMPLICIT NONE
115 !
116 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
117 !
118  INTEGER, INTENT(IN) :: IELM1,IELM2,IELMX,NPOIN,NPMAX,S,P,SIZXA
119  INTEGER, INTENT(IN) :: NDP,STOX
120  INTEGER, INTENT(INOUT) :: NPT
121  INTEGER, INTENT(IN) :: NELEM,NELMAX,LV,DIMIKM,MXPTVS,SIZGLO
122  INTEGER, INTENT(IN) :: IKLE(nelmax,*),IKLEM1(*),LIMVOI(*)
123  INTEGER, INTENT(IN) :: GLOSEG(sizglo,2)
124  CHARACTER(LEN=8), INTENT(IN) :: OP
125  CHARACTER(LEN=1),INTENT(IN) :: TYPDIA,TYPEXT
126  DOUBLE PRECISION, INTENT(INOUT) :: X(*)
127  DOUBLE PRECISION, INTENT(IN) :: Y(*),DA(*),XA(sizxa,*),C
128  DOUBLE PRECISION, INTENT(INOUT) :: W(nelmax,*)
129  LOGICAL, INTENT(IN) :: LEGO
130  TYPE(bief_mesh), INTENT(INOUT) :: MESH
131  DOUBLE PRECISION, OPTIONAL, INTENT(INOUT) :: X_ERR(*)
132  DOUBLE PRECISION, OPTIONAL, INTENT(IN) :: Y_ERR(*),DA_ERR(*)
133 !
134 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
135 !
136  INTEGER NSEG1,NSEG2,SYM,NPT2,DIM1XA
137 !
138  DOUBLE PRECISION Z(1)
139 !
140  INTRINSIC min
141 !
142 ! THESE DATA ALSO APPEAR IN MATVCT
143 !
144 ! DATA OOS/ 0 , 1 ,
145 ! * 1 , 0 ,
146 ! S=2 NOT IMPLEMENTED
147 ! * 0 , 0 ,
148 ! * 0 , 0 /
149 !
150  INTEGER :: AAS(3,3,2)
151  parameter( aas = reshape( (/
152 ! SYMMETRICAL P1-P1 EBE (S=1)
153  & 0 , 1 , 2 ,
154  & 1 , 0 , 3 ,
155  & 2 , 3 , 0 ,
156 ! SYMMETRICAL P1-P1 PRE-ASSEMBLED EBE (S=2)
157  & 0 , 1 , 3 ,
158  & 1 , 0 , 2 ,
159  & 3 , 2 , 0 /), shape=(/ 3,3,2 /) ) )
160 !
161  INTEGER :: AAQ(3,3,2)
162  parameter( aaq = reshape( (/
163 ! NONSYMMETRICAL P1-P1 EBE (S=1)
164  & 0 , 4 , 5 ,
165  & 1 , 0 , 6 ,
166  & 2 , 3 , 0 ,
167 ! NONSYMMETRICAL P1-P1 PRE-ASSEMBLED EBE (S=2)
168  & 0 , 4 , 3 ,
169  & 1 , 0 , 5 ,
170  & 6 , 2 , 0 /), shape=(/ 3,3,2 /) ) )
171 !
172  INTEGER :: BBS(4,4,2)
173  parameter( bbs = reshape( (/
174 ! SYMMETRICAL QUASI-BUBBLE QUASI-BUBBLE EBE (S=1)
175  & 0 , 1 , 2 , 3 ,
176  & 1 , 0 , 4 , 5 ,
177  & 2 , 4 , 0 , 6 ,
178  & 3 , 5 , 6 , 0 ,
179 ! SYMMETRICAL QUASI-BUBBLE QUASI-BUBBLE PRE-ASSEMBLED EBE (S=2)
180  & 0 , 4 , 6 , 1 ,
181  & 4 , 0 , 5 , 2 ,
182  & 6 , 5 , 0 , 3 ,
183  & 1 , 2 , 3 , 0 /), shape=(/ 4,4,2 /) ) )
184 !
185  INTEGER :: BBQ(4,4,2)
186  parameter( bbq = reshape( (/
187 ! NONSYMMETRICAL QUASI-BUBBLE QUASI-BUBBLE EBE (S=1)
188  & 0 , 7 , 8 , 9 ,
189  & 1 , 0 , 10 , 11 ,
190  & 2 , 4 , 0 , 12 ,
191  & 3 , 5 , 6 , 0 ,
192 ! NONSYMMETRICAL QUASI-BUBBLE QUASI-BUBBLE PRE-ASSEMBLED EBE (S=2)
193  & 0 , 10 , 6 , 7 ,
194  & 4 , 0 , 11 , 8 ,
195  & 12 , 5 , 0 , 9 ,
196  & 1 , 2 , 3 , 0 /), shape=(/ 4,4,2 /) ) )
197 !
198  INTEGER :: ABQ(3,4,2)
199  parameter( abq = reshape( (/
200 ! NONSYMMETRICAL P1 QUASI-BUBBLE EBE (S=1)
201  & 0 , 4 , 7 ,
202  & 1 , 0 , 8 ,
203  & 2 , 5 , 0 ,
204  & 3 , 6 , 9 ,
205 ! NONSYMMETRICAL P1 QUASI-BUBBLE PRE-ASSEMBLED EBE (S=2)
206  & 0 , 7 , 3 ,
207  & 1 , 0 , 8 ,
208  & 9 , 2 , 0 ,
209  & 4 , 5 , 6 /), shape=(/ 3,4,2 /) ) )
210 !
211  INTEGER :: ACQ(3,6,2)
212  parameter( acq = reshape( (/
213 ! NONSYMMETRICAL P1 P2 EBE (S=1)
214  & 0 , 6 , 11,
215  & 1 , 0 , 12,
216  & 2 , 7 , 0 ,
217  & 3 , 8 , 13,
218  & 4 , 9 , 14,
219  & 5 , 10, 15,
220 ! S=2 NOT IMPLEMENTED
221  & 0 , 0 , 0 ,
222  & 0 , 0 , 0 ,
223  & 0 , 0 , 0 ,
224  & 0 , 0 , 0 ,
225  & 0 , 0 , 0 ,
226  & 0 , 0 , 0 /), shape=(/ 3,6,2 /) ) )
227 !
228  INTEGER :: BAQ(4,3,2)
229  parameter( baq = reshape( (/
230 ! NONSYMMETRICAL QUASI-BUBBLE P1 EBE (S=1)
231  & 0 , 3 , 5 , 7 ,
232  & 1 , 0 , 6 , 8 ,
233  & 2 , 4 , 0 , 9 ,
234 ! NONSYMMETRICAL QUASI-BUBBLE P1 PRE-ASSEMBLED EBE (S=2)
235  & 0 , 7 , 3 , 4 ,
236  & 1 , 0 , 8 , 5 ,
237  & 9 , 2 , 0 , 6 /), shape=(/ 4,3,2 /) ) )
238 !
239  INTEGER :: CAQ(6,3,2)
240 ! NONSYMMETRICAL P2 P1 EBE (S=1)
241  parameter( caq = reshape( (/
242 ! NONSYMMETRICAL P2 P1 EBE (S=1)
243  & 0 , 3 , 5 , 7 , 10 , 13 ,
244  & 1 , 0 , 6 , 8 , 11 , 14 ,
245  & 2 , 4 , 0 , 9 , 12 , 15 ,
246 ! NONSYMMETRICAL P2 P1 PRE-ASSEMBLED EBE (S=2) - NOT IMPLEMENTED
247  & 0 , 0 , 0 , 0 , 0 , 0 ,
248  & 0 , 0 , 0 , 0 , 0 , 0 ,
249  & 0 , 0 , 0 , 0 , 0 , 0 /), shape=(/ 6,3,2 /) ) )
250 !
251 !-----------------------------------------------------------------------
252 !
253  IF(s.EQ.1) THEN
254 !
255 !-----------------------------------------------------------------------
256 !
257 ! TRADITIONAL EBE STORAGE AND TRADITIONAL EBE MATRIX X VECTOR PRODUCT
258 !
259 !-----------------------------------------------------------------------
260 !
261  IF(ielm1.EQ.1) THEN
262 !
263  IF(ielm2.EQ.1) THEN
264  CALL mv0202(op, x , da,typdia,
265  & xa(1,1),xa(1,2),typext, y,c,
266  & ikle(1,1),ikle(1,2),
267  & npt,nelem,w(1,1),w(1,2))
268  ELSE
269  WRITE(lu,101) ielm1,ielm2,s
270  CALL plante(1)
271  stop
272  ENDIF
273 !
274  ELSEIF(ielm1.EQ.11) THEN
275 !
276  IF(ielm2.EQ.11) THEN
277  IF(typext(1:1).EQ.'S') THEN
278  IF (modass .EQ. 1)THEN
279  CALL mv0303(op, x , da,typdia,
280  & xa(1,aas(1,2,s)),
281  & xa(1,aas(1,3,s)),
282  & xa(1,aas(2,1,s)),
283  & xa(1,aas(2,3,s)),
284  & xa(1,aas(3,1,s)),
285  & xa(1,aas(3,2,s)),
286  & typext,y,c,
287  & ikle(1,1),ikle(1,2),ikle(1,3),
288  & npt,nelem,
289  & w(1,1),w(1,2),w(1,3))
290  ELSEIF (modass .EQ. 3)THEN
291  CALL mv0303(op, x , da,typdia,
292  & xa(1,aas(1,2,s)),
293  & xa(1,aas(1,3,s)),
294  & xa(1,aas(2,1,s)),
295  & xa(1,aas(2,3,s)),
296  & xa(1,aas(3,1,s)),
297  & xa(1,aas(3,2,s)),
298  & typext,y,c,
299  & ikle(1,1),ikle(1,2),ikle(1,3),
300  & npt,nelem,
301  & w(1,1),w(1,2),w(1,3)
302  & ,x_err,y_err,da_err)
303  ENDIF
304  ELSE
305  IF (modass .EQ. 1)THEN
306  CALL mv0303(op, x , da,typdia,
307  & xa(1,aaq(1,2,s)),
308  & xa(1,aaq(1,3,s)),
309  & xa(1,aaq(2,1,s)),
310  & xa(1,aaq(2,3,s)),
311  & xa(1,aaq(3,1,s)),
312  & xa(1,aaq(3,2,s)),
313  & typext,y,c,
314  & ikle(1,1),ikle(1,2),ikle(1,3),
315  & npt,nelem,
316  & w(1,1),w(1,2),w(1,3))
317  ELSEIF (modass .EQ. 3)THEN
318  CALL mv0303(op, x , da,typdia,
319  & xa(1,aaq(1,2,s)),
320  & xa(1,aaq(1,3,s)),
321  & xa(1,aaq(2,1,s)),
322  & xa(1,aaq(2,3,s)),
323  & xa(1,aaq(3,1,s)),
324  & xa(1,aaq(3,2,s)),
325  & typext,y,c,
326  & ikle(1,1),ikle(1,2),ikle(1,3),
327  & npt,nelem,
328  & w(1,1),w(1,2),w(1,3)
329  & ,x_err,y_err,da_err)
330  ENDIF
331  ENDIF
332 !
333  ELSEIF(ielm2.EQ.12) THEN
334 !
335  CALL mv0304(op, x , da,typdia,
336  & xa(1,abq(1,2,s)),
337  & xa(1,abq(1,3,s)),
338  & xa(1,abq(1,4,s)),
339  & xa(1,abq(2,1,s)),
340  & xa(1,abq(2,3,s)),
341  & xa(1,abq(2,4,s)),
342  & xa(1,abq(3,1,s)),
343  & xa(1,abq(3,2,s)),
344  & xa(1,abq(3,4,s)),
345  & typext, y,c,
346  & ikle(1,1),ikle(1,2),ikle(1,3),ikle(1,4),
347  & npt,nelem,
348  & w(1,1),w(1,2),w(1,3),w(1,4))
349 !
350  ELSEIF(ielm2.EQ.13) THEN
351 !
352  npt2=bief_nbpts(ielm2,mesh)
353  CALL mv0306(op, x , da,typdia,
354  & xa(1,acq(1,2,s)), xa(1,acq(1,3,s)),
355  & xa(1,acq(1,4,s)), xa(1,acq(1,5,s)),
356  & xa(1,acq(1,6,s)), xa(1,acq(2,1,s)),
357  & xa(1,acq(2,3,s)), xa(1,acq(2,4,s)),
358  & xa(1,acq(2,5,s)), xa(1,acq(2,6,s)),
359  & xa(1,acq(3,1,s)), xa(1,acq(3,2,s)),
360  & xa(1,acq(3,4,s)), xa(1,acq(3,5,s)),
361  & xa(1,acq(3,6,s)),
362  & typext, y,c,
363  & ikle(1,1),ikle(1,2),ikle(1,3),
364  & ikle(1,4),ikle(1,5),ikle(1,6),
365  & npt,npt2,nelem,
366  & w(1,1),w(1,2),w(1,3),
367  & w(1,4),w(1,5),w(1,6))
368 !
369  ELSE
370  WRITE(lu,101) ielm1,ielm2,s
371  CALL plante(1)
372  stop
373  ENDIF
374 !
375  ELSEIF(ielm1.EQ.12.OR.ielm2.EQ.31.OR.ielm2.EQ.51) THEN
376 !
377  IF(ielm2.EQ.12.OR.ielm2.EQ.31.OR.ielm2.EQ.51) THEN
378  IF(stox.EQ.1) THEN
379  IF(typext(1:1).EQ.'S') THEN
380  CALL mv0404(op, x , da,typdia,
381  & xa(1,bbs(1,2,s)),
382  & xa(1,bbs(1,3,s)),
383  & xa(1,bbs(1,4,s)),
384  & xa(1,bbs(2,1,s)),
385  & xa(1,bbs(2,3,s)),
386  & xa(1,bbs(2,4,s)),
387  & xa(1,bbs(3,1,s)),
388  & xa(1,bbs(3,2,s)),
389  & xa(1,bbs(3,4,s)),
390  & xa(1,bbs(4,1,s)),
391  & xa(1,bbs(4,2,s)),
392  & xa(1,bbs(4,3,s)),
393  & typext, y,c,
394  & ikle(1,1),ikle(1,2),ikle(1,3),ikle(1,4),
395  & npt,nelem,
396  & w(1,1),w(1,2),w(1,3),w(1,4))
397  ELSE
398  CALL mv0404(op, x , da,typdia,
399  & xa(1,bbq(1,2,s)),
400  & xa(1,bbq(1,3,s)),
401  & xa(1,bbq(1,4,s)),
402  & xa(1,bbq(2,1,s)),
403  & xa(1,bbq(2,3,s)),
404  & xa(1,bbq(2,4,s)),
405  & xa(1,bbq(3,1,s)),
406  & xa(1,bbq(3,2,s)),
407  & xa(1,bbq(3,4,s)),
408  & xa(1,bbq(4,1,s)),
409  & xa(1,bbq(4,2,s)),
410  & xa(1,bbq(4,3,s)),
411  & typext, y,c,
412  & ikle(1,1),ikle(1,2),ikle(1,3),ikle(1,4),
413  & npt,nelem,
414  & w(1,1),w(1,2),w(1,3),w(1,4))
415  ENDIF
416  ELSEIF(stox.EQ.2) THEN
417  IF(typext.EQ.'Q') THEN
418  dim1xa=12
419  ELSEIF(typext.EQ.'S') THEN
420  dim1xa=6
421  ELSE
422 ! NOT USED EXCEPT FOR GIVING A DIMENSION
423  dim1xa=1
424  ENDIF
425  CALL mv0404_2(op, x , da,typdia,xa,typext, y,c,
426  & ikle(1,1),ikle(1,2),ikle(1,3),ikle(1,4),
427  & npt,nelem,w(1,1),w(1,2),w(1,3),w(1,4),dim1xa)
428  ELSE
429  WRITE(lu,*) 'MATVCT, UNKNOWN STORAGE FOR'
430  WRITE(lu,*) 'OFF-DIAGONAL TERMS:',stox
431  WRITE(lu,*) 'WHEN CALLING MV0404_2'
432  CALL plante(1)
433  stop
434  ENDIF
435  ELSEIF(ielm2.EQ.11) THEN
436  CALL mv0403(op, x , da,typdia,
437  & xa(1,baq(1,2,s)),
438  & xa(1,baq(1,3,s)),
439  & xa(1,baq(2,1,s)),
440  & xa(1,baq(2,3,s)),
441  & xa(1,baq(3,1,s)),
442  & xa(1,baq(3,2,s)),
443  & xa(1,baq(4,1,s)),
444  & xa(1,baq(4,2,s)),
445  & xa(1,baq(4,3,s)),
446  & typext, y,c,
447  & ikle(1,1),ikle(1,2),ikle(1,3),ikle(1,4),
448  & npt,nelem,
449  & w(1,1),w(1,2),w(1,3),w(1,4))
450  ELSE
451  WRITE(lu,101) ielm1,ielm2,s
452  CALL plante(1)
453  stop
454  ENDIF
455 !
456  ELSEIF(ielm1.EQ.41.OR.ielm1.EQ.13) THEN
457 !
458  IF(ielm2.EQ.41.OR.ielm2.EQ.13) THEN
459 !
460  IF(stox.EQ.1) THEN
461  CALL mv0606(op, x , da,typdia,xa,typext, y,c,
462  & ikle(1,1),ikle(1,2),ikle(1,3),
463  & ikle(1,4),ikle(1,5),ikle(1,6),
464  & npt,nelem,nelmax,
465  & w(1,1),w(1,2),w(1,3),w(1,4),w(1,5),w(1,6))
466  ELSEIF(stox.EQ.2) THEN
467  IF(typext.EQ.'Q') THEN
468  dim1xa=30
469  ELSEIF(typext.EQ.'S') THEN
470  dim1xa=15
471  ELSE
472 ! NOT USED EXCEPT FOR GIVING A DIMENSION
473  dim1xa=1
474  ENDIF
475  CALL mv0606_2(op, x , da,typdia,xa,typext, y,c,
476  & ikle(1,1),ikle(1,2),ikle(1,3),
477  & ikle(1,4),ikle(1,5),ikle(1,6),
478  & npt,nelem,
479  & w(1,1),w(1,2),w(1,3),w(1,4),w(1,5),w(1,6),
480  & dim1xa)
481  ELSE
482  WRITE(lu,*) 'MATVCT, UNKNOWN STORAGE FOR'
483  WRITE(lu,*) 'OFF-DIAGONAL TERMS:',stox
484  CALL plante(1)
485  stop
486  ENDIF
487 !
488  ELSEIF(ielm2.EQ.11) THEN
489 !
490 ! HERE IELM1=13
491  npt2=bief_nbpts(ielm1,mesh)
492  CALL mv0603(op, x , da,typdia,
493  & xa(1,caq(1,2,s)),xa(1,caq(1,3,s)),
494  & xa(1,caq(2,1,s)),xa(1,caq(2,3,s)),
495  & xa(1,caq(3,1,s)),xa(1,caq(3,2,s)),
496  & xa(1,caq(4,1,s)),xa(1,caq(4,2,s)),
497  & xa(1,caq(4,3,s)),xa(1,caq(5,1,s)),
498  & xa(1,caq(5,2,s)),xa(1,caq(5,3,s)),
499  & xa(1,caq(6,1,s)),xa(1,caq(6,2,s)),
500  & xa(1,caq(6,3,s)),
501  & typext, y,c,
502  & ikle(1,1),ikle(1,2),ikle(1,3),
503  & ikle(1,4),ikle(1,5),ikle(1,6),
504  & npt,npt2,nelem,
505  & w(1,1),w(1,2),w(1,3),
506  & w(1,4),w(1,5),w(1,6))
507 !
508  ELSE
509  WRITE(lu,101) ielm1,ielm2,s
510  CALL plante(1)
511  stop
512  ENDIF
513 !
514 ! IELM1 NOT IMPLEMENTED : ERROR
515 !
516  ELSE
517 !
518  WRITE(lu,101) ielm1,ielm2,s
519 101 FORMAT(1x,'MATVCT (BIEF) : ELEMENTS ',1i2,' AND ',1i2,/,1x,
520  & 'AND STORAGE ',1i2,' CASE NOT IMPLEMENTED')
521  CALL plante(1)
522  stop
523 !
524  ENDIF
525 !
526 ! POSSIBLE FINAL ASSEMBLY OF X
527 !
528 ! SINCE INIT = FALSE HERE, MAY NOT NEED NPT
529  npt = bief_nbpts(ielmx,mesh)
530  IF(lego) THEN
531  IF(modass .EQ. 3) THEN
532  CALL assvec(x,ikle,npt,nelem,nelmax,w,
533  & .false.,lv,.false.,z,ndp,errx=x_err)
534  ELSEIF (modass .EQ. 1) THEN
535  CALL assvec(x,ikle,npt,nelem,nelmax,w,
536  & .false.,lv,.false.,z,ndp)
537  ENDIF
538  ENDIF
539 !
540  ELSEIF(s.EQ.3.AND.p.EQ.2) THEN
541 !
542 !-----------------------------------------------------------------------
543 !
544 ! SEGMENT STORAGE AND FRONTAL MATRIX X VECTOR PRODUCT
545 !
546 !-----------------------------------------------------------------------
547 !
548  IF(ielm1.EQ.1) THEN
549 !
550  IF(ielm2.EQ.1) THEN
551 ! CALL MW0202(OP, X , DA,TYPDIA,
552 ! * XA(1,1),XA(1,2),TYPEXT, Y,C,
553 ! * IKLE(1,1),IKLE(1,2),
554 ! * NPT,NELEM,
555 ! * W(1,1),W(1,2))
556 ! ELSE
557  WRITE(lu,101) ielm1,ielm2,s
558  CALL plante(1)
559  stop
560  ENDIF
561 !
562  ELSEIF(ielm1.EQ.11) THEN
563 !
564  IF(ielm2.EQ.11) THEN
565 !
566  CALL mw0303(op, x , da,typdia,xa,typext, y,c,
567  & iklem1,dimikm,limvoi,mxptvs,npmax,npoin,w)
568 !
569  ELSEIF(ielm2.EQ.12) THEN
570 !
571 ! CALL MW0304(OP, X , DA,TYPDIA,
572 ! * XA(1,1),XA(1,2),XA(1,3),
573 ! * XA(1,4),XA(1,5),XA(1,6),
574 ! * XA(1,7),XA(1,8),XA(1,9),
575 ! * TYPEXT, Y,C,
576 ! * IKLE(1,1),IKLE(1,2),IKLE(1,3),IKLE(1,4),
577 ! * NPT,NELEM,
578 ! * W(1,1),W(1,2),W(1,3),W(1,4))
579 ! ELSE
580  WRITE(lu,101) ielm1,ielm2,s
581  CALL plante(1)
582  stop
583  ENDIF
584 !
585  ELSEIF(ielm1.EQ.12) THEN
586 !
587  IF(ielm2.EQ.12) THEN
588 ! CALL MW0404(OP, X , DA,TYPDIA,
589 ! * XA(1,1),XA(1,2),XA(1,3),XA(1,1),
590 ! * XA(1,4),XA(1,5),XA(1,2),XA(1,4),
591 ! * XA(1,6),XA(1,3),XA(1,5),XA(1,6),
592 ! * TYPEXT, Y,C,
593 ! * IKLE(1,1),IKLE(1,2),IKLE(1,3),IKLE(1,4),
594 ! * NPT,NELEM,
595 ! * W(1,1),W(1,2),W(1,3),W(1,4))
596 ! ELSEIF(IELM2.EQ.11) THEN
597 ! CALL MW0403(OP, X , DA,TYPDIA,
598 ! * XA(1,1),XA(1,2),XA(1,3),
599 ! * XA(1,4),XA(1,5),XA(1,6),
600 ! * XA(1,7),XA(1,8),XA(1,9),
601 ! * TYPEXT, Y,C,
602 ! * IKLE(1,1),IKLE(1,2),IKLE(1,3),IKLE(1,4),
603 ! * NPT,NELEM,
604 ! * W(1,1),W(1,2),W(1,3),W(1,4))
605 ! ELSE
606  WRITE(lu,101) ielm1,ielm2,s
607  CALL plante(1)
608  stop
609  ENDIF
610 !
611  ELSEIF(ielm1.EQ.41) THEN
612 !
613  IF(ielm2.EQ.41) THEN
614 !
615 ! CALL MW0606(OP, X , DA,TYPDIA,XA,TYPEXT, Y,C,
616 ! * IKLE(1,1),IKLE(1,2),IKLE(1,3),
617 ! * IKLE(1,4),IKLE(1,5),IKLE(1,6),
618 ! * NPT,NELEM,NELMAX,
619 ! * W(1,1),W(1,2),W(1,3),W(1,4),W(1,5),W(1,6))
620 ! ELSE
621  WRITE(lu,101) ielm1,ielm2,s
622  CALL plante(1)
623  stop
624  ENDIF
625 !
626 ! IELM1 NOT IMPLEMENTED : ERROR
627 !
628  ELSE
629 !
630  WRITE(lu,101) ielm1,ielm2,s
631  CALL plante(1)
632  stop
633 !
634  ENDIF
635 !
636 ! STORAGE BY SEGMENTS
637 !
638  ELSEIF(s.EQ.3.AND.p.EQ.1) THEN
639 !
640 !-----------------------------------------------------------------------
641 !
642 ! SEGMENT STORAGE AND TRADITIONAL MATRIX X VECTOR PRODUCT
643 !
644 !-----------------------------------------------------------------------
645 !
646  nseg1 = bief_nbseg(ielm1,mesh)
647  nseg2 = bief_nbseg(ielm2,mesh)
648 !
649 ! IN LINEAR-QUADRATIC RECTANGULAR MATRICES, PURELY QUADRATIC
650 ! SEGMENTS ARE NOT CONSIDERED (NUMBER 13,14 AND 15, SO 3 PER ELEMENT)
651 !
652  IF(ielm1.EQ.11.AND.ielm2.EQ.13) THEN
653  nseg2=nseg2-3*nelem
654  ELSEIF(ielm1.EQ.13.AND.ielm2.EQ.11) THEN
655  nseg1=nseg1-3*nelem
656  ENDIF
657 !
658  IF(typext(1:1).EQ.'Q') THEN
659  sym = min(nseg1,nseg2)
660  ELSE
661  sym = 0
662  ENDIF
663  CALL mvseg (op, x , da,typdia,xa(1,1),xa(sym+1,1),
664  & typext,y,c,npt,nelem,nseg1,nseg2,
665  & gloseg(1,1),gloseg(1,2),ielm1,ielm2)
666 !
667 !-----------------------------------------------------------------------
668 !
669 ! STORAGE NOT IMPLEMENTED
670 !
671 !-----------------------------------------------------------------------
672 !
673  ELSE
674  WRITE(lu,103) s,p
675 103 FORMAT(1x,'MATVCT (BIEF) : ',1i2,' AND ',1i2,/,1x,
676  & 'STORAGE AND MATRIX-VECTOR PRODUCT INCOMPATIBLE')
677  CALL plante(1)
678  stop
679  ENDIF
680 !
681 !=======================================================================
682 !
683  RETURN
684  END
685 
subroutine matvct(OP, X, DA, TYPDIA, XA, TYPEXT, Y, C, IKLE, NPT, NELEM, NELMAX, W, LEGO, IELM1, IELM2, IELMX, LV, S, P, IKLEM1, DIMIKM, LIMVOI, MXPTVS, NPMAX, NPOIN, GLOSEG, SIZGLO, SIZXA, NDP, MESH, STOX, X_ERR, Y_ERR, DA_ERR)
Definition: matvct.f:11
subroutine mv0603(OP, X, DA, TYPDIA, XA12, XA13, XA21, XA23, XA31, XA32, XA41, XA42, XA43, XA51, XA52, XA53, XA61, XA62, XA63, TYPEXT, Y, C, IKLE1, IKLE2, IKLE3, IKLE4, IKLE5, IKLE6, NPOIN, NPT2, NELEM, W1, W2, W3, W4, W5, W6)
Definition: mv0603.f:12
subroutine mv0606(OP, X, DA, TYPDIA, XA, TYPEXT, Y, C, IKLE1, IKLE2, IKLE3, IKLE4, IKLE5, IKLE6, NPOIN, NELEM, NELMAX, W1, W2, W3, W4, W5, W6)
Definition: mv0606.f:9
integer function bief_nbpts(IELM, MESH)
Definition: bief_nbpts.f:7
subroutine mv0404_2(OP, X, DA, TYPDIA, XA, TYPEXT, Y, C, IKLE1, IKLE2, IKLE3, IKLE4, NPOIN, NELEM, W1, W2, W3, W4, DIM1XA)
Definition: mv0404_2.f:8
subroutine mv0303(OP, X, DA, TYPDIA, XA12, XA13, XA21, XA23, XA31, XA32, TYPEXT, Y, C, IKLE1, IKLE2, IKLE3, NPOIN, NELEM, W1, W2, W3, X_ERR, Y_ERR, DA_ERR)
Definition: mv0303.f:9
subroutine mw0303(OP, X, DA, TYPDIA, XAS, TYPEXT, Y, C, IKLEM1, DIMIKM, LIMVOI, MXPTVS, NPMAX, NPOIN, TRAV)
Definition: mw0303.f:8
subroutine assvec(X, IKLE, NPOIN, NELEM, NELMAX, W, INIT, LV, MSK, MASKEL, NDP, ERRX)
Definition: assvec.f:7
subroutine mv0606_2(OP, X, DA, TYPDIA, XA, TYPEXT, Y, C, IKLE1, IKLE2, IKLE3, IKLE4, IKLE5, IKLE6, NPOIN, NELEM, W1, W2, W3, W4, W5, W6, DIM1XA)
Definition: mv0606_2.f:9
subroutine mv0304(OP, X, DA, TYPDIA, XA12, XA13, XA14, XA21, XA23, XA24, XA31, XA32, XA34, TYPEXT, Y, C, IKLE1, IKLE2, IKLE3, IKLE4, NPOIN, NELEM, W1, W2, W3, W4)
Definition: mv0304.f:9
subroutine mv0202(OP, X, DA, TYPDIA, XA12, XA21, TYPEXT, Y, C, IKLE1, IKLE2, NPOIN, NELEM, W1, W2)
Definition: mv0202.f:8
subroutine mv0404(OP, X, DA, TYPDIA, XA12, XA13, XA14, XA21, XA23, XA24, XA31, XA32, XA34, XA41, XA42, XA43, TYPEXT, Y, C, IKLE1, IKLE2, IKLE3, IKLE4, NPOIN, NELEM, W1, W2, W3, W4)
Definition: mv0404.f:9
integer function bief_nbseg(IELM, MESH)
Definition: bief_nbseg.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
subroutine mv0403(OP, X, DA, TYPDIA, XA12, XA13, XA21, XA23, XA31, XA32, XA41, XA42, XA43, TYPEXT, Y, C, IKLE1, IKLE2, IKLE3, IKLE4, NPOIN, NELEM, W1, W2, W3, W4)
Definition: mv0403.f:11
subroutine mv0306(OP, X, DA, TYPDIA, XA12, XA13, XA14, XA15, XA16, XA21, XA23, XA24, XA25, XA26, XA31, XA32, XA34, XA35, XA36, TYPEXT, Y, C, IKLE1, IKLE2, IKLE3, IKLE4, IKLE5, IKLE6, NPOIN, NPT2, NELEM, W1, W2, W3, W4, W5, W6)
Definition: mv0306.f:12
Definition: bief.f:3