Actual source code: zda.c
1: /*$Id: zda.c,v 1.49 2001/08/06 21:19:11 bsmith Exp $*/
3: #include src/dm/da/daimpl.h
4: #include src/fortran/custom/zpetsc.h
5: #include petscmat.h
6: #include petscda.h
8: #ifdef PETSC_HAVE_FORTRAN_CAPS
9: #define dasetblockfills_ DASETBLOCKFILLS
10: #define dasetlocalfunction_ DASETLOCALFUNCTION
11: #define dasetLocaladiforfunction_ DASETLOCALADIFORFUNCTION
12: #define dasetlocaladiformffunction_ DASETLOCALADIFORMFFUNCTION
13: #define dasetlocaljacobian_ DASETLOCALJACOBIAN
14: #define dagetlocalinfo_ DAGETLOCALINFO
15: #define dagetinterpolation_ DAGETINTERPOLATION
16: #define dacreate1d_ DACREATE1D
17: #define dacreate3d_ DACREATE3D
18: #define dacreate2d_ DACREATE2D
19: #define dadestroy_ DADESTROY
20: #define dacreateglobalvector_ DACREATEGLOBALVECTOR
21: #define dacreatenaturalvector_ DACREATENATURALVECTOR
22: #define dacreatelocalvector_ DACREATELOCALVECTOR
23: #define dagetlocalvector_ DAGETLOCALVECTOR
24: #define dagetglobalvector_ DAGETGLOBALVECTOR
25: #define darestorelocalvector_ DARESTORELOCALVECTOR
26: #define dagetscatter_ DAGETSCATTER
27: #define dagetglobalindices_ DAGETGLOBALINDICES
28: #define daview_ DAVIEW
29: #define dagetinfo_ DAGETINFO
30: #define dagetcoloring_ DAGETCOLORING
31: #define dagetmatrix_ DAGETMATRIX
32: #define dagetislocaltoglobalmapping_ DAGETISLOCALTOGLOBALMAPPING
33: #define dagetislocaltoglobalmappingblck_ DAGETISLOCALTOGLOBALMAPPINGBLCK
34: #define daload_ DALOAD
35: #define dasetfieldname_ DASETFIELDNAME
36: #define dagetfieldname_ DAGETFIELDNAME
37: #define darefine_ DAREFINE
38: #define dagetao_ DAGETAO
39: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
40: #define dasetblockfills_ dasetblockfills
41: #define dagetlocalinfo_ dagetlocalinfo
42: #define dagetlocalvector_ dagetlocalvector
43: #define dagetglobalvector_ dagetglobalvector
44: #define darestorelocalvector_ darestorelocalvector
45: #define dagetinterpolation_ dagetinterpolation
46: #define daload_ daload
47: #define dacreateglobalvector_ dacreateglobalvector
48: #define dacreatenaturalvector_ dacreatenaturalvector
49: #define dacreatelocalvector_ dacreatelocalvector
50: #define daview_ daview
51: #define dacreate1d_ dacreate1d
52: #define dacreate3d_ dacreate3d
53: #define dacreate2d_ dacreate2d
54: #define dadestroy_ dadestroy
55: #define dagetscatter_ dagetscatter
56: #define dagetglobalindices_ dagetglobalindices
57: #define dagetinfo_ dagetinfo
58: #define dagetcoloring_ dagetcoloring
59: #define dagetmatrix_ dagetmatrix
60: #define dagetislocaltoglobalmapping_ dagetislocaltoglobalmapping
61: #define dagetislocaltoglobalmappingblck_ dagetislocaltoglobalmappingblck
62: #define dasetfieldname_ dasetfieldname
63: #define dagetfieldname_ dagetfieldname
64: #define darefine_ darefine
65: #define dagetao_ dagetao
66: #define dasetlocalfunction_ dasetlocalfunction
67: #define dasetlocaladiforfunction_ dasetlocaladiforfunction
68: #define dasetlocaladiformffunction_ dasetlocaladiformffunction
69: #define dasetlocaljacobian_ dasetlocaljacobian
70: #endif
74: EXTERN_C_BEGIN
76: void PETSC_STDCALL dasetblockfills_(DA *da,int *dfill,int *ofill,int *ierr)
77: {
78: *DASetBlockFills(*da,dfill,ofill);
79: }
81: static void (PETSC_STDCALL *j1d)(DALocalInfo*,void*,void*,void*,int*);
82: static int ourlj1d(DALocalInfo *info,PetscScalar *in,Mat m,void *ptr)
83: {
84: int 0;
85: (*j1d)(info,&in[info->dof*info->gxs],&m,ptr,&ierr);
86: return 0;
87: }
89: static void (PETSC_STDCALL *j2d)(DALocalInfo*,void*,void*,void*,int*);
90: static int ourlj2d(DALocalInfo *info,PetscScalar **in,Mat m,void *ptr)
91: {
92: int 0;
93: (*j2d)(info,&in[info->gys][info->dof*info->gxs],&m,ptr,&ierr);
94: return 0;
95: }
97: static void (PETSC_STDCALL *j3d)(DALocalInfo*,void*,void*,void*,int*);
98: static int ourlj3d(DALocalInfo *info,PetscScalar ***in,Mat m,void *ptr)
99: {
100: int 0;
101: (*j3d)(info,&in[info->gzs][info->gys][info->dof*info->gxs],&m,ptr,&ierr);
102: return 0;
103: }
105: static void (PETSC_STDCALL *f1d)(DALocalInfo*,void*,void*,void*,int*);
106: static int ourlf1d(DALocalInfo *info,PetscScalar *in,PetscScalar *out,void *ptr)
107: {
108: int 0;
109: (*f1d)(info,&in[info->dof*info->gxs],&out[info->dof*info->xs],ptr,&ierr);
110: return 0;
111: }
113: static void (PETSC_STDCALL *f2d)(DALocalInfo*,void*,void*,void*,int*);
114: static int ourlf2d(DALocalInfo *info,PetscScalar **in,PetscScalar **out,void *ptr)
115: {
116: int 0;
117: (*f2d)(info,&in[info->gys][info->dof*info->gxs],&out[info->ys][info->dof*info->xs],ptr,&ierr);
118: return 0;
119: }
121: static void (PETSC_STDCALL *f3d)(DALocalInfo*,void*,void*,void*,int*);
122: static int ourlf3d(DALocalInfo *info,PetscScalar ***in,PetscScalar ***out,void *ptr)
123: {
124: int 0;
125: (*f3d)(info,&in[info->gzs][info->gys][info->dof*info->gxs],&out[info->zs][info->ys][info->dof*info->xs],ptr,&ierr);
126: return 0;
127: }
129: void PETSC_STDCALL dasetlocalfunction_(DA *da,void (PETSC_STDCALL *func)(DALocalInfo*,void*,void*,void*,int*),int *ierr)
130: {
131: int dim;
133: *DAGetInfo(*da,&dim,0,0,0,0,0,0,0,0,0,0); if (*ierr) return;
134: if (dim == 2) {
135: f2d = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,int*))func;
136: *DASetLocalFunction(*da,(DALocalFunction1)ourlf2d);
137: } else if (dim == 3) {
138: f3d = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,int*))func;
139: *DASetLocalFunction(*da,(DALocalFunction1)ourlf3d);
140: } else if (dim == 1) {
141: f1d = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,int*))func;
142: *DASetLocalFunction(*da,(DALocalFunction1)ourlf1d);
143: } else *1;
144: }
147: void PETSC_STDCALL dasetlocaladiforfunction_(DA *da,
148: void (PETSC_STDCALL *jfunc)(int*,DALocalInfo*,void*,void*,int*,void*,void*,int*,void*,int*),int *ierr)
149: {
150: (*da)->adifor_lf = (DALocalFunction1)jfunc;
151: }
153: void PETSC_STDCALL dasetlocaladiformffunction_(DA *da,
154: void (PETSC_STDCALL *jfunc)(DALocalInfo*,void*,void*,void*,void*,void*,int*),int *ierr)
155: {
156: (*da)->adiformf_lf = (DALocalFunction1)jfunc;
157: }
159: void PETSC_STDCALL dasetlocaljacobian_(DA *da,void (PETSC_STDCALL *jac)(DALocalInfo*,void*,void*,void*,int*),int *ierr)
160: {
161: int dim;
163: *DAGetInfo(*da,&dim,0,0,0,0,0,0,0,0,0,0); if (*ierr) return;
164: if (dim == 2) {
165: j2d = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,int*))jac;
166: *DASetLocalJacobian(*da,(DALocalFunction1)ourlj2d);
167: } else if (dim == 3) {
168: j3d = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,int*))jac;
169: *DASetLocalJacobian(*da,(DALocalFunction1)ourlj3d);
170: } else if (dim == 1) {
171: j1d = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,int*))jac;
172: *DASetLocalJacobian(*da,(DALocalFunction1)ourlj1d);
173: } else *1;
174: }
176: void PETSC_STDCALL dagetlocalinfo_(DA *da,DALocalInfo *ao,int *ierr)
177: {
178: *DAGetLocalInfo(*da,ao);
179: }
181: void PETSC_STDCALL dagetao_(DA *da,AO *ao,int *ierr)
182: {
183: *DAGetAO(*da,ao);
184: }
186: void PETSC_STDCALL darefine_(DA *da,MPI_Comm *comm,DA *daref, int *ierr)
187: {
188: *DARefine(*da,(MPI_Comm)PetscToPointerComm(*comm),daref);
189: }
191: void PETSC_STDCALL dagetinterpolation_(DA *dac,DA *daf,Mat *A,Vec *scale,int *ierr)
192: {
193: CHKFORTRANNULLOBJECT(scale);
194: *DAGetInterpolation(*dac,*daf,A,scale);
195: }
197: void PETSC_STDCALL dasetfieldname_(DA *da,int *nf,CHAR name PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
198: {
199: char *t;
200: FIXCHAR(name,len,t);
201: *DASetFieldName(*da,*nf,t);
202: FREECHAR(name,t);
203: }
204: void PETSC_STDCALL dagetfieldname(DA *da,int *nf,CHAR name PETSC_MIXED_LEN(len),int *ierr PETSC_END_LEN(len))
205: {
206: char *tname;
208: *DAGetFieldName(*da,*nf,&tname);
209: #if defined(PETSC_USES_CPTOFCD)
210: {
211: char *t = _fcdtocp(name); int len1 = _fcdlen(name);
212: *PetscStrncpy(t,tname,len1);
213: }
214: #else
215: *PetscStrncpy(name,tname,len);
216: #endif
217: }
219: void PETSC_STDCALL daload_(PetscViewer *viewer,int *M,int *N,int *P,DA *da,int *ierr)
220: {
221: PetscViewer v;
222: PetscPatchDefaultViewers_Fortran(viewer,v);
223: *DALoad(v,*M,*N,*P,da);
224: }
226: void PETSC_STDCALL dagetislocaltoglobalmapping_(DA *da,ISLocalToGlobalMapping *map,int *ierr)
227: {
228: *DAGetISLocalToGlobalMapping(*da,map);
229: }
231: void PETSC_STDCALL dagetislocaltoglobalmappingblck_(DA *da,ISLocalToGlobalMapping *map,int *ierr)
232: {
233: *DAGetISLocalToGlobalMappingBlck(*da,map);
234: }
236: void PETSC_STDCALL dagetcoloring_(DA *da,ISColoringType *ctype,ISColoring *coloring,int *ierr PETSC_END_LEN(len))
237: {
238: *DAGetColoring(*da,*ctype,coloring);
239: }
241: void PETSC_STDCALL dagetmatrix_(DA *da,CHAR mat_type PETSC_MIXED_LEN(len),Mat *J,int *ierr PETSC_END_LEN(len))
242: {
243: char *t;
244: FIXCHAR(mat_type,len,t);
245: *DAGetMatrix(*da,t,J);
246: FREECHAR(mat_type,t);
247: }
249: void PETSC_STDCALL daview_(DA *da,PetscViewer *vin,int *ierr)
250: {
251: PetscViewer v;
252: PetscPatchDefaultViewers_Fortran(vin,v);
253: *DAView(*da,v);
254: }
256: void PETSC_STDCALL dagetglobalindices_(DA *da,int *n,int *indices,long *ia,int *ierr)
257: {
258: int *idx;
259: *DAGetGlobalIndices(*da,n,&idx);
260: *ia = PetscIntAddressToFortran(indices,idx);
261: }
263: void PETSC_STDCALL dacreateglobalvector_(DA *da,Vec* g,int *ierr)
264: {
265: *DACreateGlobalVector(*da,g);
266: }
268: void PETSC_STDCALL dacreatenaturalvector_(DA *da,Vec* g,int *ierr)
269: {
270: *DACreateNaturalVector(*da,g);
271: }
273: void PETSC_STDCALL dacreatelocalvector_(DA *da,Vec* l,int *ierr)
274: {
275: *DACreateLocalVector(*da,l);
276: }
278: void PETSC_STDCALL dagetlocalvector_(DA *da,Vec* l,int *ierr)
279: {
280: *DAGetLocalVector(*da,l);
281: }
283: void PETSC_STDCALL dagetglobalvector_(DA *da,Vec* l,int *ierr)
284: {
285: *DAGetGlobalVector(*da,l);
286: }
288: void PETSC_STDCALL darestorelocalvector_(DA *da,Vec* l,int *ierr)
289: {
290: *DARestoreLocalVector(*da,l);
291: }
293: void PETSC_STDCALL dagetscatter_(DA *da,VecScatter *ltog,VecScatter *gtol,VecScatter *ltol,int *ierr)
294: {
295: CHKFORTRANNULLINTEGER(ltog);
296: CHKFORTRANNULLINTEGER(gtol);
297: CHKFORTRANNULLINTEGER(ltol);
298: *DAGetScatter(*da,ltog,gtol,ltol);
299: }
301: void PETSC_STDCALL dadestroy_(DA *da,int *ierr)
302: {
303: *DADestroy(*da);
304: }
306: void PETSC_STDCALL dacreate2d_(MPI_Comm *comm,DAPeriodicType *wrap,DAStencilType
307: *stencil_type,int *M,int *N,int *m,int *n,int *w,
308: int *s,int *lx,int *ly,DA *inra,int *ierr)
309: {
310: CHKFORTRANNULLINTEGER(lx);
311: CHKFORTRANNULLINTEGER(ly);
312: *DACreate2d((MPI_Comm)PetscToPointerComm(*comm),*wrap,
313: *stencil_type,*M,*N,*m,*n,*w,*s,lx,ly,inra);
314: }
316: void PETSC_STDCALL dacreate1d_(MPI_Comm *comm,DAPeriodicType *wrap,int *M,int *w,int *s,
317: int *lc,DA *inra,int *ierr)
318: {
319: CHKFORTRANNULLINTEGER(lc);
320: *DACreate1d((MPI_Comm)PetscToPointerComm(*comm),*wrap,*M,*w,*s,lc,inra);
321: }
323: void PETSC_STDCALL dacreate3d_(MPI_Comm *comm,DAPeriodicType *wrap,DAStencilType
324: *stencil_type,int *M,int *N,int *P,int *m,int *n,int *p,
325: int *w,int *s,int *lx,int *ly,int *lz,DA *inra,int *ierr)
326: {
327: CHKFORTRANNULLINTEGER(lx);
328: CHKFORTRANNULLINTEGER(ly);
329: CHKFORTRANNULLINTEGER(lz);
330: *DACreate3d((MPI_Comm)PetscToPointerComm(*comm),*wrap,*stencil_type,
331: *M,*N,*P,*m,*n,*p,*w,*s,lx,ly,lz,inra);
332: }
334: void PETSC_STDCALL dagetinfo_(DA *da,int *dim,int *M,int *N,int *P,int *m,int *n,int *p,int *w,int *s,
335: DAPeriodicType *wrap,DAStencilType *st,int *ierr)
336: {
337: CHKFORTRANNULLINTEGER(dim);
338: CHKFORTRANNULLINTEGER(M);
339: CHKFORTRANNULLINTEGER(N);
340: CHKFORTRANNULLINTEGER(P);
341: CHKFORTRANNULLINTEGER(m);
342: CHKFORTRANNULLINTEGER(n);
343: CHKFORTRANNULLINTEGER(p);
344: CHKFORTRANNULLINTEGER(w);
345: CHKFORTRANNULLINTEGER(s);
346: CHKFORTRANNULLINTEGER(wrap);
347: CHKFORTRANNULLINTEGER(st);
348: *DAGetInfo(*da,dim,M,N,P,m,n,p,w,s,wrap,st);
349: }
351: EXTERN_C_END