Actual source code: zmat.c
2: #include src/mat/impls/adj/mpi/mpiadj.h
3: #include src/fortran/custom/zpetsc.h
4: #include petscmat.h
7: #ifdef PETSC_HAVE_FORTRAN_CAPS
8: #define matpartitioningsetvertexweights_ MATPARTITIONINGSETVERTEXWEIGHTS
9: #define matsettype_ MATSETTYPE
10: #define matmpiaijgetseqaij_ MATMPIAIJGETSEQAIJ
11: #define matmpibaijgetseqbaij_ MATMPIBAIJGETSEQBAIJ
12: #define matgetrowij_ MATGETROWIJ
13: #define matrestorerowij_ MATRESTOREROWIJ
14: #define matsetfromoptions_ MATSETFROMOPTIONS
15: #define matcreateseqaijwitharrays_ MATCREATESEQAIJWITHARRAYS
16: #define matpartitioningdestroy_ MATPARTITIONINGDESTROY
17: #define matsetvalue_ MATSETVALUE
18: #define matsetvaluelocal_ MATSETVALUELOCAL
19: #define matgetrow_ MATGETROW
20: #define matrestorerow_ MATRESTOREROW
21: #define matgetordering_ MATGETORDERING
22: #define matdestroy_ MATDESTROY
23: #define matcreatempiaij_ MATCREATEMPIAIJ
24: #define matcreateseqaij_ MATCREATESEQAIJ
25: #define matcreatempibaij_ MATCREATEMPIBAIJ
26: #define matcreateseqbaij_ MATCREATESEQBAIJ
27: #define matcreatempisbaij_ MATCREATEMPISBAIJ
28: #define matcreateseqsbaij_ MATCREATESEQSBAIJ
29: #define matcreate_ MATCREATE
30: #define matcreateshell_ MATCREATESHELL
31: #define matshellsetcontext_ MATSHELLSETCONTEXT
32: #define matorderingregisterdestroy_ MATORDERINGREGISTERDESTROY
33: #define matcreatempirowbs_ MATCREATEMPIROWBS
34: #define matcreateseqbdiag_ MATCREATESEQBDIAG
35: #define matcreatempibdiag_ MATCREATEMPIBDIAG
36: #define matcreateseqdense_ MATCREATESEQDENSE
37: #define matcreatempidense_ MATCREATEMPIDENSE
38: #define matconvert_ MATCONVERT
39: #define matload_ MATLOAD
40: #define mattranspose_ MATTRANSPOSE
41: #define matgetarray_ MATGETARRAY
42: #define matrestorearray_ MATRESTOREARRAY
43: #define matgettype_ MATGETTYPE
44: #define matgetinfo_ MATGETINFO
45: #define matshellsetoperation_ MATSHELLSETOPERATION
46: #define matview_ MATVIEW
47: #define matfdcoloringcreate_ MATFDCOLORINGCREATE
48: #define matfdcoloringdestroy_ MATFDCOLORINGDESTROY
49: #define matfdcoloringsetfunctionsnes_ MATFDCOLORINGSETFUNCTIONSNES
50: #define matfdcoloringsetfunctionts_ MATFDCOLORINGSETFUNCTIONTS
51: #define matcopy_ MATCOPY
52: #define matgetsubmatrices_ MATGETSUBMATRICES
53: #define matgetcoloring_ MATGETCOLORING
54: #define matpartitioningsettype_ MATPARTITIONINGSETTYPE
55: #define matduplicate_ MATDUPLICATE
56: #define matzerorows_ MATZEROROWS
57: #define matzerorowsis_ MATZEROROWSIS
58: #define matzerorowslocal_ MATZEROROWSLOCAL
59: #define matzerorowslocalis_ MATZEROROWSLOCALIS
60: #define matpartitioningview_ MATPARTITIONINGVIEW
61: #define matpartitioningcreate_ MATPARTITIONINGCREATE
62: #define matpartitioningsetadjacency_ MATPARTITIONINGSETADJACENCY
63: #define matpartitioningapply_ MATPARTITIONINGAPPLY
64: #define matcreatempiadj_ MATCREATEMPIADJ
65: #define matsetvaluesstencil_ MATSETVALUESSTENCIL
66: #define matseqaijsetpreallocation_ MATSEQAIJSETPREALLOCATION
67: #define matmpiaijsetpreallocation_ MATMPIAIJSETPREALLOCATION
68: #define matseqbaijsetpreallocation_ MATSEQBAIJSETPREALLOCATION
69: #define matmpibaijsetpreallocation_ MATMPIBAIJSETPREALLOCATION
70: #define matseqsbaijsetpreallocation_ MATSEQSBAIJSETPREALLOCATION
71: #define matmpisbaijsetpreallocation_ MATMPISBAIJSETPREALLOCATION
72: #define matseqbdiagsetpreallocation_ MATSEQBDIAGSETPREALLOCATION
73: #define matmpibdiagsetpreallocation_ MATMPIBDIAGSETPREALLOCATION
74: #define matseqdensesetpreallocation_ MATSEQDENSESETPREALLOCATION
75: #define matmpidensesetpreallocation_ MATMPIDENSESETPREALLOCATION
76: #define matmpiadjsetpreallocation_ MATMPIADJSETPREALLOCATION
77: #define matmpirowbssetpreallocation_ MATMPIROWBSSETPREALLOCATION
78: #define matpartitioningpartysetglobal_ MATPARTITIONINGPARTYSETGLOBAL
79: #define matpartitioningpartysetlocal_ MATPARTITIONINGPARTYSETLOCAL
80: #define matpartitioningscotchsetstrategy_ MATPARTITIONINGSCOTCHSETSTRATEGY
81: #define matpartitioningscotchsetarch_ MATPARTITIONINGSCOTCHSETARCH
82: #define matpartitioningscotchsethostlist_ MATPARTITIONINGSCOTCHSETHOSTLIST
83: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
84: #define matpartitioningsetvertexweights_ matpartitioningsetvertexweights
85: #define matsettype_ matsettype
86: #define matmpiaijgetseqaij_ matmpiaijgetseqaij
87: #define matmpibaijgetseqbaij_ matmpibaijgetseqbaij
88: #define matrestorerowij_ matrestorerowij
89: #define matgetrowij_ matgetrowij
90: #define matcreateseqaijwitharrays_ matcreateseqaijwitharrays
91: #define matpartitioningdestroy_ matpartitioningdestroy
92: #define matpartitioningsettype_ matpartitioningsettype
93: #define matsetvalue_ matsetvalue
94: #define matsetvaluelocal_ matsetvaluelocal
95: #define matgetrow_ matgetrow
96: #define matrestorerow_ matrestorerow
97: #define matview_ matview
98: #define matgetinfo_ matgetinfo
99: #define matgettype_ matgettype
100: #define matdestroy_ matdestroy
101: #define matcreatempiaij_ matcreatempiaij
102: #define matcreateseqaij_ matcreateseqaij
103: #define matcreatempibaij_ matcreatempibaij
104: #define matcreateseqbaij_ matcreateseqbaij
105: #define matcreatempisbaij_ matcreatempisbaij
106: #define matcreateseqsbaij_ matcreateseqsbaij
107: #define matcreate_ matcreate
108: #define matcreateshell_ matcreateshell
109: #define matshellsetcontext_ matshellsetcontext
110: #define matorderingregisterdestroy_ matorderingregisterdestroy
111: #define matgetordering_ matgetordering
112: #define matcreatempirowbs_ matcreatempirowbs
113: #define matcreateseqbdiag_ matcreateseqbdiag
114: #define matcreatempibdiag_ matcreatempibdiag
115: #define matcreateseqdense_ matcreateseqdense
116: #define matcreatempidense_ matcreatempidense
117: #define matconvert_ matconvert
118: #define matload_ matload
119: #define mattranspose_ mattranspose
120: #define matgetarray_ matgetarray
121: #define matrestorearray_ matrestorearray
122: #define matshellsetoperation_ matshellsetoperation
123: #define matfdcoloringcreate_ matfdcoloringcreate
124: #define matfdcoloringdestroy_ matfdcoloringdestroy
125: #define matfdcoloringsetfunctionsnes_ matfdcoloringsetfunctionsnes
126: #define matfdcoloringsetfunctionts_ matfdcoloringsetfunctionts
127: #define matcopy_ matcopy
128: #define matgetsubmatrices_ matgetsubmatrices
129: #define matgetcoloring_ matgetcoloring
130: #define matduplicate_ matduplicate
131: #define matzerorows_ matzerorows
132: #define matzerorowsis_ matzerorowsis
133: #define matzerorowslocal_ matzerorowslocal
134: #define matzerorowslocalis_ matzerorowslocalis
135: #define matpartitioningview_ matpartitioningview
136: #define matpartitioningcreate_ matpartitioningcreate
137: #define matpartitioningsetadjacency_ matpartitioningsetadjacency
138: #define matpartitioningapply_ matpartitioningapply
139: #define matcreatempiadj_ matcreatempiadj
140: #define matsetfromoptions_ matsetfromoptions
141: #define matsetvaluesstencil_ matsetvaluesstencil
142: #define matseqaijsetpreallocation_ matseqaijsetpreallocation
143: #define matmpiaijsetpreallocation_ matmpiaijsetpreallocation
144: #define matseqbaijsetpreallocation_ matseqbaijsetpreallocation
145: #define matmpibaijsetpreallocation_ matmpibaijsetpreallocation
146: #define matseqsbaijsetpreallocation_ matseqsbaijsetpreallocation
147: #define matmpisbaijsetpreallocation_ matmpisbaijsetpreallocation
148: #define matseqbdiagsetpreallocation_ matseqbdiagsetpreallocation
149: #define matmpibdiagsetpreallocation_ matmpibdiagsetpreallocation
150: #define matseqdensesetpreallocation_ matseqdensesetpreallocation
151: #define matmpidensesetpreallocation_ matmpidensesetpreallocation
152: #define matmpiadjsetpreallocation_ matmpiadjsetpreallocation
153: #define matmpirowbssetpreallocation_ matmpirowbssetpreallocation
154: #define matpartitioningpartysetglobal_ matpartitioningpartysetglobal
155: #define matpartitioningpartysetlocal_ matpartitioningpartysetlocal
156: #define matpartitioningscotchsetstrategy_ matpartitioningscotchsetstrategy
157: #define matpartitioningscotchsetarch_ matpartitioningscotchsetarch
158: #define matpartitioningscotchsethostlist_ matpartitioningscotchsethostlist
159: #endif
161: #include petscts.h
164: static void (PETSC_STDCALL *f7)(TS*,double*,Vec*,Vec*,void*,PetscErrorCode*);
165: static void (PETSC_STDCALL *f8)(SNES*,Vec*,Vec*,void*,PetscErrorCode*);
169: static PetscErrorCode ourmatfdcoloringfunctionts(TS ts,double t,Vec x,Vec y,void *ctx)
170: {
171: PetscErrorCode 0;
172: (*f7)(&ts,&t,&x,&y,ctx,&ierr);
173: return ierr;
174: }
176: static PetscErrorCode ourmatfdcoloringfunctionsnes(SNES ts,Vec x,Vec y,void *ctx)
177: {
178: PetscErrorCode 0;
179: (*f8)(&ts,&x,&y,ctx,&ierr);
180: return ierr;
181: }
185: void PETSC_STDCALL matpartitioningsetvertexweights_(MatPartitioning *part,const PetscInt weights[],PetscErrorCode *ierr)
186: {
187: PetscInt len;
188: PetscInt *array;
189: *MatGetLocalSize((*part)->adj,&len,0); if (*ierr) return;
190: *PetscMalloc(len*sizeof(PetscInt),&array); if (*ierr) return;
191: *PetscMemcpy(array,weights,len*sizeof(PetscInt));if (*ierr) return;
192: *MatPartitioningSetVertexWeights(*part,array);
193: }
196: void PETSC_STDCALL matsettype_(Mat *x,CHAR type_name PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len))
197: {
198: char *t;
200: FIXCHAR(type_name,len,t);
201: *MatSetType(*x,t);
202: FREECHAR(type_name,t);
203: }
205: void PETSC_STDCALL matsetvaluesstencil_(Mat *mat,PetscInt *m,MatStencil *idxm,PetscInt *n,MatStencil *idxn,PetscScalar *v,InsertMode *addv,
206: PetscErrorCode *ierr)
207: {
208: *MatSetValuesStencil(*mat,*m,idxm,*n,idxn,v,*addv);
209: }
211: void PETSC_STDCALL matmpiaijgetseqaij_(Mat *A,Mat *Ad,Mat *Ao,PetscInt *ic,size_t *iic,PetscErrorCode *ierr)
212: {
213: PetscInt *i;
214: *MatMPIAIJGetSeqAIJ(*A,Ad,Ao,&i);if (*ierr) return;
215: *iic = PetscIntAddressToFortran(ic,i);
216: }
218: void PETSC_STDCALL matmpibaijgetseqbaij_(Mat *A,Mat *Ad,Mat *Ao,PetscInt *ic,size_t *iic,PetscErrorCode *ierr)
219: {
220: PetscInt *i;
221: *MatMPIBAIJGetSeqBAIJ(*A,Ad,Ao,&i);if (*ierr) return;
222: *iic = PetscIntAddressToFortran(ic,i);
223: }
225: void PETSC_STDCALL matgetrowij_(Mat *B,PetscInt *shift,PetscTruth *sym,PetscInt *n,PetscInt *ia,size_t *iia,PetscInt *ja,size_t *jja,
226: PetscTruth *done,PetscErrorCode *ierr)
227: {
228: PetscInt *IA,*JA;
229: *MatGetRowIJ(*B,*shift,*sym,n,&IA,&JA,done);if (*ierr) return;
230: *iia = PetscIntAddressToFortran(ia,IA);
231: *jja = PetscIntAddressToFortran(ja,JA);
232: }
234: void PETSC_STDCALL matrestorerowij_(Mat *B,PetscInt *shift,PetscTruth *sym,PetscInt *n,PetscInt *ia,size_t *iia,PetscInt *ja,size_t *jja,
235: PetscTruth *done,PetscErrorCode *ierr)
236: {
237: PetscInt *IA = PetscIntAddressFromFortran(ia,*iia),*JA = PetscIntAddressFromFortran(ja,*jja);
238: *MatRestoreRowIJ(*B,*shift,*sym,n,&IA,&JA,done);
239: }
241: void PETSC_STDCALL matsetfromoptions_(Mat *B,PetscErrorCode *ierr)
242: {
243: *MatSetFromOptions(*B);
244: }
246: void PETSC_STDCALL matcreateseqaijwitharrays_(MPI_Comm *comm,PetscInt *m,PetscInt *n,PetscInt *i,PetscInt *j,PetscScalar *a,Mat *mat,PetscErrorCode *ierr)
247: {
248: *MatCreateSeqAIJWithArrays((MPI_Comm)PetscToPointerComm(*comm),*m,*n,i,j,a,mat);
249: }
251: void PETSC_STDCALL matcreatempiadj_(MPI_Comm *comm,PetscInt *m,PetscInt *n,PetscInt *i,PetscInt *j,PetscInt *values,Mat *A,PetscErrorCode *ierr)
252: {
253: Mat_MPIAdj *adj;
255: CHKFORTRANNULLINTEGER(values);
256: *MatCreateMPIAdj((MPI_Comm)PetscToPointerComm(*comm),*m,*n,i,j,values,A);
257: adj = (Mat_MPIAdj*)(*A)->data;
258: adj->freeaij = PETSC_FALSE;
259: }
261: void PETSC_STDCALL matpartitioningdestroy_(MatPartitioning *part,PetscErrorCode *ierr)
262: {
263: *MatPartitioningDestroy(*part);
264: }
266: void PETSC_STDCALL matpartitioningcreate_(MPI_Comm *comm,MatPartitioning *part, PetscErrorCode *ierr)
267: {
268: *MatPartitioningCreate((MPI_Comm)PetscToPointerComm(*comm),part);
269: }
271: void PETSC_STDCALL matpartitioningapply_(MatPartitioning *part,IS *is,PetscErrorCode *ierr)
272: {
273: *MatPartitioningApply(*part,is);
274: }
276: void PETSC_STDCALL matpartitioningsetadjacency_(MatPartitioning *part,Mat *mat,PetscErrorCode *ierr)
277: {
278: *MatPartitioningSetAdjacency(*part,*mat);
279: }
281: void PETSC_STDCALL matpartitioningview_(MatPartitioning *part,PetscViewer *viewer, PetscErrorCode *ierr)
282: {
283: PetscViewer v;
284: PetscPatchDefaultViewers_Fortran(viewer,v);
285: *MatPartitioningView(*part,v);
286: }
288: void PETSC_STDCALL matpartitioningsettype_(MatPartitioning *part,CHAR type PETSC_MIXED_LEN(len),
289: PetscErrorCode *ierr PETSC_END_LEN(len))
290: {
291: char *t;
292: FIXCHAR(type,len,t);
293: *MatPartitioningSetType(*part,t);
294: FREECHAR(type,t);
295: }
297: void PETSC_STDCALL matgetcoloring_(Mat *mat,CHAR type PETSC_MIXED_LEN(len),ISColoring *iscoloring,
298: PetscErrorCode *ierr PETSC_END_LEN(len))
299: {
300: char *t;
301: FIXCHAR(type,len,t);
302: *MatGetColoring(*mat,t,iscoloring);
303: FREECHAR(type,t);
304: }
306: void PETSC_STDCALL matsetvalue_(Mat *mat,PetscInt *i,PetscInt *j,PetscScalar *va,InsertMode *mode,PetscErrorCode *ierr)
307: {
308: /* cannot use MatSetValue() here since that usesCHKERRQ() which has a return in it */
309: *MatSetValues(*mat,1,i,1,j,va,*mode);
310: }
312: void PETSC_STDCALL matsetvaluelocal_(Mat *mat,PetscInt *i,PetscInt *j,PetscScalar *va,InsertMode *mode,PetscErrorCode *ierr)
313: {
314: /* cannot use MatSetValueLocal() here since that usesCHKERRQ() which has a return in it */
315: *MatSetValuesLocal(*mat,1,i,1,j,va,*mode);
316: }
318: void PETSC_STDCALL matfdcoloringcreate_(Mat *mat,ISColoring *iscoloring,MatFDColoring *color,PetscErrorCode *ierr)
319: {
320: *MatFDColoringCreate(*mat,*iscoloring,color);
321: }
323: /*
324: This is a poor way of storing the column and value pointers
325: generated by MatGetRow() to be returned with MatRestoreRow()
326: but there is not natural,good place else to store them. Hence
327: Fortran programmers can only have one outstanding MatGetRows()
328: at a time.
329: */
330: static PetscErrorCode matgetrowactive = 0;
331: static const PetscInt *my_ocols = 0;
332: static const PetscScalar *my_ovals = 0;
334: void PETSC_STDCALL matgetrow_(Mat *mat,PetscInt *row,PetscInt *ncols,PetscInt *cols,PetscScalar *vals,PetscErrorCode *ierr)
335: {
336: const PetscInt **oocols = &my_ocols;
337: const PetscScalar **oovals = &my_ovals;
339: if (matgetrowactive) {
340: PetscError(__LINE__,"MatGetRow_Fortran",__FILE__,__SDIR__,1,0,
341: "Cannot have two MatGetRow() active simultaneously\n\
342: call MatRestoreRow() before calling MatGetRow() a second time");
343: *1;
344: return;
345: }
347: CHKFORTRANNULLINTEGER(cols); if (!cols) oocols = PETSC_NULL;
348: CHKFORTRANNULLSCALAR(vals); if (!vals) oovals = PETSC_NULL;
350: *MatGetRow(*mat,*row,ncols,oocols,oovals);
351: if (*ierr) return;
353: if (oocols) { *PetscMemcpy(cols,my_ocols,(*ncols)*sizeof(PetscInt)); if (*ierr) return;}
354: if (oovals) { *PetscMemcpy(vals,my_ovals,(*ncols)*sizeof(PetscScalar)); if (*ierr) return; }
355: matgetrowactive = 1;
356: }
358: void PETSC_STDCALL matrestorerow_(Mat *mat,PetscInt *row,PetscInt *ncols,PetscInt *cols,PetscScalar *vals,PetscErrorCode *ierr)
359: {
360: const PetscInt **oocols = &my_ocols;
361: const PetscScalar **oovals = &my_ovals;
362: if (!matgetrowactive) {
363: PetscError(__LINE__,"MatRestoreRow_Fortran",__FILE__,__SDIR__,1,0,
364: "Must call MatGetRow() first");
365: *1;
366: return;
367: }
368: CHKFORTRANNULLINTEGER(cols); if (!cols) oocols = PETSC_NULL;
369: CHKFORTRANNULLSCALAR(vals); if (!vals) oovals = PETSC_NULL;
371: *MatRestoreRow(*mat,*row,ncols,oocols,oovals);
372: matgetrowactive = 0;
373: }
375: void PETSC_STDCALL matview_(Mat *mat,PetscViewer *vin,PetscErrorCode *ierr)
376: {
377: PetscViewer v;
378: PetscPatchDefaultViewers_Fortran(vin,v);
379: *MatView(*mat,v);
380: }
382: void PETSC_STDCALL matcopy_(Mat *A,Mat *B,MatStructure *str,PetscErrorCode *ierr)
383: {
384: *MatCopy(*A,*B,*str);
385: }
387: void PETSC_STDCALL matgetinfo_(Mat *mat,MatInfoType *flag,double *finfo,PetscErrorCode *ierr)
388: {
389: *MatGetInfo(*mat,*flag,(MatInfo*)finfo);
390: }
392: void PETSC_STDCALL matgetarray_(Mat *mat,PetscScalar *fa,size_t *ia,PetscErrorCode *ierr)
393: {
394: PetscScalar *mm;
395: PetscInt m,n;
397: *MatGetArray(*mat,&mm); if (*ierr) return;
398: *MatGetSize(*mat,&m,&n); if (*ierr) return;
399: *PetscScalarAddressToFortran((PetscObject)*mat,fa,mm,m*n,ia); if (*ierr) return;
400: }
402: void PETSC_STDCALL matrestorearray_(Mat *mat,PetscScalar *fa,size_t *ia,PetscErrorCode *ierr)
403: {
404: PetscScalar *lx;
405: PetscInt m,n;
407: *MatGetSize(*mat,&m,&n); if (*ierr) return;
408: *PetscScalarAddressFromFortran((PetscObject)*mat,fa,*ia,m*n,&lx);if (*ierr) return;
409: *MatRestoreArray(*mat,&lx);if (*ierr) return;
410: }
412: void PETSC_STDCALL mattranspose_(Mat *mat,Mat *B,PetscErrorCode *ierr)
413: {
414: CHKFORTRANNULLOBJECT(B);
415: *MatTranspose(*mat,B);
416: }
418: void PETSC_STDCALL matload_(PetscViewer *viewer,CHAR outtype PETSC_MIXED_LEN(len),Mat *newmat,PetscErrorCode *ierr PETSC_END_LEN(len))
419: {
420: char *t;
421: PetscViewer v;
422: FIXCHAR(outtype,len,t);
423: PetscPatchDefaultViewers_Fortran(viewer,v);
424: *MatLoad(v,t,newmat);
425: FREECHAR(outtype,t);
426: }
428: void PETSC_STDCALL matconvert_(Mat *mat,CHAR outtype PETSC_MIXED_LEN(len),MatReuse *reuse,Mat *M,PetscErrorCode *ierr PETSC_END_LEN(len))
429: {
430: char *t;
431: FIXCHAR(outtype,len,t);
432: *MatConvert(*mat,t,*reuse,M);
433: FREECHAR(outtype,t);
434: }
436: void PETSC_STDCALL matcreateseqdense_(MPI_Comm *comm,PetscInt *m,PetscInt *n,PetscScalar *data,Mat *newmat,PetscErrorCode *ierr)
437: {
438: CHKFORTRANNULLSCALAR(data);
439: *MatCreateSeqDense((MPI_Comm)PetscToPointerComm(*comm),*m,*n,data,newmat);
440: }
442: void PETSC_STDCALL matcreatempidense_(MPI_Comm *comm,PetscInt *m,PetscInt *n,PetscInt *M,PetscInt *N,PetscScalar *data,Mat *newmat,
443: PetscErrorCode *ierr)
444: {
445: CHKFORTRANNULLSCALAR(data);
446: *MatCreateMPIDense((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*M,*N,data,newmat);
447: }
449: /* Fortran ignores diagv */
450: void PETSC_STDCALL matcreatempibdiag_(MPI_Comm *comm,PetscInt *m,PetscInt *M,PetscInt *N,PetscInt *nd,PetscInt *bs,
451: PetscInt *diag,PetscScalar **diagv,Mat *newmat,PetscErrorCode *ierr)
452: {
453: CHKFORTRANNULLINTEGER(diag);
454: *MatCreateMPIBDiag((MPI_Comm)PetscToPointerComm(*comm),
455: *m,*M,*N,*nd,*bs,diag,PETSC_NULL,newmat);
456: }
458: /* Fortran ignores diagv */
459: void PETSC_STDCALL matcreateseqbdiag_(MPI_Comm *comm,PetscInt *m,PetscInt *n,PetscInt *nd,PetscInt *bs,
460: PetscInt *diag,PetscScalar **diagv,Mat *newmat,PetscErrorCode *ierr)
461: {
462: CHKFORTRANNULLINTEGER(diag);
463: *MatCreateSeqBDiag((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*nd,*bs,diag,
464: PETSC_NULL,newmat);
465: }
467: #if defined(PETSC_HAVE_BLOCKSOLVE)
468: /* Fortran cannot pass in procinfo,hence ignored */
469: void PETSC_STDCALL matcreatempirowbs_(MPI_Comm *comm,PetscInt *m,PetscInt *M,PetscInt *nz,PetscInt *nnz,Mat *newmat,PetscErrorCode *ierr)
470: {
471: CHKFORTRANNULLINTEGER(nnz);
472: *MatCreateMPIRowbs((MPI_Comm)PetscToPointerComm(*comm),*m,*M,*nz,nnz,newmat);
473: }
474: #endif
476: void PETSC_STDCALL matgetordering_(Mat *mat,CHAR type PETSC_MIXED_LEN(len),IS *rperm,IS *cperm,
477: PetscErrorCode *ierr PETSC_END_LEN(len))
478: {
479: char *t;
480: FIXCHAR(type,len,t);
481: *MatGetOrdering(*mat,t,rperm,cperm);
482: FREECHAR(type,t);
483: }
485: void PETSC_STDCALL matorderingregisterdestroy_(PetscErrorCode *ierr)
486: {
487: *MatOrderingRegisterDestroy();
488: }
490: void PETSC_STDCALL matgettype_(Mat *mm,CHAR name PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len))
491: {
492: const char *tname;
494: *MatGetType(*mm,&tname);
495: #if defined(PETSC_USES_CPTOFCD)
496: {
497: char *t = _fcdtocp(name); int len1 = _fcdlen(name);
498: if (t != PETSC_NULL_CHARACTER_Fortran) {
499: *PetscStrncpy(t,tname,len1);if (*ierr) return;
500: }
501: }
502: #else
503: if (name != PETSC_NULL_CHARACTER_Fortran) {
504: *PetscStrncpy(name,tname,len);if (*ierr) return;
505: }
506: #endif
507: FIXRETURNCHAR(name,len);
509: }
511: void PETSC_STDCALL matcreate_(MPI_Comm *comm,Mat *M,PetscErrorCode *ierr)
512: {
513: *MatCreate((MPI_Comm)PetscToPointerComm(*comm),M);
514: }
516: void PETSC_STDCALL matcreateseqaij_(MPI_Comm *comm,PetscInt *m,PetscInt *n,PetscInt *nz,
517: PetscInt *nnz,Mat *newmat,PetscErrorCode *ierr)
518: {
519: CHKFORTRANNULLINTEGER(nnz);
520: *MatCreateSeqAIJ((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*nz,nnz,newmat);
521: }
523: void PETSC_STDCALL matcreateseqbaij_(MPI_Comm *comm,PetscInt *bs,PetscInt *m,PetscInt *n,PetscInt *nz,
524: PetscInt *nnz,Mat *newmat,PetscErrorCode *ierr)
525: {
526: CHKFORTRANNULLINTEGER(nnz);
527: *MatCreateSeqBAIJ((MPI_Comm)PetscToPointerComm(*comm),*bs,*m,*n,*nz,nnz,newmat);
528: }
530: void PETSC_STDCALL matcreateseqsbaij_(MPI_Comm *comm,PetscInt *bs,PetscInt *m,PetscInt *n,PetscInt *nz,
531: PetscInt *nnz,Mat *newmat,PetscErrorCode *ierr)
532: {
533: CHKFORTRANNULLINTEGER(nnz);
534: *MatCreateSeqSBAIJ((MPI_Comm)PetscToPointerComm(*comm),*bs,*m,*n,*nz,nnz,newmat);
535: }
537: void PETSC_STDCALL matfdcoloringdestroy_(MatFDColoring *mat,PetscErrorCode *ierr)
538: {
539: *MatFDColoringDestroy(*mat);
540: }
542: void PETSC_STDCALL matdestroy_(Mat *mat,PetscErrorCode *ierr)
543: {
544: *MatDestroy(*mat);
545: }
547: void PETSC_STDCALL matcreatempiaij_(MPI_Comm *comm,PetscInt *m,PetscInt *n,PetscInt *M,PetscInt *N,
548: PetscInt *d_nz,PetscInt *d_nnz,PetscInt *o_nz,PetscInt *o_nnz,Mat *newmat,PetscErrorCode *ierr)
549: {
550: CHKFORTRANNULLINTEGER(d_nnz);
551: CHKFORTRANNULLINTEGER(o_nnz);
553: *MatCreateMPIAIJ((MPI_Comm)PetscToPointerComm(*comm),
554: *m,*n,*M,*N,*d_nz,d_nnz,*o_nz,o_nnz,newmat);
555: }
556: void PETSC_STDCALL matcreatempibaij_(MPI_Comm *comm,PetscInt *bs,PetscInt *m,PetscInt *n,PetscInt *M,PetscInt *N,
557: PetscInt *d_nz,PetscInt *d_nnz,PetscInt *o_nz,PetscInt *o_nnz,Mat *newmat,PetscErrorCode *ierr)
558: {
559: CHKFORTRANNULLINTEGER(d_nnz);
560: CHKFORTRANNULLINTEGER(o_nnz);
561: *MatCreateMPIBAIJ((MPI_Comm)PetscToPointerComm(*comm),
562: *bs,*m,*n,*M,*N,*d_nz,d_nnz,*o_nz,o_nnz,newmat);
563: }
564: void PETSC_STDCALL matcreatempisbaij_(MPI_Comm *comm,PetscInt *bs,PetscInt *m,PetscInt *n,PetscInt *M,PetscInt *N,
565: PetscInt *d_nz,PetscInt *d_nnz,PetscInt *o_nz,PetscInt *o_nnz,Mat *newmat,PetscErrorCode *ierr)
566: {
567: CHKFORTRANNULLINTEGER(d_nnz);
568: CHKFORTRANNULLINTEGER(o_nnz);
569: *MatCreateMPISBAIJ((MPI_Comm)PetscToPointerComm(*comm),
570: *bs,*m,*n,*M,*N,*d_nz,d_nnz,*o_nz,o_nnz,newmat);
571: }
573: /*
574: The MatShell Matrix Vector product requires a C routine.
575: This C routine then calls the corresponding Fortran routine that was
576: set by the user.
577: */
578: void PETSC_STDCALL matcreateshell_(MPI_Comm *comm,PetscInt *m,PetscInt *n,PetscInt *M,PetscInt *N,void **ctx,Mat *mat,PetscErrorCode *ierr)
579: {
580: *MatCreateShell((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*M,*N,*ctx,mat);
581: if (*ierr) return;
582: *PetscMalloc(4*sizeof(void*),&((PetscObject)*mat)->fortran_func_pointers);
583: }
585: void PETSC_STDCALL matshellsetcontext_(Mat *mat,void**ctx, int *__ierr ){
586: *__MatShellSetContext(*mat,*ctx);
587: }
589: static PetscErrorCode ourmult(Mat mat,Vec x,Vec y)
590: {
591: PetscErrorCode 0;
592: (*(PetscErrorCode (PETSC_STDCALL *)(Mat*,Vec*,Vec*,PetscErrorCode*))(((PetscObject)mat)->fortran_func_pointers[0]))(&mat,&x,&y,&ierr);
593: return ierr;
594: }
596: static PetscErrorCode ourmulttranspose(Mat mat,Vec x,Vec y)
597: {
598: PetscErrorCode 0;
599: (*(PetscErrorCode (PETSC_STDCALL *)(Mat*,Vec*,Vec*,PetscErrorCode*))(((PetscObject)mat)->fortran_func_pointers[2]))(&mat,&x,&y,&ierr);
600: return ierr;
601: }
603: static PetscErrorCode ourmultadd(Mat mat,Vec x,Vec y,Vec z)
604: {
605: PetscErrorCode 0;
606: (*(PetscErrorCode (PETSC_STDCALL *)(Mat*,Vec*,Vec*,Vec*,PetscErrorCode*))(((PetscObject)mat)->fortran_func_pointers[1]))(&mat,&x,&y,&z,&ierr);
607: return ierr;
608: }
610: static PetscErrorCode ourmulttransposeadd(Mat mat,Vec x,Vec y,Vec z)
611: {
612: PetscErrorCode 0;
613: (*(PetscErrorCode (PETSC_STDCALL *)(Mat*,Vec*,Vec*,Vec*,PetscErrorCode*))(((PetscObject)mat)->fortran_func_pointers[3]))(&mat,&x,&y,&z,&ierr);
614: return ierr;
615: }
617: void PETSC_STDCALL matshellsetoperation_(Mat *mat,MatOperation *op,PetscErrorCode (PETSC_STDCALL *f)(Mat*,Vec*,Vec*,PetscErrorCode*),PetscErrorCode *ierr)
618: {
619: if (*op == MATOP_MULT) {
620: *MatShellSetOperation(*mat,*op,(FCNVOID)ourmult);
621: ((PetscObject)*mat)->fortran_func_pointers[0] = (FCNVOID)f;
622: } else if (*op == MATOP_MULT_TRANSPOSE) {
623: *MatShellSetOperation(*mat,*op,(FCNVOID)ourmulttranspose);
624: ((PetscObject)*mat)->fortran_func_pointers[2] = (FCNVOID)f;
625: } else if (*op == MATOP_MULT_ADD) {
626: *MatShellSetOperation(*mat,*op,(FCNVOID)ourmultadd);
627: ((PetscObject)*mat)->fortran_func_pointers[1] = (FCNVOID)f;
628: } else if (*op == MATOP_MULT_TRANSPOSE_ADD) {
629: *MatShellSetOperation(*mat,*op,(FCNVOID)ourmulttransposeadd);
630: ((PetscObject)*mat)->fortran_func_pointers[3] = (FCNVOID)f;
631: } else {
632: PetscError(__LINE__,"MatShellSetOperation_Fortran",__FILE__,__SDIR__,1,0,
633: "Cannot set that matrix operation");
634: *1;
635: }
636: }
638: /*
639: MatFDColoringSetFunction sticks the Fortran function into the fortran_func_pointers
640: this function is then accessed by ourmatfdcoloringfunction()
642: NOTE: FORTRAN USER CANNOT PUT IN A NEW J OR B currently.
644: USER CAN HAVE ONLY ONE MatFDColoring in code Because there is no place to hang f7!
645: */
648: void PETSC_STDCALL matfdcoloringsetfunctionts_(MatFDColoring *fd,void (PETSC_STDCALL *f)(TS*,double*,Vec*,Vec*,void*,PetscErrorCode*),
649: void *ctx,PetscErrorCode *ierr)
650: {
651: f7 = f;
652: *MatFDColoringSetFunction(*fd,(FCNINTVOID)ourmatfdcoloringfunctionts,ctx);
653: }
656: void PETSC_STDCALL matfdcoloringsetfunctionsnes_(MatFDColoring *fd,void (PETSC_STDCALL *f)(SNES*,Vec*,Vec*,void*,PetscErrorCode*),
657: void *ctx,PetscErrorCode *ierr)
658: {
659: f8 = f;
660: *MatFDColoringSetFunction(*fd,(FCNINTVOID)ourmatfdcoloringfunctionsnes,ctx);
661: }
663: /*
664: MatGetSubmatrices() is slightly different from C since the
665: Fortran provides the array to hold the submatrix objects,while in C that
666: array is allocated by the MatGetSubmatrices()
667: */
668: void PETSC_STDCALL matgetsubmatrices_(Mat *mat,PetscInt *n,IS *isrow,IS *iscol,MatReuse *scall,Mat *smat,PetscErrorCode *ierr)
669: {
670: Mat *lsmat;
671: PetscInt i;
673: if (*scall == MAT_INITIAL_MATRIX) {
674: *MatGetSubMatrices(*mat,*n,isrow,iscol,*scall,&lsmat);
675: for (i=0; i<*n; i++) {
676: smat[i] = lsmat[i];
677: }
678: PetscFree(lsmat);
679: } else {
680: *MatGetSubMatrices(*mat,*n,isrow,iscol,*scall,&smat);
681: }
682: }
684: void PETSC_STDCALL matduplicate_(Mat *matin,MatDuplicateOption *op,Mat *matout,PetscErrorCode *ierr)
685: {
686: *MatDuplicate(*matin,*op,matout);
687: }
689: void PETSC_STDCALL matzerorows_(Mat *mat,PetscInt *numRows,PetscInt *rows,PetscScalar *diag,PetscErrorCode *ierr)
690: {
691: CHKFORTRANNULLSCALAR(diag);
692: *MatZeroRows(*mat,*numRows,rows,*diag);
693: }
695: void PETSC_STDCALL matzerorowsis_(Mat *mat,IS *is,PetscScalar *diag,PetscErrorCode *ierr)
696: {
697: CHKFORTRANNULLSCALAR(diag);
698: *MatZeroRowsIS(*mat,*is,*diag);
699: }
701: void PETSC_STDCALL matzerorowslocal_(Mat *mat,PetscInt *numRows,PetscInt *rows,PetscScalar *diag,PetscErrorCode *ierr)
702: {
703: CHKFORTRANNULLSCALAR(diag);
704: *MatZeroRowsLocal(*mat,*numRows,rows,*diag);
705: }
707: void PETSC_STDCALL matzerorowslocalis_(Mat *mat,IS *is,PetscScalar *diag,PetscErrorCode *ierr)
708: {
709: CHKFORTRANNULLSCALAR(diag);
710: *MatZeroRowsLocalIS(*mat,*is,*diag);
711: }
713: void PETSC_STDCALL matseqaijsetpreallocation_(Mat *mat,PetscInt *nz,PetscInt *nnz,PetscErrorCode *ierr)
714: {
715: CHKFORTRANNULLINTEGER(nnz);
716: *MatSeqAIJSetPreallocation(*mat,*nz,nnz);
717: }
719: void PETSC_STDCALL matmpiaijsetpreallocation_(Mat *mat,PetscInt *d_nz,PetscInt *d_nnz,PetscInt *o_nz,PetscInt *o_nnz,PetscErrorCode *ierr)
720: {
721: CHKFORTRANNULLINTEGER(d_nnz);
722: CHKFORTRANNULLINTEGER(o_nnz);
723: *MatMPIAIJSetPreallocation(*mat,*d_nz,d_nnz,*o_nz,o_nnz);
724: }
726: void PETSC_STDCALL matseqbaijsetpreallocation_(Mat *mat,PetscInt *bs,PetscInt *nz,PetscInt *nnz,PetscErrorCode *ierr)
727: {
728: CHKFORTRANNULLINTEGER(nnz);
729: *MatSeqBAIJSetPreallocation(*mat,*bs,*nz,nnz);
730: }
732: void PETSC_STDCALL matmpibaijsetpreallocation_(Mat *mat,PetscInt *bs,PetscInt *d_nz,PetscInt *d_nnz,PetscInt *o_nz,PetscInt *o_nnz,PetscErrorCode *ierr)
733: {
734: CHKFORTRANNULLINTEGER(d_nnz);
735: CHKFORTRANNULLINTEGER(o_nnz);
736: *MatMPIBAIJSetPreallocation(*mat,*bs,*d_nz,d_nnz,*o_nz,o_nnz);
737: }
739: void PETSC_STDCALL matseqsbaijsetpreallocation_(Mat *mat,PetscInt *bs,PetscInt *nz,PetscInt *nnz,PetscErrorCode *ierr)
740: {
741: CHKFORTRANNULLINTEGER(nnz);
742: *MatSeqSBAIJSetPreallocation(*mat,*bs,*nz,nnz);
743: }
745: void PETSC_STDCALL matmpisbaijsetpreallocation_(Mat *mat,PetscInt *bs,PetscInt *d_nz,PetscInt *d_nnz,PetscInt *o_nz,PetscInt *o_nnz,PetscErrorCode *ierr)
746: {
747: CHKFORTRANNULLINTEGER(d_nnz);
748: CHKFORTRANNULLINTEGER(o_nnz);
749: *MatMPISBAIJSetPreallocation(*mat,*bs,*d_nz,d_nnz,*o_nz,o_nnz);
750: }
752: /* Fortran ignores diagv */
753: void PETSC_STDCALL matseqbdiagsetpreallocation_(Mat *mat,PetscInt *nd,PetscInt *bs,PetscInt *diag,PetscScalar **diagv,PetscErrorCode *ierr)
754: {
755: CHKFORTRANNULLINTEGER(diag);
756: *MatSeqBDiagSetPreallocation(*mat,*nd,*bs,diag,PETSC_NULL);
757: }
758: /* Fortran ignores diagv */
759: void PETSC_STDCALL matmpibdiagsetpreallocation_(Mat *mat,PetscInt *nd,PetscInt *bs,PetscInt *diag,PetscScalar **diagv,PetscErrorCode *ierr)
760: {
761: CHKFORTRANNULLINTEGER(diag);
762: *MatMPIBDiagSetPreallocation(*mat,*nd,*bs,diag,PETSC_NULL);
763: }
765: void PETSC_STDCALL matseqdensesetpreallocation_(Mat *mat,PetscScalar *data,PetscErrorCode *ierr)
766: {
767: CHKFORTRANNULLSCALAR(data);
768: *MatSeqDenseSetPreallocation(*mat,data);
769: }
770: void PETSC_STDCALL matmpidensesetpreallocation_(Mat *mat,PetscScalar *data,PetscErrorCode *ierr)
771: {
772: CHKFORTRANNULLSCALAR(data);
773: *MatMPIDenseSetPreallocation(*mat,data);
774: }
775: void PETSC_STDCALL matmpiadjsetpreallocation_(Mat *mat,PetscInt *i,PetscInt *j,PetscInt *values, PetscErrorCode *ierr)
776: {
777: CHKFORTRANNULLINTEGER(values);
778: *MatMPIAdjSetPreallocation(*mat,i,j,values);
779: }
781: #if defined(PETSC_HAVE_BLOCKSOLVE)
782: void PETSC_STDCALL matmpirowbssetpreallocation_(Mat *mat,PetscInt *nz,PetscInt *nnz,PetscErrorCode *ierr)
783: {
784: CHKFORTRANNULLINTEGER(nnz);
785: *MatMPIRowbsSetPreallocation(*mat,*nz,nnz);
786: }
787: #endif
789: #if defined(PETSC_HAVE_PARTY)
790: void PETSC_STDCALL matpartitioningpartysetglobal_(MatPartitioning *part,CHAR method PETSC_MIXED_LEN(len),
791: PetscErrorCode *ierr PETSC_END_LEN(len))
792: {
793: char *t;
794: FIXCHAR(method,len,t);
795: *MatPartitioningPartySetGlobal(*part,t);
796: FREECHAR(method,t);
797: }
799: void PETSC_STDCALL matpartitioningpartysetlocal_(MatPartitioning *part,CHAR method PETSC_MIXED_LEN(len),
800: PetscErrorCode *ierr PETSC_END_LEN(len))
801: {
802: char *t;
803: FIXCHAR(method,len,t);
804: *MatPartitioningPartySetLocal(*part,t);
805: FREECHAR(method,t);
806: }
807: #endif
809: #if defined(PETSC_HAVE_SCOTCH)
810: void PETSC_STDCALL matpartitioningscotchsetstrategy_(MatPartitioning *part,CHAR strategy PETSC_MIXED_LEN(len),
811: PetscErrorCode *ierr PETSC_END_LEN(len))
812: {
813: char *t;
814: FIXCHAR(strategy,len,t);
815: *MatPartitioningScotchSetStrategy(*part,t);
816: FREECHAR(strategy,t);
817: }
819: void PETSC_STDCALL matpartitioningscotchsetarch_(MatPartitioning *part,CHAR filename PETSC_MIXED_LEN(len),
820: PetscErrorCode *ierr PETSC_END_LEN(len))
821: {
822: char *t;
823: FIXCHAR(filename,len,t);
824: *MatPartitioningScotchSetArch(*part,t);
825: FREECHAR(filename,t);
826: }
827: #endif
829: #if defined(PETSC_HAVE_SCOTCH)
830: void PETSC_STDCALL matpartitioningscotchsethostlist_(MatPartitioning *part,CHAR filename PETSC_MIXED_LEN(len),
831: PetscErrorCode *ierr PETSC_END_LEN(len))
832: {
833: char *t;
834: FIXCHAR(filename,len,t);
835: *MatPartitioningScotchSetHostList(*part,t);
836: FREECHAR(filename,t);
837: }
838: #endif