Actual source code: zsles.c

  1: /*$Id: zsles.c,v 1.35 2001/04/02 15:10:29 bsmith Exp $*/

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

  7: #ifdef PETSC_HAVE_FORTRAN_CAPS
  8: #define slesdestroy_             SLESDESTROY
  9: #define slescreate_              SLESCREATE
 10: #define slesgetpc_               SLESGETPC
 11: #define slessetoptionsprefix_    SLESSETOPTIONSPREFIX
 12: #define slesappendoptionsprefix_ SLESAPPENDOPTIONSPREFIX
 13: #define slesgetksp_              SLESGETKSP
 14: #define slesgetoptionsprefix_    SLESGETOPTIONSPREFIX
 15: #define slesview_                SLESVIEW
 16: #define dmmgcreate_              DMMGCREATE
 17: #define dmmgdestroy_             DMMGDESTROY
 18: #define dmmgsetup_               DMMGSETUP
 19: #define dmmgsetdm_               DMMGSETDM
 20: #define dmmgview_                DMMGVIEW
 21: #define dmmgsolve_               DMMGSOLVE
 22: #define dmmgsetusematrixfree_    DMMGSETUSEMATRIXFREE
 23: #define dmmggetda_               DMMGGETDA
 24: #define dmmgsetsles_             DMMGSETSLES
 25: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 26: #define dmmgsetsles_             dmmgsetsles
 27: #define dmmgdestroy_             dmmgdestroy
 28: #define dmmgcreate_              dmmgcreate
 29: #define dmmgsetup_               dmmgsetup
 30: #define slessetoptionsprefix_    slessetoptionsprefix
 31: #define slesappendoptionsprefix_ slesappendoptionsprefix
 32: #define slesdestroy_             slesdestroy
 33: #define slescreate_              slescreate
 34: #define slesgetpc_               slesgetpc
 35: #define slesgetksp_              slesgetksp
 36: #define slesgetoptionsprefix_    slesgetoptionsprefix
 37: #define slesview_                slesview
 38: #define dmmgsetdm_               dmmgsetdm
 39: #define dmmggetda_               dmmggetda
 40: #define dmmgview_                dmmgview
 41: #define dmmgsolve_               dmmgsolve
 42: #define dmmgsetusematrixfree_    dmmgsetusematrixfree
 43: #endif

 45: EXTERN_C_BEGIN

 47: /* ----------------------------------------------------------------------------------------------------------*/
 48: static int ourrhs(DMMG dmmg,Vec vec)
 49: {
 50:   int              0;
 51:   (*(int (PETSC_STDCALL *)(DMMG*,Vec*,int*))(((PetscObject)dmmg->dm)->fortran_func_pointers[0]))(&dmmg,&vec,&ierr);
 52:   return ierr;
 53: }

 55: /*
 56:    Since DMMGSetSLES() immediately calls the matrix functions for each level we do not need to store
 57:   the mat() function inside the DMMG object
 58: */
 59: static int (PETSC_STDCALL *theirmat)(DMMG*,Mat*,int*);
 60: static int ourmat(DMMG dmmg,Mat mat)
 61: {
 62:   int              0;
 63:   (*theirmat)(&dmmg,&mat,&ierr);
 64:   return ierr;
 65: }

 67: void PETSC_STDCALL dmmgsetsles_(DMMG **dmmg,int (PETSC_STDCALL *rhs)(DMMG*,Vec*,int*),int (PETSC_STDCALL *mat)(DMMG*,Mat*,int*),int *ierr)
 68: {
 69:   int i;
 70:   theirmat = mat;
 71:   *DMMGSetSLES(*dmmg,ourrhs,ourmat);
 72:   /*
 73:     Save the fortran rhs function in the DM on each level; ourrhs() pulls it out when needed
 74:   */
 75:   for (i=0; i<(**dmmg)->nlevels; i++) {
 76:     ((PetscObject)(*dmmg)[i]->dm)->fortran_func_pointers[0] = (void (*)())rhs;
 77:   }
 78: }

 80: /* ----------------------------------------------------------------------------------------------------------*/

 82: void PETSC_STDCALL dmmggetda_(DMMG *dmmg,DA *da,int *ierr)
 83: {
 84:   *da   = (DA)(*dmmg)->dm;
 85:   *0;
 86: }

 88: void PETSC_STDCALL dmmgsetdm_(DMMG **dmmg,DM *dm,int *ierr)
 89: {
 90:   int i;
 91:   *DMMGSetDM(*dmmg,*dm);if (*ierr) return;
 92:   /* loop over the levels added a place to hang the function pointers in the DM for each level*/
 93:   for (i=0; i<(**dmmg)->nlevels; i++) {
 94:     *PetscMalloc(3*sizeof(void (*)()),&((PetscObject)(*dmmg)[i]->dm)->fortran_func_pointers);if (*ierr) return;
 95:   }
 96: }

 98: void PETSC_STDCALL dmmgview_(DMMG **dmmg,PetscViewer *viewer,int *ierr)
 99: {
100:   *DMMGView(*dmmg,*viewer);
101: }

