Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Constraints work #3850

Merged
merged 4 commits into from
May 13, 2024
Merged

Constraints work #3850

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
df6fca9
CLI option --print-constraints
roystgnr Oct 5, 2023
5937797
Correct and clarify _constraint_rows comment
roystgnr Oct 5, 2023
e2db391
More Bezier Extraction comments
roystgnr Oct 5, 2023
0a52283
Simplify IGA constraint row construction
roystgnr Oct 5, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 4 additions & 4 deletions include/mesh/mesh_base.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2000,11 +2000,11 @@ class MeshBase : public ParallelObject
// mesh, such as FE nodes whose Rational Bernstein values need to be
// constrained in terms of values on spline control nodes.
//
// _constraint_rows[constrained_node_id][i].first.first is an
// element,
// _constraint_rows[constrained_node_id][i].first.second is the
// _constraint_rows[constrained_node][i].first.first is an
// element (e.g. a NodeElem for a spline control node),
// _constraint_rows[constrained_node][i].first.second is the
// local node id of that element which is a constraining node,
// _constraint_rows[constrained_node_id][i].second is that node's
// _constraint_rows[constrained_node][i].second is that node's
// constraint coefficient.
constraint_rows_type _constraint_rows;

Expand Down
72 changes: 36 additions & 36 deletions src/mesh/exodusII_io.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -549,28 +549,42 @@ void ExodusII_IO::read (const std::string & fname)

}

// The tailing entries in each element's connectivity
// vector are indices to Exodus constraint coefficient
// rows.

// Concatenating these rows gives a matrix with
// all the constraints for the element nodes: each
// column of that matrix is the constraint coefficients
// for the node associated with that column (via the
// Exodus numbering, not the libMesh numbering).
const auto & my_constraint_rows = exio_helper->bex_cv_conn[elem_num];

for (auto elem_node_index :
make_range(elem->n_nodes()))
{
const auto & my_constraint_rows = exio_helper->bex_cv_conn[elem_num];

// New finite element node data: dot product of
// New finite element node data = dot product of
// constraint matrix columns with spline node data.
// Store that column as a key.
// Store each column's non-zero entries, along with
// the global spline node indices, as a key to
// identify shared finite element nodes.
std::vector<std::pair<dof_id_type, Real>> key;

for (auto spline_node_index :
make_range(exio_helper->bex_num_elem_cvs))
{
// Pick out a row of the element constraint matrix
const unsigned long elem_coef_vec_index =
my_constraint_rows[spline_node_index] - 1; // Exodus isn't 0-based

auto & coef_vec =
bex_constraint_vec(elem_coef_vec_index, *exio_helper);

// Get coef from this node's column intersect that row
const Real coef =
libmesh_vector_at(coef_vec, elem_node_index);

// Get the libMesh node corresponding to that row
const int gi = (elem_num)*exio_helper->bex_num_elem_cvs +
spline_node_index;
const dof_id_type libmesh_node_id =
Expand All @@ -580,47 +594,31 @@ void ExodusII_IO::read (const std::string & fname)
key.emplace_back(libmesh_node_id, coef);
}

// Have we already created this node? Connect it.
if (const auto local_node_it = local_nodes.find(key);
local_node_it != local_nodes.end())
elem->set_node(dyna_elem_defn.nodes[elem_node_index]) = local_node_it->second;
// Have we not yet created this node? Construct it,
// along with its weight and libMesh constraint row,
// then connect it.
else
{
Point p(0);
Real w = 0;
std::vector<std::pair<std::pair<const Elem *, unsigned int>, Real>> constraint_row;

for (auto spline_node_index :
make_range(exio_helper->bex_num_elem_cvs))
for (auto [libmesh_spline_node_id, coef] : key)
{
// global => libMesh index, with crazy 1-based data - see comments above
const int gi = (elem_num)*exio_helper->bex_num_elem_cvs +
spline_node_index;
const dof_id_type libmesh_node_id =
exio_helper->node_num_map[exio_helper->connect[gi] - 1] - 1;

const unsigned long elem_coef_vec_index =
my_constraint_rows[spline_node_index] - 1; // Exodus isn't 0-based

const Node & spline_node = mesh.node_ref(libmesh_node_id);

auto & coef_vec =
bex_constraint_vec(elem_coef_vec_index, *exio_helper);
const Real coef =
libmesh_vector_at(coef_vec, elem_node_index);

// We don't need to store 0 entries;
// constraint_rows is a sparse structure.
if (coef)
{
p.add_scaled(spline_node, coef);
const Real spline_w = weights_exist ?
spline_node.get_extra_datum<Real>(weight_index) : 1;
w += coef * spline_w;

const Elem * nodeelem =
libmesh_map_find(spline_nodeelem_ptrs, &spline_node);
constraint_row.emplace_back(std::make_pair(nodeelem, 0), coef);
}
const Node & spline_node = mesh.node_ref(libmesh_spline_node_id);

p.add_scaled(spline_node, coef);
const Real spline_w = weights_exist ?
spline_node.get_extra_datum<Real>(weight_index) : 1;
w += coef * spline_w;

const Elem * nodeelem =
libmesh_map_find(spline_nodeelem_ptrs, &spline_node);
constraint_row.emplace_back(std::make_pair(nodeelem, 0), coef);
}

Node *n = mesh.add_point(p, n_nodes++);
Expand All @@ -629,7 +627,9 @@ void ExodusII_IO::read (const std::string & fname)

// If we're building disconnected Bezier
// extraction elements then we don't want to
// find the new nodes to reuse later.
// find the new nodes to reuse later; each
// finite element node will connect to only one
// element.
if (!_disc_bex)
local_nodes[key] = n;
elem->set_node(dyna_elem_defn.nodes[elem_node_index]) = n;
Expand Down
2 changes: 2 additions & 0 deletions src/systems/system.C
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -485,6 +485,8 @@ void System::reinit_constraints()
get_dof_map().create_dof_constraints(_mesh, this->time);
user_constrain();
get_dof_map().process_constraints(_mesh);
if (libMesh::on_command_line ("--print-constraints"))
get_dof_map().print_dof_constraints(libMesh::out);
#endif
get_dof_map().prepare_send_list();
}
Expand Down