Actual source code: zmat.c
1: /*$Id: zmat.c,v 1.97 2001/04/07 02:33:14 bsmith Exp $*/
3: #include "src/fortran/custom/zpetsc.h"
4: #include "petscmat.h"
6: #ifdef PETSC_HAVE_FORTRAN_CAPS
7: #define matsettype_ MATSETTYPE
8: #define matmpiaijgetseqaij_ MATMPIAIJGETSEQAIJ
9: #define matmpibaijgetseqbaij_ MATMPIBAIJGETSEQBAIJ
10: #define matgetrowij_ MATGETROWIJ
11: #define matrestorerowij_ MATRESTOREROWIJ
12: #define matsetfromoptions_ MATSETFROMOPTIONS
13: #define matcreateseqaijwitharrays_ MATCREATESEQAIJWITHARRAYS
14: #define matpartitioningdestroy_ MATPARTITIONINGDESTROY
15: #define matsetvalue_ MATSETVALUE
16: #define matsetvaluelocal_ MATSETVALUELOCAL
17: #define matgetrow_ MATGETROW
18: #define matrestorerow_ MATRESTOREROW
19: #define matgetordering_ MATGETORDERING
20: #define matdestroy_ MATDESTROY
21: #define matcreatempiaij_ MATCREATEMPIAIJ
22: #define matcreateseqaij_ MATCREATESEQAIJ
23: #define matcreatempibaij_ MATCREATEMPIBAIJ
24: #define matcreateseqbaij_ MATCREATESEQBAIJ
25: #define matcreate_ MATCREATE
26: #define matcreateshell_ MATCREATESHELL
27: #define matorderingregisterdestroy_ MATORDERINGREGISTERDESTROY
28: #define matcreatempirowbs_ MATCREATEMPIROWBS
29: #define matcreateseqbdiag_ MATCREATESEQBDIAG
30: #define matcreatempibdiag_ MATCREATEMPIBDIAG
31: #define matcreateseqdense_ MATCREATESEQDENSE
32: #define matcreatempidense_ MATCREATEMPIDENSE
33: #define matconvert_ MATCONVERT
34: #define matload_ MATLOAD
35: #define mattranspose_ MATTRANSPOSE
36: #define matgetarray_ MATGETARRAY
37: #define matrestorearray_ MATRESTOREARRAY
38: #define matgettype_ MATGETTYPE
39: #define matgetinfo_ MATGETINFO
40: #define matshellsetoperation_ MATSHELLSETOPERATION
41: #define matview_ MATVIEW
42: #define matfdcoloringcreate_ MATFDCOLORINGCREATE
43: #define matfdcoloringdestroy_ MATFDCOLORINGDESTROY
44: #define matfdcoloringsetfunctionsnes_ MATFDCOLORINGSETFUNCTIONSNES
45: #define matfdcoloringsetfunctionts_ MATFDCOLORINGSETFUNCTIONTS
46: #define matcopy_ MATCOPY
47: #define matgetsubmatrices_ MATGETSUBMATRICES
48: #define matgetcoloring_ MATGETCOLORING
49: #define matpartitioningsettype_ MATPARTITIONINGSETTYPE
50: #define matduplicate_ MATDUPLICATE
51: #define matzerorows_ MATZEROROWS
52: #define matzerorowslocal_ MATZEROROWSLOCAL
53: #define matpartitioningview_ MATPARTITIONINGVIEW
54: #define matpartitioningcreate_ MATPARTITIONINGCREATE
55: #define matpartitioningsetadjacency_ MATPARTITIONINGSETADJACENCY
56: #define matpartitioningapply_ MATPARTITIONINGAPPLY
57: #define matcreatempiadj_ MATCREATEMPIADJ
58: #define matsetvaluesstencil_ MATSETVALUESSTENCIL
59: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
60: #define matsettype_ matsettype
61: #define matmpiaijgetseqaij_ matmpiaijgetseqaij
62: #define matmpibaijgetseqbaij_ matmpibaijgetseqbaij
63: #define matrestorerowij_ matrestorerowij
64: #define matgetrowij_ matgetrowij
65: #define matcreateseqaijwitharrays_ matcreateseqaijwitharrays
66: #define matpartitioningdestroy_ matpartitioningdestroy
67: #define matpartitioningsettype_ matpartitioningsettype
68: #define matsetvalue_ matsetvalue
69: #define matsetvaluelocal_ matsetvaluelocal
70: #define matgetrow_ matgetrow
71: #define matrestorerow_ matrestorerow
72: #define matview_ matview
73: #define matgetinfo_ matgetinfo
74: #define matgettype_ matgettype
75: #define matdestroy_ matdestroy
76: #define matcreatempiaij_ matcreatempiaij
77: #define matcreateseqaij_ matcreateseqaij
78: #define matcreatempibaij_ matcreatempibaij
79: #define matcreateseqbaij_ matcreateseqbaij
80: #define matcreate_ matcreate
81: #define matcreateshell_ matcreateshell
82: #define matorderingregisterdestroy_ matorderingregisterdestroy
83: #define matgetordering_ matgetordering
84: #define matcreatempirowbs_ matcreatempirowbs
85: #define matcreateseqbdiag_ matcreateseqbdiag
86: #define matcreatempibdiag_ matcreatempibdiag
87: #define matcreateseqdense_ matcreateseqdense
88: #define matcreatempidense_ matcreatempidense
89: #define matconvert_ matconvert
90: #define matload_ matload
91: #define mattranspose_ mattranspose
92: #define matgetarray_ matgetarray
93: #define matrestorearray_ matrestorearray
94: #define matshellsetoperation_ matshellsetoperation
95: #define matfdcoloringcreate_ matfdcoloringcreate
96: #define matfdcoloringdestroy_ matfdcoloringdestroy
97: #define matfdcoloringsetfunctionsnes_ matfdcoloringsetfunctionsnes
98: #define matfdcoloringsetfunctionts_ matfdcoloringsetfunctionts
99: #define matcopy_ matcopy
100: #define matgetsubmatrices_ matgetsubmatrices
101: #define matgetcoloring_ matgetcoloring
102: #define matduplicate_ matduplicate
103: #define matzerorows_ matzerorows
104: #define matzerorowslocal_ matzerorowslocal
105: #define matpartitioningview_ matpartitioningview
106: #define matpartitioningcreate_ matpartitioningcreate
107: #define matpartitioningsetadjacency_ matpartitioningsetadjacency
108: #define matpartitioningapply_ matpartitioningapply
109: #define matcreatempiadj_ matcreatempiadj
110: #define matsetfromoptions_ matsetfromoptions
111: #define matsetvaluesstencil_ matsetvaluesstencil
112: #endif
114: EXTERN_C_BEGIN
116: void PETSC_STDCALL matsettype_(Mat *x,CHAR type_name PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
117: {
118: char *t;
120: FIXCHAR(type_name,len,t);
121: *MatSetType(*x,t);
122: FREECHAR(type_name,t);
123: }
125: void PETSC_STDCALL matsetvaluesstencil_(Mat *mat,int *m,MatStencil *idxm,int *n,MatStencil *idxn,Scalar *v,InsertMode *addv,
126: int *ierr)
127: {
128: *MatSetValuesStencil(*mat,*m,idxm,*n,idxn,v,*addv);
129: }
131: void PETSC_STDCALL matmpiaijgetseqaij_(Mat *A,Mat *Ad,Mat *Ao,int *ic,long *iic,int *ierr)
132: {
133: int *i;
134: *MatMPIAIJGetSeqAIJ(*A,Ad,Ao,&i);if (*ierr) return;
135: *iic = PetscIntAddressToFortran(ic,i);
136: }
138: void PETSC_STDCALL matmpibaijgetseqbaij_(Mat *A,Mat *Ad,Mat *Ao,int *ic,long *iic,int *ierr)
139: {
140: int *i;
141: *MatMPIBAIJGetSeqBAIJ(*A,Ad,Ao,&i);if (*ierr) return;
142: *iic = PetscIntAddressToFortran(ic,i);
143: }
145: void PETSC_STDCALL matgetrowij_(Mat *B,int *shift,PetscTruth *sym,int *n,int *ia,long *iia,int *ja,long *jja,
146: PetscTruth *done,int *ierr)
147: {
148: int *IA,*JA;
149: *MatGetRowIJ(*B,*shift,*sym,n,&IA,&JA,done);if (*ierr) return;
150: *iia = PetscIntAddressToFortran(ia,IA);
151: *jja = PetscIntAddressToFortran(ja,JA);
152: }
154: void PETSC_STDCALL matrestorerowij_(Mat *B,int *shift,PetscTruth *sym,int *n,int *ia,long *iia,int *ja,long *jja,
155: PetscTruth *done,int *ierr)
156: {
157: int *IA = PetscIntAddressFromFortran(ia,*iia),*JA = PetscIntAddressFromFortran(ja,*jja);
158: *MatRestoreRowIJ(*B,*shift,*sym,n,&IA,&JA,done);
159: }
161: void PETSC_STDCALL matsetfromoptions_(Mat *B,int *ierr)
162: {
163: *MatSetFromOptions(*B);
164: }
166: void PETSC_STDCALL matcreateseqaijwitharrays_(MPI_Comm *comm,int *m,int *n,int *i,int *j,Scalar *a,Mat *mat,int *ierr)
167: {
168: *MatCreateSeqAIJWithArrays((MPI_Comm)PetscToPointerComm(*comm),*m,*n,i,j,a,mat);
169: }
171: #include "src/mat/impls/adj/mpi/mpiadj.h"
172: void PETSC_STDCALL matcreatempiadj_(MPI_Comm *comm,int *m,int *n,int *i,int *j,int *values,Mat *A,int *ierr)
173: {
174: Mat_MPIAdj *adj;
176: if (FORTRANNULLINTEGER(values)) values = PETSC_NULL;
177: *MatCreateMPIAdj((MPI_Comm)PetscToPointerComm(*comm),*m,*n,i,j,values,A);
178: adj = (Mat_MPIAdj*)(*A)->data;
179: adj->freeaij = PETSC_FALSE;
180: }
182: void PETSC_STDCALL matpartitioningdestroy_(MatPartitioning *part,int *ierr)
183: {
184: *MatPartitioningDestroy(*part);
185: }
187: void PETSC_STDCALL matpartitioningcreate_(MPI_Comm *comm,MatPartitioning *part, int *ierr)
188: {
189: *MatPartitioningCreate((MPI_Comm)PetscToPointerComm(*comm),part);
190: }
192: void PETSC_STDCALL matpartitioningapply_(MatPartitioning *part,IS *is,int *ierr)
193: {
194: *MatPartitioningApply(*part,is);
195: }
197: void PETSC_STDCALL matpartitioningsetadjacency_(MatPartitioning *part,Mat *mat,int *ierr)
198: {
199: *MatPartitioningSetAdjacency(*part,*mat);
200: }
202: void PETSC_STDCALL matpartitioningview_(MatPartitioning *part,PetscViewer *viewer, int *ierr)
203: {
204: PetscViewer v;
205: PetscPatchDefaultViewers_Fortran(viewer,v);
206: *MatPartitioningView(*part,v);
207: }
209: void PETSC_STDCALL matpartitioningsettype_(MatPartitioning *part,CHAR type PETSC_MIXED_LEN(len),
210: int *ierr PETSC_END_LEN(len))
211: {
212: char *t;
213: FIXCHAR(type,len,t);
214: *MatPartitioningSetType(*part,t);
215: FREECHAR(type,t);
216: }
218: void PETSC_STDCALL matgetcoloring_(Mat *mat,CHAR type PETSC_MIXED_LEN(len),ISColoring *iscoloring,
219: int *ierr PETSC_END_LEN(len))
220: {
221: char *t;
222: FIXCHAR(type,len,t);
223: *MatGetColoring(*mat,t,iscoloring);
224: FREECHAR(type,t);
225: }
227: void PETSC_STDCALL matsetvalue_(Mat *mat,int *i,int *j,Scalar *va,InsertMode *mode,int *ierr)
228: {
229: /* cannot use MatSetValue() here since that usesCHKERRQ() which has a return in it */
230: *MatSetValues(*mat,1,i,1,j,va,*mode);
231: }
233: void PETSC_STDCALL matsetvaluelocal_(Mat *mat,int *i,int *j,Scalar *va,InsertMode *mode,int *ierr)
234: {
235: /* cannot use MatSetValueLocal() here since that usesCHKERRQ() which has a return in it */
236: *MatSetValuesLocal(*mat,1,i,1,j,va,*mode);
237: }
239: void PETSC_STDCALL matfdcoloringcreate_(Mat *mat,ISColoring *iscoloring,MatFDColoring *color,int *ierr)
240: {
241: *MatFDColoringCreate(*mat,*iscoloring,color);
242: }
244: /*
245: This is a poor way of storing the column and value pointers
246: generated by MatGetRow() to be returned with MatRestoreRow()
247: but there is not natural,good place else to store them. Hence
248: Fortran programmers can only have one outstanding MatGetRows()
249: at a time.
250: */
251: static int matgetrowactive = 0,*my_ocols = 0;
252: static Scalar *my_ovals = 0;
254: void PETSC_STDCALL matgetrow_(Mat *mat,int *row,int *ncols,int *cols,Scalar *vals,int *ierr)
255: {
256: int **oocols = &my_ocols;
257: Scalar **oovals = &my_ovals;
259: if (matgetrowactive) {
260: PetscError(__LINE__,"MatGetRow_Fortran",__FILE__,__SDIR__,1,0,
261: "Cannot have two MatGetRow() active simultaneouslyn
262: call MatRestoreRow() before calling MatGetRow() a second time");
263: *1;
264: return;
265: }
266: if (FORTRANNULLINTEGER(cols)) oocols = PETSC_NULL;
267: if (FORTRANNULLSCALAR(vals)) oovals = PETSC_NULL;
269: *MatGetRow(*mat,*row,ncols,oocols,oovals);
270: if (*ierr) return;
272: if (oocols) { *PetscMemcpy(cols,my_ocols,(*ncols)*sizeof(int)); if (*ierr) return;}
273: if (oovals) { *PetscMemcpy(vals,my_ovals,(*ncols)*sizeof(Scalar)); if (*ierr) return; }
274: matgetrowactive = 1;
275: }
277: void PETSC_STDCALL matrestorerow_(Mat *mat,int *row,int *ncols,int *cols,Scalar *vals,int *ierr)
278: {
279: int **oocols = &my_ocols;
280: Scalar **oovals = &my_ovals;
281: if (!matgetrowactive) {
282: PetscError(__LINE__,"MatRestoreRow_Fortran",__FILE__,__SDIR__,1,0,
283: "Must call MatGetRow() first");
284: *1;
285: return;
286: }
287: if (FORTRANNULLINTEGER(cols)) oocols = PETSC_NULL;
288: if (FORTRANNULLSCALAR(vals)) oovals = PETSC_NULL;
289: *MatRestoreRow(*mat,*row,ncols,oocols,oovals);
290: matgetrowactive = 0;
291: }
293: void PETSC_STDCALL matview_(Mat *mat,PetscViewer *vin,int *ierr)
294: {
295: PetscViewer v;
296: PetscPatchDefaultViewers_Fortran(vin,v);
297: *MatView(*mat,v);
298: }
300: void PETSC_STDCALL matcopy_(Mat *A,Mat *B,MatStructure *str,int *ierr)
301: {
302: *MatCopy(*A,*B,*str);
303: }
305: void PETSC_STDCALL matgetinfo_(Mat *mat,MatInfoType *flag,double *finfo,int *ierr)
306: {
307: *MatGetInfo(*mat,*flag,(MatInfo*)finfo);
308: }
310: void PETSC_STDCALL matgetarray_(Mat *mat,Scalar *fa,long *ia,int *ierr)
311: {
312: Scalar *mm;
313: int m,n;
315: *MatGetArray(*mat,&mm); if (*ierr) return;
316: *MatGetSize(*mat,&m,&n); if (*ierr) return;
317: *PetscScalarAddressToFortran((PetscObject)*mat,fa,mm,m*n,ia); if (*ierr) return;
318: }
320: void PETSC_STDCALL matrestorearray_(Mat *mat,Scalar *fa,long *ia,int *ierr)
321: {
322: Scalar *lx;
323: int m,n;
325: *MatGetSize(*mat,&m,&n); if (*ierr) return;
326: *PetscScalarAddressFromFortran((PetscObject)*mat,fa,*ia,m*n,&lx);if (*ierr) return;
327: *MatRestoreArray(*mat,&lx);if (*ierr) return;
328: }
330: void PETSC_STDCALL mattranspose_(Mat *mat,Mat *B,int *ierr)
331: {
332: if (FORTRANNULLINTEGER(B)) B = PETSC_NULL;
333: *MatTranspose(*mat,B);
334: }
336: void PETSC_STDCALL matload_(PetscViewer *viewer,CHAR outtype PETSC_MIXED_LEN(len),Mat *newmat,int *ierr PETSC_END_LEN(len))
337: {
338: char *t;
339: PetscViewer v;
340: FIXCHAR(outtype,len,t);
341: PetscPatchDefaultViewers_Fortran(viewer,v);
342: *MatLoad(v,t,newmat);
343: FREECHAR(outtype,t);
344: }
346: void PETSC_STDCALL matconvert_(Mat *mat,CHAR outtype PETSC_MIXED_LEN(len),Mat *M,int *ierr PETSC_END_LEN(len))
347: {
348: char *t;
349: FIXCHAR(outtype,len,t);
350: *MatConvert(*mat,t,M);
351: FREECHAR(outtype,t);
352: }
354: void PETSC_STDCALL matcreateseqdense_(MPI_Comm *comm,int *m,int *n,Scalar *data,Mat *newmat,int *ierr)
355: {
356: if (FORTRANNULLSCALAR(data)) data = PETSC_NULL;
357: *MatCreateSeqDense((MPI_Comm)PetscToPointerComm(*comm),*m,*n,data,newmat);
358: }
360: void PETSC_STDCALL matcreatempidense_(MPI_Comm *comm,int *m,int *n,int *M,int *N,Scalar *data,Mat *newmat,
361: int *ierr)
362: {
363: if (FORTRANNULLSCALAR(data)) data = PETSC_NULL;
364: *MatCreateMPIDense((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*M,*N,data,newmat);
365: }
367: /* Fortran ignores diagv */
368: void PETSC_STDCALL matcreatempibdiag_(MPI_Comm *comm,int *m,int *M,int *N,int *nd,int *bs,
369: int *diag,Scalar **diagv,Mat *newmat,int *ierr)
370: {
371: *MatCreateMPIBDiag((MPI_Comm)PetscToPointerComm(*comm),
372: *m,*M,*N,*nd,*bs,diag,PETSC_NULL,newmat);
373: }
375: /* Fortran ignores diagv */
376: void PETSC_STDCALL matcreateseqbdiag_(MPI_Comm *comm,int *m,int *n,int *nd,int *bs,
377: int *diag,Scalar **diagv,Mat *newmat,int *ierr)
378: {
379: *MatCreateSeqBDiag((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*nd,*bs,diag,
380: PETSC_NULL,newmat);
381: }
383: #if defined(PETSC_HAVE_BLOCKSOLVE) && !defined(PETSC_USE_COMPLEX)
384: /* Fortran cannot pass in procinfo,hence ignored */
385: void PETSC_STDCALL matcreatempirowbs_(MPI_Comm *comm,int *m,int *M,int *nz,int *nnz,Mat *newmat,int *ierr)
386: {
387: if (FORTRANNULLINTEGER(nnz)) nnz = PETSC_NULL;
388: *MatCreateMPIRowbs((MPI_Comm)PetscToPointerComm(*comm),*m,*M,*nz,nnz,newmat);
389: }
390: #endif
392: void PETSC_STDCALL matgetordering_(Mat *mat,CHAR type PETSC_MIXED_LEN(len),IS *rperm,IS *cperm,
393: int *ierr PETSC_END_LEN(len))
394: {
395: char *t;
396: FIXCHAR(type,len,t);
397: *MatGetOrdering(*mat,t,rperm,cperm);
398: FREECHAR(type,t);
399: }
401: void PETSC_STDCALL matorderingregisterdestroy_(int *ierr)
402: {
403: *MatOrderingRegisterDestroy();
404: }
406: void PETSC_STDCALL matgettype_(Mat *mm,CHAR name PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
407: {
408: char *tname;
410: *MatGetType(*mm,&tname);
411: #if defined(PETSC_USES_CPTOFCD)
412: {
413: char *t = _fcdtocp(name); int len1 = _fcdlen(name);
414: if (t != PETSC_NULL_CHARACTER_Fortran) {
415: *PetscStrncpy(t,tname,len1);if (*ierr) return;
416: }
417: }
418: #else
419: if (name != PETSC_NULL_CHARACTER_Fortran) {
420: *PetscStrncpy(name,tname,len);if (*ierr) return;
421: }
422: #endif
423: }
425: void PETSC_STDCALL matcreate_(MPI_Comm *comm,int *m,int *n,int *M,int *N,Mat *V,int *ierr)
426: {
427: *MatCreate((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*M,*N,V);
428: }
430: void PETSC_STDCALL matcreateseqaij_(MPI_Comm *comm,int *m,int *n,int *nz,
431: int *nnz,Mat *newmat,int *ierr)
432: {
433: if (FORTRANNULLSCALAR(nnz)) {
434: PetscError(__LINE__,"MatCreateSeqAIJ_Fortran",__FILE__,__SDIR__,1,0,
435: "Cannot pass PETSC_NULL_SCALAR,use PETSC_NULL_INT");
436: *1;
437: return;
438: }
439: if (FORTRANNULLINTEGER(nnz)) nnz = PETSC_NULL;
440: *MatCreateSeqAIJ((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*nz,nnz,newmat);
441: }
443: void PETSC_STDCALL matcreateseqbaij_(MPI_Comm *comm,int *bs,int *m,int *n,int *nz,
444: int *nnz,Mat *newmat,int *ierr)
445: {
446: if (FORTRANNULLSCALAR(nnz)) {
447: PetscError(__LINE__,"MatCreateSeqBAIJ_Fortran",__FILE__,__SDIR__,1,0,
448: "Cannot pass PETSC_NULL_SCALAR,use PETSC_NULL_INT");
449: *1;
450: return;
451: }
452: if (FORTRANNULLINTEGER(nnz)) nnz = PETSC_NULL;
453: *MatCreateSeqBAIJ((MPI_Comm)PetscToPointerComm(*comm),*bs,*m,*n,*nz,nnz,newmat);
454: }
456: void PETSC_STDCALL matfdcoloringdestroy_(MatFDColoring *mat,int *ierr)
457: {
458: *MatFDColoringDestroy(*mat);
459: }
461: void PETSC_STDCALL matdestroy_(Mat *mat,int *ierr)
462: {
463: *MatDestroy(*mat);
464: }
466: void PETSC_STDCALL matcreatempiaij_(MPI_Comm *comm,int *m,int *n,int *M,int *N,
467: int *d_nz,int *d_nnz,int *o_nz,int *o_nnz,Mat *newmat,int *ierr)
468: {
469: if (FORTRANNULLSCALAR(d_nnz) || FORTRANNULLSCALAR(o_nnz)) {
470: PetscError(__LINE__,"MatCreateMPIAIJ_Fortran",__FILE__,__SDIR__,1,0,
471: "Cannot pass PETSC_NULL_SCALAR,use PETSC_NULL_INT");
472: *1;
473: return;
474: }
475: if (FORTRANNULLINTEGER(d_nnz)) d_nnz = PETSC_NULL;
476: if (FORTRANNULLINTEGER(o_nnz)) o_nnz = PETSC_NULL;
477: *MatCreateMPIAIJ((MPI_Comm)PetscToPointerComm(*comm),
478: *m,*n,*M,*N,*d_nz,d_nnz,*o_nz,o_nnz,newmat);
479: }
480: void PETSC_STDCALL matcreatempibaij_(MPI_Comm *comm,int *bs,int *m,int *n,int *M,int *N,
481: int *d_nz,int *d_nnz,int *o_nz,int *o_nnz,Mat *newmat,int *ierr)
482: {
483: if (FORTRANNULLSCALAR(d_nnz) || FORTRANNULLSCALAR(o_nnz)) {
484: PetscError(__LINE__,"MatCreateMPIBAIJ_Fortran",__FILE__,__SDIR__,1,0,
485: "Cannot pass PETSC_NULL_SCALAR,use PETSC_NULL_INT");
486: *1;
487: return;
488: }
489: if (FORTRANNULLINTEGER(d_nnz)) d_nnz = PETSC_NULL;
490: if (FORTRANNULLINTEGER(o_nnz)) o_nnz = PETSC_NULL;
491: *MatCreateMPIBAIJ((MPI_Comm)PetscToPointerComm(*comm),
492: *bs,*m,*n,*M,*N,*d_nz,d_nnz,*o_nz,o_nnz,newmat);
493: }
495: /*
496: The MatShell Matrix Vector product requires a C routine.
497: This C routine then calls the corresponding Fortran routine that was
498: set by the user.
499: */
500: void PETSC_STDCALL matcreateshell_(MPI_Comm *comm,int *m,int *n,int *M,int *N,void **ctx,Mat *mat,int *ierr)
501: {
502: *MatCreateShell((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*M,*N,*ctx,mat);
503: if (*ierr) return;
504: *PetscMalloc(sizeof(void *),&((PetscObject)*mat)->fortran_func_pointers);
505: }
507: static int ourmult(Mat mat,Vec x,Vec y)
508: {
509: int 0;
510: (*(int (PETSC_STDCALL *)(Mat*,Vec*,Vec*,int*))(((PetscObject)mat)->fortran_func_pointers[0]))(&mat,&x,&y,&ierr);
511: return ierr;
512: }
514: void PETSC_STDCALL matshellsetoperation_(Mat *mat,MatOperation *op,int (PETSC_STDCALL *f)(Mat*,Vec*,Vec*,int*),int *ierr)
515: {
516: if (*op == MATOP_MULT) {
517: *MatShellSetOperation(*mat,*op,(void(*)())ourmult);
518: ((PetscObject)*mat)->fortran_func_pointers[0] = (void(*)())f;
519: } else {
520: PetscError(__LINE__,"MatShellSetOperation_Fortran",__FILE__,__SDIR__,1,0,
521: "Cannot set that matrix operation");
522: *1;
523: }
524: }
526: #include "petscts.h"
527: /*
528: MatFDColoringSetFunction sticks the Fortran function into the fortran_func_pointers
529: this function is then accessed by ourmatfdcoloringfunction()
531: NOTE: FORTRAN USER CANNOT PUT IN A NEW J OR B currently.
533: USER CAN HAVE ONLY ONE MatFDColoring in code Because there is no place to hang f7!
534: */
536: static void (PETSC_STDCALL *f7)(TS*,double*,Vec*,Vec*,void*,int*);
538: static int ourmatfdcoloringfunctionts(TS ts,double t,Vec x,Vec y,void *ctx)
539: {
540: int 0;
541: (*f7)(&ts,&t,&x,&y,ctx,&ierr);
542: return ierr;
543: }
545: void PETSC_STDCALL matfdcoloringsetfunctionts_(MatFDColoring *fd,void (PETSC_STDCALL *f)(TS*,double*,Vec*,Vec*,void*,int*),
546: void *ctx,int *ierr)
547: {
548: f7 = f;
549: *MatFDColoringSetFunction(*fd,(int (*)(void))ourmatfdcoloringfunctionts,ctx);
550: }
552: static void (PETSC_STDCALL *f8)(SNES*,Vec*,Vec*,void*,int*);
554: static int ourmatfdcoloringfunctionsnes(SNES ts,Vec x,Vec y,void *ctx)
555: {
556: int 0;
557: (*f8)(&ts,&x,&y,ctx,&ierr);
558: return ierr;
559: }
562: void PETSC_STDCALL matfdcoloringsetfunctionsnes_(MatFDColoring *fd,void (PETSC_STDCALL *f)(SNES*,Vec*,Vec*,void*,int*),
563: void *ctx,int *ierr)
564: {
565: f8 = f;
566: *MatFDColoringSetFunction(*fd,(int (*)(void))ourmatfdcoloringfunctionsnes,ctx);
567: }
569: /*
570: MatGetSubmatrices() is slightly different from C since the
571: Fortran provides the array to hold the submatrix objects,while in C that
572: array is allocated by the MatGetSubmatrices()
573: */
574: void PETSC_STDCALL matgetsubmatrices_(Mat *mat,int *n,IS *isrow,IS *iscol,MatReuse *scall,Mat *smat,int *ierr)
575: {
576: Mat *lsmat;
577: int i;
579: if (*scall == MAT_INITIAL_MATRIX) {
580: *MatGetSubMatrices(*mat,*n,isrow,iscol,*scall,&lsmat);
581: for (i=0; i<*n; i++) {
582: smat[i] = lsmat[i];
583: }
584: PetscFree(lsmat);
585: } else {
586: *MatGetSubMatrices(*mat,*n,isrow,iscol,*scall,&smat);
587: }
588: }
590: void PETSC_STDCALL matduplicate_(Mat *matin,MatDuplicateOption *op,Mat *matout,int *ierr)
591: {
592: *MatDuplicate(*matin,*op,matout);
593: }
595: void PETSC_STDCALL matzerorows_(Mat *mat,IS *is,Scalar *diag,int *ierr)
596: {
597: if (FORTRANNULLSCALAR(diag)) diag = PETSC_NULL;
598: *MatZeroRows(*mat,*is,diag);
599: }
601: void PETSC_STDCALL matzerorowslocal_(Mat *mat,IS *is,Scalar *diag,int *ierr)
602: {
603: if (FORTRANNULLSCALAR(diag)) diag = PETSC_NULL;
604: *MatZeroRowsLocal(*mat,*is,diag);
605: }
607: EXTERN_C_END