Actual source code: zda.c

 2:  #include src/dm/da/daimpl.h
 3:  #include src/fortran/custom/zpetsc.h
 4:  #include petscmat.h
 5:  #include petscda.h

  7: #ifdef PETSC_HAVE_FORTRAN_CAPS
  8: #define dasetblockfills_             DASETBLOCKFILLS
  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 darestoreglobalvector_       DARESTOREGLOBALVECTOR
 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 darestoreglobalvector_       darestoreglobalvector
 46: #define dagetinterpolation_          dagetinterpolation
 47: #define daload_                      daload
 48: #define dacreateglobalvector_        dacreateglobalvector
 49: #define dacreatenaturalvector_       dacreatenaturalvector
 50: #define dacreatelocalvector_         dacreatelocalvector
 51: #define daview_                      daview
 52: #define dacreate1d_                  dacreate1d
 53: #define dacreate3d_                  dacreate3d
 54: #define dacreate2d_                  dacreate2d
 55: #define dadestroy_                   dadestroy
 56: #define dagetscatter_                dagetscatter
 57: #define dagetglobalindices_          dagetglobalindices
 58: #define dagetinfo_                   dagetinfo
 59: #define dagetcoloring_               dagetcoloring
 60: #define dagetmatrix_                 dagetmatrix
 61: #define dagetislocaltoglobalmapping_ dagetislocaltoglobalmapping
 62: #define dagetislocaltoglobalmappingblck_ dagetislocaltoglobalmappingblck
 63: #define dasetfieldname_              dasetfieldname
 64: #define dagetfieldname_              dagetfieldname
 65: #define darefine_                    darefine
 66: #define dagetao_                     dagetao
 67: #define dasetlocalfunction_          dasetlocalfunction
 68: #define dasetlocaladiforfunction_       dasetlocaladiforfunction
 69: #define dasetlocaladiformffunction_       dasetlocaladiformffunction
 70: #define dasetlocaljacobian_          dasetlocaljacobian
 71: #endif




 77: void PETSC_STDCALL dasetblockfills_(DA *da,PetscInt *dfill,PetscInt *ofill,PetscErrorCode *ierr)
 78: {
 79:   *DASetBlockFills(*da,dfill,ofill);
 80: }

 82: static void (PETSC_STDCALL *j1d)(DALocalInfo*,void*,void*,void*,PetscErrorCode*);
 83: static PetscErrorCode ourlj1d(DALocalInfo *info,PetscScalar *in,Mat m,void *ptr)
 84: {
 85:   PetscErrorCode 0;
 86:   (*j1d)(info,&in[info->dof*info->gxs],&m,ptr,&ierr);
 87:   return 0;
 88: }

 90: static void (PETSC_STDCALL *j2d)(DALocalInfo*,void*,void*,void*,PetscErrorCode*);
 91: static PetscErrorCode ourlj2d(DALocalInfo *info,PetscScalar **in,Mat m,void *ptr)
 92: {
 93:   PetscErrorCode 0;
 94:   (*j2d)(info,&in[info->gys][info->dof*info->gxs],&m,ptr,&ierr);
 95:   return 0;
 96: }

 98: static void (PETSC_STDCALL *j3d)(DALocalInfo*,void*,void*,void*,PetscErrorCode*);
 99: static PetscErrorCode ourlj3d(DALocalInfo *info,PetscScalar ***in,Mat m,void *ptr)
100: {
101:   PetscErrorCode 0;
102:   (*j3d)(info,&in[info->gzs][info->gys][info->dof*info->gxs],&m,ptr,&ierr);
103:   return 0;
104: }

106: static void (PETSC_STDCALL *f1d)(DALocalInfo*,void*,void*,void*,PetscErrorCode*);
107: static PetscErrorCode ourlf1d(DALocalInfo *info,PetscScalar *in,PetscScalar *out,void *ptr)
108: {
109:   PetscErrorCode 0;
110:   (*f1d)(info,&in[info->dof*info->gxs],&out[info->dof*info->xs],ptr,&ierr);
111:   return 0;
112: }

