Skip to content

Commit

Permalink
Merge pull request #10 from andrsd/dim-detect
Browse files Browse the repository at this point in the history 
Improving the procedure to detect mesh dimension
  • Loading branch information
@andrsd
andrsd authored Jan 9, 2024
2 parents 5b9ec10 + f981f40 commit 0168a69
Showing 1 changed file with 85 additions and 36 deletions.
121 changes: 85 additions & 36 deletions src/main.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -28,6 +28,25 @@ std::map<int, exodusIIcpp::NodeSet> node_sets;
// NOTE: this may be needed per block
std::map<int, int> elem_map;

struct BoundingBox {
double xmin;
double xmax;
double ymin;
double ymax;
double zmin;
double zmax;

BoundingBox() :

Check warning on line 39 in src/main.cpp

View workflow job for this annotation

GitHub Actions / c++ linter

src/main.cpp:39:5 [cppcoreguidelines-pro-type-member-init] 

constructor does not initialize these fields: xmin, xmax, ymin, ymax, zmin, zmax

Check warning on line 39 in src/main.cpp

View workflow job for this annotation

GitHub Actions / c++ linter

src/main.cpp:39:5 [modernize-use-equals-default] 

use '= default' to define a trivial default constructor
xmin(std::numeric_limits<double>::max()),
xmax(std::numeric_limits<double>::lowest()),
ymin(std::numeric_limits<double>::max()),
ymax(std::numeric_limits<double>::lowest()),
zmin(std::numeric_limits<double>::max()),
zmax(std::numeric_limits<double>::lowest())
{
}
};

namespace gmsh {
enum ElementType { EDGE2 = 1, TRI3 = 2, QUAD4 = 3, TET4 = 4, HEX8 = 5 };
}
Expand Down Expand Up @@ -92,22 +111,66 @@ read_physical_entities(const std::vector<gmshparsercpp::MshFile::PhysicalName> &
}
}

BoundingBox
compute_mesh_bounding_box()

Check warning on line 115 in src/main.cpp

View workflow job for this annotation

GitHub Actions / c++ linter

src/main.cpp:115:1 [modernize-use-trailing-return-type] 

use a trailing return type for this function
{
// this assumes that `x`, `y` and `z` contain valid coordinate data
BoundingBox bbox;
if ((x.size() == y.size()) and (y.size() == z.size())) {
auto n_nodes = x.size();
for (std::size_t i = 0; i < n_nodes; i++) {
bbox.xmin = std::min(bbox.xmin, x[i]);
bbox.xmax = std::max(bbox.xmax, x[i]);
bbox.ymin = std::min(bbox.ymin, y[i]);
bbox.ymax = std::max(bbox.ymax, y[i]);
bbox.zmin = std::min(bbox.zmin, z[i]);
bbox.zmax = std::max(bbox.zmax, z[i]);
}
return bbox;
}
else

Check warning on line 131 in src/main.cpp

View workflow job for this annotation

GitHub Actions / c++ linter

src/main.cpp:131:5 [readability-else-after-return] 

do not use 'else' after 'return'

Check warning on line 131 in src/main.cpp

View workflow job for this annotation

GitHub Actions / c++ linter

src/main.cpp:131:9 [readability-braces-around-statements] 

statement should be inside braces
throw std::logic_error(
fmt::format("Size of x ({}) must equal to size of y({}) and size of z({}).",
x.size(),
y.size(),
z.size()));
}

void
analyze_mesh(const std::vector<gmshparsercpp::MshFile::ElementBlock> & el_blks)
analyze_mesh()
{
el_blk_dim.resize(4);
for (const auto & eb : el_blks)
el_blk_dim[eb.dimension].push_back(&eb);
auto bbox = compute_mesh_bounding_box();
double x_width = bbox.xmax - bbox.xmin;
double y_width = bbox.ymax - bbox.ymin;
double z_width = bbox.zmax - bbox.zmin;

if ((std::fabs(x_width) > 0) && (std::fabs(y_width) < 1e-16) && (std::fabs(z_width) < 1e-16))

Check warning on line 147 in src/main.cpp

View workflow job for this annotation

GitHub Actions / c++ linter

src/main.cpp:147:98 [readability-braces-around-statements] 

statement should be inside braces
dim = 1;
else if ((std::fabs(x_width) > 0) && (std::fabs(y_width) > 0) && (std::fabs(z_width) < 1e-16))

Check warning on line 149 in src/main.cpp

View workflow job for this annotation

GitHub Actions / c++ linter

src/main.cpp:149:99 [readability-braces-around-statements] 

statement should be inside braces
dim = 2;
else

Check warning on line 151 in src/main.cpp

View workflow job for this annotation

GitHub Actions / c++ linter

src/main.cpp:151:9 [readability-braces-around-statements] 

statement should be inside braces
dim = 3;

// the element block with the highest dimension that has something in it will be the final mesh
// dimension
for (int i = 0; i < 4; i++)
if (!el_blk_dim[i].empty())
dim = i;
side_set_dim = dim - 1;
node_set_dim = 0;
}

