MAST
physics_discipline_base.cpp
Go to the documentation of this file.
1 /*
2  * MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
3  * Copyright (C) 2013-2019 Manav Bhatia
4  *
5  * This library is free software; you can redistribute it and/or
6  * modify it under the terms of the GNU Lesser General Public
7  * License as published by the Free Software Foundation; either
8  * version 2.1 of the License, or (at your option) any later version.
9  *
10  * This library is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13  * Lesser General Public License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with this library; if not, write to the Free Software
17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18  */
19 
20 // MAST includes
23 #include "base/parameter.h"
24 #include "base/nonlinear_system.h"
26 #include "mesh/geom_elem.h"
27 
28 // libMesh includes
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"
34 
35 
36 void
38  _side_bc_map.clear();
39  _vol_bc_map.clear();
40 }
41 
42 
43 
44 void
45 MAST::PhysicsDisciplineBase::add_side_load(libMesh::boundary_id_type bid,
47  // make sure that this boundary and load haven't already been applied
48  std::pair<MAST::SideBCMapType::const_iterator, MAST::SideBCMapType::const_iterator> it =
49  _side_bc_map.equal_range(bid);
50 
51  for ( ; it.first != it.second; it.first++)
52  libmesh_assert(it.first->second != &load);
53 
54  // displacement boundary condition needs to be hadled separately
55  _side_bc_map.insert(MAST::SideBCMapType::value_type(bid, &load));
56 }
57 
58 
59 
60 
61 void
62 MAST::PhysicsDisciplineBase::add_dirichlet_bc(libMesh::boundary_id_type bid,
64 
65  bool insert_success =
66  _dirichlet_bc_map.insert(MAST::DirichletBCMapType::value_type(bid, &load)).second;
67 
68  libmesh_assert(insert_success);
69 }
70 
71 
72 
73 void
75 constrain_subdomain_dofs_for_var(const libMesh::subdomain_id_type sid,
76  const unsigned int var) {
77 
78  std::map<libMesh::subdomain_id_type, std::vector<unsigned int> >::iterator
79  it = _subdomain_var_constraint.find(sid);
80 
81  if (it == _subdomain_var_constraint.end()) {
82  it = _subdomain_var_constraint.insert
83  (std::pair<libMesh::subdomain_id_type, std::vector<unsigned int> >
84  (sid, std::vector<unsigned int>())).first;
85  }
86 
87  it->second.push_back(var);
88 }
89 
90 
91 
92 void
93 MAST::PhysicsDisciplineBase::add_volume_load(libMesh::subdomain_id_type sid,
95  std::pair<MAST::VolumeBCMapType::iterator, MAST::VolumeBCMapType::iterator> it =
96  _vol_bc_map.equal_range(sid);
97 
98  for ( ; it.first != it.second; it.first++)
99  libmesh_assert(it.first->second != &load);
100 
101  _vol_bc_map.insert(MAST::VolumeBCMapType::value_type(sid, &load));
102 }
103 
104 
105 
106 
107 void
109 
110  libmesh_assert(!_point_loads.count(&load));
111 
112  _point_loads.insert(&load);
113 }
114 
115 
116 
117 void
118 MAST::PhysicsDisciplineBase::clear_volume_load(libMesh::subdomain_id_type bid,
120  std::pair<MAST::VolumeBCMapType::iterator, MAST::VolumeBCMapType::iterator> it =
121  _vol_bc_map.equal_range(bid);
122 
123  for ( ; it.first != it.second; it.first++)
124  if (it.first->second == &load) {
125  _vol_bc_map.erase(it.first);
126  return;
127  }
128 
129  // should not get here
130  libmesh_assert(false);
131 }
132 
133 
134 
135 void
137 set_property_for_subdomain(const libMesh::subdomain_id_type sid,
138  const MAST::ElementPropertyCardBase& prop) {
139 
140  MAST::PropertyCardMapType::const_iterator elem_p_it = _element_property.find(sid);
141  libmesh_assert(elem_p_it == _element_property.end());
142 
143  _element_property[sid] = &prop;
144 }
145 
146 
147 
149 MAST::PhysicsDisciplineBase::get_property_card(const unsigned int sid) const {
150 
151  MAST::PropertyCardMapType::const_iterator
152  elem_p_it = _element_property.find(sid);
153  libmesh_assert(elem_p_it != _element_property.end());
154 
155  return *elem_p_it->second;
156 }
157 
158 
159 
161 MAST::PhysicsDisciplineBase::get_property_card(const libMesh::Elem& elem) const {
162 
163  MAST::PropertyCardMapType::const_iterator
164  elem_p_it = _element_property.find(elem.subdomain_id());
165  libmesh_assert(elem_p_it != _element_property.end());
166 
167  return *elem_p_it->second;
168 }
169 
170 
173 
174  MAST::PropertyCardMapType::const_iterator
175  elem_p_it = _element_property.find(elem.get_reference_elem().subdomain_id());
176  libmesh_assert(elem_p_it != _element_property.end());
177 
178  return *elem_p_it->second;
179 }
180 
181 
182 
183 void
186 
187 
188  // iterate over all the dirichlet boundary conditions and add them
189  // to the system
190  MAST::DirichletBCMapType::const_iterator it = _dirichlet_bc_map.begin();
191 
192  for ( ; it != _dirichlet_bc_map.end(); it++)
193  sys.get_dof_map().add_dirichlet_boundary(it->second->dirichlet_boundary());
194 }
195 
196 
197 
198 
199 void
202 
203  // iterate over all the dirichlet boundary conditions and add them
204  // to the system
205  MAST::DirichletBCMapType::const_iterator it = _dirichlet_bc_map.begin();
206 
207  for ( ; it != _dirichlet_bc_map.end(); it++)
208  sys.get_dof_map().remove_dirichlet_boundary(it->second->dirichlet_boundary());
209 
210 }
211 
212 
213 
214 void
216 get_system_dirichlet_bc_dofs(libMesh::System& sys,
217  std::set<unsigned int>& dof_ids) const {
218 
219  dof_ids.clear();
220 
221  // first prepare a map of boundary ids and the constrained vars on that
222  // boundary
223  std::map<libMesh::boundary_id_type, std::vector<unsigned int> >
224  constrained_vars_map;
225 
226  // iterate over all the dirichlet boundary conditions and add them
227  // to the system
228  MAST::DirichletBCMapType::const_iterator it = _dirichlet_bc_map.begin();
229 
230  for ( ; it != _dirichlet_bc_map.end(); it++) {
231 
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;
235  }
236 
237 
238  //
239  // now collect the ids that correspond to the specified boundary conditions
240  //
241  // Get a constant reference to the mesh object
242  const libMesh::MeshBase& mesh = sys.get_mesh();
243 
244  // The dimension that we are running.
245  const unsigned int dim = mesh.mesh_dimension();
246 
247  // Get a constant reference to the Finite Element type
248  // for the first variable in the system.
249  libMesh::FEType fe_type = sys.get_dof_map().variable_type(0);
250 
251  const libMesh::DofMap& dof_map = sys.get_dof_map();
252 
253  libMesh::MeshBase::const_element_iterator
254  el = mesh.active_local_elements_begin(),
255  end_el = mesh.active_local_elements_end();
256 
257  for ( ; el != end_el; ++el) {
258  const libMesh::Elem* elem = *el;
259 
260  // check if the any subdomain variables have been constrained
261  std::map<libMesh::subdomain_id_type, std::vector<unsigned int> >::const_iterator
262  it = _subdomain_var_constraint.find(elem->subdomain_id());
263 
264  if (it != _subdomain_var_constraint.end()) {
265  // constrain all element dofs on this element that belong to this
266  // variable
267  for (unsigned int i=0; i<it->second.size(); i++) {
268 
269  std::vector<libMesh::dof_id_type> dof_indices;
270  dof_map.dof_indices (*el, dof_indices, it->second[i]);
271 
272  for(unsigned int ii=0; ii<dof_indices.size(); ii++)
273  dof_ids.insert(dof_indices[ii]);
274  }
275  }
276 
277 
278  // boundary condition is applied only on sides with no neighbors
279  // and if the side's boundary id has a boundary condition tag on it
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)) {
283 
284  std::vector<libMesh::boundary_id_type> bc_ids;
285  mesh.boundary_info->boundary_ids(elem, s, bc_ids);
286 
287  for (unsigned int i_bid=0; i_bid<bc_ids.size(); i_bid++)
288  if (constrained_vars_map.count(bc_ids[i_bid])) {
289 
290  const std::vector<unsigned int>&
291  vars = constrained_vars_map[bc_ids[i_bid]];
292 
293  // now iterate over each constrained variable for this boundary
294  // and collect its dofs
295  for (unsigned int i_var=0; i_var<vars.size(); i_var++) {
296 
297  std::vector<libMesh::dof_id_type> dof_indices;
298  dof_map.dof_indices (*el, dof_indices, vars[i_var]);
299 
300  // all boundary dofs are Dirichlet dofs in this case
301  std::vector<unsigned int> side_dofs;
302  libMesh::FEInterface::dofs_on_side(*el,
303  dim,
304  fe_type,
305  s,
306  side_dofs);
307 
308  for(unsigned int ii=0; ii<side_dofs.size(); ii++)
309  dof_ids.insert(dof_indices[side_dofs[ii]]);
310  }
311  } // end of boundary loop
312  } // end of side loop
313  }// end of element loop
314 
315  // also, it is likely that some of the bcs have been applied via the
316  // DofMap API by specification of row constraints. In that case,
317  // factor out the dofs that are not coupled to any other dofs
318  libMesh::DofConstraints::const_iterator
319  constraint_it = dof_map.constraint_rows_begin(),
320  constraint_end = dof_map.constraint_rows_end();
321 
322  for ( ; constraint_it != constraint_end; constraint_it++) {
323  // if the dof constraint has only one entry, then add it to the
324  // constrained set
325  if (!constraint_it->second.size())
326  dof_ids.insert(constraint_it->first);
327  }
328 }
329 
330 
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
Definition: geom_elem.cpp:54
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...
Definition: geom_elem.h:59
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