114: static void (PETSC_STDCALL *f2d)(DALocalInfo*,void*,void*,void*,PetscErrorCode*);
115: static PetscErrorCode ourlf2d(DALocalInfo *info,PetscScalar **in,PetscScalar **out,void *ptr)
116: {
117:   PetscErrorCode 0;
118:   (*f2d)(info,&in[info->gys][info->dof*info->gxs],&out[info->ys][info->dof*info->xs],ptr,&ierr);
119:   return 0;
120: }

122: static void (PETSC_STDCALL *f3d)(DALocalInfo*,void*,void*,void*,PetscErrorCode*);
123: static PetscErrorCode ourlf3d(DALocalInfo *info,PetscScalar ***in,PetscScalar ***out,void *ptr)
124: {
125:   PetscErrorCode 0;
126:   (*f3d)(info,&in[info->gzs][info->gys][info->dof*info->gxs],&out[info->zs][info->ys][info->dof*info->xs],ptr,&ierr);
127:   return 0;
128: }

130: void PETSC_STDCALL dasetlocalfunction_(DA *da,void (PETSC_STDCALL *func)(DALocalInfo*,void*,void*,void*,PetscErrorCode*),PetscErrorCode *ierr)
131: {
132:   PetscInt dim;

134:   *DAGetInfo(*da,&dim,0,0,0,0,0,0,0,0,0,0); if (*ierr) return;
135:   if (dim == 2) {
136:      f2d    = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,PetscErrorCode*))func;
137:     *DASetLocalFunction(*da,(DALocalFunction1)ourlf2d);
138:   } else if (dim == 3) {
139:      f3d    = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,PetscErrorCode*))func;
140:     *DASetLocalFunction(*da,(DALocalFunction1)ourlf3d);
141:   } else if (dim == 1) {
142:      f1d    = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,PetscErrorCode*))func;
143:     *DASetLocalFunction(*da,(DALocalFunction1)ourlf1d);
144:   } else *1;
145: }


148: void PETSC_STDCALL dasetlocaladiforfunction_(DA *da,
149: void (PETSC_STDCALL *jfunc)(PetscInt*,DALocalInfo*,void*,void*,PetscInt*,void*,void*,PetscInt*,void*,PetscErrorCode*),PetscErrorCode *ierr)
150: {
151:   (*da)->adifor_lf = (DALocalFunction1)jfunc;
152: }

154: void PETSC_STDCALL dasetlocaladiformffunction_(DA *da,
155: void (PETSC_STDCALL *jfunc)(DALocalInfo*,void*,void*,void*,void*,void*,PetscErrorCode*),PetscErrorCode *ierr)
156: {
157:   (*da)->adiformf_lf = (DALocalFunction1)jfunc;
158: }

160: void PETSC_STDCALL dasetlocaljacobian_(DA *da,void (PETSC_STDCALL *jac)(DALocalInfo*,void*,void*,void*,PetscErrorCode*),PetscErrorCode *ierr)
161: {
162:   PetscInt dim;

164:   *DAGetInfo(*da,&dim,0,0,0,0,0,0,0,0,0,0); if (*ierr) return;
165:   if (dim == 2) {
166:      j2d    = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,PetscErrorCode*))jac;
167:     *DASetLocalJacobian(*da,(DALocalFunction1)ourlj2d);
168:   } else if (dim == 3) {
169:      j3d    = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,PetscErrorCode*))jac;
170:     *DASetLocalJacobian(*da,(DALocalFunction1)ourlj3d);
171:   } else if (dim == 1) {
172:      j1d    = (void (PETSC_STDCALL *)(DALocalInfo*,void*,void*,void*,PetscErrorCode*))jac;
173:     *DASetLocalJacobian(*da,(DALocalFunction1)ourlj1d);
174:   } else *1;
175: }

177: void PETSC_STDCALL dagetlocalinfo_(DA *da,DALocalInfo *ao,PetscErrorCode *ierr)
178: {
179:   *DAGetLocalInfo(*da,ao);
180: }

182: void PETSC_STDCALL dagetao_(DA *da,AO *ao,PetscErrorCode *ierr)
183: {
184:   *DAGetAO(*da,ao);
185: }

187: void PETSC_STDCALL darefine_(DA *da,MPI_Comm *comm,DA *daref, PetscErrorCode *ierr)
188: {
189:   *DARefine(*da,(MPI_Comm)PetscToPointerComm(*comm),daref);
190: }

