Actual source code: structs.h
1: #ifndef __DECLARAT_H__
2:
4:
5: #include <stdlib.h>
6: #include <unistd.h>
7:
8: #include "macros.h"
9: #include "petscksp.h"
10: #include "petscda.h"
11: #include "petsc.h"
12:
13: #if defined(MICSC_HAVE_SILO)
14: #include "silo.h"
15: #endif
16:
17: /*E
18: St_InputBcond - Structure containing the definition of a boundary
19: condition as it is specified in the input grid file.
20:
21: Attributes:
22: + x, y, z - Coordinates of the lower corner of the region where
23: the boundary condition is defined.
24: . nx, ny, xz - Number of cells in each direction of the region where the
25: boundary condition is defined.
26: . ineighb - Neighbor index.
27: . type - Type of boundary condition (see BcondType).
28: . porosity - Porosity value for the boundary condition.
29: - next - Pointer to the next structure in a list.
30:
31: Level: advanced
32: E*/
33: typedef struct St_InputBcond {
34: PetscInt x, y, z, nx, ny, nz, ineighb;
35: BcondType type;
36: PetscReal porosity;
37: struct St_InputBcond *next;
38: } St_InputBcond;
39:
40: /*E
41: St_Input - Structure containing all the input parameters of the problem,
42: specified by the user (or set to default values).
43:
44: Attributes:
45: + log_level - Logging level
46: . log_ievery - Number of iterations performed between logs
47: . out_ievery - Number of iterations performed between
48: binary solution outputs
49: . out_idtevery - Number of time steps performed between ascii
50: solution outputs
51: . startInstantMeans - First iteration for instant means calculation
52: . startVarianceMeans - First iteration for variance means calculation
53: . initial_field - Flag to print the initial field
54: . final_field - Flag to print the final field
55: . printProps - Flag to print properties
56: . max_outer_its - Maximum number of iterations
57: . min_outer_its - Minimum number of outer iterations
58: . resnorm_max - Maximum normalized residue
59: . corrnorm_max - Maximum normalized correction vector
60: . k_restart - Value to determine the restarting ('none', 'variables',
61: 'instant_means' or 'variance_means'). See RestartType.
62: . init_it - First iteration number
63: . k_wrap - Grid periodicity
64: . n_ghost - Ghost stencil width
65: . gravity - Gravity acceleration
66: . n_order_trans - Order of the transient scheme
67: . deltat - Time step
68: . tinit - Initial time
69: . tfinal - Final time
70: . uref - Constant value por the U direction inlet
71: . mu - MU value
72: . tref - TREF value
73: . rho - RHO value
74: . cs0 - LES constant value
75: . dpdx - Use DPDX property
76: . ndim - Number of dimensions
77: . nvar - Number of variables
78: . ncoupvar - Number of coupled variables
79: . contEqMethod - Continuity equation solving method
80: . inletFunction - Inlet function
81: . useLES - Flag to determine if the LES model should be used.
82: . useKEMOD - Flag to determine if the KE model should be used.
83: . useT - Flag to determine if the temperature variable should be used.
84: . firstBcond - Pointer to the first boundary condition of the list.
85: . firstBcond - Pointer to the first boundary condition of the list.
86: - lastBcond - Pointer to the last boundary condition of the list.
87:
88: Level: advanced
89: E*/
90: typedef struct {
91: PetscInt log_level, log_ievery, out_ievery, out_idtevery;
92: PetscInt startInstantMeans, startVarianceMeans;
93: PetscTruth initial_field, final_field;
94: PetscTruth printProps;
95: PetscInt max_outer_its, min_outer_its;
96: PetscReal resnorm_max, corrnorm_max;
97: RestartType k_restart;
98: PetscInt init_it;
99: PetscInt k_wrap;
100: PetscInt n_ghost;
101: PetscReal gravity;
102: PetscInt n_order_trans;
103: PetscReal deltat, tinit, tfinal;
104: PetscReal uref, mu, tref, rho, cs0;
105: PetscTruth dpdx;
106: PetscInt ndim, nvar, ncoupvar;
107: PetscInt contEqMethod;
108: PetscErrorCode (*inletFunction)(PetscInt*, PetscReal*);
109: PetscTruth useLES, useKEMOD, useT;
110: St_InputBcond *firstBcond, *lastBcond;
111: } St_Input;
112:
113: /*E
114: St_Grid - Structure for the grid.
115:
116: Attributes:
117: + Nnode - Number of nodes (in each direction) of the global grid.
118: . Xnode - Coordinates of the nodes of the global grid.
119: - ktype - Grid type (see ktypeGrid).
120:
121: Level: advanced
122: E*/
123: typedef struct {
124: PetscInt Nnode[3];
125: PetscReal *Xnode[3];
126: ktypeGrid ktype;
127: } St_Grid;
128:
129: /*E
130: St_Trans - Structure with parameters for solving unsteady problems.
131:
132: Attributes:
133: + norder - Transient order (see TransientOrder).
134: . time - Current time.
135: - idt - Current time step.
136:
137: Level: advanced
138: E*/
139: typedef struct {
140: TransientOrder norder;
141: PetscReal time;
142: PetscInt idt;
143: } St_Trans;
144:
145: /*E
146: MeshRegion - Structure containing the attributes of a mesh region.
147:
148: Attributes:
149: + nCells - Number of cells in each direction.
150: . lowerCorner - Coordinate of the lower corner of the mesh region.
151: . upperCorner - Coordinate of the upper corner of the mesh region.
152: - totalCells - Total number of cells.
153:
154: Level: advanced
155: E*/
156: typedef struct {
157: PetscInt nCells[3];
158: PetscInt lowerCorner[3];
159: PetscInt upperCorner[3];
160: PetscInt totalCells;
161: } MeshRegion;
162:
163: /*E
164: St_Patch - Structure for the definition of a patch.
165:
166: A patch is a rectangular region in computational domain, which is
167: used, among other things, for boundary conditions and properties
168: definition.
169:
170: The patch is identified by a unique patch ID, and characterized by
171: the global mesh region occupied. Moreover, the patch is distributed
172: among the processors by function CreatePatch, and this structure
173: provides also variables for the mesh region local to each processor
174: (with and without ghost nodes).
175:
176: The structure has also a pointer to the next patch.
177:
178: Attributes:
179: + id - Patch identifier.
180: . DA - Distributed array for the patch management.
181: . x, y, z - Lower corner.
182: . nx, ny, nz - Number of elements along each coordinate.
183: . comm - Communicator for the group of processes owning the patch.
184: - Next - Pointer to the next patch in the list.
185:
186: Level: advanced
187: E*/
188: typedef struct St_Patch {
189: PetscInt id;
190: DA da;
191: PetscInt x, y, z, nx, ny, nz;
192: MPI_Comm comm;
193: struct St_Patch *Next;
194: } St_Patch;
195:
196: typedef St_Patch *Ptr_Patch;
197:
198: /*E
199: St_Prop - Structure for the definition of a patch-wise property.
200:
201: A patch-wise property is a constant or a function defined over a
202: patch. A variable-property can depend on the geometry only, or
203: on other variables or transient terms.
204:
205: When the property is variable, the structure provides a
206: pointer to an array to store the values, and to a function to
207: evaluate them.
208:
209: Attributes:
210: + id - Property identifier.
211: . Patch - Patch where the property is defined.
212: . ktype - Property type (see ktypeProp).
213: . kip - Interpolation type for the value in the faces of the cells
214: (see ktypeIp).
215: . FluxLimiter - Flux limiter.
216: . linrlx - Linear relaxation of the property. Not automatic; it has to be
217: specified when the evaluation function is set.
218: . Cons - Value of the property when constant.
219: . Funct - Function for the evaluation of the property (when not constant).
220: . array - Array where the values of the property are stored.
221: . instantMeans - Vector for the means of the property instant values.
222: . instantSums - Vector for the sum of all property instant values.
223: . varianceMeans - Vector for the means of the property variance values.
224: . values - Vector for the property instant values.
225: . numReferences - Number of references (variables or boundary conditions
226: using this property).
227: - Next - Pointer to the next structure in the list.
228:
229: Level: advanced
230: E*/
231: typedef struct St_Prop {
232: PetscInt id;
233: Ptr_Patch Patch;
234: ktypeProp ktype;
235: StorageType storage;
236: ktypeIp kip;
237: PetscErrorCode (*FluxLimiter)(PetscInt icell, PetscInt ineighb,
238: PetscReal upwdelta, PetscReal dowdelta, PetscReal *value);
239: PetscReal linrlx;
240: PetscReal value;
241: PetscErrorCode (*Funct)(struct St_Prop *Prop);
242: PetscReal *array;
243: Vec instantMeans, instantSums, varianceSums, values;
244: PetscInt numReferences;
245: struct St_Prop *Next;
246: } St_Prop;
247:
248: typedef St_Prop *Ptr_Prop;
249:
250: /*E
251: St_BcondApplication - Structure with parameters of a boundary condition
252: application.
253:
254: A boundary condition application is defined by the equation and
255: variable indexes affected by the boundary condition and the properties
256: needed for the evaluation. Thus, a boundary condition (St_Bcond) is
257: formed by one or more applications (St_BcondApplication).
258:
259: Attributes:
260: + ivar - Index of the variable to which the boundary condition applies.
261: . ieq - Index of the equation to which the boundary condition applies.
262: . CoeffA - The coefficient C in the linear source C(V - Phi)
263: . CoeffB - The CV product in the linear source C(V - Phi)
264: - next - Pointer to the next structure in the list.
265:
266: Level: advanced
267:
268: .seealso St_Bcond
269: E*/
270: typedef struct St_BcondApplication {
271: PetscInt ivar, ieq;
272: Ptr_Prop CoeffA, CoeffB;
273: struct St_BcondApplication *next;
274: } St_BcondApplication;
275:
276: /*E
277: St_Bcond - Structure with parameters of a boundary condition.
278:
279: Attributes:
280: + id - Boundary condition identifier.
281: . ineighb - Neighbor index.
282: . Patch - Patch where the boundary condition is defined.
283: . bCondType - Boundary condition type (see BcondType).
284: . applicationList - List of applications of the boundary condition.
285: . Next - Pointer to the next structure in the list.
286:
287: Level: advanced
288:
289: .seealso St_BcondApplication
290: E*/
291: typedef struct St_Bcond {
292: PetscInt id, ineighb;
293: Ptr_Patch Patch;
294: BcondType bCondType;
295: St_BcondApplication *applicationList;
296: struct St_Bcond *Next;
297: } St_Bcond;
298:
299: typedef St_Bcond *Ptr_Bcond;
300:
301: /*E
302: St_Geom - Structure with geometric properties of the mesh.
303:
304: Attributes:
305: + Surface - Surface area of the faces of the local cells.
306: . Edge - Length of the finite volume in each direction.
307: . Volume - Cell volume.
308: . wallDistance - Distance to the nearest wall;
309: - normalVec - Unit vector, normal to the direction of the
310: vector from each cell to the nearest wall cell and following
311: the flow motion (i.e: trying to maximize X component and minimize
312: Y component).
313:
314: Level: advanced
315: E*/
316: typedef struct {
317: PetscReal *Surface[3], *Edge[3], *Volume;
318: PetscReal *wallDistance, *normalVec[3];
319: } St_Geom;
320:
321: /*E
322: St_Porosity - Structure with porosity parameters.
323:
324: Attributes:
325: + id - Porosity identifier.
326: . ineighb - Face to apply the porosity.
327: . Patch - Patch where the porosity is defined.
328: . Cons - Porosity value.
329: - Next - Pointer to the next structure in the list.
330:
331: Level: advanced
332: E*/
333: typedef struct St_Porosity {
334: PetscInt id, ineighb;
335: Ptr_Patch Patch;
336: PetscReal Cons;
337: struct St_Porosity *Next;
338: } St_Porosity;
339:
340: typedef St_Porosity *Ptr_Porosity;
341:
342: /*E
343: St_Var - Structure with parameters of a dependent variable.
344:
345: Attributes:
346: + name - Name of the variable.
347: . kip - Interpolation type to the value in the faces of the cells
348: (see ktypeIp).
349: . FluxLimiter - Flux limiter.
350: . valmin, valmax - Minimum and maximum values of the variable.
351: . linrlx, fdtrlx - Relaxation values of the variable.
352: . resref, corref - Reference values of the variable (residue and
353: correction vector).
354: . Valini - Property for the initial values of the variable.
355: . Gamma - Property for the diffusion coefficient of the variable.
356: . locPhi - Local unknown vector.
357: . Array - Pointer to the array of locPhi vector object.
358: . resnorm, corrnorm - Norm values of the variable (residue and
359: correction vector).
360: - next - Pointer to the next structure in the list.
361:
362: Level: advanced
363: E*/
364: typedef struct St_Var {
365: char name[128];
366: PetscInt index;
367: ktypeIp kip;
368: PetscErrorCode (*FluxLimiter)(PetscInt icell, PetscInt ineighb,
369: PetscReal upwdelta, PetscReal dowdelta, PetscReal *value);
370: PetscReal valmin, valmax, linrlx, fdtrlx, resref, corref,
371: resnorm, corrnorm;
372: Ptr_Prop Valini, Gamma;
373: Vec locPhi;
374: PetscReal *Array;
375: struct St_Var *next;
376: } St_Var;
377:
378: typedef St_Var *Ptr_Var;
379:
380: /*E
381: St_Sys - Structure with PETSc objects for the distributed linear
382: system. Different instances of the same structure are used for the
383: coupled and segregated systems.
384:
385: Attributes:
386: + Phi - Unknown vector.
387: . B - Right side vector.
388: . A - Coefficient matrix.
389: . da - Distributed array.
390: . solver - Linear solver.
391: . ltog - Local to global indexes of the cells.
392: . ndof - Degrees of freedom (Input.ncoupvar for coupled system, 1 for
393: segregated system).
394: . ivar - Index of the first variable of the system. The last variable
395: is ivar + ndof - 1
396: . resnorm - Norm of the residue.
397: . resref - Reference value for the residue.
398: . corrnorm - Norm of the correction vector.
399: - corref - Reference value for the correction vector.
400:
401: Level: advanced
402: E*/
403: typedef struct {
404: Vec Phi, B;
405: Mat A;
406: DA da;
407: KSP solver;
408: PetscInt *ltog;
409: PetscInt ndof;
410: PetscInt ivar;
411: PetscReal resnorm, resref, corrnorm, corref;
412: } St_Sys;
413:
414: /*E
415: St_Local - Structure with local parameters of the problem.
416:
417: Attributes:
418: + rank - Processor's rank.
419: . meshRegion - Local subdomain.
420: . ghostMeshRegion - Local subdomain (including ghost nodes).
421: . icell_neighb - Local indexes of each neighbor for each cell.
422: . porosity - Porosity of each face for each cell.
423: - FirstPorosity, LastPorosity - Pointers to the first and last
424: porosities in the list.
425:
426: Level: advanced
427: E*/
428: typedef struct {
429: PetscMPIInt rank;
430: MeshRegion meshRegion, ghostMeshRegion;
431: PetscInt **icell_neighb;
432: PetscReal **porosity;
433: St_Porosity *FirstPorosity, *LastPorosity;
434: } St_Local;
435:
436: /*E
437: St_Interp - Structure with parameters to interpolate values from
438: nodes to faces, according to the convective scheme (lineal,
439: upwind). It is necessary to re-calculate them after each outer
440: iteration since they depend on the flow direction, except for
441: the linear interpolation.
442:
443: Attributes:
444: + linear - Coefficients for the linear interpolation.
445: . upwind - Coefficients for the upwind interpolation.
446: . icell_upw - Local index of the upwind cell, referring to the
447: cell (notation from the paper about Flux Limiter
448: (Waterson).
449: . icell_uupw - Local index of the up-upwind cell, referring to the
450: cell (notation from the paper about Flux Limiter
451: (Waterson).
452: . icell_dow - Local index of the downstream cell, referring to the
453: cell (notation from the paper about Flux Limiter
454: (Waterson).
455:
456: Level: advanced
457: E*/
458: typedef struct {
459: PetscReal **linear, **upwind;
460: PetscInt **icell_upw, **icell_uupw, **icell_dow;
461: } St_Interp;
462:
463: /*E
464: St_Global - Structure for the global parameters of the problem.
465:
466: Attributes:
467: + comm - MPI communicator.
468: . NPatch - Global number of patches.
469: . NBcond - Global number of boundary conditions.
470: . NPorosity - Global number of porosities.
471: . NProp - Global number of properties.
472: . FirstPatch, LastPatch - Pointer to the first and last patches
473: of the list.
474: . FirstBcond, LastBcond - Pointer to the first and last boundary
475: conditions of the list.
476: . FirstProp, LastProp - Pointer to the first and last properties
477: of the list.
478:
479: Level: advanced
480: E*/
481: typedef struct {
482: MPI_Comm comm;
483: PetscInt NPatch, NBcond, NPorosity, NProp;
484: St_Patch *FirstPatch, *LastPatch;
485: St_Bcond *FirstBcond, *LastBcond;
486: St_Prop *FirstProp, *LastProp;
487: } St_Global;
488:
489: #endif