Actual source code: matreg.c
petsc-3.14.6 2021-03-30
2: /*
3: Mechanism for register PETSc matrix types
4: */
5: #include <petsc/private/matimpl.h>
7: PetscBool MatRegisterAllCalled = PETSC_FALSE;
9: /*
10: Contains the list of registered Mat routines
11: */
12: PetscFunctionList MatList = NULL;
14: /*@C
15: MatSetType - Builds matrix object for a particular matrix type
17: Collective on Mat
19: Input Parameters:
20: + mat - the matrix object
21: - matype - matrix type
23: Options Database Key:
24: . -mat_type <method> - Sets the type; use -help for a list
25: of available methods (for instance, seqaij)
27: Notes:
28: See "${PETSC_DIR}/include/petscmat.h" for available methods
30: Level: intermediate
32: .seealso: PCSetType(), VecSetType(), MatCreate(), MatType, Mat
33: @*/
34: PetscErrorCode MatSetType(Mat mat, MatType matype)
35: {
36: PetscErrorCode ierr,(*r)(Mat);
37: PetscBool sametype,found,subclass = PETSC_FALSE;
38: MatRootName names = MatRootNameList;
43: while (names) {
44: PetscStrcmp(matype,names->rname,&found);
45: if (found) {
46: PetscMPIInt size;
47: MPI_Comm_size(PetscObjectComm((PetscObject)mat),&size);
48: if (size == 1) matype = names->sname;
49: else matype = names->mname;
50: break;
51: }
52: names = names->next;
53: }
55: PetscObjectTypeCompare((PetscObject)mat,matype,&sametype);
56: if (sametype) return(0);
58: PetscFunctionListFind(MatList,matype,&r);
59: if (!r) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown Mat type given: %s",matype);
61: if (mat->assembled && ((PetscObject)mat)->type_name) {
62: PetscStrbeginswith(matype,((PetscObject)mat)->type_name,&subclass);
63: }
64: if (subclass) {
65: MatConvert(mat,matype,MAT_INPLACE_MATRIX,&mat);
66: return(0);
67: }
68: if (mat->ops->destroy) {
69: /* free the old data structure if it existed */
70: (*mat->ops->destroy)(mat);
71: mat->ops->destroy = NULL;
73: /* should these null spaces be removed? */
74: MatNullSpaceDestroy(&mat->nullsp);
75: MatNullSpaceDestroy(&mat->nearnullsp);
76: }
77: PetscMemzero(mat->ops,sizeof(struct _MatOps));
78: mat->preallocated = PETSC_FALSE;
79: mat->assembled = PETSC_FALSE;
80: mat->was_assembled = PETSC_FALSE;
82: /*
83: Increment, rather than reset these: the object is logically the same, so its logging and
84: state is inherited. Furthermore, resetting makes it possible for the same state to be
85: obtained with a different structure, confusing the PC.
86: */
87: mat->nonzerostate++;
88: PetscObjectStateIncrease((PetscObject)mat);
90: /* create the new data structure */
91: (*r)(mat);
92: return(0);
93: }
95: /*@C
96: MatGetType - Gets the matrix type as a string from the matrix object.
98: Not Collective
100: Input Parameter:
101: . mat - the matrix
103: Output Parameter:
104: . name - name of matrix type
106: Level: intermediate
108: .seealso: MatSetType()
109: @*/
110: PetscErrorCode MatGetType(Mat mat,MatType *type)
111: {
115: *type = ((PetscObject)mat)->type_name;
116: return(0);
117: }
119: /*@C
120: MatGetVecType - Gets the vector type used by the matrix object.
122: Not Collective
124: Input Parameter:
125: . mat - the matrix
127: Output Parameter:
128: . name - name of vector type
130: Level: intermediate
132: .seealso: MatSetVecType()
133: @*/
134: PetscErrorCode MatGetVecType(Mat mat,VecType *vtype)
135: {
139: *vtype = mat->defaultvectype;
140: return(0);
141: }
143: /*@C
144: MatSetVecType - Set the vector type to be used for a matrix object
146: Collective on Mat
148: Input Parameters:
149: + mat - the matrix object
150: - vtype - vector type
152: Notes:
153: This is rarely needed in practice since each matrix object internally sets the proper vector type.
155: Level: intermediate
157: .seealso: VecSetType(), MatGetVecType()
158: @*/
159: PetscErrorCode MatSetVecType(Mat mat,VecType vtype)
160: {
165: PetscFree(mat->defaultvectype);
166: PetscStrallocpy(vtype,&mat->defaultvectype);
167: return(0);
168: }
170: /*@C
171: MatRegister - - Adds a new matrix type
173: Not Collective
175: Input Parameters:
176: + name - name of a new user-defined matrix type
177: - routine_create - routine to create method context
179: Notes:
180: MatRegister() may be called multiple times to add several user-defined solvers.
182: Sample usage:
183: .vb
184: MatRegister("my_mat",MyMatCreate);
185: .ve
187: Then, your solver can be chosen with the procedural interface via
188: $ MatSetType(Mat,"my_mat")
189: or at runtime via the option
190: $ -mat_type my_mat
192: Level: advanced
194: .seealso: MatRegisterAll()
197: Level: advanced
198: @*/
199: PetscErrorCode MatRegister(const char sname[],PetscErrorCode (*function)(Mat))
200: {
204: MatInitializePackage();
205: PetscFunctionListAdd(&MatList,sname,function);
206: return(0);
207: }
209: MatRootName MatRootNameList = NULL;
211: /*@C
212: MatRegisterRootName - Registers a name that can be used for either a sequential or its corresponding parallel matrix type. MatSetType()
213: and -mat_type will automatically use the sequential or parallel version based on the size of the MPI communicator associated with the
214: matrix.
216: Input Parameters:
217: + rname - the rootname, for example, MATAIJ
218: . sname - the name of the sequential matrix type, for example, MATSEQAIJ
219: - mname - the name of the parallel matrix type, for example, MATMPIAIJ
221: Notes: The matrix rootname should not be confused with the base type of the function PetscObjectBaseTypeCompare()
223: Developer Notes: PETSc vectors have a similar rootname that indicates PETSc should automatically select the appropriate VecType based on the
224: size of the communicator but it is implemented by simply having additional VecCreate_RootName() registerer routines that dispatch to the
225: appropriate creation routine. Why have two different ways of implementing the same functionality for different types of objects? It is
226: confusing.
228: Level: developer
230: .seealso: PetscObjectBaseTypeCompare()
232: @*/
233: PetscErrorCode MatRegisterRootName(const char rname[],const char sname[],const char mname[])
234: {
236: MatRootName names;
239: PetscNew(&names);
240: PetscStrallocpy(rname,&names->rname);
241: PetscStrallocpy(sname,&names->sname);
242: PetscStrallocpy(mname,&names->mname);
243: if (!MatRootNameList) {
244: MatRootNameList = names;
245: } else {
246: MatRootName next = MatRootNameList;
247: while (next->next) next = next->next;
248: next->next = names;
249: }
250: return(0);
251: }