Actual source code: gcreate.c
1: /*$Id: gcreate.c,v 1.127 2001/03/23 23:22:45 balay Exp $*/
3: #include petscsys.h
4: #include src/mat/matimpl.h
6: static int MatPublish_Base(PetscObject obj)
7: {
8: #if defined(PETSC_HAVE_AMS)
9: Mat mat = (Mat)obj;
11: #endif
14: #if defined(PETSC_HAVE_AMS)
15: /* if it is already published then return */
16: if (mat->amem >=0) return(0);
18: PetscObjectPublishBaseBegin(obj);
19: AMS_Memory_add_field((AMS_Memory)mat->amem,"globalsizes",&mat->M,2,AMS_INT,AMS_READ,
20: AMS_COMMON,AMS_REDUCT_UNDEF);
21: AMS_Memory_add_field((AMS_Memory)mat->amem,"localsizes",&mat->m,2,AMS_INT,AMS_READ,
22: AMS_DISTRIBUTED,AMS_REDUCT_UNDEF);
23: PetscObjectPublishBaseEnd(obj);
24: #endif
26: return(0);
27: }
30: /*@C
31: MatCreate - Creates a matrix where the type is determined
32: from the options database. Generates a parallel MPI matrix if the
33: communicator has more than one processor. The default matrix type is
34: AIJ, using the routines MatCreateSeqAIJ() and MatCreateMPIAIJ().
36: Collective on MPI_Comm
38: Input Parameters:
39: + m - number of local rows (or PETSC_DECIDE)
40: . n - number of local columns (or PETSC_DECIDE)
41: . M - number of global rows (or PETSC_DETERMINE)
42: . N - number of global columns (or PETSC_DETERMINE)
43: - comm - MPI communicator
44:
45: Output Parameter:
46: . A - the matrix
48: Options Database Keys:
49: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
50: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
51: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
52: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
53: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
54: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
55: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
56: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
57: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
59: Even More Options Database Keys:
60: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
61: for additional format-specific options.
63: Notes:
64: If PETSC_DECIDE is not used for the arguments 'm' and 'n', then the
65: user must ensure that they are chosen to be compatible with the
66: vectors. To do this, one first considers the matrix-vector product
67: 'y = A x'. The 'm' that is used in the above routine must match the
68: local size used in the vector creation routine VecCreateMPI() for 'y'.
69: Likewise, the 'n' used must match that used as the local size in
70: VecCreateMPI() for 'x'.
72: Level: beginner
74: .keywords: matrix, create
76: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
77: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
78: MatCreateSeqDense(), MatCreateMPIDense(),
79: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
80: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
81: MatConvert()
82: @*/
83: int MatCreate(MPI_Comm comm,int m,int n,int M,int N,Mat *A)
84: {
85: Mat B;
88: PetscHeaderCreate(B,_p_Mat,struct _MatOps,MAT_COOKIE,0,"Mat",comm,MatDestroy,MatView);
89: PetscLogObjectCreate(B);
91: B->m = m;
92: B->n = n;
93: B->M = M;
94: B->N = N;
96: B->preallocated = PETSC_FALSE;
97: B->bops->publish = MatPublish_Base;
98: *A = B;
99: return(0);
100: }
102: /*@C
103: MatSetFromOptions - Creates a matrix where the type is determined
104: from the options database. Generates a parallel MPI matrix if the
105: communicator has more than one processor. The default matrix type is
106: AIJ, using the routines MatSetFromOptionsSeqAIJ() and MatSetFromOptionsMPIAIJ().
108: Collective on Mat
110: Input Parameter:
111: . A - the matrix
113: Options Database Keys:
114: + -mat_type seqaij - AIJ type, uses MatCreateSeqAIJ()
115: . -mat_type mpiaij - AIJ type, uses MatCreateMPIAIJ()
116: . -mat_type seqbdiag - block diagonal type, uses MatCreateSeqBDiag()
117: . -mat_type mpibdiag - block diagonal type, uses MatCreateMPIBDiag()
118: . -mat_type mpirowbs - rowbs type, uses MatCreateMPIRowbs()
119: . -mat_type seqdense - dense type, uses MatCreateSeqDense()
120: . -mat_type mpidense - dense type, uses MatCreateMPIDense()
121: . -mat_type seqbaij - block AIJ type, uses MatCreateSeqBAIJ()
122: - -mat_type mpibaij - block AIJ type, uses MatCreateMPIBAIJ()
124: Even More Options Database Keys:
125: See the manpages for particular formats (e.g., MatCreateSeqAIJ())
126: for additional format-specific options.
128: Level: beginner
130: .keywords: matrix, create
132: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
133: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
134: MatCreateSeqDense(), MatCreateMPIDense(),
135: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
136: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
137: MatConvert()
138: @*/
139: int MatSetFromOptions(Mat B)
140: {
141: int ierr,size;
142: char mtype[256];
143: PetscTruth flg;
146: PetscOptionsGetString(B->prefix,"-mat_type",mtype,256,&flg);
147: if (flg) {
148: MatSetType(B,mtype);
149: }
150: if (!B->type_name) {
151: MPI_Comm_size(B->comm,&size);
152: if (size == 1) {
153: MatSetType(B,MATSEQAIJ);
154: } else {
155: MatSetType(B,MATMPIAIJ);
156: }
157: }
158: return(0);
159: }
161: /*@C
162: MatSetUpPreallocation
164: Collective on Mat
166: Input Parameter:
167: . A - the matrix
169: Level: beginner
171: .keywords: matrix, create
173: .seealso: MatCreateSeqAIJ((), MatCreateMPIAIJ(),
174: MatCreateSeqBDiag(),MatCreateMPIBDiag(),
175: MatCreateSeqDense(), MatCreateMPIDense(),
176: MatCreateMPIRowbs(), MatCreateSeqBAIJ(), MatCreateMPIBAIJ(),
177: MatCreateSeqSBAIJ(), MatCreateMPISBAIJ(),
178: MatConvert()
179: @*/
180: int MatSetUpPreallocation(Mat B)
181: {
182: int ierr;
185: if (B->ops->setuppreallocation) {
186: PetscLogInfo(B,"MatSetTpPreallocation: Warning not preallocating matrix storage");
187: (*B->ops->setuppreallocation)(B);
188: B->ops->setuppreallocation = 0;
189: B->preallocated = PETSC_TRUE;
190: }
191: return(0);
192: }
194: /*
195: Copies from Cs header to A
196: */
197: int MatHeaderCopy(Mat A,Mat C)
198: {
199: int ierr,refct;
200: PetscOps *Abops;
201: MatOps Aops;
202: char *mtype,*mname;
205: /* free all the interior data structures from mat */
206: (*A->ops->destroy)(A);
208: MapDestroy(A->rmap);
209: MapDestroy(A->cmap);
211: /* save the parts of A we need */
212: Abops = A->bops;
213: Aops = A->ops;
214: refct = A->refct;
215: mtype = A->type_name;
216: mname = A->name;
218: /* copy C over to A */
219: ierr = PetscMemcpy(A,C,sizeof(struct _p_Mat));
221: /* return the parts of A we saved */
222: A->bops = Abops;
223: A->ops = Aops;
224: A->qlist = 0;
225: A->refct = refct;
226: A->type_name = mtype;
227: A->name = mname;
229: PetscHeaderDestroy(C);
230: return(0);
231: }