Skip to content

Commit

Permalink
Pressure boundary condition added. Need to rotate normals.
Browse files Browse the repository at this point in the history
  • Loading branch information
@rezgarshakeri
rezgarshakeri committed Jun 27, 2022
1 parent 5c51868 commit 2a70772
Show file tree
Hide file tree
Showing 27 changed files with 405 additions and 434 deletions.
4 changes: 3 additions & 1 deletion examples/Hdiv-mixed/basis/Hdiv-hex.h
View file
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
#ifndef Hdiv_hex_h
#define Hdiv_hex_h
// Copyright (c) 2017-2018, Lawrence Livermore National Security, LLC.
// Produced at the Lawrence Livermore National Laboratory. LLNL-CODE-734707.
// All Rights reserved. See files LICENSE and NOTICE for details.
Expand Down Expand Up @@ -159,4 +161,4 @@ static void HdivBasisHex(CeedInt Q, CeedScalar *q_ref, CeedScalar *q_weights,
}
}


#endif // Hdiv_hex_h
4 changes: 3 additions & 1 deletion examples/Hdiv-mixed/basis/Hdiv-quad.h
View file
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
#ifndef Hdiv_quad_h
#define Hdiv_quad_h
// Copyright (c) 2017-2018, Lawrence Livermore National Security, LLC.
// Produced at the Lawrence Livermore National Laboratory. LLNL-CODE-734707.
// All Rights reserved. See files LICENSE and NOTICE for details.
Expand Down Expand Up @@ -84,4 +86,4 @@ static void HdivBasisQuad(CeedInt Q, CeedScalar *q_ref, CeedScalar *q_weights,
}
}


#endif // Hdiv_quad_h
4 changes: 2 additions & 2 deletions examples/Hdiv-mixed/include/cl-options.h
View file
Original file line number Diff line number Diff line change
@@ -1,9 +1,9 @@
#ifndef cloptions_h
#define cloptions_h

#include "../include/structs.h"
#include "structs.h"

// Process general command line options
PetscErrorCode ProcessCommandLineOptions(MPI_Comm comm, AppCtx app_ctx);
PetscErrorCode ProcessCommandLineOptions(AppCtx app_ctx);

#endif // cloptions_h
6 changes: 3 additions & 3 deletions examples/Hdiv-mixed/include/register-problem.h
View file
Original file line number Diff line number Diff line change
@@ -1,18 +1,18 @@
#ifndef register_problems_h
#define register_problems_h

#include "../include/structs.h"
#include "structs.h"

// Register problems to be available on the command line
PetscErrorCode RegisterProblems_Hdiv(AppCtx app_ctx);
// -----------------------------------------------------------------------------
// Set up problems function prototype
// -----------------------------------------------------------------------------
// 1) darcy2d
PetscErrorCode Hdiv_DARCY2D(ProblemData *problem_data, void *ctx);
PetscErrorCode Hdiv_DARCY2D(ProblemData problem_data, void *ctx);

// 2) darcy3d
PetscErrorCode Hdiv_DARCY3D(ProblemData *problem_data, void *ctx);
PetscErrorCode Hdiv_DARCY3D(ProblemData problem_data, void *ctx);

// 3) darcy3dprism

Expand Down
11 changes: 7 additions & 4 deletions examples/Hdiv-mixed/include/setup-boundary.h
View file
Original file line number Diff line number Diff line change
@@ -1,11 +1,11 @@
#ifndef register_boundary_h
#define register_boundary_h
#ifndef setup_boundary_h
#define setup_boundary_h

#include <petsc.h>
#include <petscdmplex.h>
#include <petscsys.h>
#include <ceed.h>
#include "../include/structs.h"
#include "structs.h"

// ---------------------------------------------------------------------------
// Create boundary label
Expand All @@ -19,4 +19,7 @@ PetscErrorCode DMAddBoundariesDirichlet(DM dm);
PetscErrorCode BoundaryDirichletMMS(PetscInt dim, PetscReal t,
const PetscReal coords[],
PetscInt num_comp_u, PetscScalar *u, void *ctx);
#endif // register_boundary_h
PetscErrorCode DMAddBoundariesPressure(Ceed ceed, CeedData ceed_data,
AppCtx app_ctx, ProblemData problem_data, DM dm,
CeedVector bc_pressure);
#endif // setup_boundary_h
4 changes: 2 additions & 2 deletions examples/Hdiv-mixed/include/setup-dm.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,12 +5,12 @@
#include <petscdmplex.h>
#include <petscsys.h>
#include <ceed.h>
#include "../include/structs.h"
#include "structs.h"

