1. petsc
  2. PETSc
  3. petsc

Commits

Toby Isaac  committed c0c1684

Move DMPlexProject*() functions to src/ksp to remove ksp dependency in dm

  • Participants
  • Parent commits 770c2a3
  • Branches knepley/fix-quadrature-order, next-oct-2014

Comments (0)

Files changed (7)

File include/petscdm.h

View file
  • Ignore whitespace
 PETSC_EXTERN PetscErrorCode DMLocalToGlobalHookAdd(DM,PetscErrorCode (*)(DM,Vec,InsertMode,Vec,void*),PetscErrorCode (*)(DM,Vec,InsertMode,Vec,void*),void*);
 PETSC_EXTERN PetscErrorCode DMGlobalToLocalBegin(DM,Vec,InsertMode,Vec);
 PETSC_EXTERN PetscErrorCode DMGlobalToLocalEnd(DM,Vec,InsertMode,Vec);
-PETSC_EXTERN PetscErrorCode DMGlobalToLocalSolve(DM,Vec,Vec);
 PETSC_EXTERN PetscErrorCode DMLocalToGlobalBegin(DM,Vec,InsertMode,Vec);
 PETSC_EXTERN PetscErrorCode DMLocalToGlobalEnd(DM,Vec,InsertMode,Vec);
 PETSC_EXTERN PetscErrorCode DMLocalToLocalBegin(DM,Vec,InsertMode,Vec);

File include/petscdmplex.h

View file
  • Ignore whitespace
 } JacActionCtx;
 
 PETSC_EXTERN PetscErrorCode DMPlexInsertBoundaryValuesFEM(DM, Vec);
-PETSC_EXTERN PetscErrorCode DMPlexProjectFunction(DM, void (**)(const PetscReal [], PetscScalar *, void *), void **, InsertMode, Vec);
 PETSC_EXTERN PetscErrorCode DMPlexProjectFunctionLocal(DM, void (**)(const PetscReal [], PetscScalar *, void *), void **, InsertMode, Vec);
-PETSC_EXTERN PetscErrorCode DMPlexProjectField(DM, Vec, void (**)(const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscReal [], PetscScalar []), InsertMode, Vec);
 PETSC_EXTERN PetscErrorCode DMPlexProjectFieldLocal(DM, Vec, void (**)(const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscReal [], PetscScalar []), InsertMode, Vec);
 PETSC_EXTERN PetscErrorCode DMPlexComputeL2Diff(DM, void (**)(const PetscReal [], PetscScalar *, void *), void **, Vec, PetscReal *);
 PETSC_EXTERN PetscErrorCode DMPlexComputeL2GradientDiff(DM, void (**)(const PetscReal [], const PetscReal [], PetscScalar *, void *), void **, Vec, const PetscReal [], PetscReal *);

File include/petscksp.h

View file
  • Ignore whitespace
 PETSC_EXTERN PetscErrorCode DMKSPSetComputeInitialGuess(DM,PetscErrorCode(*)(KSP,Vec,void*),void*);
 PETSC_EXTERN PetscErrorCode DMKSPGetComputeInitialGuess(DM,PetscErrorCode(**)(KSP,Vec,void*),void*);
 
+PETSC_EXTERN PetscErrorCode DMGlobalToLocalSolve(DM,Vec,Vec);
+PETSC_EXTERN PetscErrorCode DMPlexProjectField(DM, Vec, void (**)(const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscReal [], PetscScalar []), InsertMode, Vec);
+PETSC_EXTERN PetscErrorCode DMPlexProjectFunction(DM, void (**)(const PetscReal [], PetscScalar *, void *), void **, InsertMode, Vec);
 #endif

File src/dm/impls/plex/examples/tests/ex3.c

View file
  • Ignore whitespace
 #include <petscdmda.h>
 #include <petscfe.h>
 #include <petscds.h>
+#include <petscksp.h>
 
 typedef struct {
   PetscInt  debug;             /* The debugging level */

File src/dm/impls/plex/examples/tests/makefile

View file
  • Ignore whitespace
 	${RM} -f ex2f90.o
 
 ex3: ex3.o  chkopts
-	-${CLINKER} -o ex3 ex3.o ${PETSC_DM_LIB}
+	-${CLINKER} -o ex3 ex3.o ${PETSC_KSP_LIB}
 	${RM} -f ex3.o
 
 #--------------------------------------------------------------------------

File src/dm/impls/plex/plexfem.c

View file
  • Ignore whitespace
 }
 
 #undef __FUNCT__