192: void PETSC_STDCALL dagetinterpolation_(DA *dac,DA *daf,Mat *A,Vec *scale,PetscErrorCode *ierr)
193: {
194:   CHKFORTRANNULLOBJECT(scale);
195:   *DAGetInterpolation(*dac,*daf,A,scale);
196: }

198: void PETSC_STDCALL dasetfieldname_(DA *da,PetscInt *nf,CHAR name PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len))
199: {
200:   char *t;
201:   FIXCHAR(name,len,t);
202:   *DASetFieldName(*da,*nf,t);
203:   FREECHAR(name,t);
204: }
205: void PETSC_STDCALL dagetfieldname(DA *da,PetscInt *nf,CHAR name PETSC_MIXED_LEN(len),PetscErrorCode *ierr PETSC_END_LEN(len))
206: {
207:   char *tname;

209:   *DAGetFieldName(*da,*nf,&tname);
210: #if defined(PETSC_USES_CPTOFCD)
211:   {
212:     char *t = _fcdtocp(name); int len1 = _fcdlen(name);
213:     *PetscStrncpy(t,tname,len1);
214:   }
215: #else
216:   *PetscStrncpy(name,tname,len);
217: #endif
218: }

220: void PETSC_STDCALL daload_(PetscViewer *viewer,PetscInt *M,PetscInt *N,PetscInt *P,DA *da,PetscErrorCode *ierr)
221: {
222:   PetscViewer v;
223:   PetscPatchDefaultViewers_Fortran(viewer,v);
224:   *DALoad(v,*M,*N,*P,da);
225: }

227: void PETSC_STDCALL dagetislocaltoglobalmapping_(DA *da,ISLocalToGlobalMapping *map,PetscErrorCode *ierr)
228: {
229:   *DAGetISLocalToGlobalMapping(*da,map);
230: }

232: void PETSC_STDCALL dagetislocaltoglobalmappingblck_(DA *da,ISLocalToGlobalMapping *map,PetscErrorCode *ierr)
233: {
234:   *DAGetISLocalToGlobalMappingBlck(*da,map);
235: }

237: void PETSC_STDCALL dagetcoloring_(DA *da,ISColoringType *ctype,ISColoring *coloring,PetscErrorCode *ierr PETSC_END_LEN(len))
238: {
239:   *DAGetColoring(*da,*ctype,coloring);
240: }

242: void PETSC_STDCALL dagetmatrix_(DA *da,CHAR mat_type PETSC_MIXED_LEN(len),Mat *J,PetscErrorCode *ierr PETSC_END_LEN(len))
243: {
244:   char *t;
245:   FIXCHAR(mat_type,len,t);
246:   *DAGetMatrix(*da,t,J);
247:   FREECHAR(mat_type,t);
248: }

250: void PETSC_STDCALL daview_(DA *da,PetscViewer *vin,PetscErrorCode *ierr)
251: {
252:   PetscViewer v;
253:   PetscPatchDefaultViewers_Fortran(vin,v);
254:   *DAView(*da,v);
255: }

257: void PETSC_STDCALL dagetglobalindices_(DA *da,PetscInt *n,PetscInt *indices,size_t *ia,PetscErrorCode *ierr)
258: {
259:   PetscInt *idx;
260:   *DAGetGlobalIndices(*da,n,&idx);
261:   *ia   = PetscIntAddressToFortran(indices,idx);
262: }

264: void PETSC_STDCALL dacreateglobalvector_(DA *da,Vec* g,PetscErrorCode *ierr)
265: {
266:   *DACreateGlobalVector(*da,g);
267: }

269: void PETSC_STDCALL dacreatenaturalvector_(DA *da,Vec* g,PetscErrorCode *ierr)
270: {
271:   *DACreateNaturalVector(*da,g);
272: }

274: void PETSC_STDCALL dacreatelocalvector_(DA *da,Vec* l,PetscErrorCode *ierr)
275: {
276:   *DACreateLocalVector(*da,l);
277: }

279: void PETSC_STDCALL dagetlocalvector_(DA *da,Vec* l,PetscErrorCode *ierr)
280: {
281:   *DAGetLocalVector(*da,l);
282: }

284: void PETSC_STDCALL dagetglobalvector_(DA *da,Vec* l,PetscErrorCode *ierr)
285: {
286:   *DAGetGlobalVector(*da,l);
287: }

