Select Git revision
-
Stéphane Del Pino authored
This tool checks if chosen data are synchronized and have same values comparing two runs By now only ItemValue are checked
Stéphane Del Pino authoredThis tool checks if chosen data are synchronized and have same values comparing two runs By now only ItemValue are checked
MeshBuilderBase.cpp 8.17 KiB
#include <mesh/MeshBuilderBase.hpp>
#include <mesh/Connectivity.hpp>
#include <mesh/ConnectivityDescriptor.hpp>
#include <mesh/ConnectivityDispatcher.hpp>
#include <mesh/ConnectivityDispatcherVariant.hpp>
#include <mesh/ItemId.hpp>
#include <mesh/Mesh.hpp>
#include <mesh/MeshData.hpp>
#include <mesh/MeshDataManager.hpp>
#include <mesh/MeshVariant.hpp>
#include <utils/PugsAssert.hpp>
#include <utils/PugsMacros.hpp>
#include <vector>
template <size_t Dimension>
void
MeshBuilderBase::_dispatch()
{
if (parallel::size() == 1) {
return;
}
using ConnectivityType = Connectivity<Dimension>;
using Rd = TinyVector<Dimension>;
using MeshType = Mesh<Dimension>;
if (not m_mesh) {
ConnectivityDescriptor descriptor;
std::shared_ptr connectivity = ConnectivityType::build(descriptor);
NodeValue<Rd> xr;
m_mesh = std::make_shared<MeshVariant>(std::make_shared<const MeshType>(connectivity, xr));
}
const MeshType& mesh = *(m_mesh->get<const MeshType>());
auto p_dispatcher = std::make_shared<const ConnectivityDispatcher<Dimension>>(mesh.connectivity());
m_connectivity_dispatcher = std::make_shared<ConnectivityDispatcherVariant>(p_dispatcher);
std::shared_ptr dispatched_connectivity = p_dispatcher->dispatchedConnectivity();
NodeValue<Rd> dispatched_xr = p_dispatcher->dispatch(mesh.xr());
m_mesh = std::make_shared<MeshVariant>(std::make_shared<const MeshType>(dispatched_connectivity, dispatched_xr));
}
template void MeshBuilderBase::_dispatch<1>();
template void MeshBuilderBase::_dispatch<2>();
template void MeshBuilderBase::_dispatch<3>();
template <size_t Dimension>
void
MeshBuilderBase::_checkMesh() const
{
using MeshType = Mesh<Dimension>;
using ConnectivityType = typename MeshType::Connectivity;
if (not m_mesh) {
throw UnexpectedError("mesh is not built yet");
}
const MeshType& mesh = *(m_mesh->get<const MeshType>());
const ConnectivityType& connectivity = mesh.connectivity();
if constexpr (Dimension > 2) { // check for duplicated edges
auto edge_to_node_matrix = connectivity.edgeToNodeMatrix();
std::vector<std::vector<EdgeId>> node_edges(mesh.numberOfNodes());
for (EdgeId edge_id = 0; edge_id < mesh.numberOfEdges(); ++edge_id) {
auto node_list = edge_to_node_matrix[edge_id];
for (size_t i_node = 0; i_node < node_list.size(); ++i_node) {
node_edges[node_list[i_node]].push_back(edge_id);
}
}
for (auto&& edge_list : node_edges) {
std::sort(edge_list.begin(), edge_list.end());
}
for (EdgeId edge_id = 0; edge_id < mesh.numberOfEdges(); ++edge_id) {
auto node_list = edge_to_node_matrix[edge_id];
std::vector<EdgeId> intersection = node_edges[node_list[0]];
for (size_t i_node = 1; i_node < node_list.size(); ++i_node) {
std::vector<EdgeId> local_intersection;
std::set_intersection(intersection.begin(), intersection.end(), //
node_edges[node_list[i_node]].begin(), node_edges[node_list[i_node]].end(),
std::back_inserter(local_intersection));
std::swap(local_intersection, intersection);
if (intersection.size() < 2) {
break;
}
}
if (intersection.size() > 1) {
std::ostringstream error_msg;
error_msg << "invalid mesh.\n\tFollowing edges\n";
for (EdgeId duplicated_edge_id : intersection) {
error_msg << "\t - id=" << duplicated_edge_id << " number=" << connectivity.edgeNumber()[duplicated_edge_id]
<< '\n';
}
error_msg << "\tare duplicated";
throw NormalError(error_msg.str());
}
}
}
if constexpr (Dimension > 1) { // check for duplicated faces
auto face_to_node_matrix = connectivity.faceToNodeMatrix();
std::vector<std::vector<FaceId>> node_faces(mesh.numberOfNodes());
for (FaceId face_id = 0; face_id < mesh.numberOfFaces(); ++face_id) {
auto node_list = face_to_node_matrix[face_id];
for (size_t i_node = 0; i_node < node_list.size(); ++i_node) {
node_faces[node_list[i_node]].push_back(face_id);
}
}
for (auto&& face_list : node_faces) {
std::sort(face_list.begin(), face_list.end());
}
for (FaceId face_id = 0; face_id < mesh.numberOfFaces(); ++face_id) {
auto node_list = face_to_node_matrix[face_id];
std::vector<FaceId> intersection = node_faces[node_list[0]];
for (size_t i_node = 1; i_node < node_list.size(); ++i_node) {
std::vector<FaceId> local_intersection;
std::set_intersection(intersection.begin(), intersection.end(), //
node_faces[node_list[i_node]].begin(), node_faces[node_list[i_node]].end(),
std::back_inserter(local_intersection));
std::swap(local_intersection, intersection);
if (intersection.size() < 2) {
break;
}
}
if (intersection.size() > 1) {
std::ostringstream error_msg;
error_msg << "invalid mesh.\n\tFollowing faces\n";
for (FaceId intersection_face_id : intersection) {
error_msg << "\t - id=" << intersection_face_id
<< " number=" << connectivity.faceNumber()[intersection_face_id] << '\n';
error_msg << "\t nodes:";
for (size_t i = 0; i < face_to_node_matrix[intersection_face_id].size(); ++i) {
error_msg << ' ' << face_to_node_matrix[intersection_face_id][i];
}
error_msg << '\n';
}
error_msg << "\tare duplicated";
throw NormalError(error_msg.str());
}
}
}
auto cell_to_node_matrix = connectivity.cellToNodeMatrix();
{ // check for duplicated cells
std::vector<std::vector<CellId>> node_cells(mesh.numberOfNodes());
for (CellId cell_id = 0; cell_id < mesh.numberOfCells(); ++cell_id) {
auto node_list = cell_to_node_matrix[cell_id];
for (size_t i_node = 0; i_node < node_list.size(); ++i_node) {
node_cells[node_list[i_node]].push_back(cell_id);
}
}
for (auto&& cell_list : node_cells) {
std::sort(cell_list.begin(), cell_list.end());
}
for (CellId cell_id = 0; cell_id < mesh.numberOfCells(); ++cell_id) {
auto node_list = cell_to_node_matrix[cell_id];
std::vector<CellId> intersection = node_cells[node_list[0]];
for (size_t i_node = 1; i_node < node_list.size(); ++i_node) {
std::vector<CellId> local_intersection;
std::set_intersection(intersection.begin(), intersection.end(), //
node_cells[node_list[i_node]].begin(), node_cells[node_list[i_node]].end(),
std::back_inserter(local_intersection));
std::swap(local_intersection, intersection);
if (intersection.size() < 2) {
break;
}
}
if (intersection.size() > 1) {
std::ostringstream error_msg;
error_msg << "invalid mesh.\n\tFollowing cells\n";
for (CellId duplicated_cell_id : intersection) {
error_msg << "\t - id=" << duplicated_cell_id << " number=" << connectivity.cellNumber()[duplicated_cell_id]
<< '\n';
}
error_msg << "\tare duplicated";
throw NormalError(error_msg.str());
}
}
}
const auto& Cjr = MeshDataManager::instance().getMeshData(mesh).Cjr();
const auto& xr = mesh.xr();
for (CellId cell_id = 0; cell_id < mesh.numberOfCells(); ++cell_id) {
double cell_volume = 0;
auto cell_nodes = cell_to_node_matrix[cell_id];
for (size_t i_node = 0; i_node < cell_nodes.size(); ++i_node) {
cell_volume += dot(Cjr(cell_id, i_node), xr[cell_nodes[i_node]]);
}
if (cell_volume <= 0) {
std::ostringstream error_msg;
error_msg << "invalid mesh.\n\tThe following cell\n";
error_msg << "\t - id=" << cell_id << " number=" << connectivity.cellNumber()[cell_id] << '\n';
error_msg << "\thas non-positive volume: " << cell_volume / Dimension;
throw NormalError(error_msg.str());
}
}
}
template void MeshBuilderBase::_checkMesh<1>() const;
template void MeshBuilderBase::_checkMesh<2>() const;
template void MeshBuilderBase::_checkMesh<3>() const;