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