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