33 #include "libmesh/getpot.h"
34 #include "libmesh/elem.h"
35 #include "libmesh/fe_interface.h"
41 const std::string& physics_name )
43 _pin_value = input(
"Physics/"+physics_name+
"/pin_value", 0.0 );
45 unsigned int pin_loc_dim = input.vector_variable_size(
"Physics/"+physics_name+
"/pin_location");
48 if( pin_loc_dim > 0 && pin_loc_dim < 2 )
50 std::cerr <<
"Error: pressure pin location must be at least 2 dimensional"
55 _pin_location(0) = input(
"Physics/"+physics_name+
"/pin_location", 0.0, 0 );
56 _pin_location(1) = input(
"Physics/"+physics_name+
"/pin_location", 0.0, 1 );
58 if( pin_loc_dim == 3 )
59 _pin_location(2) = input(
"Physics/"+physics_name+
"/pin_location", 0.0, 2 );
70 const bool request_jacobian,
72 const double penalty )
75 AssemblyContext &c = libMesh::libmesh_cast_ref<AssemblyContext&>(context);
79 libMesh::DenseSubVector<libMesh::Number> &F_var = c.get_elem_residual(var);
80 libMesh::DenseSubMatrix<libMesh::Number> &K_var = c.get_elem_jacobian(var, var);
83 const unsigned int n_var_dofs = c.get_dof_indices(var).size();
85 libMesh::Number var_value = c.point_value(var,
_pin_location);
87 libMesh::FEType fe_type = c.get_element_fe(var)->get_fe_type();
89 libMesh::Point point_loc_in_masterelem =
90 libMesh::FEInterface::inverse_map(c.get_dim(), fe_type, &c.get_elem(),
_pin_location);
92 std::vector<libMesh::Real> phi(n_var_dofs);
94 for (
unsigned int i=0; i != n_var_dofs; i++)
95 phi[i] = libMesh::FEInterface::shape( c.get_dim(), fe_type, &c.get_elem(), i,
96 point_loc_in_masterelem );
98 for (
unsigned int i=0; i != n_var_dofs; i++)
100 F_var(i) += penalty*(var_value -
_pin_value)*phi[i];
103 if (request_jacobian && c.get_elem_solution_derivative())
105 libmesh_assert (c.get_elem_solution_derivative() == 1.0);
107 for (
unsigned int j=0; j != n_var_dofs; j++)
108 K_var(i,j) += penalty*phi[i]*phi[j];
unsigned int VariableIndex
More descriptive name of the type used for variable indices.
libMesh::Point _pin_location
Location we want to pin the pressure.
void pin_value(libMesh::DiffContext &context, const bool request_jacobian, const GRINS::VariableIndex var, const double penalty=1.0)
The idea here is to pin a variable to a particular value if there is a null space - e...
libMesh::Number _pin_value
Value of pressure we wish to pin.
PressurePinning(const GetPot &input, const std::string &physics_name)