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 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 dasetblockfills_             dasetblockfills
 40: #define dagetlocalinfo_              dagetlocalinfo
 41: #define dagetlocalvector_            dagetlocalvector
 42: #define dagetglobalvector_           dagetglobalvector
 43: #define darestorelocalvector_        darestorelocalvector
 44: #define dagetinterpolation_          dagetinterpolation
 45: #define daload_                      daload
 46: #define dacreateglobalvector_        dacreateglobalvector
 47: #define dacreatenaturalvector_       dacreatenaturalvector
 48: #define dacreatelocalvector_         dacreatelocalvector
 49: #define daview_                      daview
 50: #define dacreate1d_                  dacreate1d
 51: #define dacreate3d_                  dacreate3d
 52: #define dacreate2d_                  dacreate2d
 53: #define dadestroy_                   dadestroy
 54: #define dagetscatter_                dagetscatter
 55: #define dagetglobalindices_          dagetglobalindices
 56: #define dagetinfo_                   dagetinfo
 57: #define dagetcoloring_               dagetcoloring
 58: #define dagetmatrix_                 dagetmatrix
 59: #define dagetislocaltoglobalmapping_ dagetislocaltoglobalmapping
 60: #define dagetislocaltoglobalmappingblck_ dagetislocaltoglobalmappingblck
 61: #define dasetfieldname_              dasetfieldname
 62: #define dagetfieldname_              dagetfieldname
 63: #define darefine_                    darefine
 64: #define dagetao_                     dagetao
 65: #define dasetlocalfunction_          dasetlocalfunction
 66: #define dasetlocaladiforfunction_       dasetlocaladiforfunction
 67: #define dasetlocaladiformffunction_       dasetlocaladiformffunction
 68: #define dasetlocaljacobian_          dasetlocaljacobian
 69: #endif




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

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

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

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

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

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

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

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

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


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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

292: void PETSC_STDCALL dagetscatter_(DA *da,VecScatter *ltog,VecScatter *gtol,VecScatter *ltol,PetscErrorCode *ierr)
293: {
294:   CHKFORTRANNULLOBJECT(ltog);
295:   CHKFORTRANNULLOBJECT(gtol);
296:   CHKFORTRANNULLOBJECT(ltol);
297:   *DAGetScatter(*da,ltog,gtol,ltol);
298: }

300: void PETSC_STDCALL dadestroy_(DA *da,PetscErrorCode *ierr)
301: {
302:   *DADestroy(*da);
303: }

305: void PETSC_STDCALL dacreate2d_(MPI_Comm *comm,DAPeriodicType *wrap,DAStencilType
306:                   *stencil_type,PetscInt *M,PetscInt *N,PetscInt *m,PetscInt *n,PetscInt *w,
307:                   PetscInt *s,PetscInt *lx,PetscInt *ly,DA *inra,PetscErrorCode *ierr)
308: {
309:   CHKFORTRANNULLINTEGER(lx);
310:   CHKFORTRANNULLINTEGER(ly);
311:   *DACreate2d((MPI_Comm)PetscToPointerComm(*comm),*wrap,
312:                        *stencil_type,*M,*N,*m,*n,*w,*s,lx,ly,inra);
313: }

315: void PETSC_STDCALL dacreate1d_(MPI_Comm *comm,DAPeriodicType *wrap,PetscInt *M,PetscInt *w,PetscInt *s,
316:                  PetscInt *lc,DA *inra,PetscErrorCode *ierr)
317: {
318:   CHKFORTRANNULLINTEGER(lc);
319:   *DACreate1d((MPI_Comm)PetscToPointerComm(*comm),*wrap,*M,*w,*s,lc,inra);
320: }

322: void PETSC_STDCALL dacreate3d_(MPI_Comm *comm,DAPeriodicType *wrap,DAStencilType 
323:                  *stencil_type,PetscInt *M,PetscInt *N,PetscInt *P,PetscInt *m,PetscInt *n,PetscInt *p,
324:                  PetscInt *w,PetscInt *s,PetscInt *lx,PetscInt *ly,PetscInt *lz,DA *inra,PetscErrorCode *ierr)
325: {
326:   CHKFORTRANNULLINTEGER(lx);
327:   CHKFORTRANNULLINTEGER(ly);
328:   CHKFORTRANNULLINTEGER(lz);
329:   *DACreate3d((MPI_Comm)PetscToPointerComm(*comm),*wrap,*stencil_type,
330:                         *M,*N,*P,*m,*n,*p,*w,*s,lx,ly,lz,inra);
331: }

333: void PETSC_STDCALL dagetinfo_(DA *da,PetscInt *dim,PetscInt *M,PetscInt *N,PetscInt *P,PetscInt *m,PetscInt *n,PetscInt *p,PetscInt *w,PetscInt *s,
334:                 DAPeriodicType *wrap,DAStencilType *st,PetscErrorCode *ierr)
335: {
336:   CHKFORTRANNULLINTEGER(dim);
337:   CHKFORTRANNULLINTEGER(M);
338:   CHKFORTRANNULLINTEGER(N);
339:   CHKFORTRANNULLINTEGER(P);
340:   CHKFORTRANNULLINTEGER(m);
341:   CHKFORTRANNULLINTEGER(n);
342:   CHKFORTRANNULLINTEGER(p);
343:   CHKFORTRANNULLINTEGER(w);
344:   CHKFORTRANNULLINTEGER(s);
345:   CHKFORTRANNULLINTEGER(wrap);
346:   CHKFORTRANNULLINTEGER(st);
347:   *DAGetInfo(*da,dim,M,N,P,m,n,p,w,s,wrap,st);
348: }