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
8: #ifdef PETSC_HAVE_FORTRAN_CAPS
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 matorderingregisterdestroy_ MATORDERINGREGISTERDESTROY
32: #define matcreatempirowbs_ MATCREATEMPIROWBS
33: #define matcreateseqbdiag_ MATCREATESEQBDIAG
34: #define matcreatempibdiag_ MATCREATEMPIBDIAG
35: #define matcreateseqdense_ MATCREATESEQDENSE
36: #define matcreatempidense_ MATCREATEMPIDENSE
37: #define matconvert_ MATCONVERT
38: #define matload_ MATLOAD
39: #define mattranspose_ MATTRANSPOSE
40: #define matgetarray_ MATGETARRAY
41: #define matrestorearray_ MATRESTOREARRAY
42: #define matgettype_ MATGETTYPE
43: #define matgetinfo_ MATGETINFO
44: #define matshellsetoperation_ MATSHELLSETOPERATION
45: #define matview_ MATVIEW
46: #define matfdcoloringcreate_ MATFDCOLORINGCREATE
47: #define matfdcoloringdestroy_ MATFDCOLORINGDESTROY
48: #define matfdcoloringsetfunctionsnes_ MATFDCOLORINGSETFUNCTIONSNES
49: #define matfdcoloringsetfunctionts_ MATFDCOLORINGSETFUNCTIONTS
50: #define matcopy_ MATCOPY
51: #define matgetsubmatrices_ MATGETSUBMATRICES
52: #define matgetcoloring_ MATGETCOLORING
53: #define matpartitioningsettype_ MATPARTITIONINGSETTYPE
54: #define matduplicate_ MATDUPLICATE
55: #define matzerorows_ MATZEROROWS
56: #define matzerorowslocal_ MATZEROROWSLOCAL
57: #define matpartitioningview_ MATPARTITIONINGVIEW
58: #define matpartitioningcreate_ MATPARTITIONINGCREATE
59: #define matpartitioningsetadjacency_ MATPARTITIONINGSETADJACENCY
60: #define matpartitioningapply_ MATPARTITIONINGAPPLY
61: #define matcreatempiadj_ MATCREATEMPIADJ
62: #define matsetvaluesstencil_ MATSETVALUESSTENCIL
63: #define matseqaijsetpreallocation_ MATSEQAIJSETPREALLOCATION
64: #define matmpiaijsetpreallocation_ MATMPIAIJSETPREALLOCATION
65: #define matseqbaijsetpreallocation_ MATSEQBAIJSETPREALLOCATION
66: #define matmpibaijsetpreallocation_ MATMPIBAIJSETPREALLOCATION
67: #define matseqsbaijsetpreallocation_ MATSEQSBAIJSETPREALLOCATION
68: #define matmpisbaijsetpreallocation_ MATMPISBAIJSETPREALLOCATION
69: #define matseqbdiagsetpreallocation_ MATSEQBDIAGSETPREALLOCATION
70: #define matmpibdiagsetpreallocation_ MATMPIBDIAGSETPREALLOCATION
71: #define matseqdensesetpreallocation_ MATSEQDENSESETPREALLOCATION
72: #define matmpidensesetpreallocation_ MATMPIDENSESETPREALLOCATION
73: #define matmpiadjsetpreallocation_ MATMPIADJSETPREALLOCATION
74: #define matmpirowbssetpreallocation_ MATMPIROWBSSETPREALLOCATION
75: #define matpartitioningpartysetglobal_ MATPARTITIONINGPARTYSETGLOBAL
76: #define matpartitioningpartysetlocal_ MATPARTITIONINGPARTYSETLOCAL
77: #define matpartitioningscotchsetstrategy_ MATPARTITIONINGSCOTCHSETSTRATEGY
78: #define matpartitioningscotchsetarch_ MATPARTITIONINGSCOTCHSETARCH
79: #define matpartitioningscotchsethostlist_ MATPARTITIONINGSCOTCHSETHOSTLIST
80: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
81: #define matsettype_ matsettype
82: #define matmpiaijgetseqaij_ matmpiaijgetseqaij
83: #define matmpibaijgetseqbaij_ matmpibaijgetseqbaij
84: #define matrestorerowij_ matrestorerowij
85: #define matgetrowij_ matgetrowij
86: #define matcreateseqaijwitharrays_ matcreateseqaijwitharrays
87: #define matpartitioningdestroy_ matpartitioningdestroy
88: #define matpartitioningsettype_ matpartitioningsettype
89: #define matsetvalue_ matsetvalue
90: #define matsetvaluelocal_ matsetvaluelocal
91: #define matgetrow_ matgetrow
92: #define matrestorerow_ matrestorerow
93: #define matview_ matview
94: #define matgetinfo_ matgetinfo
95: #define matgettype_ matgettype
96: #define matdestroy_ matdestroy
97: #define matcreatempiaij_ matcreatempiaij
98: #define matcreateseqaij_ matcreateseqaij
99: #define matcreatempibaij_ matcreatempibaij
100: #define matcreateseqbaij_ matcreateseqbaij
101: #define matcreatempisbaij_ matcreatempisbaij
102: #define matcreateseqsbaij_ matcreateseqsbaij
103: #define matcreate_ matcreate
104: #define matcreateshell_ matcreateshell
105: #define matorderingregisterdestroy_ matorderingregisterdestroy
106: #define matgetordering_ matgetordering
107: #define matcreatempirowbs_ matcreatempirowbs
108: #define matcreateseqbdiag_ matcreateseqbdiag
109: #define matcreatempibdiag_ matcreatempibdiag
110: #define matcreateseqdense_ matcreateseqdense
111: #define matcreatempidense_ matcreatempidense
112: #define matconvert_ matconvert
113: #define matload_ matload
114: #define mattranspose_ mattranspose
115: #define matgetarray_ matgetarray
116: #define matrestorearray_ matrestorearray
117: #define matshellsetoperation_ matshellsetoperation
118: #define matfdcoloringcreate_ matfdcoloringcreate
119: #define matfdcoloringdestroy_ matfdcoloringdestroy
120: #define matfdcoloringsetfunctionsnes_ matfdcoloringsetfunctionsnes
121: #define matfdcoloringsetfunctionts_ matfdcoloringsetfunctionts
122: #define matcopy_ matcopy
123: #define matgetsubmatrices_ matgetsubmatrices
124: #define matgetcoloring_ matgetcoloring
125: #define matduplicate_ matduplicate
126: #define matzerorows_ matzerorows
127: #define matzerorowslocal_ matzerorowslocal
128: #define matpartitioningview_ matpartitioningview
129: #define matpartitioningcreate_ matpartitioningcreate
130: #define matpartitioningsetadjacency_ matpartitioningsetadjacency
131: #define matpartitioningapply_ matpartitioningapply
132: #define matcreatempiadj_ matcreatempiadj
133: #define matsetfromoptions_ matsetfromoptions
134: #define matsetvaluesstencil_ matsetvaluesstencil
135: #define matseqaijsetpreallocation_ matseqaijsetpreallocation
136: #define matmpiaijsetpreallocation_ matmpiaijsetpreallocation
137: #define matseqbaijsetpreallocation_ matseqbaijsetpreallocation
138: #define matmpibaijsetpreallocation_ matmpibaijsetpreallocation
139: #define matseqsbaijsetpreallocation_ matseqsbaijsetpreallocation
140: #define matmpisbaijsetpreallocation_ matmpisbaijsetpreallocation
141: #define matseqbdiagsetpreallocation_ matseqbdiagsetpreallocation
142: #define matmpibdiagsetpreallocation_ matmpibdiagsetpreallocation
143: #define matseqdensesetpreallocation_ matseqdensesetpreallocation
144: #define matmpidensesetpreallocation_ matmpidensesetpreallocation
145: #define matmpiadjsetpreallocation_ matmpiadjsetpreallocation
146: #define matmpirowbssetpreallocation_ matmpirowbssetpreallocation
147: #define matpartitioningpartysetglobal_ matpartitioningpartysetglobal
148: #define matpartitioningpartysetlocal_ matpartitioningpartysetlocal
149: #define matpartitioningscotchsetstrategy_ matpartitioningscotchsetstrategy
150: #define matpartitioningscotchsetarch_ matpartitioningscotchsetarch
151: #define matpartitioningscotchsethostlist_ matpartitioningscotchsethostlist
152: #endif
154: #include petscts.h
156: EXTERN_C_BEGIN
157: static void (PETSC_STDCALL *f7)(TS*,double*,Vec*,Vec*,void*,int*);
158: static void (PETSC_STDCALL *f8)(SNES*,Vec*,Vec*,void*,int*);
159: EXTERN_C_END
161: /* These are not extern C because they are passed into non-extern C user level functions */
162: static int ourmatfdcoloringfunctionts(TS ts,double t,Vec x,Vec y,void *ctx)
163: {
164: int 0;
165: (*f7)(&ts,&t,&x,&y,ctx,&ierr);
166: return ierr;
167: }
169: static int ourmatfdcoloringfunctionsnes(SNES ts,Vec x,Vec y,void *ctx)
170: {
171: int 0;
172: (*f8)(&ts,&x,&y,ctx,&ierr);
173: return ierr;
174: }
176: EXTERN_C_BEGIN
178: void PETSC_STDCALL matsettype_(Mat *x,CHAR type_name PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
179: {
180: char *t;
182: FIXCHAR(type_name,len,t);
183: *MatSetType(*x,t);
184: FREECHAR(type_name,t);
185: }
187: void PETSC_STDCALL matsetvaluesstencil_(Mat *mat,int *m,MatStencil *idxm,int *n,MatStencil *idxn,PetscScalar *v,InsertMode *addv,
188: int *ierr)
189: {
190: *MatSetValuesStencil(*mat,*m,idxm,*n,idxn,v,*addv);
191: }
193: void PETSC_STDCALL matmpiaijgetseqaij_(Mat *A,Mat *Ad,Mat *Ao,int *ic,long *iic,int *ierr)
194: {
195: int *i;
196: *MatMPIAIJGetSeqAIJ(*A,Ad,Ao,&i);if (*ierr) return;
197: *iic = PetscIntAddressToFortran(ic,i);
198: }
200: void PETSC_STDCALL matmpibaijgetseqbaij_(Mat *A,Mat *Ad,Mat *Ao,int *ic,long *iic,int *ierr)
201: {
202: int *i;
203: *MatMPIBAIJGetSeqBAIJ(*A,Ad,Ao,&i);if (*ierr) return;
204: *iic = PetscIntAddressToFortran(ic,i);
205: }
207: void PETSC_STDCALL matgetrowij_(Mat *B,int *shift,PetscTruth *sym,int *n,int *ia,long *iia,int *ja,long *jja,
208: PetscTruth *done,int *ierr)
209: {
210: int *IA,*JA;
211: *MatGetRowIJ(*B,*shift,*sym,n,&IA,&JA,done);if (*ierr) return;
212: *iia = PetscIntAddressToFortran(ia,IA);
213: *jja = PetscIntAddressToFortran(ja,JA);
214: }
216: void PETSC_STDCALL matrestorerowij_(Mat *B,int *shift,PetscTruth *sym,int *n,int *ia,long *iia,int *ja,long *jja,
217: PetscTruth *done,int *ierr)
218: {
219: int *IA = PetscIntAddressFromFortran(ia,*iia),*JA = PetscIntAddressFromFortran(ja,*jja);
220: *MatRestoreRowIJ(*B,*shift,*sym,n,&IA,&JA,done);
221: }
223: void PETSC_STDCALL matsetfromoptions_(Mat *B,int *ierr)
224: {
225: *MatSetFromOptions(*B);
226: }
228: void PETSC_STDCALL matcreateseqaijwitharrays_(MPI_Comm *comm,int *m,int *n,int *i,int *j,PetscScalar *a,Mat *mat,int *ierr)
229: {
230: *MatCreateSeqAIJWithArrays((MPI_Comm)PetscToPointerComm(*comm),*m,*n,i,j,a,mat);
231: }
233: void PETSC_STDCALL matcreatempiadj_(MPI_Comm *comm,int *m,int *n,int *i,int *j,int *values,Mat *A,int *ierr)
234: {
235: Mat_MPIAdj *adj;
237: CHKFORTRANNULLINTEGER(values);
238: *MatCreateMPIAdj((MPI_Comm)PetscToPointerComm(*comm),*m,*n,i,j,values,A);
239: adj = (Mat_MPIAdj*)(*A)->data;
240: adj->freeaij = PETSC_FALSE;
241: }
243: void PETSC_STDCALL matpartitioningdestroy_(MatPartitioning *part,int *ierr)
244: {
245: *MatPartitioningDestroy(*part);
246: }
248: void PETSC_STDCALL matpartitioningcreate_(MPI_Comm *comm,MatPartitioning *part, int *ierr)
249: {
250: *MatPartitioningCreate((MPI_Comm)PetscToPointerComm(*comm),part);
251: }
253: void PETSC_STDCALL matpartitioningapply_(MatPartitioning *part,IS *is,int *ierr)
254: {
255: *MatPartitioningApply(*part,is);
256: }
258: void PETSC_STDCALL matpartitioningsetadjacency_(MatPartitioning *part,Mat *mat,int *ierr)
259: {
260: *MatPartitioningSetAdjacency(*part,*mat);
261: }
263: void PETSC_STDCALL matpartitioningview_(MatPartitioning *part,PetscViewer *viewer, int *ierr)
264: {
265: PetscViewer v;
266: PetscPatchDefaultViewers_Fortran(viewer,v);
267: *MatPartitioningView(*part,v);
268: }
270: void PETSC_STDCALL matpartitioningsettype_(MatPartitioning *part,CHAR type PETSC_MIXED_LEN(len),
271: int *ierr PETSC_END_LEN(len))
272: {
273: char *t;
274: FIXCHAR(type,len,t);
275: *MatPartitioningSetType(*part,t);
276: FREECHAR(type,t);
277: }
279: void PETSC_STDCALL matgetcoloring_(Mat *mat,CHAR type PETSC_MIXED_LEN(len),ISColoring *iscoloring,
280: int *ierr PETSC_END_LEN(len))
281: {
282: char *t;
283: FIXCHAR(type,len,t);
284: *MatGetColoring(*mat,t,iscoloring);
285: FREECHAR(type,t);
286: }
288: void PETSC_STDCALL matsetvalue_(Mat *mat,int *i,int *j,PetscScalar *va,InsertMode *mode,int *ierr)
289: {
290: /* cannot use MatSetValue() here since that usesCHKERRQ() which has a return in it */
291: *MatSetValues(*mat,1,i,1,j,va,*mode);
292: }
294: void PETSC_STDCALL matsetvaluelocal_(Mat *mat,int *i,int *j,PetscScalar *va,InsertMode *mode,int *ierr)
295: {
296: /* cannot use MatSetValueLocal() here since that usesCHKERRQ() which has a return in it */
297: *MatSetValuesLocal(*mat,1,i,1,j,va,*mode);
298: }
300: void PETSC_STDCALL matfdcoloringcreate_(Mat *mat,ISColoring *iscoloring,MatFDColoring *color,int *ierr)
301: {
302: *MatFDColoringCreate(*mat,*iscoloring,color);
303: }
305: /*
306: This is a poor way of storing the column and value pointers
307: generated by MatGetRow() to be returned with MatRestoreRow()
308: but there is not natural,good place else to store them. Hence
309: Fortran programmers can only have one outstanding MatGetRows()
310: at a time.
311: */
312: static int matgetrowactive = 0,*my_ocols = 0;
313: static PetscScalar *my_ovals = 0;
315: void PETSC_STDCALL matgetrow_(Mat *mat,int *row,int *ncols,int *cols,PetscScalar *vals,int *ierr)
316: {
317: int **oocols = &my_ocols;
318: PetscScalar **oovals = &my_ovals;
320: if (matgetrowactive) {
321: PetscError(__LINE__,"MatGetRow_Fortran",__FILE__,__SDIR__,1,0,
322: "Cannot have two MatGetRow() active simultaneously\n\
323: call MatRestoreRow() before calling MatGetRow() a second time");
324: *1;
325: return;
326: }
328: CHKFORTRANNULLINTEGER(cols); if (!cols) oocols = PETSC_NULL;
329: CHKFORTRANNULLSCALAR(vals); if (!vals) oovals = PETSC_NULL;
331: *MatGetRow(*mat,*row,ncols,oocols,oovals);
332: if (*ierr) return;
334: if (oocols) { *PetscMemcpy(cols,my_ocols,(*ncols)*sizeof(int)); if (*ierr) return;}
335: if (oovals) { *PetscMemcpy(vals,my_ovals,(*ncols)*sizeof(PetscScalar)); if (*ierr) return; }
336: matgetrowactive = 1;
337: }
339: void PETSC_STDCALL matrestorerow_(Mat *mat,int *row,int *ncols,int *cols,PetscScalar *vals,int *ierr)
340: {
341: int **oocols = &my_ocols;
342: PetscScalar **oovals = &my_ovals;
343: if (!matgetrowactive) {
344: PetscError(__LINE__,"MatRestoreRow_Fortran",__FILE__,__SDIR__,1,0,
345: "Must call MatGetRow() first");
346: *1;
347: return;
348: }
349: CHKFORTRANNULLINTEGER(cols); if (!cols) oocols = PETSC_NULL;
350: CHKFORTRANNULLSCALAR(vals); if (!vals) oovals = PETSC_NULL;
352: *MatRestoreRow(*mat,*row,ncols,oocols,oovals);
353: matgetrowactive = 0;
354: }
356: void PETSC_STDCALL matview_(Mat *mat,PetscViewer *vin,int *ierr)
357: {
358: PetscViewer v;
359: PetscPatchDefaultViewers_Fortran(vin,v);
360: *MatView(*mat,v);
361: }
363: void PETSC_STDCALL matcopy_(Mat *A,Mat *B,MatStructure *str,int *ierr)
364: {
365: *MatCopy(*A,*B,*str);
366: }
368: void PETSC_STDCALL matgetinfo_(Mat *mat,MatInfoType *flag,double *finfo,int *ierr)
369: {
370: *MatGetInfo(*mat,*flag,(MatInfo*)finfo);
371: }
373: void PETSC_STDCALL matgetarray_(Mat *mat,PetscScalar *fa,long *ia,int *ierr)
374: {
375: PetscScalar *mm;
376: int m,n;
378: *MatGetArray(*mat,&mm); if (*ierr) return;
379: *MatGetSize(*mat,&m,&n); if (*ierr) return;
380: *PetscScalarAddressToFortran((PetscObject)*mat,fa,mm,m*n,ia); if (*ierr) return;
381: }
383: void PETSC_STDCALL matrestorearray_(Mat *mat,PetscScalar *fa,long *ia,int *ierr)
384: {
385: PetscScalar *lx;
386: int m,n;
388: *MatGetSize(*mat,&m,&n); if (*ierr) return;
389: *PetscScalarAddressFromFortran((PetscObject)*mat,fa,*ia,m*n,&lx);if (*ierr) return;
390: *MatRestoreArray(*mat,&lx);if (*ierr) return;
391: }
393: void PETSC_STDCALL mattranspose_(Mat *mat,Mat *B,int *ierr)
394: {
395: CHKFORTRANNULLINTEGER(B);
396: *MatTranspose(*mat,B);
397: }
399: void PETSC_STDCALL matload_(PetscViewer *viewer,CHAR outtype PETSC_MIXED_LEN(len),Mat *newmat,int *ierr PETSC_END_LEN(len))
400: {
401: char *t;
402: PetscViewer v;
403: FIXCHAR(outtype,len,t);
404: PetscPatchDefaultViewers_Fortran(viewer,v);
405: *MatLoad(v,t,newmat);
406: FREECHAR(outtype,t);
407: }
409: void PETSC_STDCALL matconvert_(Mat *mat,CHAR outtype PETSC_MIXED_LEN(len),Mat *M,int *ierr PETSC_END_LEN(len))
410: {
411: char *t;
412: FIXCHAR(outtype,len,t);
413: *MatConvert(*mat,t,M);
414: FREECHAR(outtype,t);
415: }
417: void PETSC_STDCALL matcreateseqdense_(MPI_Comm *comm,int *m,int *n,PetscScalar *data,Mat *newmat,int *ierr)
418: {
419: CHKFORTRANNULLSCALAR(data);
420: *MatCreateSeqDense((MPI_Comm)PetscToPointerComm(*comm),*m,*n,data,newmat);
421: }
423: void PETSC_STDCALL matcreatempidense_(MPI_Comm *comm,int *m,int *n,int *M,int *N,PetscScalar *data,Mat *newmat,
424: int *ierr)
425: {
426: CHKFORTRANNULLSCALAR(data);
427: *MatCreateMPIDense((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*M,*N,data,newmat);
428: }
430: /* Fortran ignores diagv */
431: void PETSC_STDCALL matcreatempibdiag_(MPI_Comm *comm,int *m,int *M,int *N,int *nd,int *bs,
432: int *diag,PetscScalar **diagv,Mat *newmat,int *ierr)
433: {
434: CHKFORTRANNULLINTEGER(diag);
435: *MatCreateMPIBDiag((MPI_Comm)PetscToPointerComm(*comm),
436: *m,*M,*N,*nd,*bs,diag,PETSC_NULL,newmat);
437: }
439: /* Fortran ignores diagv */
440: void PETSC_STDCALL matcreateseqbdiag_(MPI_Comm *comm,int *m,int *n,int *nd,int *bs,
441: int *diag,PetscScalar **diagv,Mat *newmat,int *ierr)
442: {
443: CHKFORTRANNULLINTEGER(diag);
444: *MatCreateSeqBDiag((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*nd,*bs,diag,
445: PETSC_NULL,newmat);
446: }
448: #if defined(PETSC_HAVE_BLOCKSOLVE) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
449: /* Fortran cannot pass in procinfo,hence ignored */
450: void PETSC_STDCALL matcreatempirowbs_(MPI_Comm *comm,int *m,int *M,int *nz,int *nnz,Mat *newmat,int *ierr)
451: {
452: CHKFORTRANNULLINTEGER(nnz);
453: *MatCreateMPIRowbs((MPI_Comm)PetscToPointerComm(*comm),*m,*M,*nz,nnz,newmat);
454: }
455: #endif
457: void PETSC_STDCALL matgetordering_(Mat *mat,CHAR type PETSC_MIXED_LEN(len),IS *rperm,IS *cperm,
458: int *ierr PETSC_END_LEN(len))
459: {
460: char *t;
461: FIXCHAR(type,len,t);
462: *MatGetOrdering(*mat,t,rperm,cperm);
463: FREECHAR(type,t);
464: }
466: void PETSC_STDCALL matorderingregisterdestroy_(int *ierr)
467: {
468: *MatOrderingRegisterDestroy();
469: }
471: void PETSC_STDCALL matgettype_(Mat *mm,CHAR name PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
472: {
473: char *tname;
475: *MatGetType(*mm,&tname);
476: #if defined(PETSC_USES_CPTOFCD)
477: {
478: char *t = _fcdtocp(name); int len1 = _fcdlen(name);
479: if (t != PETSC_NULL_CHARACTER_Fortran) {
480: *PetscStrncpy(t,tname,len1);if (*ierr) return;
481: }
482: }
483: #else
484: if (name != PETSC_NULL_CHARACTER_Fortran) {
485: *PetscStrncpy(name,tname,len);if (*ierr) return;
486: }
487: #endif
488: }
490: void PETSC_STDCALL matcreate_(MPI_Comm *comm,int *m,int *n,int *M,int *N,Mat *V,int *ierr)
491: {
492: *MatCreate((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*M,*N,V);
493: }
495: void PETSC_STDCALL matcreateseqaij_(MPI_Comm *comm,int *m,int *n,int *nz,
496: int *nnz,Mat *newmat,int *ierr)
497: {
498: CHKFORTRANNULLINTEGER(nnz);
499: *MatCreateSeqAIJ((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*nz,nnz,newmat);
500: }
502: void PETSC_STDCALL matcreateseqbaij_(MPI_Comm *comm,int *bs,int *m,int *n,int *nz,
503: int *nnz,Mat *newmat,int *ierr)
504: {
505: CHKFORTRANNULLINTEGER(nnz);
506: *MatCreateSeqBAIJ((MPI_Comm)PetscToPointerComm(*comm),*bs,*m,*n,*nz,nnz,newmat);
507: }
509: void PETSC_STDCALL matcreateseqsbaij_(MPI_Comm *comm,int *bs,int *m,int *n,int *nz,
510: int *nnz,Mat *newmat,int *ierr)
511: {
512: CHKFORTRANNULLINTEGER(nnz);
513: *MatCreateSeqSBAIJ((MPI_Comm)PetscToPointerComm(*comm),*bs,*m,*n,*nz,nnz,newmat);
514: }
516: void PETSC_STDCALL matfdcoloringdestroy_(MatFDColoring *mat,int *ierr)
517: {
518: *MatFDColoringDestroy(*mat);
519: }
521: void PETSC_STDCALL matdestroy_(Mat *mat,int *ierr)
522: {
523: *MatDestroy(*mat);
524: }
526: void PETSC_STDCALL matcreatempiaij_(MPI_Comm *comm,int *m,int *n,int *M,int *N,
527: int *d_nz,int *d_nnz,int *o_nz,int *o_nnz,Mat *newmat,int *ierr)
528: {
529: CHKFORTRANNULLINTEGER(d_nnz);
530: CHKFORTRANNULLINTEGER(o_nnz);
532: *MatCreateMPIAIJ((MPI_Comm)PetscToPointerComm(*comm),
533: *m,*n,*M,*N,*d_nz,d_nnz,*o_nz,o_nnz,newmat);
534: }
535: void PETSC_STDCALL matcreatempibaij_(MPI_Comm *comm,int *bs,int *m,int *n,int *M,int *N,
536: int *d_nz,int *d_nnz,int *o_nz,int *o_nnz,Mat *newmat,int *ierr)
537: {
538: CHKFORTRANNULLINTEGER(d_nnz);
539: CHKFORTRANNULLINTEGER(o_nnz);
540: *MatCreateMPIBAIJ((MPI_Comm)PetscToPointerComm(*comm),
541: *bs,*m,*n,*M,*N,*d_nz,d_nnz,*o_nz,o_nnz,newmat);
542: }
543: void PETSC_STDCALL matcreatempisbaij_(MPI_Comm *comm,int *bs,int *m,int *n,int *M,int *N,
544: int *d_nz,int *d_nnz,int *o_nz,int *o_nnz,Mat *newmat,int *ierr)
545: {
546: CHKFORTRANNULLINTEGER(d_nnz);
547: CHKFORTRANNULLINTEGER(o_nnz);
548: *MatCreateMPISBAIJ((MPI_Comm)PetscToPointerComm(*comm),
549: *bs,*m,*n,*M,*N,*d_nz,d_nnz,*o_nz,o_nnz,newmat);
550: }
552: /*
553: The MatShell Matrix Vector product requires a C routine.
554: This C routine then calls the corresponding Fortran routine that was
555: set by the user.
556: */
557: void PETSC_STDCALL matcreateshell_(MPI_Comm *comm,int *m,int *n,int *M,int *N,void **ctx,Mat *mat,int *ierr)
558: {
559: *MatCreateShell((MPI_Comm)PetscToPointerComm(*comm),*m,*n,*M,*N,*ctx,mat);
560: if (*ierr) return;
561: *PetscMalloc(4*sizeof(void *),&((PetscObject)*mat)->fortran_func_pointers);
562: }
564: static int ourmult(Mat mat,Vec x,Vec y)
565: {
566: int 0;
567: (*(int (PETSC_STDCALL *)(Mat*,Vec*,Vec*,int*))(((PetscObject)mat)->fortran_func_pointers[0]))(&mat,&x,&y,&ierr);
568: return ierr;
569: }
571: static int ourmulttranspose(Mat mat,Vec x,Vec y)
572: {
573: int 0;
574: (*(int (PETSC_STDCALL *)(Mat*,Vec*,Vec*,int*))(((PetscObject)mat)->fortran_func_pointers[2]))(&mat,&x,&y,&ierr);
575: return ierr;
576: }
578: static int ourmultadd(Mat mat,Vec x,Vec y,Vec z)
579: {
580: int 0;
581: (*(int (PETSC_STDCALL *)(Mat*,Vec*,Vec*,Vec*,int*))(((PetscObject)mat)->fortran_func_pointers[1]))(&mat,&x,&y,&z,&ierr);
582: return ierr;
583: }
585: static int ourmulttransposeadd(Mat mat,Vec x,Vec y,Vec z)
586: {
587: int 0;
588: (*(int (PETSC_STDCALL *)(Mat*,Vec*,Vec*,Vec*,int*))(((PetscObject)mat)->fortran_func_pointers[3]))(&mat,&x,&y,&z,&ierr);
589: return ierr;
590: }
592: void PETSC_STDCALL matshellsetoperation_(Mat *mat,MatOperation *op,int (PETSC_STDCALL *f)(Mat*,Vec*,Vec*,int*),int *ierr)
593: {
594: if (*op == MATOP_MULT) {
595: *MatShellSetOperation(*mat,*op,(FCNVOID)ourmult);
596: ((PetscObject)*mat)->fortran_func_pointers[0] = (FCNVOID)f;
597: } else if (*op == MATOP_MULT_TRANSPOSE) {
598: *MatShellSetOperation(*mat,*op,(FCNVOID)ourmulttranspose);
599: ((PetscObject)*mat)->fortran_func_pointers[2] = (FCNVOID)f;
600: } else if (*op == MATOP_MULT_ADD) {
601: *MatShellSetOperation(*mat,*op,(FCNVOID)ourmultadd);
602: ((PetscObject)*mat)->fortran_func_pointers[1] = (FCNVOID)f;
603: } else if (*op == MATOP_MULT_TRANSPOSE_ADD) {
604: *MatShellSetOperation(*mat,*op,(FCNVOID)ourmulttransposeadd);
605: ((PetscObject)*mat)->fortran_func_pointers[3] = (FCNVOID)f;
606: } else {
607: PetscError(__LINE__,"MatShellSetOperation_Fortran",__FILE__,__SDIR__,1,0,
608: "Cannot set that matrix operation");
609: *1;
610: }
611: }
613: /*
614: MatFDColoringSetFunction sticks the Fortran function into the fortran_func_pointers
615: this function is then accessed by ourmatfdcoloringfunction()
617: NOTE: FORTRAN USER CANNOT PUT IN A NEW J OR B currently.
619: USER CAN HAVE ONLY ONE MatFDColoring in code Because there is no place to hang f7!
620: */
623: void PETSC_STDCALL matfdcoloringsetfunctionts_(MatFDColoring *fd,void (PETSC_STDCALL *f)(TS*,double*,Vec*,Vec*,void*,int*),
624: void *ctx,int *ierr)
625: {
626: f7 = f;
627: *MatFDColoringSetFunction(*fd,(FCNINTVOID)ourmatfdcoloringfunctionts,ctx);
628: }
631: void PETSC_STDCALL matfdcoloringsetfunctionsnes_(MatFDColoring *fd,void (PETSC_STDCALL *f)(SNES*,Vec*,Vec*,void*,int*),
632: void *ctx,int *ierr)
633: {
634: f8 = f;
635: *MatFDColoringSetFunction(*fd,(FCNINTVOID)ourmatfdcoloringfunctionsnes,ctx);
636: }
638: /*
639: MatGetSubmatrices() is slightly different from C since the
640: Fortran provides the array to hold the submatrix objects,while in C that
641: array is allocated by the MatGetSubmatrices()
642: */
643: void PETSC_STDCALL matgetsubmatrices_(Mat *mat,int *n,IS *isrow,IS *iscol,MatReuse *scall,Mat *smat,int *ierr)
644: {
645: Mat *lsmat;
646: int i;
648: if (*scall == MAT_INITIAL_MATRIX) {
649: *MatGetSubMatrices(*mat,*n,isrow,iscol,*scall,&lsmat);
650: for (i=0; i<*n; i++) {
651: smat[i] = lsmat[i];
652: }
653: PetscFree(lsmat);
654: } else {
655: *MatGetSubMatrices(*mat,*n,isrow,iscol,*scall,&smat);
656: }
657: }
659: void PETSC_STDCALL matduplicate_(Mat *matin,MatDuplicateOption *op,Mat *matout,int *ierr)
660: {
661: *MatDuplicate(*matin,*op,matout);
662: }
664: void PETSC_STDCALL matzerorows_(Mat *mat,IS *is,PetscScalar *diag,int *ierr)
665: {
666: CHKFORTRANNULLSCALAR(diag);
667: *MatZeroRows(*mat,*is,diag);
668: }
670: void PETSC_STDCALL matzerorowslocal_(Mat *mat,IS *is,PetscScalar *diag,int *ierr)
671: {
672: CHKFORTRANNULLSCALAR(diag);
673: *MatZeroRowsLocal(*mat,*is,diag);
674: }
676: void PETSC_STDCALL matseqaijsetpreallocation_(Mat *mat,int *nz,int *nnz,int *ierr)
677: {
678: CHKFORTRANNULLINTEGER(nnz);
679: CHKFORTRANNULLINTEGER(nnz);
680: *MatSeqAIJSetPreallocation(*mat,*nz,nnz);
681: }
683: void PETSC_STDCALL matmpiaijsetpreallocation_(Mat *mat,int *d_nz,int *d_nnz,int *o_nz,int *o_nnz,int *ierr)
684: {
685: CHKFORTRANNULLINTEGER(d_nnz);
686: CHKFORTRANNULLINTEGER(o_nnz);
687: *MatMPIAIJSetPreallocation(*mat,*d_nz,d_nnz,*o_nz,o_nnz);
688: }
690: void PETSC_STDCALL matseqbaijsetpreallocation_(Mat *mat,int *bs,int *nz,int *nnz,int *ierr)
691: {
692: CHKFORTRANNULLINTEGER(nnz);
693: CHKFORTRANNULLINTEGER(nnz);
694: *MatSeqBAIJSetPreallocation(*mat,*bs,*nz,nnz);
695: }
697: void PETSC_STDCALL matmpibaijsetpreallocation_(Mat *mat,int *bs,int *d_nz,int *d_nnz,int *o_nz,int *o_nnz,int *ierr)
698: {
699: CHKFORTRANNULLINTEGER(d_nnz);
700: CHKFORTRANNULLINTEGER(o_nnz);
701: *MatMPIBAIJSetPreallocation(*mat,*bs,*d_nz,d_nnz,*o_nz,o_nnz);
702: }
704: void PETSC_STDCALL matseqsbaijsetpreallocation_(Mat *mat,int *bs,int *nz,int *nnz,int *ierr)
705: {
706: CHKFORTRANNULLINTEGER(nnz);
707: CHKFORTRANNULLINTEGER(nnz);
708: *MatSeqSBAIJSetPreallocation(*mat,*bs,*nz,nnz);
709: }
711: void PETSC_STDCALL matmpisbaijsetpreallocation_(Mat *mat,int *bs,int *d_nz,int *d_nnz,int *o_nz,int *o_nnz,int *ierr)
712: {
713: CHKFORTRANNULLINTEGER(d_nnz);
714: CHKFORTRANNULLINTEGER(o_nnz);
715: *MatMPISBAIJSetPreallocation(*mat,*bs,*d_nz,d_nnz,*o_nz,o_nnz);
716: }
718: /* Fortran ignores diagv */
719: void PETSC_STDCALL matseqbdiagsetpreallocation_(Mat *mat,int *nd,int *bs,int *diag,PetscScalar **diagv,int *ierr)
720: {
721: CHKFORTRANNULLINTEGER(diag);
722: *MatSeqBDiagSetPreallocation(*mat,*nd,*bs,diag,PETSC_NULL);
723: }
724: /* Fortran ignores diagv */
725: void PETSC_STDCALL matmpibdiagsetpreallocation_(Mat *mat,int *nd,int *bs,int *diag,PetscScalar **diagv,int *ierr)
726: {
727: CHKFORTRANNULLINTEGER(diag);
728: *MatMPIBDiagSetPreallocation(*mat,*nd,*bs,diag,PETSC_NULL);
729: }
731: void PETSC_STDCALL matseqdensesetpreallocation_(Mat *mat,PetscScalar *data,int *ierr)
732: {
733: CHKFORTRANNULLSCALAR(data);
734: *MatSeqDenseSetPreallocation(*mat,data);
735: }
736: void PETSC_STDCALL matmpidensesetpreallocation_(Mat *mat,PetscScalar *data,int *ierr)
737: {
738: CHKFORTRANNULLSCALAR(data);
739: *MatMPIDenseSetPreallocation(*mat,data);
740: }
741: void PETSC_STDCALL matmpiadjsetpreallocation_(Mat *mat,int *i,int *j,int *values, int *ierr)
742: {
743: CHKFORTRANNULLINTEGER(values);
744: *MatMPIAdjSetPreallocation(*mat,i,j,values);
745: }
747: #if defined(PETSC_HAVE_BLOCKSOLVE) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
748: void PETSC_STDCALL matmpirowbssetpreallocation_(Mat *mat,int *nz,int *nnz,int *ierr)
749: {
750: CHKFORTRANNULLINTEGER(nnz);
751: *MatMPIRowbsSetPreallocation(*mat,*nz,nnz);
752: }
753: #endif
755: #if defined(PETSC_HAVE_PARTY)
756: void PETSC_STDCALL matpartitioningpartysetglobal_(MatPartitioning *part,CHAR method PETSC_MIXED_LEN(len),
757: int *ierr PETSC_END_LEN(len))
758: {
759: char *t;
760: FIXCHAR(method,len,t);
761: *MatPartitioningPartySetGlobal(*part,t);
762: FREECHAR(method,t);
763: }
765: void PETSC_STDCALL matpartitioningpartysetlocal_(MatPartitioning *part,CHAR method PETSC_MIXED_LEN(len),
766: int *ierr PETSC_END_LEN(len))
767: {
768: char *t;
769: FIXCHAR(method,len,t);
770: *MatPartitioningPartySetLocal(*part,t);
771: FREECHAR(method,t);
772: }
773: #endif
775: #if defined(PETSC_HAVE_SCOTCH)
776: void PETSC_STDCALL matpartitioningscotchsetstrategy_(MatPartitioning *part,CHAR strategy PETSC_MIXED_LEN(len),
777: int *ierr PETSC_END_LEN(len))
778: {
779: char *t;
780: FIXCHAR(strategy,len,t);
781: *MatPartitioningScotchSetStrategy(*part,t);
782: FREECHAR(strategy,t);
783: }
785: void PETSC_STDCALL matpartitioningscotchsetarch_(MatPartitioning *part,CHAR filename PETSC_MIXED_LEN(len),
786: int *ierr PETSC_END_LEN(len))
787: {
788: char *t;
789: FIXCHAR(filename,len,t);
790: *MatPartitioningScotchSetArch(*part,t);
791: FREECHAR(filename,t);
792: }
793: #endif
795: #if defined(PETSC_HAVE_SCOTCH)
796: void PETSC_STDCALL matpartitioningscotchsethostlist_(MatPartitioning *part,CHAR filename PETSC_MIXED_LEN(len),
797: int *ierr PETSC_END_LEN(len))
798: {
799: char *t;
800: FIXCHAR(filename,len,t);
801: *MatPartitioningScotchSetHostList(*part,t);
802: FREECHAR(filename,t);
803: }
804: #endif
806: EXTERN_C_END