Actual source code: dmstagimpl.h
1: #pragma once
3: #include <petscdmstag.h>
4: #include <petsc/private/dmimpl.h>
6: #define DMSTAG_MAX_DIM 3
7: #define DMSTAG_MAX_STRATA DMSTAG_MAX_DIM + 1
9: /* This value is 1 + 3^DMSTAG_MAX_DIM */
10: #define DMSTAG_NUMBER_LOCATIONS 28
12: typedef struct {
13: /* Fields which may require being set before DMSetUp() is called, set by DMStagInitialize().
14: Some may be adjusted by DMSetUp() */
15: PetscInt N[DMSTAG_MAX_DIM]; /* Global dimensions (elements) */
16: PetscInt n[DMSTAG_MAX_DIM]; /* Local dimensions (elements) */
17: PetscInt *l[DMSTAG_MAX_DIM]; /* Elements/rank in each direction */
18: PetscInt dof[DMSTAG_MAX_STRATA]; /* Dof per point for each stratum */
19: DMStagStencilType stencilType; /* Elementwise stencil type */
20: PetscInt stencilWidth; /* Elementwise ghost width */
21: DMBoundaryType boundaryType[DMSTAG_MAX_DIM]; /* Physical domain ghosting type */
22: PetscMPIInt nRanks[DMSTAG_MAX_DIM]; /* Ranks in each direction */
24: /* Fields unrelated to setup */
25: DMType coordinateDMType; /* DM type to create for coordinates */
26: PetscInt refineFactor[DMSTAG_MAX_DIM]; /* Ratio used in refining and coarsening */
28: /* Data above is copied by DMStagDuplicateWithoutSetup(), while data below is not */
30: /* Fields populated by DMSetUp() */
31: PetscInt nGhost[DMSTAG_MAX_DIM]; /* Local dimensions (w/ ghosts) */
32: PetscInt start[DMSTAG_MAX_DIM]; /* First element number */
33: PetscInt startGhost[DMSTAG_MAX_DIM]; /* First element number (w/ ghosts) */
34: PetscMPIInt rank[DMSTAG_MAX_DIM]; /* Location in grid of ranks */
35: PetscMPIInt *neighbors; /* dim^3 local ranks */
36: VecScatter gtol; /* Global --> Local */
37: VecScatter ltog_injective; /* Local --> Global, injective */
38: VecScatter ltol; /* Local --> Local */
39: PetscInt *locationOffsets; /* Offsets for points in loc. rep. */
41: /* Additional convenience fields populated by DMSetUp() (easily computed from the above) */
42: PetscInt entriesPerElement; /* Entries stored with each element */
43: PetscInt entries; /* Local number of entries */
44: PetscInt entriesGhost; /* Local numbers of entries w/ ghosts */
45: PetscBool firstRank[DMSTAG_MAX_DIM]; /* First rank in this dim? */
46: PetscBool lastRank[DMSTAG_MAX_DIM]; /* Last rank in this dim? */
47: } DM_Stag;
49: PETSC_INTERN PetscErrorCode DMCreateMatrix_Stag_1D_AIJ_Assemble(DM, Mat);
50: PETSC_INTERN PetscErrorCode DMCreateMatrix_Stag_2D_AIJ_Assemble(DM, Mat);
51: PETSC_INTERN PetscErrorCode DMCreateMatrix_Stag_3D_AIJ_Assemble(DM, Mat);
52: PETSC_INTERN PetscErrorCode DMStagDuplicateWithoutSetup(DM, MPI_Comm, DM *);
53: PETSC_INTERN PetscErrorCode DMStagInitialize(DMBoundaryType, DMBoundaryType, DMBoundaryType, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, DMStagStencilType, PetscInt, const PetscInt[], const PetscInt[], const PetscInt[], DM);
54: PETSC_INTERN PetscErrorCode DMSetUp_Stag_1d(DM);
55: PETSC_INTERN PetscErrorCode DMSetUp_Stag_2d(DM);
56: PETSC_INTERN PetscErrorCode DMSetUp_Stag_3d(DM);
57: PETSC_INTERN PetscErrorCode DMStagRestrictSimple_1d(DM, Vec, DM, Vec);
58: PETSC_INTERN PetscErrorCode DMStagRestrictSimple_2d(DM, Vec, DM, Vec);
59: PETSC_INTERN PetscErrorCode DMStagRestrictSimple_3d(DM, Vec, DM, Vec);
60: PETSC_INTERN PetscErrorCode DMStagPopulateInterpolation1d_Internal(DM, DM, Mat);
61: PETSC_INTERN PetscErrorCode DMStagPopulateInterpolation2d_Internal(DM, DM, Mat);
62: PETSC_INTERN PetscErrorCode DMStagPopulateInterpolation3d_Internal(DM, DM, Mat);
63: PETSC_INTERN PetscErrorCode DMStagPopulateLocalToGlobalInjective_1d(DM);
64: PETSC_INTERN PetscErrorCode DMStagPopulateLocalToGlobalInjective_2d(DM);
65: PETSC_INTERN PetscErrorCode DMStagPopulateLocalToGlobalInjective_3d(DM);
66: PETSC_INTERN PetscErrorCode DMStagPopulateLocalToLocal1d_Internal(DM);
67: PETSC_INTERN PetscErrorCode DMStagPopulateLocalToLocal2d_Internal(DM);
68: PETSC_INTERN PetscErrorCode DMStagPopulateLocalToLocal3d_Internal(DM);
69: PETSC_INTERN PetscErrorCode DMStagPopulateRestriction1d_Internal(DM, DM, Mat);
70: PETSC_INTERN PetscErrorCode DMStagPopulateRestriction2d_Internal(DM, DM, Mat);
71: PETSC_INTERN PetscErrorCode DMStagPopulateRestriction3d_Internal(DM, DM, Mat);
72: PETSC_INTERN PetscErrorCode DMStagSetUniformCoordinatesExplicit_1d(DM, PetscReal, PetscReal);
73: PETSC_INTERN PetscErrorCode DMStagSetUniformCoordinatesExplicit_2d(DM, PetscReal, PetscReal, PetscReal, PetscReal);
74: PETSC_INTERN PetscErrorCode DMStagSetUniformCoordinatesExplicit_3d(DM, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal);
75: PETSC_INTERN PetscErrorCode DMStagStencilLocationCanonicalize(DMStagStencilLocation, DMStagStencilLocation *);