// ---------------------------------------------------------------------------
// Set-up DM
// ---------------------------------------------------------------------------
PetscErrorCode CreateDistributedDM(MPI_Comm comm, ProblemData *problem_data,
PetscErrorCode CreateDistributedDM(MPI_Comm comm, ProblemData problem_data,
DM *dm);

#endif // setupdm_h
4 changes: 2 additions & 2 deletions examples/Hdiv-mixed/include/setup-libceed.h
View file
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
#ifndef setuplibceed_h
#define setuplibceed_h

#include "../include/structs.h"
#include "structs.h"

// Convert PETSc MemType to libCEED MemType
CeedMemType MemTypeP2C(PetscMemType mtype);
Expand All @@ -18,7 +18,7 @@ PetscErrorCode CreateRestrictionFromPlexOriented(Ceed ceed, DM dm, CeedInt P,
CeedElemRestriction *elem_restr_oriented, CeedElemRestriction *elem_restr);
// Set up libCEED for a given degree
PetscErrorCode SetupLibceed(DM dm, Ceed ceed, AppCtx app_ctx,
ProblemData *problem_data,
ProblemData problem_data,
PetscInt rhs_loc_size, CeedData ceed_data,
CeedVector rhs_ceed, CeedVector *target);
#endif // setuplibceed_h
4 changes: 2 additions & 2 deletions examples/Hdiv-mixed/include/setup-solvers.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@

#include "structs.h"

PetscErrorCode SetupCommonCtx(MPI_Comm comm, DM dm, Ceed ceed,
PetscErrorCode SetupCommonCtx(DM dm, Ceed ceed,
CeedData ceed_data,
OperatorApplyContext op_apply_ctx);
PetscErrorCode SetupJacobianOperatorCtx(CeedData ceed_data,
Expand All @@ -24,7 +24,7 @@ PetscErrorCode PDESolver(CeedData ceed_data, VecType vec_type, SNES snes,
PetscErrorCode ComputeL2Error(CeedData ceed_data, Vec U, CeedVector target,
CeedScalar *l2_error_u, CeedScalar *l2_error_p,
OperatorApplyContext op_apply_ctx);
PetscErrorCode PrintOutput(MPI_Comm comm, Ceed ceed,
PetscErrorCode PrintOutput(Ceed ceed,
CeedMemType mem_type_backend,
SNES snes, KSP ksp,
Vec U, CeedScalar l2_error_u,
Expand Down
46 changes: 32 additions & 14 deletions examples/Hdiv-mixed/include/structs.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -7,13 +7,14 @@
// Application context from user command line options
typedef struct AppCtx_ *AppCtx;
struct AppCtx_ {
MPI_Comm comm;
// Degree of polynomial (1 only), extra quadrature pts
PetscInt degree;
PetscInt q_extra;
// For applying traction BCs
PetscInt bc_traction_count;
PetscInt bc_traction_faces[16];
PetscScalar bc_traction_vector[16][3];
PetscInt bc_pressure_count;
PetscInt bc_faces[16]; // face ID
PetscScalar bc_pressure_value[16];
// Problem type arguments
PetscFunctionList problems;
char problem_name[PETSC_MAX_PATH_LEN];
Expand All @@ -23,13 +24,12 @@ struct AppCtx_ {
// libCEED data struct
typedef struct CeedData_ *CeedData;
struct CeedData_ {
CeedBasis basis_x, basis_u, basis_p;
CeedBasis basis_x, basis_u, basis_p, basis_u_face;
CeedElemRestriction elem_restr_x, elem_restr_u, elem_restr_U_i,
elem_restr_p;
CeedQFunction qf_residual, qf_jacobian, qf_error;
CeedOperator op_residual, op_jacobian, op_error;
CeedVector x_ceed, y_ceed;
CeedQFunctionContext pq2d_context;
CeedVector x_ceed, y_ceed, x_coord;
};

// 1) darcy2d
Expand All @@ -50,13 +50,31 @@ struct DARCY3DContext_ {
};
#endif

// 4) richard
// 3) richard2d
#ifndef PHYSICS_RICHARD2D_STRUCT
#define PHYSICS_RICHARD2D_STRUCT
typedef struct RICHARD2DContext_ *RICHARD2DContext;
struct RICHARD2DContext_ {
CeedScalar kappa;
};
#endif

// 4) richard3d
#ifndef PHYSICS_RICHARD3D_STRUCT
#define PHYSICS_RICHARD3D_STRUCT
typedef struct RICHARD3DContext_ *RICHARD3DContext;
struct RICHARD3DContext_ {
CeedScalar kappa;
};
#endif

// Struct that contains all enums and structs used for the physics of all problems
typedef struct Physics_ *Physics;
struct Physics_ {
DARCY2DContext darcy2d_ctx;
DARCY3DContext darcy3d_ctx;
RICHARD2DContext richard2d_ctx;
RICHARD3DContext richard3d_ctx;
};

// PETSc operator contexts
Expand All @@ -71,15 +89,15 @@ struct OperatorApplyContext_ {
};

// Problem specific data
typedef struct {
CeedQFunctionUser setup_rhs, residual, jacobian, setup_error,
setup_true, setup_face_geo;
const char *setup_rhs_loc, *residual_loc, *jacobian_loc,
*setup_error_loc, *setup_true_loc, *setup_face_geo_loc;
typedef struct ProblemData_ *ProblemData;
struct ProblemData_ {
CeedQFunctionUser force, residual, jacobian, error,
setup_true, bc_pressure;
const char *force_loc, *residual_loc, *jacobian_loc,
*error_loc, *setup_true_loc, *bc_pressure_loc;
CeedQuadMode quadrature_mode;
CeedInt elem_node, dim, q_data_size_face;
PetscErrorCode (*setup_ctx)(Ceed, CeedData, Physics);

} ProblemData;
};

#endif // structs_h
32 changes: 20 additions & 12 deletions examples/Hdiv-mixed/main.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,9 +23,11 @@
//
// Build with: make
// Run with:
// ./main -pc_type svd -problem darcy2d -dm_plex_dim 2 -dm_plex_box_faces 4,4
// ./main -pc_type svd -problem darcy3d -dm_plex_dim 3 -dm_plex_box_faces 4,4,4
// ./main -pc_type svd -problem darcy3d -dm_plex_filename /path to the mesh file
// ./main -pc_type svd -problem darcy2d -dm_plex_dim 2 -dm_plex_box_faces 4,4
// ./main -pc_type svd -problem darcy3d -dm_plex_dim 3 -dm_plex_box_faces 4,4,4
// ./main -pc_type svd -problem darcy3d -dm_plex_filename /path to the mesh file
// ./main -pc_type svd -problem darcy2d -dm_plex_dim 2 -dm_plex_box_faces 4,4 -bc_pressure 1
// ./main -pc_type svd -problem darcy2d -dm_plex_dim 2 -dm_plex_box_faces 4,4 -bc_pressure 1,2,3,4
const char help[] = "Solve H(div)-mixed problem using PETSc and libCEED\n";

#include "main.h"
Expand All @@ -42,7 +44,7 @@ int main(int argc, char **argv) {
AppCtx app_ctx;
PetscCall( PetscCalloc1(1, &app_ctx) );

ProblemData *problem_data = NULL;
ProblemData problem_data = NULL;
PetscCall( PetscCalloc1(1, &problem_data) );

CeedData ceed_data;
Expand All @@ -62,13 +64,14 @@ int main(int argc, char **argv) {

// -- Process general command line options
MPI_Comm comm = PETSC_COMM_WORLD;
PetscCall( ProcessCommandLineOptions(comm, app_ctx) );
app_ctx->comm = comm;
PetscCall( ProcessCommandLineOptions(app_ctx) );

// ---------------------------------------------------------------------------
// Choose the problem from the list of registered problems
// ---------------------------------------------------------------------------
{
PetscErrorCode (*p)(ProblemData *, void *);
PetscErrorCode (*p)(ProblemData, void *);
PetscCall( PetscFunctionListFind(app_ctx->problems, app_ctx->problem_name,
&p) );
if (!p) SETERRQ(PETSC_COMM_SELF, 1, "Problem '%s' not found",
Expand Down Expand Up @@ -107,12 +110,12 @@ int main(int argc, char **argv) {


// ---------------------------------------------------------------------------
// Create local RHS vector
// Create local Force vector
// ---------------------------------------------------------------------------
Vec F_loc;
PetscInt F_loc_size;
CeedScalar *f;
CeedVector force_ceed, target;
CeedVector force_ceed, target, bc_pressure;
PetscMemType force_mem_type;
PetscCall( DMCreateLocalVector(dm, &F_loc) );
// Local size for libCEED
Expand All @@ -121,15 +124,18 @@ int main(int argc, char **argv) {
PetscCall( VecGetArrayAndMemType(F_loc, &f, &force_mem_type) );
CeedVectorCreate(ceed, F_loc_size, &force_ceed);
CeedVectorSetArray(force_ceed, MemTypeP2C(force_mem_type), CEED_USE_POINTER, f);

CeedVectorCreate(ceed, F_loc_size, &bc_pressure);
CeedVectorSetArray(bc_pressure, MemTypeP2C(force_mem_type), CEED_USE_POINTER,
f);
// ---------------------------------------------------------------------------
// Setup libCEED - Compute local F and true solution (target)
// ---------------------------------------------------------------------------
// -- Set up libCEED objects
PetscCall( SetupLibceed(dm, ceed, app_ctx, problem_data,
F_loc_size, ceed_data, force_ceed, &target) );
//CeedVectorView(force_ceed, "%12.8f", stdout);

PetscCall( DMAddBoundariesPressure(ceed, ceed_data, app_ctx, problem_data, dm,
bc_pressure) );
// ---------------------------------------------------------------------------
// Create global F
// ---------------------------------------------------------------------------
Expand All @@ -147,9 +153,10 @@ int main(int argc, char **argv) {
SNES snes;
KSP ksp;
Vec U;
op_apply_ctx->comm = comm;
PetscCall( SNESCreate(comm, &snes) );
PetscCall( SNESGetKSP(snes, &ksp) );
PetscCall( SetupCommonCtx(comm, dm, ceed, ceed_data, op_apply_ctx) );
PetscCall( SetupCommonCtx(dm, ceed, ceed_data, op_apply_ctx) );
PetscCall( PDESolver(ceed_data, vec_type, snes, ksp, F, &U, op_apply_ctx) );
//VecView(U, PETSC_VIEWER_STDOUT_WORLD);

Expand All @@ -163,7 +170,7 @@ int main(int argc, char **argv) {
// ---------------------------------------------------------------------------
// Print output results
// ---------------------------------------------------------------------------
PetscCall( PrintOutput(comm, ceed, mem_type_backend,
PetscCall( PrintOutput(ceed, mem_type_backend,
snes, ksp, U, l2_error_u, l2_error_p, app_ctx) );

// ---------------------------------------------------------------------------
Expand Down Expand Up @@ -201,6 +208,7 @@ int main(int argc, char **argv) {

// Free libCEED objects
CeedVectorDestroy(&force_ceed);
CeedVectorDestroy(&bc_pressure);
CeedVectorDestroy(&target);
PetscCall( CeedDataDestroy(ceed_data) );
CeedDestroy(&ceed);
Expand Down
27 changes: 13 additions & 14 deletions examples/Hdiv-mixed/problems/darcy2d.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,15 +17,14 @@
/// @file
/// Utility functions for setting up Darcy problem in 2D

#include "../include/setup-libceed.h"
#include "../include/register-problem.h"
#include "../qfunctions/darcy-rhs2d.h"
#include "../qfunctions/darcy-force2d.h"
#include "../qfunctions/darcy-mass2d.h"
#include "../qfunctions/darcy-error2d.h"
#include "../qfunctions/face-geo2d.h"
#include "../qfunctions/pressure-boundary2d.h"

// Hdiv_DARCY2D is registered in cl-option.c
PetscErrorCode Hdiv_DARCY2D(ProblemData *problem_data, void *ctx) {
PetscErrorCode Hdiv_DARCY2D(ProblemData problem_data, void *ctx) {
Physics phys = *(Physics *)ctx;
MPI_Comm comm = PETSC_COMM_WORLD;
PetscFunctionBeginUser;
Expand All @@ -39,16 +38,16 @@ PetscErrorCode Hdiv_DARCY2D(ProblemData *problem_data, void *ctx) {
problem_data->elem_node = 4;
problem_data->q_data_size_face = 3;
problem_data->quadrature_mode = CEED_GAUSS;
problem_data->setup_rhs = SetupDarcyRhs2D;
problem_data->setup_rhs_loc = SetupDarcyRhs2D_loc;
problem_data->residual = SetupDarcyMass2D;
problem_data->residual_loc = SetupDarcyMass2D_loc;
problem_data->jacobian = SetupJacobianDarcyMass2D;
problem_data->jacobian_loc = SetupJacobianDarcyMass2D_loc;
problem_data->setup_error = SetupDarcyError2D;
problem_data->setup_error_loc = SetupDarcyError2D_loc;
problem_data->setup_face_geo = SetupFaceGeo2D;
problem_data->setup_face_geo_loc = SetupFaceGeo2D_loc;
problem_data->force = DarcyForce2D;
problem_data->force_loc = DarcyForce2D_loc;
problem_data->residual = DarcyMass2D;
problem_data->residual_loc = DarcyMass2D_loc;
problem_data->jacobian = JacobianDarcyMass2D;
problem_data->jacobian_loc = JacobianDarcyMass2D_loc;
problem_data->error = DarcyError2D;
problem_data->error_loc = DarcyError2D_loc;
problem_data->bc_pressure = BCPressure2D;
problem_data->bc_pressure_loc = BCPressure2D_loc;

// ------------------------------------------------------
// Command line Options
Expand Down
Loading

0 comments on commit 2a70772

Please sign in to comment.