Actual source code: user.h
2: #include <petscsys.h>
3: #include <petsc/private/fortranimpl.h>
5: #define max_colors 200
6: #define max_nbtran 20
8: #define REAL double
10: typedef struct gxy { /* GRID STRUCTURE */
11: int nnodes; /* Global number of nodes */
12: int ncell; /* Global number of cells */
13: int nedge; /* Global number of edges */
14: int ncolor; /* Number of colors */
15: int nccolor; /* Number of colors for cells */
16: int nncolor; /* Number of colors for nodes */
17: int ncount[max_colors]; /* No. of faces in color */
18: int nccount[max_colors]; /* No. of cells in color */
19: int nncount[max_colors]; /* No. of nodes in color */
20: int iup; /* if 1, upward int coef reqd */
21: int idown; /* if 1, dwnwrd int coef reqd */
23: int nsface; /* Total # of solid faces */
24: int nvface; /* Total # of viscous faces */
25: int nfface; /* Total # of far field faces */
26: int nsnode; /* Total # of solid nodes */
27: int nvnode; /* Total # of viscous nodes */
28: int nfnode; /* Total # of far field nodes */
29: int jvisc; /* 0 = Euler */
30: /* 1 = laminar no visc LHS */
31: /* 2 = laminar w/ visc LHS */
32: /* 3 = turb B-B no visc LHS */
33: /* 4 = turb B-B w/ visc LHS */
34: /* 5 = turb Splrt no visc LHS */
35: /* 6 = turb Splrt w/ visc LHS */
36: int ileast; /* 1 = Lst square gradient */
37: int nsets; /* No of levels for scheduling*/
38: int *eptr; /* edge pointers */
39: int *isface; /* Face # of solid faces */
40: int *ifface; /* Face # of far field faces */
41: int *ivface; /* Face # of viscous faces */
42: int *isford; /* Copies of isface, ifface, */
43: int *ifford; /* and ivface used for */
44: int *ivford; /* ordering */
45: int *isnode; /* Node # of solid nodes */
46: int *ivnode; /* Node # of viscous nodes */
47: int *ifnode; /* Node # of far field nodes */
48: int *nflag; /* Node flag */
49: int *nnext; /* Next node */
50: int *nneigh; /* Neighbor of a node */
51: int *c2n; /* Cell-to-node pointers */
52: int *c2e; /* Cell-to-edge pointers */
53: int *c2c; /* Cell-to-cell pointers */
54: int *ctag; /* Cell tags */
55: int *csearch; /* Cell search list */
56: int *cenc; /* Enclosing cell for node */
57: int *clist; /* Colored list of cells */
58: int *iupdate; /* Tells whether to update */
59: int *sface; /* Nodes for solid faces */
60: int *vface; /* Nodes for viscous faces */
61: int *fface; /* Nodes for far field faces */
62: int *icount; /* # of surrounding nodes */
63: int *isetc; /* Nodes in each level */
64: int *iset; /* Actual nodes for levels */
65: /* Forward substitution for ILU */
66: int *nlcol; /* No of edge colors for sets */
67: int *nlcount; /* How many edges in each colr*/
68: int *lvface; /* Edges that influence a set */
69: /* Back substitution for ILU */
70: int *nbcol; /* No of edge colors for sets */
71: int *nbcount; /* How many edges in each colr*/
72: int *lbface; /* Edges that influence a set */
73: REAL *xyz; /* Node Coordinates */
74: REAL *area; /* Area of control volumes */
75: /*REAL *gradx, *grady, *gradz;*/ /* Gradients */
76: REAL *cdt; /* Local time step */
77: REAL *qcp, *rcp; /* Two work arrays */
78: REAL *ff; /* MG forcing function */
79: REAL *dfp, *dfm; /* Flux Jacobians */
80: REAL *dft1, *dft2; /* Turb mod linearization */
81: REAL *slen; /* Generalized distance */
82: REAL *turbre; /* nu x turb Reynolds # */
83: REAL *amut; /* Turbulent mu (viscosity) */
84: REAL *turbres; /* Turbulent residual */
85: REAL *turbff; /* Turbulent forcing function */
86: REAL *turbold; /* Turbulent unknown (for MG) */
87: REAL *sxn, *syn, *szn, *sa; /* Normals at solid nodes */
88: REAL *vxn, *vyn, *vzn, *va; /* Normals at viscous nodes */
89: REAL *fxn, *fyn, *fzn, *fa; /* Normals at far field nodes */
90: REAL *xyzn; /* Normal to faces and length */
91: REAL *us, *vs, *ws, *as; /* For linearizing viscous */
92: REAL *phi; /* Flux limiter */
93: REAL *rxy; /* */
95: int *icoefup; /* Surrounding nodes */
96: REAL *rcoefup; /* Coefficients */
97: int *icoefdn; /* Surrounding nodes */
98: REAL *rcoefdn; /* Coefficients */
99: REAL *AP; /* Array for GMRES */
100: REAL *Fgm; /* Big array for GMRES */
101: REAL *Xgm; /* Another GMRES array */
102: REAL *temr; /* Temporary array */
103: REAL *ALU; /* Big array for ILU(0) */
104: int *ia, *iau, *ja, *fhelp; /* Stuff for ILU(0) */
106: /*
107: * stuff to read in daves grid file
108: */
109: int nnbound,nvbound,nfbound,nnfacet,nvfacet,nffacet,ntte;
110: int *ncolorn,*countn,*ncolorv,*countv,*ncolorf,*countf;
111: int *nntet,*nnpts,*nvtet,*nvpts,*nftet,*nfpts;
112: int *f2ntn,*f2ntv,*f2ntf;
114: /* PETSc datastructures and other related info */
115: Vec qnode; /* Global distributed solution vector*/
116: Vec qnodeLoc; /* Local sequential solution vector*/
117: Vec dq; /* Delta Q */
118: Vec qold; /* Global distributed solution vector*/
119: Vec res; /* Residual */
120: Vec grad; /* Gradient Vector */
121: Vec gradLoc; /* Local Gradient Vector */
122: Vec B; /* Right hand side */
123: Mat A; /* Left hand side */
124: VecScatter scatter, gradScatter; /* Scatter between local and global vectors */
125: int *loc2pet; /* local to PETSc mapping */
126: int *loc2glo; /* local to global mapping */
127: int *v2p; /* Vertex to processor mapping */
128: int *sface_bit, *vface_bit;
129: int nnodesLoc, nedgeLoc, nvertices; /* nnodesLoc=number of owned nodes, nedgeLoc=number of edges touching owned nodes, nvertices=includes ghost nodes */
130: int nsnodeLoc, nvnodeLoc, nfnodeLoc;
131: int nnfacetLoc, nvfacetLoc, nffacetLoc;
133: /* global arrays needed for post processing */
134: /*int *indGlo, *isnodeGlo, *ivnodeGlo, *f2ntnGlo, *f2ntvGlo;
135: REAL *xGlo, *yGlo, *zGlo;
136: Vec qGlo;
137: VecScatter scatterGlo;*/
139: #if defined(_OPENMP)
140: #if defined(HAVE_REDUNDANT_WORK)
141: REAL *resd;
142: #else
143: int nedgeAllThr;
144: int *part_thr,*nedge_thr,*edge_thr;
145: REAL *xyzn_thr;
146: #endif
147: #endif
148: } GRID; /* Grids */
149: /*============================*/
151: /*============================*/
152: typedef struct { /* GENERAL INFORMATION */
153: REAL title[20]; /* Title line */
154: REAL beta; /* Artificial Compress. Param */
155: REAL alpha; /* Angle of attack */
156: REAL Re; /* Reynolds number */
157: REAL dt; /* Input cfl */
158: REAL tot; /* total computer time */
159: REAL res0; /* Begining residual */
160: REAL resc; /* Current residual */
161: int ntt; /* A counter */
162: int mseq; /* Mesh sequencing */
163: int ivisc; /* 0 = Euler */
164: /* 1 = laminar no visc LHS */
165: /* 2 = laminar w/ visc LHS */
166: /* 3 = turb BB no visc LHS */
167: /* 4 = turb BB w/ visc LHS */
168: /* 5 = turb SA w/ visc LHS */
169: /* 6 = turb SA w/ visc LHS */
170: int irest; /* for restarts irest = 1 */
171: int icyc; /* iterations completed */
172: int ihane; /* ihane = 0 for van leer fds */
173: /* = 1 for hanel flux */
174: /* = 2 for Roe's fds */
175: int ntturb; /* Counter for turbulence */
176: } CINFO; /* COMMON INFO */
177: /*============================*/
179: /*============================*/
180: typedef struct { /* FLOW SOLVER VARIABLES */
181: REAL cfl1; /* Starting CFL number */
182: REAL cfl2; /* Ending CFL number */
183: int nsmoth; /* How many its for Res smooth*/
184: int iflim; /* 1=use limiter 0=no limiter */
185: int itran; /* 1=transition (spalart only)*/
186: int nbtran; /* No. of transition points */
187: int jupdate; /* For freezing Jacobians */
188: int nstage; /* Number of subiterations */
189: int ncyct; /* Subiterations for turb mod */
190: int iramp; /* Ramp CFL over iramp iters */
191: int nitfo; /* Iterations first order */
192: } CRUNGE; /* COMMON RUNGE */
193: /*============================*/
195: typedef struct { /*============================*/
196: REAL gtol; /* linear system tolerence */
197: int icycle; /* Number of GMRES iterations */
198: int nsrch; /* Dimension of Krylov */
199: int ilu0; /* 1 for ILU(0) */
200: int ifcn; /* 0=fcn2 1=fcneval(nwt Krlv) */
201: } CGMCOM; /* COMMON GMCOM */
202: /*============================*/
204: extern int set_up_grid(GRID*);
205: extern int write_fine_grid(GRID*);
207: /* ================= Fortran routines called from C ======================= */
208: #if defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
209: # define f77name(ucase,lcase,lcbar) lcbar
210: #elif defined(PETSC_HAVE_FORTRAN_CAPS)
211: # define f77name(ucase,lcase,lcbar) ucase
212: #else
213: # define f77name(ucase,lcase,lcbar) lcase
214: #endif
215: #define f77INFO f77name(INFO,info,info_)
216: #define f77RUNGE f77name(RUNGE,runge,runge_)
217: #define f77GMCOM f77name(GMCOM,gmcom,gmcom_)
218: #define f77FORLINK f77name(FORLINK,forlink,forlink_)
219: #define f77OPENM f77name(OPENM,openm,openm_)
220: #define f77READR1 f77name(READR1,readr1,readr1_)
221: #define f77READR2 f77name(READR2,readr2,readr2_)
222: #define f77READR3 f77name(READR3,readr3,readr3_)
223: #define f77RDGPAR f77name(RDGPAR,rdgpar,rdgpar_)
224: #define f77README f77name(README,readme,readme_)
225: #define f77COLORCJ f77name(COLORCJ,colorcj,colorcj_)
226: #define f77COLORCGS f77name(COLORCGS,colorcgs,colorcgs_)
227: #define f77BNDORD f77name(BNDORD,bndord,bndord_)
228: #define f77FINDIN f77name(FINDIN,findin,findin_)
229: #define f77ELMORD f77name(ELMORD,elmord,elmord_)
230: #define f77BNSHFT f77name(BNSHFT,bnshft,bnshft_)
231: #define f77VNSHFT f77name(VNSHFT,vnshft,vnshft_)
232: #define f77NSHIFT f77name(NSHIFT,nshift,nshift_)
233: #define f77NEIGHBR f77name(NEIGHBR,neighbr,neighbr_)
234: #define f77NSTACK f77name(NSTACK,nstack,nstack_)
235: #define f77GTCPTR f77name(GTCPTR,gtcptr,gtcptr_)
236: #define f77GTENCC f77name(GTENCC,gtencc,gtencc_)
237: #define f77INCOEF f77name(INCOEF,incoef,incoef_)
238: #define f77INTERP1 f77name(INTERP1,interp1,interp1_)
239: #define f77INTERP4 f77name(INTERP4,interp4,interp4_)
240: #define f77RCOLL1 f77name(RCOLL1,rcoll1,rcoll1_)
241: #define f77RCOLL f77name(RCOLL,rcoll,rcoll_)
242: #define f77INIT f77name(INIT,init,init_)
243: #define f77SUMGS f77name(SUMGS,sumgs,sumgs_)
244: #define f77GETAREA f77name(GETAREA,getarea,getarea_)
245: #define f77INFOTRN f77name(INFOTRN,infotrn,infotrn_)
246: #define f77GETRES f77name(GETRES,getres,getres_)
247: #define f77L2NORM f77name(L2NORM,l2norm,l2norm_)
248: #define f77FORCE f77name(FORCE,force,force_)
249: #define f77UPDATE f77name(UPDATE,update,update_)
250: #define f77WREST f77name(WREST,wrest,wrest_)
251: #define f77RREST f77name(RREST,rrest,rrest_)
252: #define f77PLLAN f77name(PLLAN,pllan,pllan_)
253: #define f77FLLAN f77name(FLLAN,fllan,fllan_)
254: #define f77TECFLO f77name(TECFLO,tecflo,tecflo_)
255: #define f77FASFLO f77name(FASFLO,fasflo,fasflo_)
256: #define f77BC f77name(BC,bc,bc_)
257: #define f77CLINK f77name(CLINK,clink,clink_)
258: #define f77SLENGTH f77name(SLENGTH,slength,slength_)
259: #define f77GETNDEX f77name(GETNDEX,getndex,getndex_)
260: #define f77CHANGEV f77name(CHANGEV,changev,changev_)
261: #define f77CHANGEP f77name(CHANGEP,changep,changep_)
262: #define f77TURBER f77name(TURBER,turber,turber_)
263: #define f77TURBRES f77name(TURBRES,turbres,turbres_)
264: #define f77SPALART f77name(SPALART,spalart,spalart_)
265: #define f77SPALRES f77name(SPALRES,spalres,spalres_)
266: #define f77PLOTURB f77name(PLOTURB,ploturb,ploturb_)
267: #define f77GETSKIN f77name(GETSKIN,getskin,getskin_)
268: #define f77GETC2N f77name(GETC2N,getc2n,getc2n_)
269: #define f77VWEIGHT f77name(VWEIGHT,vweight,vweight_)
270: #define f77PLOTCP f77name(PLOTCP,plotcp,plotcp_)
271: #define f77CORRSM f77name(CORRSM,corrsm,corrsm_)
272: #define f77CORRSM1 f77name(CORRSM1,corrsm1,corrsm1_)
274: #define f77GETIA f77name(GETIA,getia,getia_)
275: #define f77GETJA f77name(GETJA,getja,getja_)
276: #define f77SORTER f77name(SORTER,sorter,sorter_)
277: #define f77BLKILU f77name(BLKILU,blkilu,blkilu_)
278: #define f77BLKSOL f77name(BLKSOL,blksol,blksol_)
279: #define f77GETLEVEL f77name(GETLEVEL,getlevel,getlevel_)
280: #define f77LVCOLOR f77name(LVCOLOR,lvcolor,lvcolor_)
281: #define f77LBCOLOR f77name(LBCOLOR,lbcolor,lbcolor_)
283: #define f77FILLA f77name(FILLA,filla,filla_)
284: #define f77LSTGS f77name(LSTGS,lstgs,lstgs_)
285: #define f77IREAD f77name(IREAD,iread,iread_)
286: #define f77RREAD f77name(RREAD,rread,rread_)
288: EXTERN_C_BEGIN
289: extern void f77FORLINK(void);
290: extern void f77OPENM(int*);
291: extern void f77FILLA(int*,int*,int*,int*,int*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,
292: PetscScalar*,int*,int*,int*,int*,int*,int*,PetscScalar*,Mat*,PetscScalar*,PetscScalar*,PetscScalar*,
293: PetscScalar*,int*,int*);
294: extern void f77READR1(int*,int*);
295: extern void f77SUMGS(int*,int*,int*,PetscScalar*,PetscScalar*,int*,int*);
296: extern void f77INIT(int*,PetscScalar*,PetscScalar*,PetscScalar*,int*,int*,int*);
297: extern void f77LSTGS(int*,int*,int*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,int*,int*);
298: extern void f77GETRES(int*,int*,int*,int*,int*,int*,int*,int*,int*,int*,int*,int*,int*,
299: int*,int*,int*,int*,int*,int*,int*,int*,int*,int*,int*,int*,int*,
300: int*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,
301: PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,
302: PetscScalar*,PetscScalar*,PetscScalar*,int*,int*,
303: PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,int*,
304: #if defined(_OPENMP)
305: int*,
306: #if defined(HAVE_EDGE_COLORING)
307: int*, int*,
308: #elif defined(HAVE_REDUNDANT_WORK)
309: PetscScalar*,
310: #else
311: int*,
312: int*,int*,int*,PetscScalar*,
313: #endif
314: #endif
315: int*,int*,int*);
316: extern void f77FORCE(int*,int*,int*,int*,int*,int*,int*,int*,int*,int*,int*,PetscScalar*,PetscScalar*,int*,int*,PetscScalar*,PetscScalar*,PetscScalar*,int*,int*);
317: extern void f77GETIA(int*,int*,int*,int*,int*,int*);
318: extern void f77GETJA(int*,int*,int*,int*,int*,int*,int*);
319: extern void f77TECFLO(int* nnodes,int* nvbound,int* nvfacet,int* nvnode,
320: PetscScalar* x,PetscScalar* y,PetscScalar* z,
321: PetscScalar* qnode, int* nvpts, int* nvtet,
322: int* f2ntv, int* ivnode,
323: int* timeStep, int* rank, int* openFile, int* closeFile,
324: int* boundaryType,PetscScalar* title);
325: EXTERN_C_END