-#define __FUNCT__ "DMPlexProjectFunction"
-/*@C
-  DMPlexProjectFunction - This projects the given function into the function space provided.
-
-  Input Parameters:
-+ dm      - The DM
-. funcs   - The coordinate functions to evaluate, one per field
-. ctxs    - Optional array of contexts to pass to each coordinate function.  ctxs itself may be null.
-- mode    - The insertion mode for values
-
-  Output Parameter:
-. X - vector
-
-  Level: developer
-
-.seealso: DMPlexComputeL2Diff()
-@*/
-PetscErrorCode DMPlexProjectFunction(DM dm, void (**funcs)(const PetscReal [], PetscScalar *, void *), void **ctxs, InsertMode mode, Vec X)
-{
-  Vec            localX;
-  PetscErrorCode ierr;
-
-  PetscFunctionBegin;
-  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
-  ierr = DMGetLocalVector(dm, &localX);CHKERRQ(ierr);
-  ierr = DMPlexProjectFunctionLocal(dm, funcs, ctxs, mode, localX);CHKERRQ(ierr);
-  ierr = DMLocalToGlobalBegin(dm, localX, mode, X);CHKERRQ(ierr);
-  ierr = DMLocalToGlobalEnd(dm, localX, mode, X);CHKERRQ(ierr);
-  if (mode == INSERT_VALUES || mode == INSERT_ALL_VALUES || mode == INSERT_BC_VALUES) {
-    Mat cMat;
-
-    ierr = DMPlexGetConstraintMatrix(dm, &cMat);CHKERRQ(ierr);
-    if (cMat) {
-      ierr = DMGlobalToLocalSolve(dm, localX, X);CHKERRQ(ierr);
-    }
-  }
-  ierr = DMRestoreLocalVector(dm, &localX);CHKERRQ(ierr);
-  PetscFunctionReturn(0);
-}
-
-#undef __FUNCT__
 #define __FUNCT__ "DMPlexProjectFieldLocal"
 PetscErrorCode DMPlexProjectFieldLocal(DM dm, Vec localU, void (**funcs)(const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscReal [], PetscScalar []), InsertMode mode, Vec localX)
 {
 }
 
 #undef __FUNCT__
-#define __FUNCT__ "DMPlexProjectField"
-/*@C
-  DMPlexProjectField - This projects the given function of the fields into the function space provided.
-
-  Input Parameters:
-+ dm      - The DM
-. U       - The input field vector
-. funcs   - The functions to evaluate, one per field
-- mode    - The insertion mode for values
-
-  Output Parameter:
-. X       - The output vector
-
-  Level: developer
-
-.seealso: DMPlexProjectFunction(), DMPlexComputeL2Diff()
-@*/
-PetscErrorCode DMPlexProjectField(DM dm, Vec U, void (**funcs)(const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscReal [], PetscScalar []), InsertMode mode, Vec X)
-{
-  Vec            localX, localU;
-  PetscErrorCode ierr;
-
-  PetscFunctionBegin;
-  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
-  ierr = DMGetLocalVector(dm, &localX);CHKERRQ(ierr);
-  ierr = DMGetLocalVector(dm, &localU);CHKERRQ(ierr);
-  ierr = DMGlobalToLocalBegin(dm, U, INSERT_VALUES, localU);CHKERRQ(ierr);
-  ierr = DMGlobalToLocalEnd(dm, U, INSERT_VALUES, localU);CHKERRQ(ierr);
-  ierr = DMPlexProjectFieldLocal(dm, localU, funcs, mode, localX);CHKERRQ(ierr);
-  ierr = DMLocalToGlobalBegin(dm, localX, mode, X);CHKERRQ(ierr);
-  ierr = DMLocalToGlobalEnd(dm, localX, mode, X);CHKERRQ(ierr);
-  if (mode == INSERT_VALUES || mode == INSERT_ALL_VALUES || mode == INSERT_BC_VALUES) {
-    Mat cMat;
-
-    ierr = DMPlexGetConstraintMatrix(dm, &cMat);CHKERRQ(ierr);
-    if (cMat) {
-      ierr = DMGlobalToLocalSolve(dm, localX, X);CHKERRQ(ierr);
-    }
-  }
-  ierr = DMRestoreLocalVector(dm, &localX);CHKERRQ(ierr);
-  ierr = DMRestoreLocalVector(dm, &localU);CHKERRQ(ierr);
-  PetscFunctionReturn(0);
-}
-
-#undef __FUNCT__
 #define __FUNCT__ "DMPlexInsertBoundaryValuesFEM"
 PetscErrorCode DMPlexInsertBoundaryValuesFEM(DM dm, Vec localX)
 {

File src/ksp/ksp/utils/dmproject.c

View file
  • Ignore whitespace
 
+#include <petsc-private/petscimpl.h>
 #include <petscdm.h>     /*I "petscdm.h" I*/
 #include <petscdmplex.h> /*I "petscdmplex.h" I*/
 #include <petscksp.h>    /*I "petscksp.h" I*/
   ierr = MatDestroy(&CtC);CHKERRQ(ierr);
   PetscFunctionReturn(0);
 }
