Actual source code: zda.c

  1: /*$Id: zda.c,v 1.42 2001/04/07 15:21:07 bsmith Exp $*/

  3: #include "src/fortran/custom/zpetsc.h"
  4: #include "petscmat.h"
  5: #include "petscda.h"

  7: #ifdef PETSC_HAVE_FORTRAN_CAPS
  8: #define dagetinterpolation_          DAGETINTERPOLATION
  9: #define dacreate1d_                  DACREATE1D
 10: #define dacreate3d_                  DACREATE3D
 11: #define dacreate2d_                  DACREATE2D
 12: #define dadestroy_                   DADESTROY
 13: #define dacreateglobalvector_        DACREATEGLOBALVECTOR
 14: #define dacreatelocalvector_         DACREATELOCALVECTOR
 15: #define dagetlocalvector_            DAGETLOCALVECTOR
 16: #define darestorelocalvector_        DARESTORELOCALVECTOR
 17: #define dagetscatter_                DAGETSCATTER
 18: #define dagetglobalindices_          DAGETGLOBALINDICES
 19: #define daview_                      DAVIEW
 20: #define dagetinfo_                   DAGETINFO
 21: #define dagetcoloring_               DAGETCOLORING
 22: #define dagetislocaltoglobalmapping_ DAGETISLOCALTOGLOBALMAPPING
 23: #define dagetislocaltoglobalmappingblck_ DAGETISLOCALTOGLOBALMAPPINGBLCK
 24: #define daload_                      DALOAD
 25: #define dasetfieldname_              DASETFIELDNAME
 26: #define dagetfieldname_              DAGETFIELDNAME
 27: #define darefine_                    DAREFINE
 28: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 29: #define dagetlocalvector_            dagetlocalvector
 30: #define darestorelocalvector_        darestorelocalvector
 31: #define dagetinterpolation_          dagetinterpolation
 32: #define daload_                      daload
 33: #define dacreateglobalvector_        dacreateglobalvector
 34: #define dacreatelocalvector_         dacreatelocalvector
 35: #define daview_                      daview
 36: #define dacreate1d_                  dacreate1d
 37: #define dacreate3d_                  dacreate3d
 38: #define dacreate2d_                  dacreate2d
 39: #define dadestroy_                   dadestroy
 40: #define dagetscatter_                dagetscatter
 41: #define dagetglobalindices_          dagetglobalindices
 42: #define dagetinfo_                   dagetinfo
 43: #define dagetcoloring_               dagetcoloring
 44: #define dagetislocaltoglobalmapping_ dagetislocaltoglobalmapping
 45: #define dagetislocaltoglobalmappingblck_ dagetislocaltoglobalmappingblck
 46: #define dasetfieldname_              dasetfieldname
 47: #define dagetfieldname_              dagetfieldname
 48: #define darefine_                    darefine
 49: #endif

 51: EXTERN_C_BEGIN

 53: void PETSC_STDCALL darefine_(DA *da,MPI_Comm *comm,DA *daref, int *ierr )
 54: {
 55:   *DARefine(*da,(MPI_Comm)PetscToPointerComm(*comm),daref);
 56: }

 58: void PETSC_STDCALL dagetinterpolation_(DA *dac,DA *daf,Mat *A,Vec *scale,int *ierr)
 59: {
 60:   if (FORTRANNULLOBJECT(scale)) scale = PETSC_NULL;
 61:   *DAGetInterpolation(*dac,*daf,A,scale);
 62: }

 64: void PETSC_STDCALL dasetfieldname_(DA *da,int *nf,CHAR name PETSC_MIXED_LEN(len),
 65:                                    int *ierr PETSC_END_LEN(len))
 66: {
 67:   char *t;
 68:   FIXCHAR(name,len,t);
 69:   *DASetFieldName(*da,*nf,t);
 70:   FREECHAR(name,t);
 71: }
 72: void PETSC_STDCALL dagetfieldname(DA *da,int *nf,CHAR name PETSC_MIXED_LEN(len),
 73:                                   int *ierr PETSC_END_LEN(len))
 74: {
 75:   char *tname;

 77:   *DAGetFieldName(*da,*nf,&tname);
 78: #if defined(PETSC_USES_CPTOFCD)
 79:   {
 80:     char *t = _fcdtocp(name); int len1 = _fcdlen(name);
 81:     *PetscStrncpy(t,tname,len1);
 82:   }
 83: #else
 84:   *PetscStrncpy(name,tname,len);
 85: #endif
 86: }

 88: void PETSC_STDCALL daload_(PetscViewer *viewer,int *M,int *N,int *P,DA *da,int *ierr)
 89: {
 90:   PetscViewer v;
 91:   PetscPatchDefaultViewers_Fortran(viewer,v);
 92:   *DALoad(v,*M,*N,*P,da);
 93: }

 95: void PETSC_STDCALL dagetislocaltoglobalmapping_(DA *da,ISLocalToGlobalMapping *map,int *ierr)
 96: {
 97:   *DAGetISLocalToGlobalMapping(*da,map);
 98: }

100: void PETSC_STDCALL dagetislocaltoglobalmappingblck_(DA *da,ISLocalToGlobalMapping *map,int *ierr)
101: {
102:   *DAGetISLocalToGlobalMappingBlck(*da,map);
103: }

