Actual source code: hypre.c
1: /*
2: Provides an interface to the LLNL package hypre
3: */
4: #include <petsc/private/petscimpl.h>
5: #include <petscpc.h>
7: /*@
8: PCHYPRESetDiscreteGradient - Set the discrete gradient matrix for `PCHYPRE` type of AMS or ADS
10: Collective
12: Input Parameters:
13: + pc - the preconditioning context
14: - G - the discrete gradient
16: Level: intermediate
18: Notes:
19: `G` should have as many rows as the number of edges and as many columns as the number of vertices in the mesh
21: Each row of `G` has 2 nonzeros, with column indexes being the global indexes of edge's endpoints: matrix entries are +1 and -1 depending on edge orientation
23: Developer Note:
24: This automatically converts the matrix to `MATHYPRE` if it is not already of that type
26: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteCurl()`
27: @*/
28: PetscErrorCode PCHYPRESetDiscreteGradient(PC pc, Mat G)
29: {
30: PetscFunctionBegin;
33: PetscCheckSameComm(pc, 1, G, 2);
34: PetscTryMethod(pc, "PCHYPRESetDiscreteGradient_C", (PC, Mat), (pc, G));
35: PetscFunctionReturn(PETSC_SUCCESS);
36: }
38: /*@
39: PCHYPRESetDiscreteCurl - Set the discrete curl matrix for `PCHYPRE` type of ADS
41: Collective
43: Input Parameters:
44: + pc - the preconditioning context
45: - C - the discrete curl
47: Level: intermediate
49: Notes:
50: `C` should have as many rows as the number of faces and as many columns as the number of edges in the mesh
52: Each row of `C` has as many nonzeros as the number of edges of a face, with column indexes being the global indexes of the corresponding edge.
53: Matrix entries are +1 and -1 depending on edge orientation with respect to the face orientation
55: Developer Notes:
56: This automatically converts the matrix to `MATHYPRE` if it is not already of that type
58: If this is only for `PCHYPRE` type of ADS it should be called `PCHYPREADSSetDiscreteCurl()`
60: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteGradient()`
61: @*/
62: PetscErrorCode PCHYPRESetDiscreteCurl(PC pc, Mat C)
63: {
64: PetscFunctionBegin;
67: PetscCheckSameComm(pc, 1, C, 2);
68: PetscTryMethod(pc, "PCHYPRESetDiscreteCurl_C", (PC, Mat), (pc, C));
69: PetscFunctionReturn(PETSC_SUCCESS);
70: }
72: /*@
73: PCHYPRESetInterpolations - Set the interpolation matrices for `PCHYPRE` type of AMS or ADS
75: Collective
77: Input Parameters:
78: + pc - the preconditioning context
79: . dim - the dimension of the problem, only used in AMS
80: . RT_PiFull - Raviart-Thomas interpolation matrix
81: . RT_Pi - x/y/z component of Raviart-Thomas interpolation matrix
82: . ND_PiFull - Nedelec interpolation matrix
83: - ND_Pi - x/y/z component of Nedelec interpolation matrix
85: Level: intermediate
87: Notes:
88: For AMS, only Nedelec interpolation matrices are needed, the Raviart-Thomas interpolation matrices can be set to `NULL`.
90: For ADS, both type of interpolation matrices are needed.
92: Developer Note:
93: This automatically converts the matrix to `MATHYPRE` if it is not already of that type
95: .seealso: [](ch_ksp), `PCHYPRE`
96: @*/
97: PetscErrorCode PCHYPRESetInterpolations(PC pc, PetscInt dim, Mat RT_PiFull, Mat RT_Pi[], Mat ND_PiFull, Mat ND_Pi[])
98: {
99: PetscInt i;
101: PetscFunctionBegin;
103: if (RT_PiFull) {
105: PetscCheckSameComm(pc, 1, RT_PiFull, 3);
106: }
107: if (RT_Pi) {
108: PetscAssertPointer(RT_Pi, 4);
109: for (i = 0; i < dim; ++i) {
110: if (RT_Pi[i]) {
112: PetscCheckSameComm(pc, 1, RT_Pi[i], 4);
113: }
114: }
115: }
116: if (ND_PiFull) {
118: PetscCheckSameComm(pc, 1, ND_PiFull, 5);
119: }
120: if (ND_Pi) {
121: PetscAssertPointer(ND_Pi, 6);
122: for (i = 0; i < dim; ++i) {
123: if (ND_Pi[i]) {
125: PetscCheckSameComm(pc, 1, ND_Pi[i], 6);
126: }
127: }
128: }
129: PetscTryMethod(pc, "PCHYPRESetInterpolations_C", (PC, PetscInt, Mat, Mat[], Mat, Mat[]), (pc, dim, RT_PiFull, RT_Pi, ND_PiFull, ND_Pi));
130: PetscFunctionReturn(PETSC_SUCCESS);
131: }
133: /*@
134: PCHYPRESetAlphaPoissonMatrix - Set the vector Poisson matrix for `PCHYPRE` of type AMS
136: Collective
138: Input Parameters:
139: + pc - the preconditioning context
140: - A - the matrix
142: Level: intermediate
144: Note:
145: `A` should be obtained by discretizing the vector valued Poisson problem with linear finite elements
147: Developer Notes:
148: This automatically converts the matrix to `MATHYPRE` if it is not already of that type
150: If this is only for `PCHYPRE` type of AMS it should be called `PCHYPREAMSSetAlphaPoissonMatrix()`
152: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetBetaPoissonMatrix()`
153: @*/
154: PetscErrorCode PCHYPRESetAlphaPoissonMatrix(PC pc, Mat A)
155: {
156: PetscFunctionBegin;
159: PetscCheckSameComm(pc, 1, A, 2);
160: PetscTryMethod(pc, "PCHYPRESetPoissonMatrix_C", (PC, Mat, PetscBool), (pc, A, PETSC_TRUE));
161: PetscFunctionReturn(PETSC_SUCCESS);
162: }
164: /*@
165: PCHYPRESetBetaPoissonMatrix - Set the Poisson matrix for `PCHYPRE` of type AMS
167: Collective
169: Input Parameters:
170: + pc - the preconditioning context
171: - A - the matrix, or `NULL` to turn it off
173: Level: intermediate
175: Note:
176: `A` should be obtained by discretizing the Poisson problem with linear finite elements.
178: Developer Notes:
179: This automatically converts the matrix to `MATHYPRE` if it is not already of that type
181: If this is only for `PCHYPRE` type of AMS it should be called `PCHYPREAMSPCHYPRESetBetaPoissonMatrix()`
183: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetAlphaPoissonMatrix()`
184: @*/
185: PetscErrorCode PCHYPRESetBetaPoissonMatrix(PC pc, Mat A)
186: {
187: PetscFunctionBegin;
189: if (A) {
191: PetscCheckSameComm(pc, 1, A, 2);
192: }
193: PetscTryMethod(pc, "PCHYPRESetPoissonMatrix_C", (PC, Mat, PetscBool), (pc, A, PETSC_FALSE));
194: PetscFunctionReturn(PETSC_SUCCESS);
195: }
197: /*@
198: PCHYPRESetEdgeConstantVectors - Set the representation of the constant vector fields in the edge element basis for `PCHYPRE` of type AMS
200: Collective
202: Input Parameters:
203: + pc - the preconditioning context
204: . ozz - vector representing (1,0,0) (or (1,0) in 2D)
205: . zoz - vector representing (0,1,0) (or (0,1) in 2D)
206: - zzo - vector representing (0,0,1) (use NULL in 2D)
208: Level: intermediate
210: Developer Note:
211: If this is only for `PCHYPRE` type of AMS it should be called `PCHYPREAMSSetEdgeConstantVectors()`
213: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetAlphaPoissonMatrix()`
214: @*/
215: PetscErrorCode PCHYPRESetEdgeConstantVectors(PC pc, Vec ozz, Vec zoz, Vec zzo)
216: {
217: PetscFunctionBegin;
222: PetscCheckSameComm(pc, 1, ozz, 2);
223: PetscCheckSameComm(pc, 1, zoz, 3);
224: if (zzo) PetscCheckSameComm(pc, 1, zzo, 4);
225: PetscTryMethod(pc, "PCHYPRESetEdgeConstantVectors_C", (PC, Vec, Vec, Vec), (pc, ozz, zoz, zzo));
226: PetscFunctionReturn(PETSC_SUCCESS);
227: }
229: /*@
230: PCHYPREAMSSetInteriorNodes - Set the list of interior nodes to a zero-conductivity region for `PCHYPRE` of type AMS
232: Collective
234: Input Parameters:
235: + pc - the preconditioning context
236: - interior - vector. node is interior if its entry in the array is 1.0.
238: Level: intermediate
240: Note:
241: This calls `HYPRE_AMSSetInteriorNodes()`
243: .seealso: [](ch_ksp), `PCHYPRE`, `PCHYPRESetDiscreteGradient()`, `PCHYPRESetDiscreteCurl()`, `PCHYPRESetAlphaPoissonMatrix()`
244: @*/
245: PetscErrorCode PCHYPREAMSSetInteriorNodes(PC pc, Vec interior)
246: {
247: PetscFunctionBegin;
250: PetscCheckSameComm(pc, 1, interior, 2);
251: PetscTryMethod(pc, "PCHYPREAMSSetInteriorNodes_C", (PC, Vec), (pc, interior));
252: PetscFunctionReturn(PETSC_SUCCESS);
253: }
255: /*@
256: PCHYPRESetType - Sets which hypre preconditioner you wish to use
258: Input Parameters:
259: + pc - the preconditioner context
260: - name - either euclid, ilu, pilut, parasails, boomeramg, ams, or ads
262: Options Database Key:
263: . pc_hypre_type - One of euclid, ilu, pilut, parasails, boomeramg, ams, or ads
265: Level: intermediate
267: .seealso: [](ch_ksp), `PCCreate()`, `PCSetType()`, `PCType`, `PC`, `PCHYPRE`
268: @*/
269: PetscErrorCode PCHYPRESetType(PC pc, const char name[])
270: {
271: PetscFunctionBegin;
273: PetscAssertPointer(name, 2);
274: PetscTryMethod(pc, "PCHYPRESetType_C", (PC, const char[]), (pc, name));
275: PetscFunctionReturn(PETSC_SUCCESS);
276: }
278: /*@C
279: PCHYPREGetCFMarkers - Gets CF marker arrays for all levels (except the finest level)
281: Logically Collective
283: Input Parameter:
284: . pc - the preconditioner context
286: Output Parameters:
287: + n_per_level - the number of nodes per level (size of `num_levels`)
288: - CFMarkers - the Coarse/Fine Boolean arrays (size of `num_levels` - 1)
290: Level: advanced
292: Note:
293: Caller is responsible for memory management of `n_per_level` and `CFMarkers` pointers. That is they should free them with `PetscFree()` when no longer needed.
295: .seealso: [](ch_ksp), `PC`, `PCMG`, `PCMGGetRestriction()`, `PCMGSetInterpolation()`, `PCMGGetRScale()`, `PCMGGetInterpolation()`, `PCGetInterpolations()`
296: @*/
297: PetscErrorCode PCHYPREGetCFMarkers(PC pc, PetscInt *n_per_level[], PetscBT *CFMarkers[])
298: {
299: PetscFunctionBegin;
301: PetscAssertPointer(n_per_level, 2);
302: PetscAssertPointer(CFMarkers, 3);
303: PetscUseMethod(pc, "PCHYPREGetCFMarkers_C", (PC, PetscInt *[], PetscBT *[]), (pc, n_per_level, CFMarkers));
304: PetscFunctionReturn(PETSC_SUCCESS);
305: }
307: /*@
308: PCHYPREGetType - Gets which hypre preconditioner you are using
310: Input Parameter:
311: . pc - the preconditioner context
313: Output Parameter:
314: . name - either euclid, ilu, pilut, parasails, boomeramg, ams, or ads
316: Level: intermediate
318: .seealso: [](ch_ksp), `PCCreate()`, `PCHYPRESetType()`, `PCType`, `PC`, `PCHYPRE`
319: @*/
320: PetscErrorCode PCHYPREGetType(PC pc, const char *name[])
321: {
322: PetscFunctionBegin;
324: PetscAssertPointer(name, 2);
325: PetscTryMethod(pc, "PCHYPREGetType_C", (PC, const char *[]), (pc, name));
326: PetscFunctionReturn(PETSC_SUCCESS);
327: }
329: /*@
330: PCMGGalerkinSetMatProductAlgorithm - Set type of sparse matrix-matrix product for hypre's BoomerAMG to use on GPUs
332: Logically Collective
334: Input Parameters:
335: + pc - the hypre context
336: - name - one of 'cusparse', 'hypre'
338: Options Database Key:
339: . -pc_mg_galerkin_mat_product_algorithm <cusparse,hypre> - Type of sparse matrix-matrix product to use in hypre
341: Level: intermediate
343: Developer Note:
344: How the name starts with `PCMG`, should it not be `PCHYPREBoomerAMG`?
346: .seealso: [](ch_ksp), `PCHYPRE`, `PCMGGalerkinGetMatProductAlgorithm()`
347: @*/
348: PetscErrorCode PCMGGalerkinSetMatProductAlgorithm(PC pc, const char name[])
349: {
350: PetscFunctionBegin;
352: PetscTryMethod(pc, "PCMGGalerkinSetMatProductAlgorithm_C", (PC, const char[]), (pc, name));
353: PetscFunctionReturn(PETSC_SUCCESS);
354: }
356: /*@
357: PCMGGalerkinGetMatProductAlgorithm - Get type of sparse matrix-matrix product for hypre's BoomerAMG to use on GPUs
359: Not Collective
361: Input Parameter:
362: . pc - the multigrid context
364: Output Parameter:
365: . name - one of 'cusparse', 'hypre'
367: Level: intermediate
369: .seealso: [](ch_ksp), `PCHYPRE`, `PCMGGalerkinSetMatProductAlgorithm()`
370: @*/
371: PetscErrorCode PCMGGalerkinGetMatProductAlgorithm(PC pc, const char *name[])
372: {
373: PetscFunctionBegin;
375: PetscTryMethod(pc, "PCMGGalerkinGetMatProductAlgorithm_C", (PC, const char *[]), (pc, name));
376: PetscFunctionReturn(PETSC_SUCCESS);
377: }