const std::vector<gmshparsercpp::MshFile::MultiDEntity> *
get_entities_by_dim(const gmshparsercpp::MshFile & f, int dim)

Check warning on line 159 in src/main.cpp

View workflow job for this annotation

GitHub Actions / c++ linter

src/main.cpp:159:1 [modernize-use-trailing-return-type] 

use a trailing return type for this function

Check warning on line 159 in src/main.cpp

View workflow job for this annotation

GitHub Actions / c++ linter

src/main.cpp:159:21 [bugprone-easily-swappable-parameters] 

2 adjacent parameters of 'get_entities_by_dim' of similar type are easily swapped by mistake
{
const std::vector<gmshparsercpp::MshFile::MultiDEntity> * ents = nullptr;
switch (dim) {
case 1:
return &f.get_curve_entities();
case 2:
return &f.get_surface_entities();
case 3:
return &f.get_volume_entities();
default:
throw std::runtime_error(fmt::format("Unsupported spatial dimension {}.", dim));
}
}

void
build_coordinates(const std::vector<gmshparsercpp::MshFile::Node> & nodes)
{
Expand All @@ -127,6 +190,14 @@ build_coordinates(const std::vector<gmshparsercpp::MshFile::Node> & nodes)
}
}

void
build_element_block_dim(const std::vector<gmshparsercpp::MshFile::ElementBlock> & el_blks)
{
el_blk_dim.resize(4);
for (const auto & eb : el_blks)
el_blk_dim[eb.dimension].push_back(&eb);
}

void
build_element_blocks(const std::vector<const gmshparsercpp::MshFile::ElementBlock *> & el_blks,
const std::vector<gmshparsercpp::MshFile::MultiDEntity> & entities)
Expand Down Expand Up @@ -278,37 +349,15 @@ read_gmsh_file(const std::string & file_name)
read_physical_entities(phys_entities);

const std::vector<gmshparsercpp::MshFile::Node> & nodes = f.get_nodes();
const std::vector<gmshparsercpp::MshFile::ElementBlock> & el_blks = f.get_element_blocks();
analyze_mesh(el_blks);
build_coordinates(nodes);
const std::vector<gmshparsercpp::MshFile::ElementBlock> & el_blks = f.get_element_blocks();
analyze_mesh();
build_element_block_dim(el_blks);

const std::vector<gmshparsercpp::MshFile::MultiDEntity> * ents = nullptr;
switch (dim) {
case 1:
ents = &f.get_curve_entities();
break;
case 2:
ents = &f.get_surface_entities();
break;
case 3:
ents = &f.get_volume_entities();
break;
default:
throw std::runtime_error(fmt::format("Unsupported spatial dimension {}.", dim));
}
const std::vector<gmshparsercpp::MshFile::MultiDEntity> * ents = get_entities_by_dim(f, dim);
build_element_blocks(el_blk_dim[dim], *ents);

const std::vector<gmshparsercpp::MshFile::MultiDEntity> * sideset_ents = nullptr;
switch (side_set_dim) {
case 1:
sideset_ents = &f.get_curve_entities();
break;
case 2:
sideset_ents = &f.get_surface_entities();
break;
default:
throw std::runtime_error(fmt::format("Unsupported side set dimension {}.", side_set_dim));
}
const std::vector<gmshparsercpp::MshFile::MultiDEntity> * sideset_ents = get_entities_by_dim(f, side_set_dim);
build_side_sets(el_blk_dim[side_set_dim], *sideset_ents);
}

Expand Down

0 comments on commit 0168a69

Please sign in to comment.