+
+#undef __FUNCT__
+#define __FUNCT__ "DMPlexProjectFunction"
+/*@C
+  DMPlexProjectFunction - This projects the given function into the function space provided.
+
+  Input Parameters:
++ dm      - The DM
+. funcs   - The coordinate functions to evaluate, one per field
+. ctxs    - Optional array of contexts to pass to each coordinate function.  ctxs itself may be null.
+- mode    - The insertion mode for values
+
+  Output Parameter:
+. X - vector
+
+  Level: developer
+
+.seealso: DMPlexComputeL2Diff()
+@*/
+PetscErrorCode DMPlexProjectFunction(DM dm, void (**funcs)(const PetscReal [], PetscScalar *, void *), void **ctxs, InsertMode mode, Vec X)
+{
+  Vec            localX;
+  PetscErrorCode ierr;
+
+  PetscFunctionBegin;
+  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
+  ierr = DMGetLocalVector(dm, &localX);CHKERRQ(ierr);
+  ierr = DMPlexProjectFunctionLocal(dm, funcs, ctxs, mode, localX);CHKERRQ(ierr);
+  ierr = DMLocalToGlobalBegin(dm, localX, mode, X);CHKERRQ(ierr);
+  ierr = DMLocalToGlobalEnd(dm, localX, mode, X);CHKERRQ(ierr);
+  if (mode == INSERT_VALUES || mode == INSERT_ALL_VALUES || mode == INSERT_BC_VALUES) {
+    Mat cMat;
+
+    ierr = DMPlexGetConstraintMatrix(dm, &cMat);CHKERRQ(ierr);
+    if (cMat) {
+      ierr = DMGlobalToLocalSolve(dm, localX, X);CHKERRQ(ierr);
+    }
+  }
+  ierr = DMRestoreLocalVector(dm, &localX);CHKERRQ(ierr);
+  PetscFunctionReturn(0);
+}
+
+#undef __FUNCT__
+#define __FUNCT__ "DMPlexProjectField"
+/*@C
+  DMPlexProjectField - This projects the given function of the fields into the function space provided.
+
+  Input Parameters:
++ dm      - The DM
+. U       - The input field vector
+. funcs   - The functions to evaluate, one per field
+- mode    - The insertion mode for values
+
+  Output Parameter:
+. X       - The output vector
+
+  Level: developer
+
+.seealso: DMPlexProjectFunction(), DMPlexComputeL2Diff()
+@*/
+PetscErrorCode DMPlexProjectField(DM dm, Vec U, void (**funcs)(const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscReal [], PetscScalar []), InsertMode mode, Vec X)
+{
+  Vec            localX, localU;
+  PetscErrorCode ierr;
+
+  PetscFunctionBegin;
+  PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
+  ierr = DMGetLocalVector(dm, &localX);CHKERRQ(ierr);
+  ierr = DMGetLocalVector(dm, &localU);CHKERRQ(ierr);
+  ierr = DMGlobalToLocalBegin(dm, U, INSERT_VALUES, localU);CHKERRQ(ierr);
+  ierr = DMGlobalToLocalEnd(dm, U, INSERT_VALUES, localU);CHKERRQ(ierr);
+  ierr = DMPlexProjectFieldLocal(dm, localU, funcs, mode, localX);CHKERRQ(ierr);
+  ierr = DMLocalToGlobalBegin(dm, localX, mode, X);CHKERRQ(ierr);
+  ierr = DMLocalToGlobalEnd(dm, localX, mode, X);CHKERRQ(ierr);
+  if (mode == INSERT_VALUES || mode == INSERT_ALL_VALUES || mode == INSERT_BC_VALUES) {
+    Mat cMat;
+
+    ierr = DMPlexGetConstraintMatrix(dm, &cMat);CHKERRQ(ierr);
+    if (cMat) {
+      ierr = DMGlobalToLocalSolve(dm, localX, X);CHKERRQ(ierr);
+    }
+  }
+  ierr = DMRestoreLocalVector(dm, &localX);CHKERRQ(ierr);
+  ierr = DMRestoreLocalVector(dm, &localU);CHKERRQ(ierr);
+  PetscFunctionReturn(0);
+}
+