29 #include "libmesh/dof_map.h" 30 #include "libmesh/dirichlet_boundaries.h" 31 #include "libmesh/fe_interface.h" 32 #include "libmesh/dirichlet_boundaries.h" 33 #include "libmesh/elem.h" 48 std::pair<MAST::SideBCMapType::const_iterator, MAST::SideBCMapType::const_iterator> it =
51 for ( ; it.first != it.second; it.first++)
52 libmesh_assert(it.first->second != &load);
55 _side_bc_map.insert(MAST::SideBCMapType::value_type(bid, &load));
66 _dirichlet_bc_map.insert(MAST::DirichletBCMapType::value_type(bid, &load)).second;
68 libmesh_assert(insert_success);
76 const unsigned int var) {
78 std::map<libMesh::subdomain_id_type, std::vector<unsigned int> >::iterator
83 (std::pair<libMesh::subdomain_id_type, std::vector<unsigned int> >
84 (sid, std::vector<unsigned int>())).first;
87 it->second.push_back(var);
95 std::pair<MAST::VolumeBCMapType::iterator, MAST::VolumeBCMapType::iterator> it =
98 for ( ; it.first != it.second; it.first++)
99 libmesh_assert(it.first->second != &load);
101 _vol_bc_map.insert(MAST::VolumeBCMapType::value_type(sid, &load));
120 std::pair<MAST::VolumeBCMapType::iterator, MAST::VolumeBCMapType::iterator> it =
123 for ( ; it.first != it.second; it.first++)
124 if (it.first->second == &load) {
130 libmesh_assert(
false);
140 MAST::PropertyCardMapType::const_iterator elem_p_it =
_element_property.find(sid);
151 MAST::PropertyCardMapType::const_iterator
155 return *elem_p_it->second;
163 MAST::PropertyCardMapType::const_iterator
167 return *elem_p_it->second;
174 MAST::PropertyCardMapType::const_iterator
178 return *elem_p_it->second;
193 sys.get_dof_map().add_dirichlet_boundary(it->second->dirichlet_boundary());
208 sys.get_dof_map().remove_dirichlet_boundary(it->second->dirichlet_boundary());
217 std::set<unsigned int>& dof_ids)
const {
223 std::map<libMesh::boundary_id_type, std::vector<unsigned int> >
224 constrained_vars_map;
232 libMesh::DirichletBoundary& dirichlet_b = it->second->dirichlet_boundary();
233 sys.get_dof_map().add_dirichlet_boundary(dirichlet_b);
234 constrained_vars_map[it->first] = dirichlet_b.variables;
242 const libMesh::MeshBase& mesh = sys.get_mesh();
245 const unsigned int dim = mesh.mesh_dimension();
249 libMesh::FEType fe_type = sys.get_dof_map().variable_type(0);
251 const libMesh::DofMap& dof_map = sys.get_dof_map();
253 libMesh::MeshBase::const_element_iterator
254 el = mesh.active_local_elements_begin(),
255 end_el = mesh.active_local_elements_end();
257 for ( ; el != end_el; ++el) {
258 const libMesh::Elem* elem = *el;
261 std::map<libMesh::subdomain_id_type, std::vector<unsigned int> >::const_iterator
267 for (
unsigned int i=0; i<it->second.size(); i++) {
269 std::vector<libMesh::dof_id_type> dof_indices;
270 dof_map.dof_indices (*el, dof_indices, it->second[i]);
272 for(
unsigned int ii=0; ii<dof_indices.size(); ii++)
273 dof_ids.insert(dof_indices[ii]);
280 for (
unsigned int s=0; s<elem->n_sides(); s++)
281 if ((*el)->neighbor_ptr(s) ==
nullptr &&
282 mesh.boundary_info->n_boundary_ids(elem, s)) {
284 std::vector<libMesh::boundary_id_type> bc_ids;
285 mesh.boundary_info->boundary_ids(elem, s, bc_ids);
287 for (
unsigned int i_bid=0; i_bid<bc_ids.size(); i_bid++)
288 if (constrained_vars_map.count(bc_ids[i_bid])) {
290 const std::vector<unsigned int>&
291 vars = constrained_vars_map[bc_ids[i_bid]];
295 for (
unsigned int i_var=0; i_var<vars.size(); i_var++) {
297 std::vector<libMesh::dof_id_type> dof_indices;
298 dof_map.dof_indices (*el, dof_indices, vars[i_var]);
301 std::vector<unsigned int> side_dofs;
302 libMesh::FEInterface::dofs_on_side(*el,
308 for(
unsigned int ii=0; ii<side_dofs.size(); ii++)
309 dof_ids.insert(dof_indices[side_dofs[ii]]);
318 libMesh::DofConstraints::const_iterator
319 constraint_it = dof_map.constraint_rows_begin(),
320 constraint_end = dof_map.constraint_rows_end();
322 for ( ; constraint_it != constraint_end; constraint_it++) {
325 if (!constraint_it->second.size())
326 dof_ids.insert(constraint_it->first);
void add_volume_load(libMesh::subdomain_id_type bid, MAST::BoundaryConditionBase &load)
adds the specified volume loads for the elements with subdomain tag s_id
void clear_volume_load(libMesh::subdomain_id_type sid, MAST::BoundaryConditionBase &load)
clear the specified volume load from the applied loads
This class implements a system for solution of nonlinear systems.
void get_system_dirichlet_bc_dofs(libMesh::System &sys, std::set< unsigned int > &dof_ids) const
Prepares a list of the constrained dofs for system sys and returns in dof_ids.
MAST::VolumeBCMapType _vol_bc_map
volume boundary condition map of boundary id and load
void init_system_dirichlet_bc(MAST::NonlinearSystem &sys) const
initializes the system for dirichlet boundary conditions
void add_dirichlet_bc(libMesh::boundary_id_type bid, MAST::DirichletBoundaryCondition &load)
adds the specified Dirichlet boundary condition for the boundary with tag b_id
void add_side_load(libMesh::boundary_id_type bid, MAST::BoundaryConditionBase &load)
adds the specified side loads for the boudnary with tag b_id
virtual const libMesh::Elem & get_reference_elem() const
void clear_loads()
clear the loads and pointer to static solution system for this structural model
std::map< libMesh::subdomain_id_type, std::vector< unsigned int > > _subdomain_var_constraint
variables constrained on subdomain
This class allows for the specification of load associated with specified nodes in a user-provided se...
void constrain_subdomain_dofs_for_var(const libMesh::subdomain_id_type sid, const unsigned int var)
constrain dofs on a subdomain to zero
const MAST::ElementPropertyCardBase & get_property_card(const libMesh::Elem &elem) const
get property card for the specified element
This class acts as a wrapper around libMesh::Elem for the purpose of providing a uniform interface fo...
void add_point_load(MAST::PointLoadCondition &load)
adds the specified point load
MAST::DirichletBCMapType _dirichlet_bc_map
Dirichlet boundary condition map of boundary id and load.
void set_property_for_subdomain(const libMesh::subdomain_id_type sid, const MAST::ElementPropertyCardBase &prop)
sets the same property for all elements in the specified subdomain
MAST::PointLoadSetType _point_loads
point loads
MAST::PropertyCardMapType _element_property
map of element property cards for each element
MAST::SideBCMapType _side_bc_map
side boundary condition map of boundary id and load
void clear_system_dirichlet_bc(MAST::NonlinearSystem &sys) const
clears the system dirichlet boundary conditions