105: void PETSC_STDCALL dagetcoloring_(DA *da,ISColoringType *ctype,MatType *mtype,ISColoring *coloring,Mat *J,int *ierr)
106: {
107:   if (FORTRANNULLOBJECT(coloring)) coloring = PETSC_NULL;
108:   if (FORTRANNULLOBJECT(J))        J        = PETSC_NULL;
109:   *DAGetColoring(*da,*ctype,*mtype,coloring,J);
110: }

112: void PETSC_STDCALL daview_(DA *da,PetscViewer *vin,int *ierr)
113: {
114:   PetscViewer v;
115:   PetscPatchDefaultViewers_Fortran(vin,v);
116:   *DAView(*da,v);
117: }

119: void PETSC_STDCALL dagetglobalindices_(DA *da,int *n,int *indices,long *ia,int *ierr)
120: {
121:   int *idx;
122:   *DAGetGlobalIndices(*da,n,&idx);
123:   *ia   = PetscIntAddressToFortran(indices,idx);
124: }

126: void PETSC_STDCALL dacreateglobalvector_(DA *da,Vec* g,int *ierr)
127: {
128:   *DACreateGlobalVector(*da,g);
129: }

131: void PETSC_STDCALL dacreatelocalvector_(DA *da,Vec* l,int *ierr)
132: {
133:   *DACreateLocalVector(*da,l);
134: }

136: void PETSC_STDCALL dagetlocalvector_(DA *da,Vec* l,int *ierr)
137: {
138:   *DAGetLocalVector(*da,l);
139: }

141: void PETSC_STDCALL darestorelocalvector_(DA *da,Vec* l,int *ierr)
142: {
143:   *DARestoreLocalVector(*da,l);
144: }

146: void PETSC_STDCALL dagetscatter_(DA *da,VecScatter *ltog,VecScatter *gtol,VecScatter *ltol,int *ierr)
147: {
148:   if (!FORTRANNULLINTEGER(ltog)) ltog = PETSC_NULL;
149:   if (!FORTRANNULLINTEGER(gtol)) gtol = PETSC_NULL;
150:   if (!FORTRANNULLINTEGER(ltol)) ltol = PETSC_NULL;
151:   *DAGetScatter(*da,ltog,gtol,ltol);
152: }

154: void PETSC_STDCALL dadestroy_(DA *da,int *ierr)
155: {
156:   *DADestroy(*da);
157: }

159: void PETSC_STDCALL dacreate2d_(MPI_Comm *comm,DAPeriodicType *wrap,DAStencilType
160:                   *stencil_type,int *M,int *N,int *m,int *n,int *w,
161:                   int *s,int *lx,int *ly,DA *inra,int *ierr)
162: {
163:   if (FORTRANNULLINTEGER(lx)) lx = PETSC_NULL;
164:   if (FORTRANNULLINTEGER(ly)) ly = PETSC_NULL;
165:   *DACreate2d((MPI_Comm)PetscToPointerComm(*comm),*wrap,
166:                        *stencil_type,*M,*N,*m,*n,*w,*s,lx,ly,inra);
167: }

169: void PETSC_STDCALL dacreate1d_(MPI_Comm *comm,DAPeriodicType *wrap,int *M,int *w,int *s,
170:                  int *lc,DA *inra,int *ierr)
171: {
172:   if (FORTRANNULLINTEGER(lc)) lc = PETSC_NULL;
173:   *DACreate1d((MPI_Comm)PetscToPointerComm(*comm),*wrap,*M,*w,*s,lc,inra);
174: }

176: void PETSC_STDCALL dacreate3d_(MPI_Comm *comm,DAPeriodicType *wrap,DAStencilType 
177:                  *stencil_type,int *M,int *N,int *P,int *m,int *n,int *p,
178:                  int *w,int *s,int *lx,int *ly,int *lz,DA *inra,int *ierr)
179: {
180:   if (FORTRANNULLINTEGER(lx)) lx = PETSC_NULL;
181:   if (FORTRANNULLINTEGER(ly)) ly = PETSC_NULL;
182:   if (FORTRANNULLINTEGER(lz)) lz = PETSC_NULL;
183:   *DACreate3d((MPI_Comm)PetscToPointerComm(*comm),*wrap,*stencil_type,
184:                         *M,*N,*P,*m,*n,*p,*w,*s,lx,ly,lz,inra);
185: }

187: void PETSC_STDCALL dagetinfo_(DA *da,int *dim,int *M,int *N,int *P,int *m,int *n,int *p,int *w,int *s,
188:                 DAPeriodicType *wrap,DAStencilType *st,int *ierr)
189: {
190:   if (FORTRANNULLINTEGER(dim)) dim  = PETSC_NULL;
191:   if (FORTRANNULLINTEGER(M))   M    = PETSC_NULL;
192:   if (FORTRANNULLINTEGER(N))   N    = PETSC_NULL;
193:   if (FORTRANNULLINTEGER(P))   P    = PETSC_NULL;
194:   if (FORTRANNULLINTEGER(m))   m    = PETSC_NULL;
195:   if (FORTRANNULLINTEGER(n))   n    = PETSC_NULL;
196:   if (FORTRANNULLINTEGER(p))   p    = PETSC_NULL;
197:   if (FORTRANNULLINTEGER(w))   w    = PETSC_NULL;
198:   if (FORTRANNULLINTEGER(s))   s    = PETSC_NULL;
199:   if (FORTRANNULLINTEGER(wrap))wrap = PETSC_NULL;
200:   if (FORTRANNULLINTEGER(st))  st   = PETSC_NULL;
201:   *DAGetInfo(*da,dim,M,N,P,m,n,p,w,s,wrap,st);
202: }

204: EXTERN_C_END