Actual source code: zdmdasnesf.c

  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

 12: static struct {
 13:   PetscFortranCallbackId lf1d;
 14:   PetscFortranCallbackId lf2d;
 15:   PetscFortranCallbackId lf3d;
 16:   PetscFortranCallbackId lj1d;
 17:   PetscFortranCallbackId lj2d;
 18:   PetscFortranCallbackId lj3d;
 19: } _cb;

 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;

 27:   PetscFunctionBegin;
 28:   PetscCall(DMGetDMSNES(info->da, &sdm));
 29:   PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, _cb.lj1d, (PetscVoidFn **)&func, &ctx));
 30:   PetscCallFortranVoidFunction((*func)(info, &in[info->dof * info->gxs], &A, &m, ctx, &ierr));
 31:   PetscFunctionReturn(PETSC_SUCCESS);
 32: }

 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;

 39:   PetscFunctionBegin;
 40:   PetscCall(DMGetDMSNES(info->da, &sdm));
 41:   PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, _cb.lj2d, (PetscVoidFn **)&func, &ctx));
 42:   PetscCallFortranVoidFunction((*func)(info, &in[info->gys][info->dof * info->gxs], &A, &m, ctx, &ierr));
 43:   PetscFunctionReturn(PETSC_SUCCESS);
 44: }

 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;

 51:   PetscFunctionBegin;
 52:   PetscCall(DMGetDMSNES(info->da, &sdm));
 53:   PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, _cb.lj2d, (PetscVoidFn **)&func, &ctx));
 54:   PetscCallFortranVoidFunction((*func)(info, &in[info->gzs][info->gys][info->dof * info->gxs], &A, &m, ctx, &ierr));
 55:   PetscFunctionReturn(PETSC_SUCCESS);
 56: }

 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;

 63:   *ierr = DMGetDMSNESWrite(*da, &sdm);
 64:   if (*ierr) return;
 65:   *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
 66:   if (*ierr) return;
 67:   if (dim == 2) {
 68:     *ierr = PetscObjectSetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.lj2d, (PetscVoidFn *)jac, ctx);
 69:     if (*ierr) return;
 70:     *ierr = DMDASNESSetJacobianLocal(*da, (PetscErrorCode (*)(DMDALocalInfo *, void *, Mat, Mat, void *))sourlj2d, NULL);
 71:   } else if (dim == 3) {
 72:     *ierr = PetscObjectSetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.lj3d, (PetscVoidFn *)jac, ctx);
 73:     if (*ierr) return;
 74:     *ierr = DMDASNESSetJacobianLocal(*da, (PetscErrorCode (*)(DMDALocalInfo *, void *, Mat, Mat, void *))sourlj3d, NULL);
 75:   } else if (dim == 1) {
 76:     *ierr = PetscObjectSetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.lj1d, (PetscVoidFn *)jac, ctx);
 77:     if (*ierr) return;
 78:     *ierr = DMDASNESSetJacobianLocal(*da, (PetscErrorCode (*)(DMDALocalInfo *, void *, Mat, Mat, void *))sourlj1d, NULL);
 79:   } else *ierr = PETSC_ERR_ARG_OUTOFRANGE;
 80: }

 82: /************************************************/

 84: static PetscErrorCode sourlf1d(DMDALocalInfo *info, PetscScalar *in, PetscScalar *out, void *ptr)
 85: {
 86:   void (*func)(DMDALocalInfo *, PetscScalar *, PetscScalar *, void *, PetscErrorCode *), *ctx;
 87:   DMSNES sdm;

 89:   PetscFunctionBegin;
 90:   PetscCall(DMGetDMSNES(info->da, &sdm));
 91:   PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, _cb.lf1d, (PetscVoidFn **)&func, &ctx));
 92:   PetscCallFortranVoidFunction((*func)(info, &in[info->dof * info->gxs], &out[info->dof * info->xs], ctx, &ierr));
 93:   PetscFunctionReturn(PETSC_SUCCESS);
 94: }

 96: static PetscErrorCode sourlf2d(DMDALocalInfo *info, PetscScalar **in, PetscScalar **out, void *ptr)
 97: {
 98:   void (*func)(DMDALocalInfo *, PetscScalar *, PetscScalar *, void *, PetscErrorCode *), *ctx;
 99:   DMSNES sdm;

101:   PetscFunctionBegin;
102:   PetscCall(DMGetDMSNES(info->da, &sdm));
103:   PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, _cb.lf2d, (PetscVoidFn **)&func, &ctx));
104:   PetscCallFortranVoidFunction((*func)(info, &in[info->gys][info->dof * info->gxs], &out[info->ys][info->dof * info->xs], ctx, &ierr));
105:   PetscFunctionReturn(PETSC_SUCCESS);
106: }

108: static PetscErrorCode sourlf3d(DMDALocalInfo *info, PetscScalar ***in, PetscScalar ***out, void *ptr)
109: {
110:   void (*func)(DMDALocalInfo *, PetscScalar *, PetscScalar *, void *, PetscErrorCode *), *ctx;
111:   DMSNES sdm;

113:   PetscFunctionBegin;
114:   PetscCall(DMGetDMSNES(info->da, &sdm));
115:   PetscCall(PetscObjectGetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, _cb.lf3d, (PetscVoidFn **)&func, &ctx));
116:   PetscCallFortranVoidFunction((*func)(info, &in[info->gzs][info->gys][info->dof * info->gxs], &out[info->zs][info->ys][info->dof * info->xs], ctx, &ierr));
117:   PetscFunctionReturn(PETSC_SUCCESS);
118: }

120: PETSC_EXTERN void dmdasnessetfunctionlocal_(DM *da, InsertMode *mode, void (*func)(DMDALocalInfo *, void *, void *, void *, PetscErrorCode *), void *ctx, PetscErrorCode *ierr)
121: {
122:   DMSNES   sdm;
123:   PetscInt dim;

125:   *ierr = DMGetDMSNESWrite(*da, &sdm);
126:   if (*ierr) return;
127:   *ierr = DMDAGetInfo(*da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
128:   if (*ierr) return;
129:   if (dim == 2) {
130:     *ierr = PetscObjectSetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.lf2d, (PetscVoidFn *)func, ctx);
131:     if (*ierr) return;
132:     *ierr = DMDASNESSetFunctionLocal(*da, *mode, (PetscErrorCode (*)(DMDALocalInfo *, void *, void *, void *))sourlf2d, NULL);
133:   } else if (dim == 3) {
134:     *ierr = PetscObjectSetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.lf3d, (PetscVoidFn *)func, ctx);
135:     if (*ierr) return;
136:     *ierr = DMDASNESSetFunctionLocal(*da, *mode, (PetscErrorCode (*)(DMDALocalInfo *, void *, void *, void *))sourlf3d, NULL);
137:   } else if (dim == 1) {
138:     *ierr = PetscObjectSetFortranCallback((PetscObject)sdm, PETSC_FORTRAN_CALLBACK_SUBTYPE, &_cb.lf1d, (PetscVoidFn *)func, ctx);
139:     if (*ierr) return;
140:     *ierr = DMDASNESSetFunctionLocal(*da, *mode, (PetscErrorCode (*)(DMDALocalInfo *, void *, void *, void *))sourlf1d, NULL);
141:   } else *ierr = PETSC_ERR_ARG_OUTOFRANGE;
142: }