289: void PETSC_STDCALL darestorelocalvector_(DA *da,Vec* l,PetscErrorCode *ierr)
290: {
291:   *DARestoreLocalVector(*da,l);
292: }

294: void PETSC_STDCALL darestoreglobalvector_(DA *da,Vec* g,PetscErrorCode *ierr)
295: {
296:   *DARestoreGlobalVector(*da,g);
297: }

299: void PETSC_STDCALL dagetscatter_(DA *da,VecScatter *ltog,VecScatter *gtol,VecScatter *ltol,PetscErrorCode *ierr)
300: {
301:   CHKFORTRANNULLOBJECT(ltog);
302:   CHKFORTRANNULLOBJECT(gtol);
303:   CHKFORTRANNULLOBJECT(ltol);
304:   *DAGetScatter(*da,ltog,gtol,ltol);
305: }

307: void PETSC_STDCALL dadestroy_(DA *da,PetscErrorCode *ierr)
308: {
309:   *DADestroy(*da);
310: }

312: void PETSC_STDCALL dacreate2d_(MPI_Comm *comm,DAPeriodicType *wrap,DAStencilType
313:                   *stencil_type,PetscInt *M,PetscInt *N,PetscInt *m,PetscInt *n,PetscInt *w,
314:                   PetscInt *s,PetscInt *lx,PetscInt *ly,DA *inra,PetscErrorCode *ierr)
315: {
316:   CHKFORTRANNULLINTEGER(lx);
317:   CHKFORTRANNULLINTEGER(ly);
318:   *DACreate2d((MPI_Comm)PetscToPointerComm(*comm),*wrap,
319:                        *stencil_type,*M,*N,*m,*n,*w,*s,lx,ly,inra);
320: }

322: void PETSC_STDCALL dacreate1d_(MPI_Comm *comm,DAPeriodicType *wrap,PetscInt *M,PetscInt *w,PetscInt *s,
323:                  PetscInt *lc,DA *inra,PetscErrorCode *ierr)
324: {
325:   CHKFORTRANNULLINTEGER(lc);
326:   *DACreate1d((MPI_Comm)PetscToPointerComm(*comm),*wrap,*M,*w,*s,lc,inra);
327: }

329: void PETSC_STDCALL dacreate3d_(MPI_Comm *comm,DAPeriodicType *wrap,DAStencilType 
330:                  *stencil_type,PetscInt *M,PetscInt *N,PetscInt *P,PetscInt *m,PetscInt *n,PetscInt *p,
331:                  PetscInt *w,PetscInt *s,PetscInt *lx,PetscInt *ly,PetscInt *lz,DA *inra,PetscErrorCode *ierr)
332: {
333:   CHKFORTRANNULLINTEGER(lx);
334:   CHKFORTRANNULLINTEGER(ly);
335:   CHKFORTRANNULLINTEGER(lz);
336:   *DACreate3d((MPI_Comm)PetscToPointerComm(*comm),*wrap,*stencil_type,
337:                         *M,*N,*P,*m,*n,*p,*w,*s,lx,ly,lz,inra);
338: }

340: void PETSC_STDCALL dagetinfo_(DA *da,PetscInt *dim,PetscInt *M,PetscInt *N,PetscInt *P,PetscInt *m,PetscInt *n,PetscInt *p,PetscInt *w,PetscInt *s,
341:                 DAPeriodicType *wrap,DAStencilType *st,PetscErrorCode *ierr)
342: {
343:   CHKFORTRANNULLINTEGER(dim);
344:   CHKFORTRANNULLINTEGER(M);
345:   CHKFORTRANNULLINTEGER(N);
346:   CHKFORTRANNULLINTEGER(P);
347:   CHKFORTRANNULLINTEGER(m);
348:   CHKFORTRANNULLINTEGER(n);
349:   CHKFORTRANNULLINTEGER(p);
350:   CHKFORTRANNULLINTEGER(w);
351:   CHKFORTRANNULLINTEGER(s);
352:   CHKFORTRANNULLINTEGER(wrap);
353:   CHKFORTRANNULLINTEGER(st);
354:   *DAGetInfo(*da,dim,M,N,P,m,n,p,w,s,wrap,st);
355: }