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