KTH framework for Nek5000 toolboxes; testing version  0.0.1
subs2.f
Go to the documentation of this file.
1  SUBROUTINE setaxw1 (IFAXWG)
2 C
3  include 'SIZE'
4  include 'WZ'
5 C
6  LOGICAL IFAXWG
7 C
8  IF (ifaxwg) THEN
9  CALL copy (w3m1,w2am1,lx1*ly1)
10  ELSE
11  CALL copy (w3m1,w2cm1,lx1*ly1)
12  ENDIF
13 C
14  RETURN
15  END
16  SUBROUTINE setaxw2 (IFAXWG)
17 C-----------------------------------------------------------------------
18  include 'SIZE'
19  include 'WZ'
20 C
21  LOGICAL IFAXWG
22 C
23  IF (ifaxwg) THEN
24  CALL copy (w3m2,w2am2,lx2*ly2)
25  ELSE
26  CALL copy (w3m2,w2cm2,lx2*ly2)
27  ENDIF
28 C
29  RETURN
30  END
31 C-----------------------------------------------------------------------
32  SUBROUTINE stnrinv
33 C
34 C Calculate 2nd and 3rd strain-rate invariants
35 C
36  include 'SIZE'
37  include 'SOLN'
38  include 'TSTEP'
39  common /screv/ ei2(lx1,ly1,lz1,lelt)
40  $ , ei3(lx1,ly1,lz1,lelt)
41  common /ctmp1/ exx(lx1,ly1,lz1,lelt)
42  $ , exy(lx1,ly1,lz1,lelt)
43  $ , eyy(lx1,ly1,lz1,lelt)
44  $ , ezz(lx1,ly1,lz1,lelt)
45  common /ctmp0/ exz(lx1,ly1,lz1,lelt)
46  $ , eyz(lx1,ly1,lz1,lelt)
47 c
48  ntot1 = lx1*ly1*lz1*nelv
49  CALL rzero (ei2,ntot1)
50  CALL rzero (ei3,ntot1)
51  IF (istep.EQ.0) RETURN
52 C
53  matmod = 0
54  CALL stnrate (vx,vy,vz,nelv,matmod)
55 C
56  IF (ldim.EQ.2) THEN
57  CALL col3 (ei2,exx,eyy,ntot1)
58  CALL subcol3 (ei2,exy,exy,ntot1)
59  CALL rzero (ei3,ntot1)
60  ELSE
61  const = 2.0
62  CALL col4 (ei3,exx,eyy,ezz,ntot1)
63  CALL col4 (ei2,exy,exz,eyz,ntot1)
64  CALL add2s2 (ei3,ei2,const,ntot1)
65  CALL subcol4 (ei3,exx,eyz,eyz,ntot1)
66  CALL subcol4 (ei3,eyy,exz,exz,ntot1)
67  CALL subcol4 (ei3,ezz,exy,exy,ntot1)
68  CALL col3 (ei2,exx,eyy,ntot1)
69  CALL addcol3 (ei2,exx,ezz,ntot1)
70  CALL addcol3 (ei2,eyy,ezz,ntot1)
71  CALL subcol3 (ei2,exy,exy,ntot1)
72  CALL subcol3 (ei2,exz,exz,ntot1)
73  CALL subcol3 (ei2,eyz,eyz,ntot1)
74  ENDIF
75 C
76  RETURN
77  END
78 C-----------------------------------------------------------------------
79  SUBROUTINE opdot (DP,A1,A2,A3,B1,B2,B3,N)
80 C
81  include 'SIZE'
82 C
83  dimension dp(lx1,ly1,lz1,1)
84  $ , a1(lx1,ly1,lz1,1)
85  $ , a2(lx1,ly1,lz1,1)
86  $ , a3(lx1,ly1,lz1,1)
87  $ , b1(lx1,ly1,lz1,1)
88  $ , b2(lx1,ly1,lz1,1)
89  $ , b3(lx1,ly1,lz1,1)
90 C
91  IF (ldim.EQ.2) THEN
92  CALL vdot2 (dp,a1,a2,b1,b2,n)
93  ELSE
94  CALL vdot3 (dp,a1,a2,a3,b1,b2,b3,n)
95  ENDIF
96 C
97  RETURN
98  END
99 C-----------------------------------------------------------------------
100  SUBROUTINE opadds (A1,A2,A3,B1,B2,B3,CONST,N,ISC)
101 C
102  include 'SIZE'
103 C
104  dimension a1(lx1,ly1,lz1,1)
105  $ , a2(lx1,ly1,lz1,1)
106  $ , a3(lx1,ly1,lz1,1)
107  $ , b1(lx1,ly1,lz1,1)
108  $ , b2(lx1,ly1,lz1,1)
109  $ , b3(lx1,ly1,lz1,1)
110 C
111  IF (isc.EQ.1) THEN
112  CALL add2s1 (a1,b1,const,n)
113  CALL add2s1 (a2,b2,const,n)
114  IF (ldim.EQ.3) CALL add2s1 (a3,b3,const,n)
115  ELSEIF (isc.EQ.2) THEN
116  CALL add2s2 (a1,b1,const,n)
117  CALL add2s2 (a2,b2,const,n)
118  IF (ldim.EQ.3) CALL add2s2 (a3,b3,const,n)
119  ENDIF
120 C
121  RETURN
122  END
123 C-----------------------------------------------------------------------
124  SUBROUTINE facexs (A,B,IFACE1,IOP)
125 C
126 C IOP = 0
127 C Extract scalar A from B on face IFACE1.
128 C
129 C IOP = 1
130 C Extract scalar B from A on face IFACE1.
131 C
132 C A has the (NX,NY,NFACE) data structure
133 C B has the (NX,NY,NZ) data structure
134 C IFACE1 is in the preprocessor notation
135 C IFACE is the dssum notation.
136 C
137  include 'SIZE'
138  include 'TOPOL'
139 C
140  dimension a(lx1,ly1),b(lx1,ly1,lz1)
141 C
142  CALL dsset(lx1,ly1,lz1)
143  iface = eface1(iface1)
144  js1 = skpdat(1,iface)
145  jf1 = skpdat(2,iface)
146  jskip1 = skpdat(3,iface)
147  js2 = skpdat(4,iface)
148  jf2 = skpdat(5,iface)
149  jskip2 = skpdat(6,iface)
150 C
151  i = 0
152  IF (iop.EQ.0) THEN
153  DO 100 j2=js2,jf2,jskip2
154  DO 100 j1=js1,jf1,jskip1
155  i = i+1
156  a(i,1) = b(j1,j2,1)
157  100 CONTINUE
158  ELSE
159  DO 150 j2=js2,jf2,jskip2
160  DO 150 j1=js1,jf1,jskip1
161  i = i+1
162  b(j1,j2,1) = a(i,1)
163  150 CONTINUE
164  ENDIF
165 C
166  RETURN
167  END
168 C-----------------------------------------------------------------------
169  SUBROUTINE facexv (A1,A2,A3,B1,B2,B3,IFACE1,IOP)
170 C
171 C IOP = 0
172 C Extract vector (A1,A2,A3) from (B1,B2,B3) on face IFACE1.
173 C
174 C IOP = 1
175 C Extract vector (B1,B2,B3) from (A1,A2,A3) on face IFACE1.
176 C
177 C A1, A2, A3 have the (NX,NY,NFACE) data structure
178 C B1, B2, B3 have the (NX,NY,NZ) data structure
179 C IFACE1 is in the preprocessor notation
180 C IFACE is the dssum notation.
181 C
182  include 'SIZE'
183  include 'TOPOL'
184 C
185  dimension a1(lx1,ly1),a2(lx1,ly1),a3(lx1,ly1),
186  $ b1(lx1,ly1,lz1),b2(lx1,ly1,lz1),b3(lx1,ly1,lz1)
187 C
188  CALL dsset(lx1,ly1,lz1)
189  iface = eface1(iface1)
190  js1 = skpdat(1,iface)
191  jf1 = skpdat(2,iface)
192  jskip1 = skpdat(3,iface)
193  js2 = skpdat(4,iface)
194  jf2 = skpdat(5,iface)
195  jskip2 = skpdat(6,iface)
196  i = 0
197 C
198  IF (iop.EQ.0) THEN
199  DO 100 j2=js2,jf2,jskip2
200  DO 100 j1=js1,jf1,jskip1
201  i = i+1
202  a1(i,1) = b1(j1,j2,1)
203  a2(i,1) = b2(j1,j2,1)
204  a3(i,1) = b3(j1,j2,1)
205  100 CONTINUE
206  ELSE
207  DO 150 j2=js2,jf2,jskip2
208  DO 150 j1=js1,jf1,jskip1
209  i = i+1
210  b1(j1,j2,1) = a1(i,1)
211  b2(j1,j2,1) = a2(i,1)
212  b3(j1,j2,1) = a3(i,1)
213  150 CONTINUE
214  ENDIF
215 C
216  RETURN
217  END
218 C-----------------------------------------------------------------------
219  SUBROUTINE facsub2 (A1,A2,A3,B1,B2,B3,IFACE1)
220 C
221 C Subtract B1,B2,B3 from A1,A2,A3 on surface IFACE1 of element IE.
222 C
223 C A1, A2, A3 have the (NX,NY,NZ) data structure
224 C B1, B2, B3 have the (NX,NY,NFACE) data structure
225 C IFACE1 is in the preprocessor notation
226 C IFACE is the dssum notation.
227 C
228  include 'SIZE'
229  include 'TOPOL'
230 C
231  dimension a1(lx1,ly1,lz1),a2(lx1,ly1,lz1),a3(lx1,ly1,lz1),
232  $ b1(lx1,ly1),b2(lx1,ly1),b3(lx1,ly1)
233 C
234  CALL dsset(lx1,ly1,lz1)
235  iface = eface1(iface1)
236  js1 = skpdat(1,iface)
237  jf1 = skpdat(2,iface)
238  jskip1 = skpdat(3,iface)
239  js2 = skpdat(4,iface)
240  jf2 = skpdat(5,iface)
241  jskip2 = skpdat(6,iface)
242 C
243  i = 0
244  DO 100 j2=js2,jf2,jskip2
245  DO 100 j1=js1,jf1,jskip1
246  i = i+1
247  a1(j1,j2,1) = a1(j1,j2,1) - b1(i,1)
248  a2(j1,j2,1) = a2(j1,j2,1) - b2(i,1)
249  a3(j1,j2,1) = a3(j1,j2,1) - b3(i,1)
250  100 CONTINUE
251 C
252  RETURN
253  END
254 C-----------------------------------------------------------------------
255  SUBROUTINE gammasf (H1,H2)
256 C-----------------------------------------------------------------------
257 C
258 C Compute lagrest eigenvalue of coupled Helmholtz operator
259 C
260 C-----------------------------------------------------------------------
261  include 'SIZE'
262  include 'EIGEN'
263  include 'INPUT'
264  include 'MASS'
265  include 'MVGEOM'
266  include 'SOLN'
267  include 'TSTEP'
268  include 'WZ'
269  common /scrmg/ ae1(lx1,ly1,lz1,lelv)
270  $ , ae2(lx1,ly1,lz1,lelv)
271  $ , ae3(lx1,ly1,lz1,lelv)
272  common /scruz/ e1(lx1,ly1,lz1,lelv)
273  $ , e2(lx1,ly1,lz1,lelv)
274  $ , e3(lx1,ly1,lz1,lelv)
275 C
276  dimension h1(lx1,ly1,lz1,1),h2(lx1,ly1,lz1,1)
277 C
278  ntot1 = lx1*ly1*lz1*nelv
279  imesh = 1
280  matmod = 0
281 C
282  IF (istep.EQ.0) THEN
283  eigga = 0.0
284  CALL stx1sf
285  ELSE
286  CALL copy (e1,ev1,ntot1)
287  CALL copy (e2,ev2,ntot1)
288  IF (ldim.EQ.3) CALL copy (e3,ev3,ntot1)
289  ENDIF
290 C
291  evnew = eigga
292 C
293  DO 1000 iter=1,nmxe
294 C
295  CALL axhmsf (ae1,ae2,ae3,e1,e2,e3,h1,h2,matmod)
296  CALL rmask (ae1,ae2,ae3,nelv)
297  CALL opdssum (ae1,ae2,ae3)
298 C
299  evold = evnew
300  evnew = glsc3(e1,ae1,vmult,ntot1) + glsc3(e2,ae2,vmult,ntot1)
301  IF (ldim.EQ.3) evnew = evnew + glsc3(e3,ae3,vmult,ntot1)
302  crit = abs( (evnew - evold)/evnew )
303  IF ( crit .LT. tolev ) GOTO 2000
304 C
305  CALL col3 (e1,binvm1,ae1,ntot1)
306  CALL col3 (e2,binvm1,ae2,ntot1)
307  IF (ldim.EQ.3) CALL col3 (e3,binvm1,ae3,ntot1)
308  xx = glsc3(e1,ae1,vmult,ntot1) + glsc3(e2,ae2,vmult,ntot1)
309  IF (ldim.EQ.3) xx = xx + glsc3(e3,ae3,vmult,ntot1)
310  IF (xx .LT. 0.0) GO TO 9000
311 C
312  xnorm=1./sqrt( xx )
313  CALL cmult (e1,xnorm,ntot1)
314  CALL cmult (e2,xnorm,ntot1)
315  IF (ldim.EQ.3) CALL cmult (e3,xnorm,ntot1)
316 C
317  1000 CONTINUE
318 C
319 C Save eigenvalue for all cases.
320 C Save eigenvectors for deforming geometries.
321 C
322  2000 eigga = evnew
323  IF (ifmvbd) THEN
324  CALL copy (ev1,e1,ntot1)
325  CALL copy (ev2,e2,ntot1)
326  IF (ldim.EQ.3) CALL copy (ev3,e3,ntot1)
327  ENDIF
328 C
329  RETURN
330 C
331  9000 CONTINUE
332  IF (nid.EQ.0)
333  $WRITE ( 6,*) ' Non +ve def. A-operator detected during eigenvalue
334  $computation : tran(x)Ax =',xx
335  CALL emerxit
336  call exitt
337 C
338  END
339 C-----------------------------------------------------------------------
340  SUBROUTINE cmult2 (A,B,CONST,N)
341  dimension a(1),b(1)
342  DO 100 i=1,n
343  a(i)=b(i)*const
344  100 CONTINUE
345  RETURN
346  END
347 C-----------------------------------------------------------------------
348  SUBROUTINE add3s (A,B,C,CONST,N)
349  dimension a(1),b(1),c(1)
350  DO 100 i=1,n
351  a(i)=b(i)+const*c(i)
352  100 CONTINUE
353  RETURN
354  END
355 C-----------------------------------------------------------------------
356  SUBROUTINE emerxit
357 C
358  include 'SIZE'
359  include 'INPUT'
360  include 'TSTEP'
361  include 'PARALLEL'
362 C
363 C Try to hang in there on the first few time steps (pff 8/92)
364  IF (iftran.AND.istep.LT.9) RETURN
365 C
366  lastep = 1
367  CALL prepost(.true.,' ')
368 C
369  IF (np.EQ.1) THEN
370  WRITE (6,*) ' '
371  WRITE (6,*)
372  $ ' Emergency exit:',istep,' time =',time
373  WRITE (6,*)
374  $ ' Latest solution and data are dumped for post-processing.'
375  WRITE (6,*) ' *** STOP ***'
376  ELSE
377  WRITE (6,*) ' '
378  WRITE (6,*) nid,
379  $ ' Emergency exit:',istep,' time =',time
380  WRITE (6,*)
381  $ ' Latest solution and data are dumped for post-processing.'
382  WRITE (6,*) ' *** STOP ***'
383  ENDIF
384 C
385  call runstat
386 c
387  call exitt
388  RETURN
389  END
390 C-----------------------------------------------------------------------
391  SUBROUTINE faccvs (A1,A2,A3,B,IFACE1)
392 C
393 C Collocate scalar B with vector A, components A1,A2,A3,
394 C on the surface IFACE1 of an element.
395 C
396 C A1,A2,A3 have the (NX,NY,NZ) data structure
397 C B has the (NX,NY,IFACE) data structure
398 C IFACE1 is in the preprocessor notation
399 C IFACE is the dssum notation.
400 C
401  include 'SIZE'
402  include 'TOPOL'
403  dimension a1(lx1,ly1,lz1),a2(lx1,ly1,lz1),a3(lx1,ly1,lz1),
404  $ b(lx1,ly1)
405 C
406 C Set up counters
407 C
408  CALL dsset(lx1,ly1,lz1)
409  iface = eface1(iface1)
410  js1 = skpdat(1,iface)
411  jf1 = skpdat(2,iface)
412  jskip1 = skpdat(3,iface)
413  js2 = skpdat(4,iface)
414  jf2 = skpdat(5,iface)
415  jskip2 = skpdat(6,iface)
416  i = 0
417 C
418  IF (ldim.EQ.2) THEN
419  DO 100 j2=js2,jf2,jskip2
420  DO 100 j1=js1,jf1,jskip1
421  i = i+1
422  a1(j1,j2,1) = a1(j1,j2,1)*b(i,1)
423  a2(j1,j2,1) = a2(j1,j2,1)*b(i,1)
424  100 CONTINUE
425  ELSE
426  DO 200 j2=js2,jf2,jskip2
427  DO 200 j1=js1,jf1,jskip1
428  i = i+1
429  a1(j1,j2,1) = a1(j1,j2,1)*b(i,1)
430  a2(j1,j2,1) = a2(j1,j2,1)*b(i,1)
431  a3(j1,j2,1) = a3(j1,j2,1)*b(i,1)
432  200 CONTINUE
433  ENDIF
434 C
435  RETURN
436  END
437 C-----------------------------------------------------------------------
438  SUBROUTINE stx1sf
439 C------------------------------------------------------------------
440 C
441 C Compute startvector for finding an eigenvalue on mesh 1.
442 C Normalization: XT*B*X = 1
443 C
444 C------------------------------------------------------------------
445  include 'SIZE'
446  include 'MASS'
447  include 'SOLN'
448  common /scrmg/ ae1(lx1,ly1,lz1,lelv)
449  $ , ae2(lx1,ly1,lz1,lelv)
450  $ , ae3(lx1,ly1,lz1,lelv)
451  common /scruz/ e1(lx1,ly1,lz1,lelv)
452  $ , e2(lx1,ly1,lz1,lelv)
453  $ , e3(lx1,ly1,lz1,lelv)
454 C
455  ntot1 = lx1*ly1*lz1*nelv
456  CALL rzero3 (e1 ,e2 ,e3 ,ntot1)
457  CALL rzero3 (ae1,ae2,ae3,ntot1)
458 C
459  CALL copy (e1,bm1,ntot1)
460  CALL copy (e2,bm1,ntot1)
461  IF (ldim.EQ.3) CALL copy (e3,bm1,ntot1)
462 C
463  CALL rmask (e1,e2,e3,nelv)
464  CALL col3 (ae1,bm1,e1,ntot1)
465  CALL col3 (ae2,bm1,e2,ntot1)
466  IF (ldim.EQ.3) CALL col3 (ae3,bm1,e3,ntot1)
467 C
468  CALL opdssum (ae1,ae2,ae3)
469 C
470  xx = glsc3(e1,ae1,vmult,ntot1) + glsc3(e2,ae2,vmult,ntot1)
471  IF (ldim.EQ.3) xx = xx + glsc3(e3,ae3,vmult,ntot1)
472  xnorm = 1./sqrt(xx)
473  CALL cmult (e1,xnorm,ntot1)
474  CALL cmult (e2,xnorm,ntot1)
475  IF (ldim.EQ.3) CALL cmult (e3,xnorm,ntot1)
476 
477 c call exitti ('quit in stx1sf$,',nel)
478 
479 C
480  RETURN
481  END
482 C-----------------------------------------------------------------------
483  SUBROUTINE solvel
484 C
485  include 'SIZE'
486  include 'GEOM'
487  include 'SOLN'
488  include 'TSTEP'
489 C
490  DO 100 iel=1,nelv
491  DO 100 k=1,lz1
492  DO 100 j=1,ly1
493  DO 100 i=1,lx1
494  CALL vsoln (vx(i,j,k,iel),vy(i,j,k,iel),vz(i,j,k,iel),
495  $ xm1(i,j,k,iel),ym1(i,j,k,iel),zm1(i,j,k,iel),pi)
496  100 CONTINUE
497 C
498  RETURN
499  END
500 C-----------------------------------------------------------------------
501  SUBROUTINE vsoln (UX,UY,UZ,X,Y,Z,PI)
502 C
503 C URR=(1.-0.75/SQRT(X**2+Y**2)+0.0625/(SQRT(X**2+Y**2)**3))*(X/
504 C $ SQRT(X**2+Y**2))
505 C UTETA=-(1.-0.375/SQRT(X**2+Y**2)-0.03125/(SQRT(X**2+Y**2)**3))*
506 C $ (Y/SQRT(X**2+Y**2))
507 C UX=URR*(X/SQRT(X**2+Y**2))-UTETA*(Y/SQRT(X**2+Y**2))
508 C UY=URR*(Y/SQRT(X**2+Y**2))+UTETA*(X/SQRT(X**2+Y**2))
509 C
510 C UX = 2.*COS( PI*X )
511 C UY = PI*Y*SIN( PI*X )
512 C
513  ux = 0.0
514  uy = 0.0
515 C
516  RETURN
517  END
518 C-----------------------------------------------------------------------
519  SUBROUTINE solpres
520 C
521  include 'SIZE'
522  include 'GEOM'
523  include 'SOLN'
524  include 'TSTEP'
525 C
526  DO 100 iel=1,nelv
527  DO 100 k=1,lz2
528  DO 100 j=1,ly2
529  DO 100 i=1,lx2
530  CALL prsoln (pr(i,j,k,iel),xm2(i,j,k,iel),ym2(i,j,k,iel),
531  * zm2(i,j,k,iel),pi)
532  100 CONTINUE
533 C
534  RETURN
535  END
536 C-----------------------------------------------------------------------
537  SUBROUTINE prsoln (P,X,Y,Z,PI)
538 C
539 C R = SQRT( X**2 + Y**2 )
540 C CS = X/R
541 C P = -0.75 * CS / R**2
542 C
543 C P = -SIN( PI*X )*COS( PI*Y )
544 C
545  p = 0.0
546 C
547  RETURN
548  END
549 C-----------------------------------------------------------------------
550  SUBROUTINE printel (TA,A,IEL)
551 C
552  include 'SIZE'
553  dimension ta(lx1,ly1,lz1,lelt)
554  CHARACTER A*10
555 C
556  lz1i = 1
557  lz1j = lz1
558  lz1inc = 1
559 C
560  WRITE (21,*) 'ELEMENT NUMBER ',iel
561  DO 101 ipl=lz1i,lz1j,lz1inc
562  CALL outm1 (ta,a,lz1,iel,ipl)
563  101 CONTINUE
564 C
565  RETURN
566  END
567 C-----------------------------------------------------------------------
568  SUBROUTINE printv (TA,A,NEL)
569 C-----------------------------------------------------------------------
570 C
571 C Store the results
572 C
573 C-----------------------------------------------------------------------
574  include 'SIZE'
575  dimension ta(lx1,ly1,lz1,lelt)
576  CHARACTER A*10
577 C
578  lz1i = 1
579  lz1j = lz1
580  lz1inc = 1
581 C
582  DO 9001 iel = 1,nel
583  WRITE (21,*) 'ELEMENT NUMBER ',iel
584  DO 101 ipl=lz1i,lz1j,lz1inc
585  CALL outm1 (ta,a,lz1,iel,ipl)
586  101 CONTINUE
587  9001 CONTINUE
588 C
589  RETURN
590  END
591 C-----------------------------------------------------------------------
592  SUBROUTINE outf1 (X,TXT,IEL,IFC)
593  include 'SIZE'
594  dimension x(lx1,lz1,6,lelt)
595  CHARACTER*10 TXT
596 C
597  nface = 2*ldim
598  nzi = lz1
599  nzj = 1
600  nzinc = -1
601  nxi = 1
602  nxj = lx1
603  nxinc = 1
604 C
605  WRITE(21,106) txt,ifc,nface
606  DO 100 j=nzi,nzj,nzinc
607  WRITE(21,105) (x(i,j,ifc,iel),i=nxi,nxj,nxinc)
608  100 CONTINUE
609 C
610  105 FORMAT(5e15.6)
611  106 FORMAT(///,5x,' ^ ',/,
612  $ 5x,' S | ',/,
613  $ 5x,' | ',a10,/,
614  $ 5x,' +----> ','Plane = ',i2,'/',i2,/,
615  $ 5x,' R ',/)
616 C
617  RETURN
618  END
619 C-----------------------------------------------------------------------
620  SUBROUTINE outm1 (X,TXT,NP,IEL,IP)
621  include 'SIZE'
622  dimension x(lx1,ly1,lz1,lelt)
623  CHARACTER*10 TXT
624 C
625  nyi = ly1
626  nyj = 1
627  nyinc = -1
628  nxi = 1
629  nxj = lx1
630  nxinc = 1
631 C
632  WRITE(6,106) txt,ip,np
633  DO 100 j=nyi,nyj,nyinc
634  WRITE(6,105) (x(i,j,ip,iel),i=nxi,nxj,nxinc)
635  100 CONTINUE
636 C
637 c 105 FORMAT(1p8e10.3)
638  105 FORMAT(8f10.3)
639  106 FORMAT(///,5x,' ^ ',/,
640  $ 5x,' Y | ',/,
641  $ 5x,' | ',a10,/,
642  $ 5x,' +----> ','Plane = ',i2,'/',i2,/,
643  $ 5x,' X ',/)
644 C
645  RETURN
646  END
647 C-----------------------------------------------------------------------
648  SUBROUTINE outm2 (X,TXT,NP,IEL,IP)
649  include 'SIZE'
650  dimension x(lx2,ly2,lz2,lelv)
651  CHARACTER*10 TXT
652 C
653  nyi = ly2
654  nyj = 1
655  nyinc = -1
656  nxi = 1
657  nxj = lx2
658  nxinc = 1
659 C
660  WRITE(21,106) txt,ip,np
661  DO 100 j=nyi,nyj,nyinc
662  WRITE(21,105) (x(i,j,ip,iel),i=nxi,nxj,nxinc)
663  100 CONTINUE
664 C
665  105 FORMAT(5e15.6)
666  106 FORMAT(///,5x,' ^ ',/,
667  $ 5x,' Y | ',/,
668  $ 5x,' | ',a10,/,
669  $ 5x,' +----> ','Plane = ',i2,'/',i2,/,
670  $ 5x,' X ',/)
671 C
672  RETURN
673  END
674 C-----------------------------------------------------------------------
675  SUBROUTINE stsmask (C1MASK,C2MASK,C3MASK)
676 C
677  include 'SIZE'
678  include 'GEOM'
679  include 'TSTEP'
680  include 'INPUT'
681  common /screv/ hfmask(lx1,lz1,6,lelt)
682  $ , hvmask(lx1,ly1,lz1,lelt)
683 C
684  dimension c1mask(lx1,ly1,lz1,1)
685  $ , c2mask(lx1,ly1,lz1,1)
686  $ , c3mask(lx1,ly1,lz1,1)
687  INTEGER IMDATA
688  SAVE imdata
689  DATA imdata /0/
690 C
691  ifld = ifield
692  nel = nelfld(ifield)
693 C
694  IF (imdata.EQ.0) THEN
695  CALL setcdat
696  imdata=1
697  ENDIF
698 C
699  IF (ifld.EQ.1) CALL skipcnr (nel)
700  CALL sethmsk (hvmask,hfmask,ifld,nel)
701  CALL setmlog (hvmask,hfmask,ifld,nel)
702  CALL setmask (c1mask,c2mask,c3mask,hvmask,nel)
703  IF (iflmsf(ifld)) CALL setcsys (hvmask,hfmask,nel)
704  IF (ifld.EQ.0) CALL fixwmsk (c2mask,c3mask,hvmask,hfmask,nel)
705 
706  if (ifaxis.and.ifld.eq.1) call fixmska (c1mask,c2mask,c3mask)
707 C
708  RETURN
709  END
710 C-----------------------------------------------------------------------
711  SUBROUTINE updmsys (IFLD)
712 C
713  include 'SIZE'
714  include 'GEOM'
715  include 'TSTEP'
716  common /screv/ hfmask(lx1,lz1,6,lelt)
717  $ , hvmask(lx1,ly1,lz1,lelt)
718 C
719  IF (.NOT.iflmsf(ifld)) RETURN
720 C
721  nel = nelfld(ifld)
722  CALL sethmsk (hvmask,hfmask,ifld,nel)
723  CALL setcsys (hvmask,hfmask,nel)
724 C
725  RETURN
726  END
727 C-----------------------------------------------------------------------
728  SUBROUTINE sethmsk (HVMASK,HFMASK,IFLD,NEL)
729 C
730  include 'SIZE'
731  include 'INPUT'
732  include 'TSTEP'
733 C
734  dimension hvmask(lx1,ly1,lz1,1)
735  $ , hfmask(lx1,lz1,6,1)
736  CHARACTER CB*3
737 C
738  ntot1 = lx1*ly1*lz1*nel
739  nxz1 = lx1*lz1
740  ntotf = nxz1*6*nel
741  nface = 2*ldim
742  const = 5.0
743  CALL cfill (hvmask,const,ntot1)
744  CALL cfill (hfmask,const,ntotf)
745 C
746  IF (ifld.EQ.1) THEN
747 C
748  DO 110 iel=1,nel
749  DO 110 ifc=1,nface
750  cb=cbc(ifc,iel,ifld)
751  IF (cb.EQ.'ON ' .OR. cb.EQ.'on ' .or.
752  $ cb.EQ.'MM ' .OR. cb.EQ.'mm ' ) THEN
753  CALL facev (hvmask,iel,ifc,3.0,lx1,ly1,lz1)
754  CALL cfill (hfmask(1,1,ifc,iel),3.0,nxz1)
755  ENDIF
756  110 CONTINUE
757 C
758  DO 120 iel=1,nel
759  DO 120 ifc=1,nface
760  cb=cbc(ifc,iel,ifld)
761  IF (cb.EQ.'SYM' .OR. cb.EQ.'A ' .OR. cb.EQ.'WS ' .OR.
762  $ cb.EQ.'ws ' .OR. cb.EQ.'WSL' .OR. cb.EQ.'wsl' .OR.
763  $ cb.EQ.'SH ' .OR. cb.EQ.'sh ' .OR. cb.EQ.'SHL' .OR.
764  $ cb.EQ.'shl') THEN
765  CALL facev (hvmask,iel,ifc,2.0,lx1,ly1,lz1)
766  CALL cfill (hfmask(1,1,ifc,iel),2.0,nxz1)
767  ENDIF
768  120 CONTINUE
769 C
770  DO 130 iel=1,nel
771  DO 130 ifc=1,nface
772  cb=cbc(ifc,iel,ifld)
773  IF (cb.EQ.'MF ' .OR. cb.EQ.'V ' .OR. cb.EQ.'v ' .OR.
774  $ cb.EQ.'VL ' .OR. cb.EQ.'vl ' .OR. cb(1:2).EQ.'mv') THEN
775  CALL facev (hvmask,iel,ifc,1.0,lx1,ly1,lz1)
776  CALL cfill (hfmask(1,1,ifc,iel),1.0,nxz1)
777  ENDIF
778  130 CONTINUE
779 C
780  DO 140 iel=1,nel
781  DO 140 ifc=1,nface
782  cb=cbc(ifc,iel,ifld)
783  IF (cb.EQ.'W ') THEN
784  CALL facev (hvmask,iel,ifc,0.0,lx1,ly1,lz1)
785  CALL cfill (hfmask(1,1,ifc,iel),0.0,nxz1)
786  ENDIF
787  140 CONTINUE
788 C
789  ELSE
790 C
791  DO 210 iel=1,nel
792  DO 210 ifc=1,nface
793  cb=cbc(ifc,iel,ifld)
794  IF (cb.EQ.'SYM') THEN
795  CALL facev (hvmask,iel,ifc,2.0,lx1,ly1,lz1)
796  CALL cfill (hfmask(1,1,ifc,iel),2.0,nxz1)
797  ENDIF
798  210 CONTINUE
799 C
800 c write(6,*) 'MASK this is ifield:',ifield
801  DO 220 iel=1,nel
802  DO 220 ifc=1,nface
803  cb=cbc(ifc,iel,ifld)
804  IF (cb(1:1).EQ.'M' .OR. cb(1:1).EQ.'m') THEN
805 c CALL FACEV (HVMASK,IEL,IFC,1.0,lx1,ly1,lz1)
806 c CALL CFILL (HFMASK(1,1,IFC,IEL),1.0,NXZ1)
807  CALL facev (hvmask,iel,ifc,2.0,lx1,ly1,lz1)
808  CALL cfill (hfmask(1,1,ifc,iel),2.0,nxz1)
809  ENDIF
810  220 CONTINUE
811 C
812  DO 230 iel=1,nel
813  DO 230 ifc=1,nface
814  cb=cbc(ifc,iel,ifld)
815  IF (cb.EQ.'FIX') THEN
816  CALL facev (hvmask,iel,ifc,0.0,lx1,ly1,lz1)
817  CALL cfill (hfmask(1,1,ifc,iel),0.0,nxz1)
818  ENDIF
819  230 CONTINUE
820 C
821  ENDIF
822 C
823  CALL dsop (hvmask,'MNA',lx1,ly1,lz1)
824 C
825  RETURN
826  END
827 C-----------------------------------------------------------------------
828  SUBROUTINE skipcnr (NEL)
829 C
830  include 'SIZE'
831  include 'GEOM'
832  include 'INPUT'
833  common /indxfc/ mcrfc(4,6)
834  $ , mfccr(3,8)
835  $ , megcr(3,8)
836  $ , mfceg(2,12)
837  $ , mcreg(2,12)
838  $ , mcrrst(3,8)
839  $ , midrst(3,12)
840  $ , mcrind(8)
841  $ , medind(2,4)
842  $ , ntefc(2,12)
843  $ , ntcrf(2,3)
844 C
845  nface = 2*ldim
846  ncrfc = nface - 2
847  nmxcr = 8*nel
848  CALL lfalse (ifnskp,nmxcr)
849 C
850  DO 100 iel=1,nel
851  DO 100 ifc=1,nface
852  IF (cdof(ifc,iel).EQ.'1') THEN
853  icr=mcrfc(1,ifc)
854  ifnskp(icr,iel)=.true.
855  ELSEIF (cdof(ifc,iel).EQ.'2') THEN
856  icr=mcrfc(2,ifc)
857  ifnskp(icr,iel)=.true.
858  ELSEIF (cdof(ifc,iel).EQ.'3') THEN
859  icr=mcrfc(3,ifc)
860  ifnskp(icr,iel)=.true.
861  ELSEIF (cdof(ifc,iel).EQ.'4') THEN
862  icr=mcrfc(4,ifc)
863  ifnskp(icr,iel)=.true.
864  ENDIF
865  IF (cdof(ifc,iel).EQ.'*') THEN
866  DO 160 icrfc=1,ncrfc
867  icr=mcrfc(icrfc,ifc)
868  ifnskp(icr,iel)=.true.
869  160 CONTINUE
870  ENDIF
871  100 CONTINUE
872 C
873  RETURN
874  END
875 C-----------------------------------------------------------------------
876  SUBROUTINE setmask (C1MASK,C2MASK,C3MASK,HVMASK,NEL)
877 C
878  include 'SIZE'
879  include 'TOTAL' !XXXXX
880 C
881  dimension hvmask(lx1,ly1,lz1,1)
882  $ , c1mask(lx1,ly1,lz1,1)
883  $ , c2mask(lx1,ly1,lz1,1)
884  $ , c3mask(lx1,ly1,lz1,1)
885 C
886  ntot1 = lx1*ly1*lz1*nel
887  CALL rzero3 (c1mask,c2mask,c3mask,ntot1)
888 C
889  DO 100 iel=1,nel
890  DO 100 iz=1,lz1
891  DO 100 iy=1,ly1
892  DO 100 ix=1,lx1
893  hmv=abs( hvmask(ix,iy,iz,iel) )
894  IF (hmv .GT. 2.9) THEN
895  c1mask(ix,iy,iz,iel) = 1.0
896  ENDIF
897  IF ((hmv.GT.1.9 .AND. hmv.LT.2.1) .OR. hmv.GT.4.9) THEN
898  c2mask(ix,iy,iz,iel) = 1.0
899  ENDIF
900  100 CONTINUE
901 C
902  IF (ldim.EQ.3) CALL copy (c3mask,c2mask,ntot1)
903 C
904  RETURN
905  END
906 C-----------------------------------------------------------------------
907  SUBROUTINE setmlog (HVMASK,HFMASK,IFLD,NEL)
908 C
909  include 'SIZE'
910  include 'GEOM'
911  common /indxfc/ mcrfc(4,6)
912  $ , mfccr(3,8)
913  $ , megcr(3,8)
914  $ , mfceg(2,12)
915  $ , mcreg(2,12)
916  $ , mcrrst(3,8)
917  $ , midrst(3,12)
918  $ , mcrind(8)
919  $ , medind(2,4)
920  $ , ntefc(2,12)
921  $ , ntcrf(2,3)
922 C
923  dimension hvmask(lx1,ly1,lz1,1)
924  $ , hfmask(lx1,lz1,6,1)
925 C
926  nface = 2*ldim
927  nedge = 12
928  ncrnr = 2**ldim
929  ntotf = nface*nel
930  ntots = nedge*nel
931  ntotc = ncrnr*nel
932  epsa = 1.e-6
933 C
934  iflmsf(ifld) = .false.
935  iflmse(ifld) = .false.
936  iflmsc(ifld) = .false.
937  CALL lfalse (ifmsfc(1,1,ifld),ntotf)
938  CALL lfalse (ifmseg(1,1,ifld),ntots)
939  CALL lfalse (ifmscr(1,1,ifld),ntotc)
940 C
941  DO 100 iel=1,nel
942  DO 100 ifc=1,nface
943  hmf = abs( hfmask(1,1,ifc,iel) )
944  IF (hmf .GT. 1.9 .AND. hmf .LT. 3.1 ) THEN
945  iflmsf(ifld) = .true.
946  ifmsfc(ifc,iel,ifld) = .true.
947  ENDIF
948  100 CONTINUE
949  CALL gllog(iflmsf(ifld),.true.)
950 C
951  IF (ldim.EQ.3) THEN
952  DO 200 iel=1,nel
953  DO 200 isd=1,nedge
954  ix = midrst(1,isd)
955  iy = midrst(2,isd)
956  iz = midrst(3,isd)
957  hmv = abs( hvmask(ix,iy,iz,iel) )
958  IF (hmv .LT. 1.9 .OR. hmv .GT. 3.1) GOTO 200
959  idiff = 0
960  DO 220 ii=1,2
961  ifc = mfceg(ii,isd)
962  hmf = abs( hfmask(1,1,ifc,iel) )
963  IF (abs(hmv - hmf) .GT. epsa) idiff=idiff + 1
964  220 CONTINUE
965  IF (idiff.EQ.2) THEN
966  iflmse(ifld) = .true.
967  ifmseg(isd,iel,ifld) = .true.
968  ENDIF
969  200 CONTINUE
970  CALL gllog(iflmse(ifld),.true.)
971  ENDIF
972 C
973  DO 300 iel=1,nel
974  DO 300 icr=1,ncrnr
975  ix = mcrrst(1,icr)
976  iy = mcrrst(2,icr)
977  iz = mcrrst(3,icr)
978  hmv = abs( hvmask(ix,iy,iz,iel) )
979  IF (hmv .LT. 1.9 .OR. hmv .GT. 3.1) GOTO 300
980  idiff = 0
981  DO 330 ii=1,ldim
982  ifc = mfccr(ii,icr)
983  hmf = abs( hfmask(1,1,ifc,iel) )
984  IF (abs(hmv - hmf) .GT. epsa) idiff=idiff + 1
985  330 CONTINUE
986  IF (ldim.EQ.3) THEN
987  DO 360 ii=1,ldim
988  isd = megcr(ii,icr)
989  ixs = midrst(1,isd)
990  iys = midrst(2,isd)
991  izs = midrst(3,isd)
992  hms = abs( hvmask(ixs,iys,izs,iel) )
993  IF (abs(hmv - hms) .GT. epsa) idiff=idiff + 1
994  360 CONTINUE
995  ENDIF
996  IF ( (ldim.EQ.2 .AND. idiff.EQ.2) .OR.
997  $ (ldim.EQ.3 .AND. idiff.EQ.6) ) THEN
998  iflmsc(ifld) = .true.
999  ifmscr(icr,iel,ifld) = .true.
1000  ENDIF
1001  300 CONTINUE
1002  CALL gllog(iflmsc(ifld),.true.)
1003 C
1004  RETURN
1005  END
1006 C-----------------------------------------------------------------------
1007  SUBROUTINE setcsys (HVMASK,HFMASK,NEL)
1008 C
1009  include 'SIZE'
1010  include 'GEOM'
1011  include 'TSTEP'
1012  dimension hvmask(lx1,ly1,lz1,1)
1013  $ , hfmask(lx1,lz1,6,1)
1014 C
1015  nface = 2*ldim
1016  ntot1 = lx1*ly1*lz1*nel
1017 C
1018  CALL rzero3 (vnx,vny,vnz,ntot1)
1019  CALL rzero3 (v1x,v1y,v1z,ntot1)
1020  CALL rzero3 (v2x,v2y,v2z,ntot1)
1021 C
1022  DO 10 iel=1,nel
1023  DO 10 ifc=1,nface
1024  hmf = abs( hfmask(1,1,ifc,iel) )
1025  IF (hmf .GT. 1.9 .AND. hmf .LT. 3.1)
1026  $ CALL facexv(unx(1,1,ifc,iel),uny(1,1,ifc,iel),unz(1,1,ifc,iel),
1027  $ vnx(1,1,1,iel),vny(1,1,1,iel),vnz(1,1,1,iel),ifc,1)
1028  10 CONTINUE
1029 C
1030  IF (ldim.EQ.2) THEN
1031  CALL comavn2 (hvmask,hfmask,nel)
1032  ELSE
1033  CALL comavn3 (hvmask,hfmask,nel)
1034  ENDIF
1035 C
1036  RETURN
1037  END
1038 C-----------------------------------------------------------------------
1039  SUBROUTINE comavn2 (HVMASK,HFMASK,NEL)
1040 C
1041  include 'SIZE'
1042  include 'GEOM'
1043  common /indxfc/ mcrfc(4,6)
1044  $ , mfccr(3,8)
1045  $ , megcr(3,8)
1046  $ , mfceg(2,12)
1047  $ , mcreg(2,12)
1048  $ , mcrrst(3,8)
1049  $ , midrst(3,12)
1050  $ , mcrind(8)
1051  $ , medind(2,4)
1052  $ , ntefc(2,12)
1053  $ , ntcrf(2,3)
1054 C
1055  dimension hvmask(lx1,ly1,lz1,1)
1056  $ , hfmask(lx1,lz1,6,1)
1057 C
1058  ntot1 = lx1*ly1*lz1*nel
1059  nface = 2*ldim
1060  ncrnr = 2**ldim
1061  epsa = 1.0e-06
1062  CALL rzero (vnz,ntot1)
1063 C
1064  iz = 1
1065  DO 100 iel=1,nel
1066  DO 100 icr=1,ncrnr
1067  ix = mcrrst(1,icr)
1068  iy = mcrrst(2,icr)
1069  hmv = abs( hvmask(ix,iy,iz,iel) )
1070  IF (hmv .LT. 1.9 .OR. hmv .GT. 3.1) GOTO 100
1071  vnx(ix,iy,iz,iel) = 0.0
1072  vny(ix,iy,iz,iel) = 0.0
1073  100 CONTINUE
1074 C
1075  DO 200 iel=1,nel
1076  DO 200 icr=1,ncrnr
1077  IF (ifnskp(icr,iel)) GOTO 200
1078  ix = mcrrst(1,icr)
1079  iy = mcrrst(2,icr)
1080  hmv = abs( hvmask(ix,iy,iz,iel) )
1081  IF (hmv .LT. 1.9 .OR. hmv .GT. 3.1) GOTO 200
1082  DO 220 ii=1,2
1083  ifc = mfccr(ii,icr)
1084  hmf = abs( hfmask(1,1,ifc,iel) )
1085  IF (abs(hmv - hmf) .LT. epsa) THEN
1086  ir = mcrrst(ii,icr)
1087  vnx(ix,iy,iz,iel)=vnx(ix,iy,iz,iel) + unx(ir,iz,ifc,iel)
1088  vny(ix,iy,iz,iel)=vny(ix,iy,iz,iel) + uny(ir,iz,ifc,iel)
1089  ENDIF
1090  220 CONTINUE
1091  200 CONTINUE
1092 C
1093  CALL dssum (vnx,lx1,ly1,lz1)
1094  CALL dssum (vny,lx1,ly1,lz1)
1095  CALL unitvec (vnx,vny,vnz,ntot1)
1096 C
1097  CALL copy (v1y,vnx,ntot1)
1098  CALL copy (v1x,vny,ntot1)
1099  CALL chsign (v1x,ntot1)
1100 C
1101  RETURN
1102  END
1103 C-----------------------------------------------------------------------
1104  SUBROUTINE comavn3 (HVMASK,HFMASK,NEL)
1105 C
1106  include 'SIZE'
1107  include 'GEOM'
1108  common /scrcg/ vnmag(lx1,ly1,lz1,lelt)
1109  common /indxfc/ mcrfc(4,6)
1110  $ , mfccr(3,8)
1111  $ , megcr(3,8)
1112  $ , mfceg(2,12)
1113  $ , mcreg(2,12)
1114  $ , mcrrst(3,8)
1115  $ , midrst(3,12)
1116  $ , mcrind(8)
1117  $ , medind(2,4)
1118  $ , ntefc(2,12)
1119  $ , ntcrf(2,3)
1120 C
1121  dimension hvmask(lx1,ly1,lz1,1)
1122  $ , hfmask(lx1,lz1,6,1)
1123 C
1124  ntot1 = lx1*ly1*lz1*nel
1125  nface = 2*ldim
1126  ncrnr = 2**ldim
1127  nedge = 12
1128  nmid =(lx1 + 1)/2
1129  nxm1 = lx1 - 1
1130  epsa = 1.0e-06
1131  epsn = 1.0e-03
1132 C
1133  DO 100 iel=1,nel
1134  DO 100 isd=1,nedge
1135  ix = midrst(1,isd)
1136  iy = midrst(2,isd)
1137  iz = midrst(3,isd)
1138  hmv = abs( hvmask(ix,iy,iz,iel) )
1139  IF (hmv .LT. 1.9 .OR. hmv .GT. 3.1) GOTO 100
1140  CALL edgindv (lv1,lv2,lvskip,isd)
1141  DO 120 i=2,nxm1
1142  lv = lv1 + (i-1)*lvskip
1143  vnx(lv,1,1,iel) = 0.0
1144  vny(lv,1,1,iel) = 0.0
1145  vnz(lv,1,1,iel) = 0.0
1146  120 CONTINUE
1147  100 CONTINUE
1148 C
1149  DO 150 iel=1,nel
1150  DO 150 icr=1,ncrnr
1151  ix = mcrrst(1,icr)
1152  iy = mcrrst(2,icr)
1153  iz = mcrrst(3,icr)
1154  hmv = abs( hvmask(ix,iy,iz,iel) )
1155  IF (hmv .LT. 1.9 .OR. hmv .GT. 3.1) GOTO 150
1156  vnx(ix,iy,iz,iel) = 0.0
1157  vny(ix,iy,iz,iel) = 0.0
1158  vnz(ix,iy,iz,iel) = 0.0
1159  150 CONTINUE
1160 C
1161 C (1) All Edges
1162 C
1163  DO 200 iel=1,nel
1164  DO 200 isd=1,nedge
1165  ix = midrst(1,isd)
1166  iy = midrst(2,isd)
1167  iz = midrst(3,isd)
1168  hmv = abs( hvmask(ix,iy,iz,iel) )
1169  IF (hmv .LT. 1.9 .OR. hmv .GT. 3.1) GOTO 200
1170  DO 220 ii=1,2
1171  ifc = mfceg(ii,isd)
1172  hmf = abs( hfmask(1,1,ifc,iel) )
1173  IF (abs(hmv - hmf) .LT. epsa) THEN
1174  CALL edgindv (lv1,lv2,lvskip,isd)
1175  CALL edgindf (lf1,lf2,lfskip,isd,ii)
1176  DO 240 i=2,nxm1
1177  lv = lv1 + (i-1)*lvskip
1178  lf = lf1 + (i-1)*lfskip
1179  vnx(lv,1,1,iel)=vnx(lv,1,1,iel)+unx(lf,1,ifc,iel)
1180  vny(lv,1,1,iel)=vny(lv,1,1,iel)+uny(lf,1,ifc,iel)
1181  vnz(lv,1,1,iel)=vnz(lv,1,1,iel)+unz(lf,1,ifc,iel)
1182  240 CONTINUE
1183  ENDIF
1184  220 CONTINUE
1185  200 CONTINUE
1186 C
1187 C (2) All Corners
1188 C
1189  DO 300 iel=1,nel
1190  DO 300 icr=1,ncrnr
1191  ix = mcrrst(1,icr)
1192  iy = mcrrst(2,icr)
1193  iz = mcrrst(3,icr)
1194  hmv = abs( hvmask(ix,iy,iz,iel) )
1195  IF (hmv .LT. 1.9 .OR. hmv .GT. 3.1) GOTO 300
1196  DO 320 ii=1,3
1197  ifc = mfccr(ii,icr)
1198  hmf = abs( hfmask(1,1,ifc,iel) )
1199  IF (abs(hmv - hmf) .LT. epsa) THEN
1200  ira = ntcrf(1,ii)
1201  isa = ntcrf(2,ii)
1202  ir = mcrrst(ira,icr)
1203  is = mcrrst(isa,icr)
1204  vnx(ix,iy,iz,iel)=vnx(ix,iy,iz,iel)+unx(ir,is,ifc,iel)
1205  vny(ix,iy,iz,iel)=vny(ix,iy,iz,iel)+uny(ir,is,ifc,iel)
1206  vnz(ix,iy,iz,iel)=vnz(ix,iy,iz,iel)+unz(ir,is,ifc,iel)
1207  ENDIF
1208  320 CONTINUE
1209  300 CONTINUE
1210 C
1211  CALL dssum (vnx,lx1,ly1,lz1)
1212  CALL dssum (vny,lx1,ly1,lz1)
1213  CALL dssum (vnz,lx1,ly1,lz1)
1214  CALL unitvec (vnx,vny,vnz,ntot1)
1215  CALL vdot3 (vnmag,vnx,vny,vnz,vnx,vny,vnz,ntot1)
1216 C
1217  DO 500 iel=1,nel
1218  DO 500 ifc=1,nface
1219  hmf = abs( hfmask(1,1,ifc,iel) )
1220  IF (hmf .LT. 1.9 .OR. hmf .GT. 3.1) GOTO 500
1221  CALL facind2 (js1,jf1,jskip1,js2,jf2,jskip2,ifc)
1222  DO 520 j2=js2,jf2,jskip2
1223  DO 520 j1=js1,jf1,jskip1
1224  IF (vnmag(j1,j2,1,iel) .LT. epsa) GOTO 520
1225  vlzdif = abs(vnz(j1,j2,1,iel)) - 1.0
1226  IF (abs(vlzdif) .LT. epsn) THEN
1227  v1x(j1,j2,1,iel) = 1.0
1228  v1y(j1,j2,1,iel) = 0.0
1229  v1z(j1,j2,1,iel) = 0.0
1230  ELSE
1231  ssn = sqrt(vnx(j1,j2,1,iel)**2 + vny(j1,j2,1,iel)**2)
1232  v1x(j1,j2,1,iel) = -vny(j1,j2,1,iel) / ssn
1233  v1y(j1,j2,1,iel) = vnx(j1,j2,1,iel) / ssn
1234  v1z(j1,j2,1,iel) = 0.0
1235  ENDIF
1236  520 CONTINUE
1237  500 CONTINUE
1238 C
1239  CALL dssum (v1x,lx1,ly1,lz1)
1240  CALL dssum (v1y,lx1,ly1,lz1)
1241  CALL dssum (v1z,lx1,ly1,lz1)
1242  CALL unitvec (v1x,v1y,v1z,ntot1)
1243 C
1244  CALL vcross (v2x,v2y,v2z,vnx,vny,vnz,v1x,v1y,v1z,ntot1)
1245 C
1246  RETURN
1247  END
1248 C-----------------------------------------------------------------------
1249  SUBROUTINE fixwmsk (W2MASK,W3MASK,HVMASK,HFMASK,NEL)
1250 C
1251  include 'SIZE'
1252  dimension hvmask(lx1,ly1,lz1,1)
1253  $ , hfmask(lx1,lz1,6,1)
1254 C
1255  dimension w2mask(lx1,ly1,lz1,1)
1256  $ , w3mask(lx1,ly1,lz1,1)
1257 C
1258  IF (ldim.EQ.2) THEN
1259  CALL fxwms2 (w2mask,hvmask,hfmask,nel)
1260  ELSE
1261  CALL fxwms3 (w2mask,w3mask,hvmask,hfmask,nel)
1262  ENDIF
1263 C
1264  CALL dsop(w2mask,'MUL',lx1,ly1,lz1)
1265  IF (ldim.EQ.3) CALL dsop(w3mask,'MUL',lx1,ly1,lz1)
1266 C
1267  RETURN
1268  END
1269 C-----------------------------------------------------------------------
1270  SUBROUTINE fxwms2 (W2MASK,HVMASK,HFMASK,NEL)
1271 C
1272  include 'SIZE'
1273  include 'GEOM'
1274  common /indxfc/ mcrfc(4,6)
1275  $ , mfccr(3,8)
1276  $ , megcr(3,8)
1277  $ , mfceg(2,12)
1278  $ , mcreg(2,12)
1279  $ , mcrrst(3,8)
1280  $ , midrst(3,12)
1281  $ , mcrind(8)
1282  $ , medind(2,4)
1283  $ , ntefc(2,12)
1284  $ , ntcrf(2,3)
1285 C
1286  dimension w2mask(lx1,ly1,lz1,1)
1287  $ , hvmask(lx1,ly1,lz1,1)
1288  $ , hfmask(lx1,lz1,6,1)
1289 C
1290  ncrnr = 2**ldim
1291  epsa = 1.0e-06
1292 C
1293  iz = 1
1294  DO 100 iel=1,nel
1295  DO 100 icr=1,ncrnr
1296  ix = mcrrst(1,icr)
1297  iy = mcrrst(2,icr)
1298  hmv = abs( hvmask(ix,iy,iz,iel) )
1299  IF (hmv .LT. 1.9 .OR. hmv .GT. 2.1) GOTO 100
1300  DO 120 ii=1,2
1301  ifc = mfccr(ii,icr)
1302  hmf = abs( hfmask(1,1,ifc,iel) )
1303  IF (abs(hmv - hmf) .LT. epsa) THEN
1304  ir = mcrrst(ii,icr)
1305  dot = vnx(ix,iy,iz,iel)*unx(ir,iz,ifc,iel) +
1306  $ vny(ix,iy,iz,iel)*uny(ir,iz,ifc,iel)
1307  IF (dot .LT. 0.99) THEN
1308  w2mask(ix,iy,iz,iel) = 0.0
1309  GOTO 100
1310  ENDIF
1311  ENDIF
1312  120 CONTINUE
1313  100 CONTINUE
1314 C
1315  RETURN
1316  END
1317 C-----------------------------------------------------------------------
1318  SUBROUTINE fxwms3 (W2MASK,W3MASK,HVMASK,HFMASK,NEL)
1319 C
1320  include 'SIZE'
1321  include 'GEOM'
1322  common /indxfc/ mcrfc(4,6)
1323  $ , mfccr(3,8)
1324  $ , megcr(3,8)
1325  $ , mfceg(2,12)
1326  $ , mcreg(2,12)
1327  $ , mcrrst(3,8)
1328  $ , midrst(3,12)
1329  $ , mcrind(8)
1330  $ , medind(2,4)
1331  $ , ntefc(2,12)
1332  $ , ntcrf(2,3)
1333 C
1334  dimension w2mask(lx1,ly1,lz1,1)
1335  $ , w3mask(lx1,ly1,lz1,1)
1336  $ , hvmask(lx1,ly1,lz1,1)
1337  $ , hfmask(lx1,lz1,6,1)
1338 C
1339  ncrnr = 2**ldim
1340  nedge = 12
1341  nmid = (lx1 + 1)/2
1342  epsa = 1.0e-06
1343 C
1344  DO 100 iel=1,nel
1345  DO 100 isd=1,12
1346  ix = midrst(1,isd)
1347  iy = midrst(2,isd)
1348  iz = midrst(3,isd)
1349  hmv = abs( hvmask(ix,iy,iz,iel) )
1350  IF (hmv .LT. 1.9 .OR. hmv .GT. 2.1) GOTO 100
1351  DO 120 ii=1,2
1352  ifc = mfceg(ii,isd)
1353  hmf = abs( hfmask(1,1,ifc,iel) )
1354  IF (abs(hmv - hmf) .LT. epsa) THEN
1355  CALL edgindf (lf1,lf2,lfskip,isd,ii)
1356  lf = lf1 + (nmid-1)*lfskip
1357  dot = vnx(ix,iy,iz,iel)*unx(lf,1,ifc,iel) +
1358  $ vny(ix,iy,iz,iel)*uny(lf,1,ifc,iel) +
1359  $ vnz(ix,iy,iz,iel)*unz(lf,1,ifc,iel)
1360  IF (dot .LT. 0.99) THEN
1361  CALL edgindv (lv1,lv2,lvskip,isd)
1362  DO 140 lv=lv1,lv2,lvskip
1363  w3mask(lv,1,1,iel) = 0.0
1364  140 CONTINUE
1365  ENDIF
1366  ENDIF
1367  120 CONTINUE
1368  100 CONTINUE
1369 C
1370 C All Corners
1371 C
1372  DO 300 iel=1,nel
1373  DO 300 icr=1,ncrnr
1374  ix = mcrrst(1,icr)
1375  iy = mcrrst(2,icr)
1376  iz = mcrrst(3,icr)
1377  hmv = abs( hvmask(ix,iy,iz,iel) )
1378  IF (hmv .LT. 1.9 .OR. hmv .GT. 2.1) GOTO 300
1379  DO 320 ii=1,3
1380  ifc = mfccr(ii,icr)
1381  hmf = abs( hfmask(1,1,ifc,iel) )
1382  IF (abs(hmv - hmf) .LT. epsa) THEN
1383  ira = ntcrf(1,ii)
1384  isa = ntcrf(2,ii)
1385  ir = mcrrst(ira,icr)
1386  is = mcrrst(isa,icr)
1387  dot = vnx(ix,iy,iz,iel)*unx(ir,is,ifc,iel) +
1388  $ vny(ix,iy,iz,iel)*uny(ir,is,ifc,iel) +
1389  $ vnz(ix,iy,iz,iel)*unz(ir,is,ifc,iel)
1390  ENDIF
1391  IF (dot .LT. 0.99) THEN
1392  w2mask(ix,iy,iz,iel) = 0.0
1393  GOTO 300
1394  ENDIF
1395  320 CONTINUE
1396  300 CONTINUE
1397 C
1398  RETURN
1399  END
1400 C-----------------------------------------------------------------------
1401  SUBROUTINE setcdat
1402 C
1403  include 'SIZE'
1404  common /indxfc/ mcrfc(4,6)
1405  $ , mfccr(3,8)
1406  $ , megcr(3,8)
1407  $ , mfceg(2,12)
1408  $ , mcreg(2,12)
1409  $ , mcrrst(3,8)
1410  $ , midrst(3,12)
1411  $ , mcrind(8)
1412  $ , medind(2,4)
1413  $ , ntefc(2,12)
1414  $ , ntcrf(2,3)
1415 C
1416  nmid = (lx1 +1)/2
1417 C
1418 C Corners on faces
1419 C
1420  mcrfc(1,1) = 1
1421  mcrfc(2,1) = 2
1422  mcrfc(3,1) = 6
1423  mcrfc(4,1) = 5
1424  mcrfc(1,2) = 2
1425  mcrfc(2,2) = 3
1426  mcrfc(3,2) = 7
1427  mcrfc(4,2) = 6
1428  mcrfc(1,3) = 3
1429  mcrfc(2,3) = 4
1430  mcrfc(3,3) = 8
1431  mcrfc(4,3) = 7
1432  mcrfc(1,4) = 4
1433  mcrfc(2,4) = 1
1434  mcrfc(3,4) = 5
1435  mcrfc(4,4) = 8
1436  mcrfc(1,5) = 4
1437  mcrfc(2,5) = 3
1438  mcrfc(3,5) = 2
1439  mcrfc(4,5) = 1
1440  mcrfc(1,6) = 5
1441  mcrfc(2,6) = 6
1442  mcrfc(3,6) = 7
1443  mcrfc(4,6) = 8
1444 C
1445 C Faces at corners
1446 C
1447  mfccr(1,1) = 4
1448  mfccr(2,1) = 1
1449  mfccr(3,1) = 5
1450  mfccr(1,2) = 1
1451  mfccr(2,2) = 2
1452  mfccr(3,2) = 5
1453  mfccr(1,3) = 2
1454  mfccr(2,3) = 3
1455  mfccr(3,3) = 5
1456  mfccr(1,4) = 3
1457  mfccr(2,4) = 4
1458  mfccr(3,4) = 5
1459  mfccr(1,5) = 4
1460  mfccr(2,5) = 1
1461  mfccr(3,5) = 6
1462  mfccr(1,6) = 1
1463  mfccr(2,6) = 2
1464  mfccr(3,6) = 6
1465  mfccr(1,7) = 2
1466  mfccr(2,7) = 3
1467  mfccr(3,7) = 6
1468  mfccr(1,8) = 3
1469  mfccr(2,8) = 4
1470  mfccr(3,8) = 6
1471 C
1472 C Edges at corners
1473 C
1474  megcr(1,1) = 4
1475  megcr(2,1) = 1
1476  megcr(3,1) = 9
1477  megcr(1,2) = 1
1478  megcr(2,2) = 2
1479  megcr(3,2) = 10
1480  megcr(1,3) = 2
1481  megcr(2,3) = 3
1482  megcr(3,3) = 11
1483  megcr(1,4) = 3
1484  megcr(2,4) = 4
1485  megcr(3,4) = 12
1486  megcr(1,5) = 8
1487  megcr(2,5) = 5
1488  megcr(3,5) = 9
1489  megcr(1,6) = 5
1490  megcr(2,6) = 6
1491  megcr(3,6) = 10
1492  megcr(1,7) = 6
1493  megcr(2,7) = 7
1494  megcr(3,7) = 11
1495  megcr(1,8) = 7
1496  megcr(2,8) = 8
1497  megcr(3,8) = 12
1498 C
1499 C Faces on edges
1500 C
1501  mfceg(1,1) = 1
1502  mfceg(2,1) = 5
1503  mfceg(1,2) = 2
1504  mfceg(2,2) = 5
1505  mfceg(1,3) = 3
1506  mfceg(2,3) = 5
1507  mfceg(1,4) = 4
1508  mfceg(2,4) = 5
1509  mfceg(1,5) = 1
1510  mfceg(2,5) = 6
1511  mfceg(1,6) = 2
1512  mfceg(2,6) = 6
1513  mfceg(1,7) = 3
1514  mfceg(2,7) = 6
1515  mfceg(1,8) = 4
1516  mfceg(2,8) = 6
1517  mfceg(1,9) = 4
1518  mfceg(2,9) = 1
1519  mfceg(1,10) = 1
1520  mfceg(2,10) = 2
1521  mfceg(1,11) = 2
1522  mfceg(2,11) = 3
1523  mfceg(1,12) = 3
1524  mfceg(2,12) = 4
1525 C
1526 C Corners at edges
1527 C
1528  mcreg(1,1) = 1
1529  mcreg(2,1) = 2
1530  mcreg(1,2) = 2
1531  mcreg(2,2) = 3
1532  mcreg(1,3) = 4
1533  mcreg(2,3) = 3
1534  mcreg(1,4) = 1
1535  mcreg(2,4) = 4
1536  mcreg(1,5) = 5
1537  mcreg(2,5) = 6
1538  mcreg(1,6) = 6
1539  mcreg(2,6) = 7
1540  mcreg(1,7) = 8
1541  mcreg(2,7) = 7
1542  mcreg(1,8) = 5
1543  mcreg(2,8) = 8
1544  mcreg(1,9) = 1
1545  mcreg(2,9) = 5
1546  mcreg(1,10) = 2
1547  mcreg(2,10) = 6
1548  mcreg(1,11) = 3
1549  mcreg(2,11) = 7
1550  mcreg(1,12) = 4
1551  mcreg(2,12) = 8
1552 C
1553 C Corner indices (Vol array)
1554 C
1555  mcrrst(1,1) = 1
1556  mcrrst(2,1) = 1
1557  mcrrst(3,1) = 1
1558  mcrrst(1,2) = lx1
1559  mcrrst(2,2) = 1
1560  mcrrst(3,2) = 1
1561  mcrrst(1,3) = lx1
1562  mcrrst(2,3) = lx1
1563  mcrrst(3,3) = 1
1564  mcrrst(1,4) = 1
1565  mcrrst(2,4) = lx1
1566  mcrrst(3,4) = 1
1567  mcrrst(1,5) = 1
1568  mcrrst(2,5) = 1
1569  mcrrst(3,5) = lx1
1570  mcrrst(1,6) = lx1
1571  mcrrst(2,6) = 1
1572  mcrrst(3,6) = lx1
1573  mcrrst(1,7) = lx1
1574  mcrrst(2,7) = lx1
1575  mcrrst(3,7) = lx1
1576  mcrrst(1,8) = 1
1577  mcrrst(2,8) = lx1
1578  mcrrst(3,8) = lx1
1579 C
1580 C Mid-edge indcies (Vol array)
1581 C
1582  midrst(1,1) = nmid
1583  midrst(1,2) = lx1
1584  midrst(1,3) = nmid
1585  midrst(1,4) = 1
1586  midrst(1,5) = nmid
1587  midrst(1,6) = lx1
1588  midrst(1,7) = nmid
1589  midrst(1,8) = 1
1590  midrst(1,9) = 1
1591  midrst(1,10) = lx1
1592  midrst(1,11) = lx1
1593  midrst(1,12) = 1
1594  midrst(2,1) = 1
1595  midrst(2,2) = nmid
1596  midrst(2,3) = lx1
1597  midrst(2,4) = nmid
1598  midrst(2,5) = 1
1599  midrst(2,6) = nmid
1600  midrst(2,7) = lx1
1601  midrst(2,8) = nmid
1602  midrst(2,9) = 1
1603  midrst(2,10) = 1
1604  midrst(2,11) = lx1
1605  midrst(2,12) = lx1
1606  midrst(3,1) = 1
1607  midrst(3,2) = 1
1608  midrst(3,3) = 1
1609  midrst(3,4) = 1
1610  midrst(3,5) = lx1
1611  midrst(3,6) = lx1
1612  midrst(3,7) = lx1
1613  midrst(3,8) = lx1
1614  midrst(3,9) = nmid
1615  midrst(3,10) = nmid
1616  midrst(3,11) = nmid
1617  midrst(3,12) = nmid
1618 C
1619 C 1-D corners indices (Vol array)
1620 C
1621  mcrind(1) = 1
1622  mcrind(2) = lx1
1623  mcrind(3) = lx1**2
1624  mcrind(7) = lx1**3
1625  mcrind(4) = mcrind(3) - lx1 + 1
1626  mcrind(5) = mcrind(7) - mcrind(3) + 1
1627  mcrind(6) = mcrind(5) + lx1 - 1
1628  mcrind(8) = mcrind(7) - lx1 + 1
1629 C
1630 C 1-D edge indices (Face array)
1631 C
1632  medind(1,1) = 1
1633  medind(2,1) = lx1
1634  medind(1,2) = lx1**2 - lx1 + 1
1635  medind(2,2) = lx1**2
1636  medind(1,3) = 1
1637  medind(2,3) = medind(1,2)
1638  medind(1,4) = lx1
1639  medind(2,4) = lx1**2
1640 C
1641 C 1-D edge index type (Face array)
1642 C
1643  ntefc(1,1) = 1
1644  ntefc(2,1) = 1
1645  ntefc(1,2) = 1
1646  ntefc(2,2) = 4
1647  ntefc(1,3) = 1
1648  ntefc(2,3) = 2
1649  ntefc(1,4) = 1
1650  ntefc(2,4) = 3
1651  ntefc(1,5) = 2
1652  ntefc(2,5) = 1
1653  ntefc(1,6) = 2
1654  ntefc(2,6) = 4
1655  ntefc(1,7) = 2
1656  ntefc(2,7) = 2
1657  ntefc(1,8) = 2
1658  ntefc(2,8) = 3
1659  ntefc(1,9) = 3
1660  ntefc(2,9) = 3
1661  ntefc(1,10) = 4
1662  ntefc(2,10) = 3
1663  ntefc(1,11) = 4
1664  ntefc(2,11) = 4
1665  ntefc(1,12) = 3
1666  ntefc(2,12) = 4
1667 C
1668 C Corner index address on face in MCRRST
1669 C
1670  ntcrf(1,1) = 1
1671  ntcrf(2,1) = 3
1672  ntcrf(1,2) = 2
1673  ntcrf(2,2) = 3
1674  ntcrf(1,3) = 1
1675  ntcrf(2,3) = 2
1676 C
1677  RETURN
1678  END
1679 C-----------------------------------------------------------------------
1680  SUBROUTINE edgindf (LF1,LF2,LFSKIP,ISD,IFCN)
1681 C
1682  include 'SIZE'
1683  common /indxfc/ mcrfc(4,6)
1684  $ , mfccr(3,8)
1685  $ , megcr(3,8)
1686  $ , mfceg(2,12)
1687  $ , mcreg(2,12)
1688  $ , mcrrst(3,8)
1689  $ , midrst(3,12)
1690  $ , mcrind(8)
1691  $ , medind(2,4)
1692  $ , ntefc(2,12)
1693  $ , ntcrf(2,3)
1694 C
1695  ityp = ntefc(ifcn,isd)
1696 C
1697  lf1 = medind(1,ityp)
1698  lf2 = medind(2,ityp)
1699 C
1700  lfskip = 1
1701  IF (ityp .GE. 3) lfskip = lx1
1702 C
1703  RETURN
1704  END
1705 C-----------------------------------------------------------------------
1706  SUBROUTINE edgindv (LV1,LV2,LVSKIP,ISD)
1707 C
1708  include 'SIZE'
1709  common /indxfc/ mcrfc(4,6)
1710  $ , mfccr(3,8)
1711  $ , megcr(3,8)
1712  $ , mfceg(2,12)
1713  $ , mcreg(2,12)
1714  $ , mcrrst(3,8)
1715  $ , midrst(3,12)
1716  $ , mcrind(8)
1717  $ , medind(2,4)
1718  $ , ntefc(2,12)
1719  $ , ntcrf(2,3)
1720 C
1721  iodd = isd - isd/2*2
1722  icr1 = mcreg(1,isd)
1723  icr2 = mcreg(2,isd)
1724 C
1725  lv1 = mcrind(icr1)
1726  lv2 = mcrind(icr2)
1727 C
1728  IF (isd .GE. 9) THEN
1729  lvskip = lx1**2
1730  ELSE
1731  IF (iodd.EQ.0) THEN
1732  lvskip = lx1
1733  ELSE
1734  lvskip = 1
1735  ENDIF
1736  ENDIF
1737 C
1738  RETURN
1739  END
1740 C-----------------------------------------------------------------------
1741  SUBROUTINE setcdof
1742 C
1743  include 'SIZE'
1744  include 'INPUT'
1745 C
1746  nface = 2*ldim
1747 C
1748  DO 100 iel=1,nelt
1749  DO 100 ifc=1,nface
1750  cdof(ifc,iel)=cbc(ifc,iel,0)(1:1)
1751  100 CONTINUE
1752 C
1753  RETURN
1754  END
1755 C-----------------------------------------------------------------------
1756  SUBROUTINE amask (VB1,VB2,VB3,V1,V2,V3,NEL)
1757 C
1758  include 'SIZE'
1759  include 'GEOM'
1760  include 'INPUT'
1761  include 'SOLN'
1762  include 'TSTEP'
1763  common /scrsf/ a1mask(lx1,ly1,lz1,lelt)
1764  $ , a2mask(lx1,ly1,lz1,lelt)
1765  $ , a3mask(lx1,ly1,lz1,lelt)
1766  common /ctmp0/ wa(lx1,ly1,lz1,lelt)
1767 C
1768  dimension vb1(lx1,ly1,lz1,1)
1769  $ , vb2(lx1,ly1,lz1,1)
1770  $ , vb3(lx1,ly1,lz1,1)
1771  $ , v1(lx1,ly1,lz1,1)
1772  $ , v2(lx1,ly1,lz1,1)
1773  $ , v3(lx1,ly1,lz1,1)
1774 C
1775  ntot1 = lx1*ly1*lz1*nel
1776  CALL rone (wa,ntot1)
1777  CALL copy (vb1,v1,ntot1)
1778  CALL copy (vb2,v2,ntot1)
1779  IF (ldim.EQ.3) CALL copy (vb3,v3,ntot1)
1780 C
1781  IF (ifield.EQ.1) THEN
1782  CALL sub3 (a1mask,wa,v1mask,ntot1)
1783  CALL sub3 (a2mask,wa,v2mask,ntot1)
1784  IF (ldim.EQ.3) CALL sub3 (a3mask,wa,v3mask,ntot1)
1785  ELSEIF (ifield.EQ.ifldmhd) THEN
1786  CALL sub3 (a1mask,wa,b1mask,ntot1)
1787  CALL sub3 (a2mask,wa,b2mask,ntot1)
1788  IF (ldim.EQ.3) CALL sub3 (a3mask,wa,b3mask,ntot1)
1789  ELSE
1790  CALL sub3 (a1mask,wa,w1mask,ntot1)
1791  CALL sub3 (a2mask,wa,w2mask,ntot1)
1792  IF (ldim.EQ.3) CALL sub3 (a3mask,wa,w3mask,ntot1)
1793  ENDIF
1794 C
1795  CALL qmask (vb1,vb2,vb3,a1mask,a2mask,a3mask,nel)
1796 C
1797  RETURN
1798  END
1799 C-----------------------------------------------------------------------
1800  SUBROUTINE rmask (R1,R2,R3,NEL)
1801 C
1802  include 'SIZE'
1803  include 'INPUT'
1804  include 'SOLN'
1805  include 'TSTEP'
1806  include 'MVGEOM'
1807 C
1808  dimension r1(lx1,ly1,lz1,1)
1809  $ , r2(lx1,ly1,lz1,1)
1810  $ , r3(lx1,ly1,lz1,1)
1811 C
1812 c if (ifsplit) then
1813 c call opmask(r1,r2,r3)
1814 c return
1815 c endif
1816 
1817 c call outfldro (v1mask,'v1mask rmk',0)
1818 c call outfldro (v2mask,'v2mask rmk',1)
1819 
1820  IF (ifield.EQ.1) THEN
1821  CALL qmask (r1,r2,r3,v1mask,v2mask,v3mask,nel)
1822  ELSEIF (ifield.eq.ifldmhd) then
1823  CALL qmask (r1,r2,r3,b1mask,b2mask,b3mask,nel)
1824  ELSE
1825  CALL qmask (r1,r2,r3,w1mask,w2mask,w3mask,nel)
1826  ENDIF
1827 
1828 c call outfldro (r1,'r1 rmask',0)
1829 c call outfldro (r2,'r2 rmask',1)
1830 c call exitti ('quit in rmask$,',nel)
1831 
1832  RETURN
1833  END
1834 C-----------------------------------------------------------------------
1835  SUBROUTINE qmask (R1,R2,R3,R1MASK,R2MASK,R3MASK,NEL)
1836 C
1837  include 'SIZE'
1838  include 'GEOM'
1839  include 'TSTEP'
1840  common /ctmp1/ s1(lx1,ly1,lz1,lelt)
1841  $ , s2(lx1,ly1,lz1,lelt)
1842  $ , s3(lx1,ly1,lz1,lelt)
1843 C
1844  dimension r1(lx1,ly1,lz1,1)
1845  $ , r2(lx1,ly1,lz1,1)
1846  $ , r3(lx1,ly1,lz1,1)
1847  $ , r1mask(lx1,ly1,lz1,1)
1848  $ , r2mask(lx1,ly1,lz1,1)
1849  $ , r3mask(lx1,ly1,lz1,1)
1850 C
1851  ntot1 = lx1*ly1*lz1*nel
1852 C
1853 C (0) Collocate Volume Mask
1854 C
1855  CALL copy (s1,r1,ntot1)
1856  CALL copy (s2,r2,ntot1)
1857  CALL col2 (r1,r1mask,ntot1)
1858  CALL col2 (r2,r2mask,ntot1)
1859  IF (ldim.EQ.3) THEN
1860  CALL copy (s3,r3,ntot1)
1861  CALL col2 (r3,r3mask,ntot1)
1862  ENDIF
1863 C
1864 C (1) Face Mask
1865 C
1866  IF (iflmsf(ifield)) THEN
1867  IF (ldim.EQ.2) THEN
1868  CALL fcmsk2 (r1,r2,s1,s2,r1mask,r2mask,nel)
1869  ELSE
1870  CALL fcmsk3 (r1,r2,r3,s1,s2,s3,r1mask,r2mask,r3mask,nel)
1871  ENDIF
1872  ENDIF
1873 C
1874 C (2) Edge Mask (3-D only)
1875 C
1876  IF (ldim.EQ.3 .AND. iflmse(ifield))
1877  $ CALL egmask (r1,r2,r3,s1,s2,s3,r1mask,r2mask,r3mask,nel)
1878 C
1879 C (3) Corner Mask
1880 C
1881  IF (iflmsc(ifield)) THEN
1882  IF (ldim.EQ.2) THEN
1883  CALL crmsk2 (r1,r2,s1,s2,r1mask,r2mask,nel)
1884  ELSE
1885  CALL crmsk3 (r1,r2,r3,s1,s2,s3,r1mask,r2mask,r3mask,nel)
1886  ENDIF
1887  ENDIF
1888 C
1889  RETURN
1890  END
1891 C-----------------------------------------------------------------------
1892  SUBROUTINE fcmsk2 (R1,R2,S1,S2,R1MASK,R2MASK,NEL)
1893 C
1894  include 'SIZE'
1895  include 'GEOM'
1896  include 'TSTEP'
1897  dimension r1(lx1,ly1,lz1,1)
1898  $ , r2(lx1,ly1,lz1,1)
1899  $ , s1(lx1,ly1,lz1,1)
1900  $ , s2(lx1,ly1,lz1,1)
1901  $ , r1mask(lx1,ly1,lz1,1)
1902  $ , r2mask(lx1,ly1,lz1,1)
1903 C
1904  nface = 2*ldim
1905 C
1906  DO 100 iel=1,nel
1907  DO 100 ifc=1,nface
1908  IF (.NOT.ifmsfc(ifc,iel,ifield)) GO TO 100
1909  CALL facind2 (js1,jf1,jskip1,js2,jf2,jskip2,ifc)
1910  DO 120 j2=js2,jf2,jskip2
1911  DO 120 j1=js1,jf1,jskip1
1912  rnor = ( s1(j1,j2,1,iel)*vnx(j1,j2,1,iel) +
1913  $ s2(j1,j2,1,iel)*vny(j1,j2,1,iel) ) *
1914  $ r1mask(j1,j2,1,iel)
1915  rtn1 = ( s1(j1,j2,1,iel)*v1x(j1,j2,1,iel) +
1916  $ s2(j1,j2,1,iel)*v1y(j1,j2,1,iel) ) *
1917  $ r2mask(j1,j2,1,iel)
1918  r1(j1,j2,1,iel) = rnor*vnx(j1,j2,1,iel) +
1919  $ rtn1*v1x(j1,j2,1,iel)
1920  r2(j1,j2,1,iel) = rnor*vny(j1,j2,1,iel) +
1921  $ rtn1*v1y(j1,j2,1,iel)
1922  120 CONTINUE
1923  100 CONTINUE
1924 C
1925  RETURN
1926  END
1927 C-----------------------------------------------------------------------
1928  SUBROUTINE fcmsk3 (R1,R2,R3,S1,S2,S3,R1MASK,R2MASK,R3MASK,NEL)
1929 C
1930  include 'SIZE'
1931  include 'GEOM'
1932  include 'TSTEP'
1933  dimension r1(lx1,ly1,lz1,1)
1934  $ , r2(lx1,ly1,lz1,1)
1935  $ , r3(lx1,ly1,lz1,1)
1936  $ , s1(lx1,ly1,lz1,1)
1937  $ , s2(lx1,ly1,lz1,1)
1938  $ , s3(lx1,ly1,lz1,1)
1939  $ , r1mask(lx1,ly1,lz1,1)
1940  $ , r2mask(lx1,ly1,lz1,1)
1941  $ , r3mask(lx1,ly1,lz1,1)
1942 C
1943  nface = 2*ldim
1944 C
1945  DO 100 iel=1,nel
1946  DO 100 ifc=1,nface
1947  IF (.NOT.ifmsfc(ifc,iel,ifield)) GO TO 100
1948  CALL facind2 (js1,jf1,jskip1,js2,jf2,jskip2,ifc)
1949  DO 120 j2=js2,jf2,jskip2
1950  DO 120 j1=js1,jf1,jskip1
1951  rnor = ( s1(j1,j2,1,iel)*vnx(j1,j2,1,iel) +
1952  $ s2(j1,j2,1,iel)*vny(j1,j2,1,iel) +
1953  $ s3(j1,j2,1,iel)*vnz(j1,j2,1,iel) ) *
1954  $ r1mask(j1,j2,1,iel)
1955  rtn1 = ( s1(j1,j2,1,iel)*v1x(j1,j2,1,iel) +
1956  $ s2(j1,j2,1,iel)*v1y(j1,j2,1,iel) +
1957  $ s3(j1,j2,1,iel)*v1z(j1,j2,1,iel) ) *
1958  $ r2mask(j1,j2,1,iel)
1959  rtn2 = ( s1(j1,j2,1,iel)*v2x(j1,j2,1,iel) +
1960  $ s2(j1,j2,1,iel)*v2y(j1,j2,1,iel) +
1961  $ s3(j1,j2,1,iel)*v2z(j1,j2,1,iel) ) *
1962  $ r3mask(j1,j2,1,iel)
1963  r1(j1,j2,1,iel) = rnor*vnx(j1,j2,1,iel) +
1964  $ rtn1*v1x(j1,j2,1,iel) +
1965  $ rtn2*v2x(j1,j2,1,iel)
1966  r2(j1,j2,1,iel) = rnor*vny(j1,j2,1,iel) +
1967  $ rtn1*v1y(j1,j2,1,iel) +
1968  $ rtn2*v2y(j1,j2,1,iel)
1969  r3(j1,j2,1,iel) = rnor*vnz(j1,j2,1,iel) +
1970  $ rtn1*v1z(j1,j2,1,iel) +
1971  $ rtn2*v2z(j1,j2,1,iel)
1972  120 CONTINUE
1973  100 CONTINUE
1974 C
1975  RETURN
1976  END
1977 C-----------------------------------------------------------------------
1978  SUBROUTINE egmask (R1,R2,R3,S1,S2,S3,R1MASK,R2MASK,R3MASK,NEL)
1979 C
1980  include 'SIZE'
1981  include 'GEOM'
1982  include 'TSTEP'
1983  dimension r1(lx1,ly1,lz1,1)
1984  $ , r2(lx1,ly1,lz1,1)
1985  $ , r3(lx1,ly1,lz1,1)
1986  $ , s1(lx1,ly1,lz1,1)
1987  $ , s2(lx1,ly1,lz1,1)
1988  $ , s3(lx1,ly1,lz1,1)
1989  $ , r1mask(lx1,ly1,lz1,1)
1990  $ , r2mask(lx1,ly1,lz1,1)
1991  $ , r3mask(lx1,ly1,lz1,1)
1992 C
1993  nedge = 12
1994 C
1995  DO 100 iel=1,nel
1996  DO 100 isd=1,nedge
1997  IF (.NOT.ifmseg(isd,iel,ifield)) GOTO 100
1998  CALL edgindv (lv1,lv2,lvskip,isd)
1999  DO 120 lv=lv1,lv2,lvskip
2000  rnor = ( s1(lv,1,1,iel)*vnx(lv,1,1,iel) +
2001  $ s2(lv,1,1,iel)*vny(lv,1,1,iel) +
2002  $ s3(lv,1,1,iel)*vnz(lv,1,1,iel) ) *
2003  $ r1mask(lv,1,1,iel)
2004  rtn1 = ( s1(lv,1,1,iel)*v1x(lv,1,1,iel) +
2005  $ s2(lv,1,1,iel)*v1y(lv,1,1,iel) +
2006  $ s3(lv,1,1,iel)*v1z(lv,1,1,iel) ) *
2007  $ r2mask(lv,1,1,iel)
2008  rtn2 = ( s1(lv,1,1,iel)*v2x(lv,1,1,iel) +
2009  $ s2(lv,1,1,iel)*v2y(lv,1,1,iel) +
2010  $ s3(lv,1,1,iel)*v2z(lv,1,1,iel) ) *
2011  $ r3mask(lv,1,1,iel)
2012  r1(lv,1,1,iel) = rnor*vnx(lv,1,1,iel) +
2013  $ rtn1*v1x(lv,1,1,iel) +
2014  $ rtn2*v2x(lv,1,1,iel)
2015  r2(lv,1,1,iel) = rnor*vny(lv,1,1,iel) +
2016  $ rtn1*v1y(lv,1,1,iel) +
2017  $ rtn2*v2y(lv,1,1,iel)
2018  r3(lv,1,1,iel) = rnor*vnz(lv,1,1,iel) +
2019  $ rtn1*v1z(lv,1,1,iel) +
2020  $ rtn2*v2z(lv,1,1,iel)
2021  120 CONTINUE
2022  100 CONTINUE
2023 C
2024  RETURN
2025  END
2026 C-----------------------------------------------------------------------
2027  SUBROUTINE crmsk2 (R1,R2,S1,S2,R1MASK,R2MASK,NEL)
2028 C
2029  include 'SIZE'
2030  include 'GEOM'
2031  include 'TSTEP'
2032  common /indxfc/ mcrfc(4,6)
2033  $ , mfccr(3,8)
2034  $ , megcr(3,8)
2035  $ , mfceg(2,12)
2036  $ , mcreg(2,12)
2037  $ , mcrrst(3,8)
2038  $ , midrst(3,12)
2039  $ , mcrind(8)
2040  $ , medind(2,4)
2041  $ , ntefc(2,12)
2042  $ , ntcrf(2,3)
2043  dimension r1(lx1,ly1,lz1,1)
2044  $ , r2(lx1,ly1,lz1,1)
2045  $ , s1(lx1,ly1,lz1,1)
2046  $ , s2(lx1,ly1,lz1,1)
2047  $ , r1mask(lx1,ly1,lz1,1)
2048  $ , r2mask(lx1,ly1,lz1,1)
2049 C
2050  ncrnr = 2**ldim
2051 C
2052  DO 100 iel=1,nel
2053  DO 100 icr=1,ncrnr
2054  IF (.NOT.ifmscr(icr,iel,ifield)) GO TO 100
2055  ix = mcrrst(1,icr)
2056  iy = mcrrst(2,icr)
2057  iz = 1
2058  rnor = ( s1(ix,iy,iz,iel)*vnx(ix,iy,iz,iel) +
2059  $ s2(ix,iy,iz,iel)*vny(ix,iy,iz,iel) ) *
2060  $ r1mask(ix,iy,iz,iel)
2061  rtn1 = ( s1(ix,iy,iz,iel)*v1x(ix,iy,iz,iel) +
2062  $ s2(ix,iy,iz,iel)*v1y(ix,iy,iz,iel) ) *
2063  $ r2mask(ix,iy,iz,iel)
2064  r1(ix,iy,iz,iel) = rnor*vnx(ix,iy,iz,iel) +
2065  $ rtn1*v1x(ix,iy,iz,iel)
2066  r2(ix,iy,iz,iel) = rnor*vny(ix,iy,iz,iel) +
2067  $ rtn1*v1y(ix,iy,iz,iel)
2068  100 CONTINUE
2069 C
2070  RETURN
2071  END
2072 C-----------------------------------------------------------------------
2073  SUBROUTINE crmsk3 (R1,R2,R3,S1,S2,S3,R1MASK,R2MASK,R3MASK,NEL)
2074 C
2075  include 'SIZE'
2076  include 'GEOM'
2077  include 'TSTEP'
2078  common /indxfc/ mcrfc(4,6)
2079  $ , mfccr(3,8)
2080  $ , megcr(3,8)
2081  $ , mfceg(2,12)
2082  $ , mcreg(2,12)
2083  $ , mcrrst(3,8)
2084  $ , midrst(3,12)
2085  $ , mcrind(8)
2086  $ , medind(2,4)
2087  $ , ntefc(2,12)
2088  $ , ntcrf(2,3)
2089  dimension r1(lx1,ly1,lz1,1)
2090  $ , r2(lx1,ly1,lz1,1)
2091  $ , r3(lx1,ly1,lz1,1)
2092  $ , s1(lx1,ly1,lz1,1)
2093  $ , s2(lx1,ly1,lz1,1)
2094  $ , s3(lx1,ly1,lz1,1)
2095  $ , r1mask(lx1,ly1,lz1,1)
2096  $ , r2mask(lx1,ly1,lz1,1)
2097  $ , r3mask(lx1,ly1,lz1,1)
2098 C
2099  ncrnr = 2**ldim
2100 C
2101  DO 100 iel=1,nel
2102  DO 100 icr=1,ncrnr
2103  IF (.NOT.ifmscr(icr,iel,ifield)) GO TO 100
2104  ix = mcrrst(1,icr)
2105  iy = mcrrst(2,icr)
2106  iz = mcrrst(3,icr)
2107  rnor = ( s1(ix,iy,iz,iel)*vnx(ix,iy,iz,iel) +
2108  $ s2(ix,iy,iz,iel)*vny(ix,iy,iz,iel) +
2109  $ s3(ix,iy,iz,iel)*vnz(ix,iy,iz,iel) ) *
2110  $ r1mask(ix,iy,iz,iel)
2111  rtn1 = ( s1(ix,iy,iz,iel)*v1x(ix,iy,iz,iel) +
2112  $ s2(ix,iy,iz,iel)*v1y(ix,iy,iz,iel) +
2113  $ s3(ix,iy,iz,iel)*v1z(ix,iy,iz,iel) ) *
2114  $ r2mask(ix,iy,iz,iel)
2115  rtn2 = ( s1(ix,iy,iz,iel)*v2x(ix,iy,iz,iel) +
2116  $ s2(ix,iy,iz,iel)*v2y(ix,iy,iz,iel) +
2117  $ s3(ix,iy,iz,iel)*v2z(ix,iy,iz,iel) ) *
2118  $ r3mask(ix,iy,iz,iel)
2119  r1(ix,iy,iz,iel) = rnor*vnx(ix,iy,iz,iel) +
2120  $ rtn1*v1x(ix,iy,iz,iel) +
2121  $ rtn2*v2x(ix,iy,iz,iel)
2122  r2(ix,iy,iz,iel) = rnor*vny(ix,iy,iz,iel) +
2123  $ rtn1*v1y(ix,iy,iz,iel) +
2124  $ rtn2*v2y(ix,iy,iz,iel)
2125  r3(ix,iy,iz,iel) = rnor*vnz(ix,iy,iz,iel) +
2126  $ rtn1*v1z(ix,iy,iz,iel) +
2127  $ rtn2*v2z(ix,iy,iz,iel)
2128  100 CONTINUE
2129 C
2130  RETURN
2131  END
2132 C-----------------------------------------------------------------------
2133  subroutine getsnormal(sn,ix,iy,iz,iside,e)
2134 
2135 c calculate surface normal
2136 
2137  include 'SIZE'
2138  include 'GEOM'
2139  include 'TOPOL'
2140 
2141  real sn(3)
2142  integer e,f
2143 
2144  f = eface1(iside)
2145 
2146  if (1.le.f.and.f.le.2) then ! "r face"
2147  sn(1) = unx(iy,iz,iside,e)
2148  sn(2) = uny(iy,iz,iside,e)
2149  sn(3) = unz(iy,iz,iside,e)
2150  elseif (3.le.f.and.f.le.4) then ! "s face"
2151  sn(1) = unx(ix,iz,iside,e)
2152  sn(2) = uny(ix,iz,iside,e)
2153  sn(3) = unz(ix,iz,iside,e)
2154  elseif (5.le.f.and.f.le.6) then ! "t face"
2155  sn(1) = unx(ix,iy,iside,e)
2156  sn(2) = uny(ix,iy,iside,e)
2157  sn(3) = unz(ix,iy,iside,e)
2158  endif
2159 
2160  return
2161  end
2162 
2163  subroutine fixmska (c1mask,c2mask,c3mask)
2164 
2165 c fixes masks for A/SYM face corners
2166 
2167  include 'SIZE'
2168  include 'INPUT'
2169 
2170  real c1mask(lx1,ly1,lz1,1)
2171  $ ,c2mask(lx1,ly1,lz1,1)
2172  $ ,c3mask(lx1,ly1,lz1,1)
2173 
2174  common /ctmp0/ im1(lx1,ly1,lz1),im2(lx1,ly1,lz1)
2175  integer e,f,val,im1,im2
2176 
2177  character*3 cb
2178 
2179  n = lx1*ly1*lz1
2180 
2181  nface = 2*ldim
2182 
2183  do e=1,nelv
2184  call izero (im1,n)
2185  call izero (im2,n)
2186  do f=1,nface
2187  cb = cbc(f,e,1)
2188  if (cb.eq.'SYM') call iface_e(im1,f,1,lx1,ly1,lz1)
2189  if (cb.eq.'A ') call iface_e(im2,f,2,lx1,ly1,lz1)
2190  enddo
2191  call icol2(im2,im1,n)
2192 
2193  k = 1
2194  do j=1,ly1,ly1-1
2195  do i=1,lx1,lx1-1
2196  if ( im2(i,j,k) .eq. 2) then ! corner of SYM & 'A ' faces
2197  c1mask(i,j,k,e) = 0.
2198  c2mask(i,j,k,e) = 0.
2199  endif
2200  enddo
2201  enddo
2202  enddo
2203 
2204  return
2205  end
2206 c-----------------------------------------------------------------------
2207  subroutine icol2(a,b,n)
2208  integer a(1),b(1)
2209 
2210  do i=1,n
2211  a(i)=a(i)*b(i)
2212  enddo
2213 
2214  return
2215  end
2216 c-----------------------------------------------------------------------
2217  subroutine iface_e(a,iface,val,nx,ny,nz)
2218 
2219 C Assign the value VAL to face(IFACE,IE) of array A.
2220 C IFACE is the input in the pre-processor ordering scheme.
2221 
2222  include 'SIZE'
2223  integer a(nx,ny,nz),val
2224  call facind (kx1,kx2,ky1,ky2,kz1,kz2,nx,ny,nz,iface)
2225  do 100 iz=kz1,kz2
2226  do 100 iy=ky1,ky2
2227  do 100 ix=kx1,kx2
2228  a(ix,iy,iz)=val
2229  100 continue
2230  return
2231  end
2232 c-----------------------------------------------------------------------
2233  function op_vlsc2_wt(b1,b2,b3,x1,x2,x3,wt)
2234  include 'SIZE'
2235  include 'INPUT'
2236  include 'TSTEP'
2237  real b1(1),b2(1),b3(1),x1(1),x2(1),x3(1),wt(1)
2238 
2239  nel = nelfld(ifield)
2240  n = lx1*ly1*lz1*nel
2241 
2242  s = 0
2243  if (if3d) then
2244  do i=1,n
2245  s=s+wt(i)*(b1(i)*x1(i)+b2(i)*x2(i)+b3(i)*x3(i))
2246  enddo
2247  else
2248  do i=1,n
2249  s=s+wt(i)*(b1(i)*x1(i)+b2(i)*x2(i))
2250  enddo
2251  endif
2252  op_vlsc2_wt = s
2253 
2254  return
2255  end
2256 c-----------------------------------------------------------------------
2257  function op_glsc2_wt(b1,b2,b3,x1,x2,x3,wt)
2258  include 'SIZE'
2259  include 'INPUT'
2260  include 'TSTEP'
2261  real b1(1),b2(1),b3(1),x1(1),x2(1),x3(1),wt(1)
2262 
2263  nel = nelfld(ifield)
2264  n = lx1*ly1*lz1*nel
2265 
2266  s = 0
2267  if (if3d) then
2268  do i=1,n
2269  s=s+wt(i)*(b1(i)*x1(i)+b2(i)*x2(i)+b3(i)*x3(i))
2270  enddo
2271  else
2272  do i=1,n
2273  s=s+wt(i)*(b1(i)*x1(i)+b2(i)*x2(i))
2274  enddo
2275  endif
2276  op_glsc2_wt = glsum(s,1)
2277 
2278  return
2279  end
2280 c-----------------------------------------------------------------------
unitvec
subroutine unitvec(X, Y, Z, N)
Definition: bdry.f:1831
lfalse
subroutine lfalse(IFA, N)
Definition: bdry.f:1813
rzero3
subroutine rzero3(A, B, C, N)
Definition: bdry.f:1821
facind2
subroutine facind2(JS1, JF1, JSKIP1, JS2, JF2, JSKIP2, IFC)
Definition: bdry.f:2098
exitt
void exitt()
Definition: comm_mpi.f:604
facev
subroutine facev(a, ie, iface, val, nx, ny, nz)
Definition: connect1.f:1077
facind
subroutine facind(kx1, kx2, ky1, ky2, kz1, kz2, nx, ny, nz, iface)
Definition: connect1.f:1028
dsset
subroutine dsset(nx, ny, nz)
Definition: connect1.f:553
runstat
subroutine runstat
Definition: drive2.f:915
dsop
subroutine dsop(u, op, nx, ny, nz)
Definition: dssum.f:101
dssum
subroutine dssum(u, nx, ny, nz)
Definition: dssum.f:34
dot
real function dot(V1, V2, N)
Definition: genxyz.f:885
col3
subroutine col3(a, b, c, n)
Definition: math.f:598
addcol3
subroutine addcol3(a, b, c, n)
Definition: math.f:654
glsc3
function glsc3(a, b, mult, n)
Definition: math.f:776
col2
subroutine col2(a, b, n)
Definition: math.f:564
add2s2
subroutine add2s2(a, b, c1, n)
Definition: math.f:690
copy
subroutine copy(a, b, n)
Definition: math.f:260
vdot2
subroutine vdot2(dot, u1, u2, v1, v2, n)
Definition: math.f:447
col4
subroutine col4(a, b, c, d, n)
Definition: math.f:120
izero
subroutine izero(a, n)
Definition: math.f:215
subcol3
subroutine subcol3(a, b, c, n)
Definition: math.f:186
gllog
subroutine gllog(la, lb)
Definition: math.f:986
glsum
function glsum(x, n)
Definition: math.f:861
subcol4
subroutine subcol4(a, b, c, d, n)
Definition: math.f:197
sub3
subroutine sub3(a, b, c, n)
Definition: math.f:175
cmult
subroutine cmult(a, const, n)
Definition: math.f:315
cfill
subroutine cfill(a, b, n)
Definition: math.f:244
rone
subroutine rone(a, n)
Definition: math.f:230
add2s1
subroutine add2s1(a, b, c1, n)
Definition: math.f:678
vdot3
subroutine vdot3(dot, u1, u2, u3, v1, v2, v3, n)
Definition: math.f:462
vcross
subroutine vcross(u1, u2, u3, v1, v2, v3, w1, w2, w3, n)
Definition: math.f:430
rzero
subroutine rzero(a, n)
Definition: math.f:208
chsign
subroutine chsign(a, n)
Definition: math.f:305
opdssum
subroutine opdssum(a, b, c)
Definition: navier1.f:2582
prepost
subroutine prepost(ifdoin, prefin)
Definition: prepost.f:91
axhmsf
subroutine axhmsf(au1, au2, au3, u1, u2, u3, h1, h2, matmod)
Definition: subs1.f:1241
stnrate
subroutine stnrate(u1, u2, u3, nel, matmod)
Definition: subs1.f:1319
amask
subroutine amask(VB1, VB2, VB3, V1, V2, V3, NEL)
Definition: subs2.f:1757
edgindf
subroutine edgindf(LF1, LF2, LFSKIP, ISD, IFCN)
Definition: subs2.f:1681
solvel
subroutine solvel
Definition: subs2.f:484
vsoln
subroutine vsoln(UX, UY, UZ, X, Y, Z, PI)
Definition: subs2.f:502
printel
subroutine printel(TA, A, IEL)
Definition: subs2.f:551
crmsk3
subroutine crmsk3(R1, R2, R3, S1, S2, S3, R1MASK, R2MASK, R3MASK, NEL)
Definition: subs2.f:2074
fxwms3
subroutine fxwms3(W2MASK, W3MASK, HVMASK, HFMASK, NEL)
Definition: subs2.f:1319
comavn2
subroutine comavn2(HVMASK, HFMASK, NEL)
Definition: subs2.f:1040
comavn3
subroutine comavn3(HVMASK, HFMASK, NEL)
Definition: subs2.f:1105
solpres
subroutine solpres
Definition: subs2.f:520
facsub2
subroutine facsub2(A1, A2, A3, B1, B2, B3, IFACE1)
Definition: subs2.f:220
setaxw2
subroutine setaxw2(IFAXWG)
Definition: subs2.f:17
egmask
subroutine egmask(R1, R2, R3, S1, S2, S3, R1MASK, R2MASK, R3MASK, NEL)
Definition: subs2.f:1979
skipcnr
subroutine skipcnr(NEL)
Definition: subs2.f:829
setcdof
subroutine setcdof
Definition: subs2.f:1742
cmult2
subroutine cmult2(A, B, CONST, N)
Definition: subs2.f:341
crmsk2
subroutine crmsk2(R1, R2, S1, S2, R1MASK, R2MASK, NEL)
Definition: subs2.f:2028
emerxit
subroutine emerxit
Definition: subs2.f:357
opadds
subroutine opadds(A1, A2, A3, B1, B2, B3, CONST, N, ISC)
Definition: subs2.f:101
op_glsc2_wt
function op_glsc2_wt(b1, b2, b3, x1, x2, x3, wt)
Definition: subs2.f:2258
stx1sf
subroutine stx1sf
Definition: subs2.f:439
stnrinv
subroutine stnrinv
Definition: subs2.f:33
outf1
subroutine outf1(X, TXT, IEL, IFC)
Definition: subs2.f:593
setmlog
subroutine setmlog(HVMASK, HFMASK, IFLD, NEL)
Definition: subs2.f:908
getsnormal
subroutine getsnormal(sn, ix, iy, iz, iside, e)
Definition: subs2.f:2134
outm1
subroutine outm1(X, TXT, NP, IEL, IP)
Definition: subs2.f:621
qmask
subroutine qmask(R1, R2, R3, R1MASK, R2MASK, R3MASK, NEL)
Definition: subs2.f:1836
facexv
subroutine facexv(A1, A2, A3, B1, B2, B3, IFACE1, IOP)
Definition: subs2.f:170
opdot
subroutine opdot(DP, A1, A2, A3, B1, B2, B3, N)
Definition: subs2.f:80
setcdat
subroutine setcdat
Definition: subs2.f:1402
setmask
subroutine setmask(C1MASK, C2MASK, C3MASK, HVMASK, NEL)
Definition: subs2.f:877
fcmsk2
subroutine fcmsk2(R1, R2, S1, S2, R1MASK, R2MASK, NEL)
Definition: subs2.f:1893
fcmsk3
subroutine fcmsk3(R1, R2, R3, S1, S2, S3, R1MASK, R2MASK, R3MASK, NEL)
Definition: subs2.f:1929
stsmask
subroutine stsmask(C1MASK, C2MASK, C3MASK)
Definition: subs2.f:676
prsoln
subroutine prsoln(P, X, Y, Z, PI)
Definition: subs2.f:538
edgindv
subroutine edgindv(LV1, LV2, LVSKIP, ISD)
Definition: subs2.f:1707
icol2
subroutine icol2(a, b, n)
Definition: subs2.f:2208
add3s
subroutine add3s(A, B, C, CONST, N)
Definition: subs2.f:349
faccvs
subroutine faccvs(A1, A2, A3, B, IFACE1)
Definition: subs2.f:392
outm2
subroutine outm2(X, TXT, NP, IEL, IP)
Definition: subs2.f:649
updmsys
subroutine updmsys(IFLD)
Definition: subs2.f:712
printv
subroutine printv(TA, A, NEL)
Definition: subs2.f:569
op_vlsc2_wt
function op_vlsc2_wt(b1, b2, b3, x1, x2, x3, wt)
Definition: subs2.f:2234
setaxw1
subroutine setaxw1(IFAXWG)
Definition: subs2.f:2
facexs
subroutine facexs(A, B, IFACE1, IOP)
Definition: subs2.f:125
setcsys
subroutine setcsys(HVMASK, HFMASK, NEL)
Definition: subs2.f:1008
fixmska
subroutine fixmska(c1mask, c2mask, c3mask)
Definition: subs2.f:2164
fxwms2
subroutine fxwms2(W2MASK, HVMASK, HFMASK, NEL)
Definition: subs2.f:1271
gammasf
subroutine gammasf(H1, H2)
Definition: subs2.f:256
rmask
subroutine rmask(R1, R2, R3, NEL)
Definition: subs2.f:1801
fixwmsk
subroutine fixwmsk(W2MASK, W3MASK, HVMASK, HFMASK, NEL)
Definition: subs2.f:1250
sethmsk
subroutine sethmsk(HVMASK, HFMASK, IFLD, NEL)
Definition: subs2.f:729
iface_e
subroutine iface_e(a, iface, val, nx, ny, nz)
Definition: subs2.f:2218