103: void PETSC_STDCALL dmmgsolve_(DMMG **dmmg,int *ierr)
104: {
105:   *DMMGSolve(*dmmg);
106: }

108: void PETSC_STDCALL dmmgsetusematrixfree_(DMMG **dmmg,int *ierr)
109: {
110:   *DMMGSetUseMatrixFree(*dmmg);
111: }

113: void PETSC_STDCALL dmmgcreate_(MPI_Comm *comm,int *nlevels,void *user,DMMG **dmmg,int *ierr)
114: {
115:   *DMMGCreate((MPI_Comm)PetscToPointerComm(*comm),*nlevels,user,dmmg);
116: }

118: void PETSC_STDCALL dmmgdestroy_(DMMG **dmmg,int *ierr)
119: {
120:   *DMMGDestroy(*dmmg);
121: }

123: void PETSC_STDCALL dmmgsetup_(DMMG **dmmg,int *ierr)
124: {
125:   *DMMGSetUp(*dmmg);
126: }

128: void PETSC_STDCALL slesview_(SLES *sles,PetscViewer *viewer, int *ierr)
129: {
130:   PetscViewer v;
131:   PetscPatchDefaultViewers_Fortran(viewer,v);
132:   *SLESView(*sles,v);
133: }

135: void PETSC_STDCALL slessetoptionsprefix_(SLES *sles,CHAR prefix PETSC_MIXED_LEN(len),
136:                                          int *ierr PETSC_END_LEN(len))
137: {
138:   char *t;

140:   FIXCHAR(prefix,len,t);
141:   *SLESSetOptionsPrefix(*sles,t);
142:   FREECHAR(prefix,t);
143: }

145: void PETSC_STDCALL slesappendoptionsprefix_(SLES *sles,CHAR prefix PETSC_MIXED_LEN(len),
146:                                             int *ierr PETSC_END_LEN(len))
147: {
148:   char *t;

150:   FIXCHAR(prefix,len,t);
151:   *SLESAppendOptionsPrefix(*sles,t);
152:   FREECHAR(prefix,t);
153: }

155: void PETSC_STDCALL slesgetksp_(SLES *sles,KSP *ksp,int *ierr)
156: {
157:   *SLESGetKSP(*sles,ksp);
158: }

160: void PETSC_STDCALL slesgetpc_(SLES *sles,PC *pc,int *ierr)
161: {
162:   *SLESGetPC(*sles,pc);
163: }

165: void PETSC_STDCALL slesdestroy_(SLES *sles,int *ierr)
166: {
167:   *SLESDestroy(*sles);
168: }

170: void PETSC_STDCALL slescreate_(MPI_Comm *comm,SLES *outsles,int *ierr)
171: {
172:   *SLESCreate((MPI_Comm)PetscToPointerComm(*comm),outsles);

174: }

176: void PETSC_STDCALL slesgetoptionsprefix_(SLES *sles,CHAR prefix PETSC_MIXED_LEN(len),
177:                                          int *ierr PETSC_END_LEN(len))
178: {
179:   char *tname;

181:   *SLESGetOptionsPrefix(*sles,&tname);
182: #if defined(PETSC_USES_CPTOFCD)
183:   {
184:     char *t = _fcdtocp(prefix); int len1 = _fcdlen(prefix);
185:     *PetscStrncpy(t,tname,len1);
186:   }
187: #else
188:   *PetscStrncpy(prefix,tname,len);
189: #endif
190: }

192: EXTERN_C_END