42 #include "libmesh/exodusII_io.h" 43 #include "libmesh/libmesh.h" 44 #include "libmesh/mesh.h" 45 #include "libmesh/mesh_generation.h" 46 #include "libmesh/linear_implicit_system.h" 47 #include "libmesh/equation_systems.h" 48 #include "libmesh/enum_xdr_mode.h" 49 #include "libmesh/getpot.h" 50 #include "libmesh/mesh_refinement.h" 53 #include "libmesh/fe.h" 54 #include "libmesh/inf_fe.h" 55 #include "libmesh/inf_elem_builder.h" 58 #include "libmesh/quadrature_gauss.h" 62 #include "libmesh/sparse_matrix.h" 63 #include "libmesh/numeric_vector.h" 64 #include "libmesh/dense_matrix.h" 65 #include "libmesh/dense_vector.h" 69 #include "libmesh/dof_map.h" 72 #include "libmesh/node.h" 75 #include "libmesh/elem.h" 83 const std::string & system_name);
86 int main (
int argc,
char ** argv)
92 #ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS 93 libmesh_example_requires(
false,
"--enable-ifem");
97 libmesh_example_requires(3 <= LIBMESH_DIM,
"3D support");
104 GetPot input(argc, argv);
106 const unsigned int nx = input(
"nx", 4),
123 #ifdef LIBMESH_HAVE_EXODUS_API 147 for (
auto & elem :
mesh.element_ptr_range())
148 if (elem->infinite())
149 elem->subdomain_id() = 1;
156 #ifdef LIBMESH_HAVE_EXODUS_API 187 equation_systems.
get_system(
"Wave").add_variable(
"p", fe_type);
200 equation_systems.
init();
202 #ifdef LIBMESH_ENABLE_AMR 204 const unsigned int nr = input(
"nr", 0);
209 equation_systems.
reinit();
227 equation_systems.
write (
"eqn_sys.dat",
WRITE);
234 #endif // else part of ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS 240 const std::string & libmesh_dbg_var(system_name))
244 libmesh_assert_equal_to (system_name,
"Wave");
250 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 271 const FEType & fe_type = dof_map.variable_type(0);
285 fe->attach_quadrature_rule (&qrule);
292 inf_fe->attach_quadrature_rule (&qrule);
311 std::vector<dof_id_type> dof_indices;
319 for (
const auto & elem :
mesh.active_local_element_ptr_range())
325 dof_map.dof_indices (elem, dof_indices);
327 const unsigned int n_dofs =
328 cast_int<unsigned int>(dof_indices.size());
346 if (elem->infinite())
390 const std::vector<Real> & JxW = cfe->get_JxWxdecay_sq();
393 const std::vector<std::vector<Real>> & phi = cfe->get_phi_over_decayxR();
398 const std::vector<std::vector<RealGradient>> & dphi = cfe->get_dphi_over_decayxR();
409 const std::vector<RealGradient> & dphase = cfe->get_dphase();
410 const std::vector<Real> &
weight = cfe->get_Sobolev_weightxR_sq();
411 const std::vector<RealGradient> & dweight = cfe->get_Sobolev_dweightxR_sq();
422 Ke.
resize (n_dofs, n_dofs);
423 Ce.
resize (n_dofs, n_dofs);
424 Me.
resize (n_dofs, n_dofs);
432 unsigned int max_qp = cfe->n_quadrature_points();
435 for (
unsigned int qp=0; qp<max_qp; qp++)
441 const unsigned int n_sf =
459 for (
unsigned int i=0; i<n_sf; i++)
460 for (
unsigned int j=0; j<n_sf; j++)
465 (dweight[qp] * phi[i][qp]
474 (dphase[qp] * dphi[j][qp])
477 (dweight[qp] * dphase[qp])
478 * phi[i][qp] * phi[j][qp]
480 (dphi[i][qp] * dphase[qp])
487 (1. - (dphase[qp] * dphase[qp]))
488 * phi[i][qp] * phi[j][qp] *
weight[qp]
497 Ke.
add(1./speed , Ce);
498 Ke.
add(1./(speed*speed), Me);
502 dof_map.constrain_element_matrix(Ke, dof_indices);
509 for (
const auto & node :
mesh.local_node_ptr_range())
515 unsigned int dn = node->dof_number(0,0,0);
525 #endif //ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
This is the EquationSystems class.
int main(int argc, char **argv)
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
void write(std::string_view name, const XdrMODE, const unsigned int write_flags=(WRITE_DATA), bool partition_agnostic=true) const
Write the systems to disk using the XDR data format.
static constexpr Real TOLERANCE
The ExodusII_IO class implements reading meshes in the ExodusII file format from Sandia National Labs...
Manages consistently variables, degrees of freedom, coefficient vectors, matrices and linear solvers ...
void resize(const unsigned int n)
Resize the vector.
virtual void add_vector(const T *v, const std::vector< numeric_index_type > &dof_indices)
Computes , where v is a pointer and each dof_indices[i] specifies where to add value v[i]...
void print_info(std::ostream &os=libMesh::out) const
Prints information about the equation systems, by default to libMesh::out.
NumericVector< Number > * rhs
The system matrix.
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
The libMesh namespace provides an interface to certain functionality in the library.
const T_sys & get_system(std::string_view name) const
This is the MeshBase class.
This class handles the numbering of degrees of freedom on a mesh.
void libmesh_ignore(const Args &...)
virtual void add_matrix(const DenseMatrix< T > &dm, const std::vector< numeric_index_type > &rows, const std::vector< numeric_index_type > &cols)=0
Add the full matrix dm to the SparseMatrix.
virtual void find_neighbors(const bool reset_remote_elements=false, const bool reset_current_list=true)=0
Locate element face (edge in 2D) neighbors.
const T & get(std::string_view) const
Implements (adaptive) mesh refinement algorithms for a MeshBase.
void print_info(std::ostream &os=libMesh::out, const unsigned int verbosity=0, const bool global=true) const
Prints relevant information about the mesh.
virtual void reinit()
Handle any mesh changes and reinitialize all the systems on the updated mesh.
boostcopy::enable_if_c< ScalarTraits< T2 >::value, void >::type add(const T2 factor, const DenseMatrix< T3 > &mat)
Adds factor times mat to this matrix.
void init(triangulateio &t)
Initializes the fields of t to nullptr/0 as necessary.
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
This class is used to build infinite elements on top of an existing mesh.
void assemble_wave(EquationSystems &es, const std::string &system_name)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Point build_inf_elem(const bool be_verbose=false)
Build infinite elements atop a volume-based mesh, determine origin automatically. ...
virtual void write(const std::string &fname) override
This method implements writing a mesh to a specified file.
T & set(const std::string &)
const MeshBase & get_mesh() const
void resize(const unsigned int new_m, const unsigned int new_n)
Resizes the matrix to the specified size and calls zero().
This class implements specific orders of Gauss quadrature.
unsigned int mesh_dimension() const
Parameters parameters
Data structure holding arbitrary parameters.
static std::unique_ptr< FEGenericBase > build_InfFE(const unsigned int dim, const FEType &type)
Builds a specific infinite element type.
virtual void init()
Initialize all the systems.
virtual void add(const numeric_index_type i, const T value)=0
Adds value to the vector entry specified by i.
virtual System & add_system(std::string_view system_type, std::string_view name)
Add the system of type system_type named name to the systems array.
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
const DofMap & get_dof_map() const
const SparseMatrix< Number > & get_system_matrix() const
This class forms the foundation from which generic finite elements may be derived.
void uniformly_refine(unsigned int n=1)
Uniformly refines the mesh n times.