1 subroutine Au1MFG (ypre, y, ac, x, 2 & rmes, res, uBrg, 3 & BDiag, iBC, BC, 4 & iper, ilwork, shp, 5 & shgl, shpb, 6 & shglb) 7c 8c---------------------------------------------------------------------- 9c 10c This routine performs a matrix-vector product for the Matrix-Free 11c Implicit/Iterative solver using a one-sided scheme. 12c 13c input: 14c y (nshg,ndof) : Y-variables 15c ypre (nshg,nflow) : preconditioned Y-variables 16c (perturbed, no-scalars) 17c x (numnp,nsd) : node coordinates 18c rmes (nshg,nflow) : modified residual 19c res (nshg,nflow) : residual 20c uBrg (nshg,nflow) : Krylov space vector 21c BDiag (nshg,nflow,nflow) : block-diagonal preconditioner 22c iBC (nshg) : BC codes 23c BC (nshg,ndofBC) : BC constraint parameters 24c engBC (nshg) : energy for BC on density or pressure 25c shp(b) (nshape,ngauss) : element shape functions (boundary) 26c shgl(b)(nsd,nshape,ngauss) : local gradients of shape functions 27c 28c output: 29c uBrg (nshg,nflow) : Krylov space vectors 30c 31c 32c Zdenek Johan, Winter 1991. (Fortran 90) 33c---------------------------------------------------------------------- 34c 35 include "common.h" 36 include "mpif.h" 37 include "auxmpi.h" 38c 39 dimension y(nshg,ndof), ypre(nshg,nflow), 40 & x(numnp,nsd), ac(nshg,ndof), 41 & rmes(nshg,nflow), ytmp(nshg,nflow), 42 & res(nshg,nflow), uBrg(nshg,nflow), 43 & BDiag(nshg,nflow,nflow), iBC(nshg), 44 & BC(nshg,ndofBC), iper(nshg) 45c 46 dimension uBtmp(nshg,nflow), tmpBC(nshg), 47 & ilwork(nlwork) 48c 49 50 dimension shp(MAXTOP,maxsh,MAXQPT), 51 & shgl(MAXTOP,nsd,maxsh,MAXQPT), 52 & shpb(MAXTOP,maxsh,MAXQPT), 53 & shglb(MAXTOP,nsd,maxsh,MAXQPT) 54 55c$$$ dimension shp(nshl,ngauss), shgl(nsd,nshl,ngauss), 56c$$$ & shpb(nshl,ngaussb), 57c$$$ & shglb(nsd,nshl,ngaussb) 58c 59c.... calculate Rmod(y + eps u_i) 60c 61 uBtmp = zero 62c 63c call yshuffle(ypre, 'new2old ') 64 uBrg = ypre + eGMRES * uBrg 65c 66 call i3LU (BDiag, uBrg, 'backward') 67c 68 call yshuffle(uBrg, 'old2new ') 69c 70 call itrBC (uBrg, uBrg, iBC, BC, iper, ilwork) 71c 72 call itrRes (uBrg, y, 73 & x, shp, 74 & shgl, iBC, 75 & BC, shpb, 76 & shglb, uBtmp, 77 & iper, ilwork, ac) 78c 79 call i3LU (BDiag, uBtmp, 'forward ') 80c 81c.... calculate ( Rmod(y + eps u_i) - Rmod(y) ) / eps 82c 83 uBrg = ( uBtmp - rmes ) / eGMRES 84c ... before returning lets put ypre back in the new format 85c call yshuffle(ypre, 'old2new ') 86c 87c.... flop count 88c 89 ! flops = flops + 4*nflow*nshg 90c 91c.... end 92c 93 return 94 end 95 96 97 98 99