1 #include <petsc/private/fortranimpl.h> 2 #include <petsc/private/dmdaimpl.h> 3 #include <petsc/private/snesimpl.h> 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define dmdasnessetjacobianlocal_ DMDASNESSETJACOBIANLOCAL 6 #define dmdasnessetfunctionlocal_ DMDASNESSETFUNCTIONLOCAL 7 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 8 #define dmdasnessetjacobianlocal_ dmdasnessetjacobianlocal 9 #define dmdasnessetfunctionlocal_ dmdasnessetfunctionlocal 10 #endif 11 12 static struct { 13 PetscFortranCallbackId lf1d; 14 PetscFortranCallbackId lf2d; 15 PetscFortranCallbackId lf3d; 16 PetscFortranCallbackId lj1d; 17 PetscFortranCallbackId lj2d; 18 PetscFortranCallbackId lj3d; 19 } _cb; 20 21 /************************************************/ 22 static PetscErrorCode sourlj1d(DMDALocalInfo *info,PetscScalar *in,Mat A,Mat m,void *ptr) 23 { 24 void (*func)(DMDALocalInfo*,PetscScalar*,Mat*,Mat*,void*,PetscErrorCode*),*ctx; 25 DMSNES sdm; 26 27 PetscFunctionBegin; 28 PetscCall(DMGetDMSNES(info->da,&sdm)); 29 PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,_cb.lj1d,(PetscVoidFunction*)&func,&ctx)); 30 PetscCallFortranVoidFunction((*func)(info,&in[info->dof*info->gxs],&A,&m,ctx,&ierr)); 31 PetscFunctionReturn(0); 32 } 33 34 static PetscErrorCode sourlj2d(DMDALocalInfo *info,PetscScalar **in,Mat A,Mat m,void *ptr) 35 { 36 void (*func)(DMDALocalInfo*,PetscScalar*,Mat*,Mat*,void*,PetscErrorCode*),*ctx; 37 DMSNES sdm; 38 39 PetscFunctionBegin; 40 PetscCall(DMGetDMSNES(info->da,&sdm)); 41 PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,_cb.lj2d,(PetscVoidFunction*)&func,&ctx)); 42 PetscCallFortranVoidFunction((*func)(info,&in[info->gys][info->dof*info->gxs],&A,&m,ctx,&ierr)); 43 PetscFunctionReturn(0); 44 } 45 46 static PetscErrorCode sourlj3d(DMDALocalInfo *info,PetscScalar ***in,Mat A,Mat m,void *ptr) 47 { 48 void (*func)(DMDALocalInfo*,PetscScalar*,Mat*,Mat*,void*,PetscErrorCode*),*ctx; 49 DMSNES sdm; 50 51 PetscFunctionBegin; 52 PetscCall(DMGetDMSNES(info->da,&sdm)); 53 PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,_cb.lj2d,(PetscVoidFunction*)&func,&ctx)); 54 PetscCallFortranVoidFunction((*func)(info,&in[info->gzs][info->gys][info->dof*info->gxs],&A,&m,ctx,&ierr)); 55 PetscFunctionReturn(0); 56 } 57 58 PETSC_EXTERN void dmdasnessetjacobianlocal_(DM *da,void (*jac)(DMDALocalInfo*,void*,void*,void*,void*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr) 59 { 60 DMSNES sdm; 61 PetscInt dim; 62 63 *ierr = DMGetDMSNESWrite(*da,&sdm); if (*ierr) return; 64 *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,0,0,0,0,0,0); if (*ierr) return; 65 if (dim == 2) { 66 *ierr = PetscObjectSetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,&_cb.lj2d,(PetscVoidFunction)jac,ctx); if (*ierr) return; 67 *ierr = DMDASNESSetJacobianLocal(*da,(PetscErrorCode (*)(DMDALocalInfo*,void*,Mat,Mat,void*))sourlj2d,NULL); 68 } else if (dim == 3) { 69 *ierr = PetscObjectSetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,&_cb.lj3d,(PetscVoidFunction)jac,ctx); if (*ierr) return; 70 *ierr = DMDASNESSetJacobianLocal(*da,(PetscErrorCode (*)(DMDALocalInfo*,void*,Mat,Mat,void*))sourlj3d,NULL); 71 } else if (dim == 1) { 72 *ierr = PetscObjectSetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,&_cb.lj1d,(PetscVoidFunction)jac,ctx); if (*ierr) return; 73 *ierr = DMDASNESSetJacobianLocal(*da,(PetscErrorCode (*)(DMDALocalInfo*,void*,Mat,Mat,void*))sourlj1d,NULL); 74 } else *ierr = 1; 75 } 76 77 /************************************************/ 78 79 static PetscErrorCode sourlf1d(DMDALocalInfo *info,PetscScalar *in,PetscScalar *out,void *ptr) 80 { 81 void (*func)(DMDALocalInfo*,PetscScalar*,PetscScalar*,void*,PetscErrorCode*),*ctx; 82 DMSNES sdm; 83 84 PetscFunctionBegin; 85 PetscCall(DMGetDMSNES(info->da,&sdm)); 86 PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,_cb.lf1d,(PetscVoidFunction*)&func,&ctx)); 87 PetscCallFortranVoidFunction((*func)(info,&in[info->dof*info->gxs],&out[info->dof*info->xs],ctx,&ierr)); 88 PetscFunctionReturn(0); 89 } 90 91 static PetscErrorCode sourlf2d(DMDALocalInfo *info,PetscScalar **in,PetscScalar **out,void *ptr) 92 { 93 void (*func)(DMDALocalInfo*,PetscScalar*,PetscScalar*,void*,PetscErrorCode*),*ctx; 94 DMSNES sdm; 95 96 PetscFunctionBegin; 97 PetscCall(DMGetDMSNES(info->da,&sdm)); 98 PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,_cb.lf2d,(PetscVoidFunction*)&func,&ctx)); 99 PetscCallFortranVoidFunction((*func)(info,&in[info->gys][info->dof*info->gxs],&out[info->ys][info->dof*info->xs],ctx,&ierr)); 100 PetscFunctionReturn(0); 101 } 102 103 static PetscErrorCode sourlf3d(DMDALocalInfo *info,PetscScalar ***in,PetscScalar ***out,void *ptr) 104 { 105 void (*func)(DMDALocalInfo*,PetscScalar*,PetscScalar*,void*,PetscErrorCode*),*ctx; 106 DMSNES sdm; 107 108 PetscFunctionBegin; 109 PetscCall(DMGetDMSNES(info->da,&sdm)); 110 PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,_cb.lf3d,(PetscVoidFunction*)&func,&ctx)); 111 PetscCallFortranVoidFunction((*func)(info,&in[info->gzs][info->gys][info->dof*info->gxs],&out[info->zs][info->ys][info->dof*info->xs],ctx,&ierr)); 112 PetscFunctionReturn(0); 113 } 114 115 PETSC_EXTERN void dmdasnessetfunctionlocal_(DM *da,InsertMode *mode,void (*func)(DMDALocalInfo*,void*,void*,void*,PetscErrorCode*),void *ctx,PetscErrorCode *ierr) 116 { 117 DMSNES sdm; 118 PetscInt dim; 119 120 *ierr = DMGetDMSNESWrite(*da,&sdm); if (*ierr) return; 121 *ierr = DMDAGetInfo(*da,&dim,0,0,0,0,0,0,0,0,0,0,0,0); if (*ierr) return; 122 if (dim == 2) { 123 *ierr = PetscObjectSetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,&_cb.lf2d,(PetscVoidFunction)func,ctx); if (*ierr) return; 124 *ierr = DMDASNESSetFunctionLocal(*da,*mode,(PetscErrorCode (*)(DMDALocalInfo*,void*,void*,void*))sourlf2d,NULL); 125 } else if (dim == 3) { 126 *ierr = PetscObjectSetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,&_cb.lf3d,(PetscVoidFunction)func,ctx); if (*ierr) return; 127 *ierr = DMDASNESSetFunctionLocal(*da,*mode,(PetscErrorCode (*)(DMDALocalInfo*,void*,void*,void*))sourlf3d,NULL); 128 } else if (dim == 1) { 129 *ierr = PetscObjectSetFortranCallback((PetscObject)sdm,PETSC_FORTRAN_CALLBACK_SUBTYPE,&_cb.lf1d,(PetscVoidFunction)func,ctx); if (*ierr) return; 130 *ierr = DMDASNESSetFunctionLocal(*da,*mode,(PetscErrorCode (*)(DMDALocalInfo*,void*,void*,void*))sourlf1d,NULL); 131 } else *ierr = 1; 132 } 133