Skip to content

Commit

Permalink
whack-a-mole
Browse files Browse the repository at this point in the history
  • Loading branch information
@landinjm
landinjm committed Nov 26, 2024
1 parent 981e8dc commit 54c08c3
Show file tree
Hide file tree
Showing 5 changed files with 145 additions and 128 deletions.
3 changes: 2 additions & 1 deletion include/matrixFreePDE.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@
#define MATRIXFREEPDE_H

// general headers
#include <filesystem>
#include <fstream>
#include <iterator>
#include <sstream>
Expand Down Expand Up @@ -481,7 +482,7 @@ class MatrixFreePDE : public Subscriptor
void
computeIntegral(double &integratedField,
int index,
std::vector<vectorType *> postProcessedSet);
std::vector<vectorType *> variableSet);

// variables for time dependent problems
/*Flag used to see if invM, time stepping in run(), etc are necessary*/
Expand Down
207 changes: 107 additions & 100 deletions src/FloodFiller/FloodFiller.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,14 +4,15 @@

template <int dim, int degree>
void
FloodFiller<dim, degree>::calcGrainSets(dealii::FESystem<dim> &finite_element,
dealii::DoFHandler<dim> &dof_handler,
vectorType *solution_field,
double threshold_lower,
double threshold_upper,
int min_id,
unsigned int order_parameter_index,
std::vector<GrainSet<dim>> &grain_sets)
FloodFiller<dim, degree>::calcGrainSets(
[[maybe_unused]] dealii::FESystem<dim> &finite_element,
dealii::DoFHandler<dim> &dof_handler,
vectorType *solution_field,
double threshold_lower,
double threshold_upper,
int min_id,
unsigned int order_parameter_index,
std::vector<GrainSet<dim>> &grain_sets)
{
unsigned int grain_index = 0;

Expand Down Expand Up @@ -76,6 +77,7 @@ FloodFiller<dim, degree>::calcGrainSets(dealii::FESystem<dim> &finite_eleme
mergeSplitGrains(grain_sets);
}

// NOLINTBEGIN(misc-no-recursion)
template <int dim, int degree>
template <typename T>
void
Expand All @@ -89,97 +91,102 @@ FloodFiller<dim, degree>::recursiveFloodFill(T cell,
std::vector<GrainSet<dim>> &grain_sets,
bool &grain_assigned)
{
if (cell != cell_end)
if (cell == cell_end)
{
// Check if the cell has been marked yet
bool cellMarked = cell->user_flag_set();
return;
}

// Check if the cell has been marked yet
if (cell->user_flag_set())
{
return;
}

if (!cellMarked)
if (cell->has_children())
{
// Call recursiveFloodFill on the element's children
for (unsigned int n_child = 0; n_child < cell->n_children(); n_child++)
{
if (cell->has_children())
{
// Call recursiveFloodFill on the element's children
for (unsigned int n_child = 0; n_child < cell->n_children(); n_child++)
{
recursiveFloodFill<T>(cell->child(n_child),
cell_end,
solution_field,
threshold_lower,
threshold_upper,
min_id,
grain_index,
grain_sets,
grain_assigned);
}
}
else
{
if (cell->is_locally_owned())
{
cell->set_user_flag();
recursiveFloodFill<T>(cell->child(n_child),
cell_end,
solution_field,
threshold_lower,
threshold_upper,
min_id,
grain_index,
grain_sets,
grain_assigned);
}
}
else
{
if (!cell->is_locally_owned())
{
return;
}

dealii::FEValues<dim> fe_values(*fe, quadrature, dealii::update_values);
cell->set_user_flag();

std::vector<double> var_values(num_quad_points);
std::vector<dealii::Point<dim>> q_point_list(num_quad_points);
dealii::FEValues<dim> fe_values(*fe, quadrature, dealii::update_values);

// Get the most common value for the element
fe_values.reinit(cell);
fe_values.get_function_values(*solution_field, var_values);
std::vector<double> var_values(num_quad_points);
std::vector<dealii::Point<dim>> q_point_list(num_quad_points);

double ele_val;
std::map<double, int> quadratureValues;
int maxNumberSeen = 0;
double mostCommonQPointValue = -1;
for (unsigned int q_point = 0; q_point < num_quad_points; ++q_point)
{
// Add the number of times that var_values[q_point] has
// been seen
if (var_values[q_point] > min_id)
{
++quadratureValues[var_values[q_point]];
}
if (quadratureValues[var_values[q_point]] > maxNumberSeen)
{
maxNumberSeen = quadratureValues[var_values[q_point]];
mostCommonQPointValue = var_values[q_point];
}
}
ele_val = mostCommonQPointValue;
// Get the most common value for the element
fe_values.reinit(cell);
fe_values.get_function_values(*solution_field, var_values);

if (ele_val > threshold_lower && ele_val < threshold_upper)
{
grain_assigned = true;
std::map<double, int> quadratureValues;
int maxNumberSeen = 0;
double mostCommonQPointValue = -1;
for (unsigned int q_point = 0; q_point < num_quad_points; ++q_point)
{
// Add the number of times that var_values[q_point] has
// been seen
if (var_values[q_point] > min_id)
{
++quadratureValues[var_values[q_point]];
}
if (quadratureValues[var_values[q_point]] > maxNumberSeen)
{
maxNumberSeen = quadratureValues[var_values[q_point]];
mostCommonQPointValue = var_values[q_point];
}
}
double ele_val = mostCommonQPointValue;

std::vector<dealii::Point<dim>> vertex_list;
for (unsigned int vertex_index = 0;
vertex_index < dealii::Utilities::fixed_power<dim>(2.0);
vertex_index++)
{
vertex_list.push_back(cell->vertex(vertex_index));
}
grain_sets.back().addVertexList(vertex_list);
if (ele_val > threshold_lower && ele_val < threshold_upper)
{
grain_assigned = true;

// Call recursiveFloodFill on the element's neighbors
for (unsigned int n_child = 0; n_child < 2 * dim; n_child++)
{
recursiveFloodFill<T>(cell->neighbor(n_child),
cell_end,
solution_field,
threshold_lower,
threshold_upper,
min_id,
grain_index,
grain_sets,
grain_assigned);
}
}
}
std::vector<dealii::Point<dim>> vertex_list;
for (unsigned int vertex_index = 0;
vertex_index < dealii::Utilities::fixed_power<dim>(2.0);
vertex_index++)
{
vertex_list.push_back(cell->vertex(vertex_index));
}
grain_sets.back().addVertexList(vertex_list);

// Call recursiveFloodFill on the element's neighbors
for (unsigned int n_child = 0; n_child < 2 * dim; n_child++)
{
recursiveFloodFill<T>(cell->neighbor(n_child),
cell_end,
solution_field,
threshold_lower,
threshold_upper,
min_id,
grain_index,
grain_sets,
grain_assigned);
}
}
}
}

// NOLINTEND(misc-no-recursion)

// =================================================================================
// All-to-all communication of the grain sets
// =================================================================================
Expand All @@ -188,7 +195,7 @@ void
FloodFiller<dim, degree>::createGlobalGrainSetList(
std::vector<GrainSet<dim>> &grain_sets) const
{
int numProcs = dealii::Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);
unsigned int numProcs = dealii::Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD);

unsigned int num_grains_local = grain_sets.size();

Expand Down Expand Up @@ -229,7 +236,7 @@ FloodFiller<dim, degree>::createGlobalGrainSetList(
MPI_Allgather(&num_grains_local,
1,
MPI_INT,
&num_grains_per_core[0],
num_grains_per_core.data(),
1,
MPI_INT,
MPI_COMM_WORLD);
Expand All @@ -246,24 +253,24 @@ FloodFiller<dim, degree>::createGlobalGrainSetList(

std::vector<unsigned int> order_parameters_global(num_grains_global, 0);

MPI_Allgatherv(&order_parameters[0],
MPI_Allgatherv(order_parameters.data(),
num_grains_local,
MPI_UNSIGNED,
&order_parameters_global[0],
&num_grains_per_core[0],
&offset[0],
order_parameters_global.data(),
num_grains_per_core.data(),
offset.data(),
MPI_UNSIGNED,
MPI_COMM_WORLD);

// Communicate the number of elements
std::vector<unsigned int> num_elements_global(num_grains_global, 0);

MPI_Allgatherv(&num_elements[0],
MPI_Allgatherv(num_elements.data(),
num_grains_local,
MPI_UNSIGNED,
&num_elements_global[0],
&num_grains_per_core[0],
&offset[0],
num_elements_global.data(),
num_grains_per_core.data(),
offset.data(),
MPI_UNSIGNED,
MPI_COMM_WORLD);

Expand All @@ -277,7 +284,7 @@ FloodFiller<dim, degree>::createGlobalGrainSetList(
std::vector<int> num_vertices_per_core;

unsigned int g = 0;
for (int i = 0; i < numProcs; i++)
for (unsigned int i = 0; i < numProcs; i++)
{
int num_vert_single_core = 0;
for (int j = 0; j < num_grains_per_core.at(i); j++)
Expand All @@ -291,17 +298,17 @@ FloodFiller<dim, degree>::createGlobalGrainSetList(
}

offset.at(0) = 0;
for (int n = 1; n < numProcs; n++)
for (unsigned int n = 1; n < numProcs; n++)
{
offset[n] = offset[n - 1] + num_vertices_per_core[n - 1];
}

MPI_Allgatherv(&vertices[0],
MPI_Allgatherv(vertices.data(),
num_vertices,
MPI_DOUBLE,
&vertices_global[0],
&num_vertices_per_core[0],
&offset[0],
vertices_global.data(),
num_vertices_per_core.data(),
offset.data(),
MPI_DOUBLE,
MPI_COMM_WORLD);

Expand Down
4 changes: 2 additions & 2 deletions src/matrixfree/AdaptiveRefinement.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -90,12 +90,12 @@ AdaptiveRefinement<dim, degree>::adaptive_refinement_criterion()
for (const auto &criterion : userInputs.refinement_criteria)
{
// Get the values and/or gradients
if (update_values & update_flags)
if (static_cast<bool>(update_values & update_flags))
{
fe_values.get_function_values(*solutionSet[criterion.variable_index],
values);
}
if (update_gradients & update_flags)
if (static_cast<bool>(update_gradients & update_flags))
{
fe_values.get_function_gradients(*solutionSet[criterion.variable_index],
gradients);
Expand Down
8 changes: 4 additions & 4 deletions src/matrixfree/boundaryConditions.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -44,9 +44,9 @@ MatrixFreePDE<dim, degree>::applyNeumannBCs()
FEFaceValues<dim> fe_face_values(*fe,
face_quadrature_formula,
update_values | update_JxW_values);
const unsigned int n_face_q_points = face_quadrature_formula.size(),
dofs_per_cell = fe->dofs_per_cell;
Vector<double> cell_rhs(dofs_per_cell);
const unsigned int n_face_q_points = face_quadrature_formula.size();
const unsigned int dofs_per_cell = fe->dofs_per_cell;
Vector<double> cell_rhs(dofs_per_cell);
std::vector<types::global_dof_index> local_dof_indices(dofs_per_cell);

// Loop over each face on a boundary
Expand Down Expand Up @@ -225,7 +225,7 @@ MatrixFreePDE<dim, degree>::setPeriodicity()
periodic_pair = true;
}
}
if (periodic_pair == true)
if (periodic_pair)
{
GridTools::collect_periodic_faces(triangulation,
/*b_id1*/ 2 * i,
Expand Down
Loading

0 comments on commit 54c08c3

Please sign in to comment.