From 6768d75c7a34045941d50a847017a5363c245664 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Wed, 11 Sep 2024 13:31:59 -0600 Subject: [PATCH 01/28] Basic outline for checkNonMatchingEdges -Adds checkNonMatchingEdges function to MeshDiagnosticsGenerator.C/.h -Includes basic outline of algorithm -Includes nested for-loops for looping over all mesh elements refs #26596 Adding in sudo-code for the algorithm that finds edge overlap -Still need a method to get (x,y,z) coordinates from a node -Also added the function to the header file Updated method for finding intersecting edges -Implemented method to find closest line connecting two edges -If the line length is < tol it counts as intersecting -Only checks is edges don't share any nodes -Needs improvement to ensure it doesn't double check edges -Needs improvement to ensure it loops over all edges even if they're in different blocks/volumes -Need to improve tolerance so it's not hard-coded -Need to add tests Fixed issue where no nodes were intersecting -created new issue where every node is now intersecting --- .../meshgenerators/MeshDiagnosticsGenerator.h | 4 + .../include/utils/MeshBaseDiagnosticsUtils.h | 4 + .../meshgenerators/MeshDiagnosticsGenerator.C | 85 +++++++++++++++++ .../src/utils/MeshBaseDiagnosticsUtils.C | 95 +++++++++++++++++++ 4 files changed, 188 insertions(+) diff --git a/framework/include/meshgenerators/MeshDiagnosticsGenerator.h b/framework/include/meshgenerators/MeshDiagnosticsGenerator.h index 7a7ef2cab0ee..9bc214768562 100644 --- a/framework/include/meshgenerators/MeshDiagnosticsGenerator.h +++ b/framework/include/meshgenerators/MeshDiagnosticsGenerator.h @@ -49,6 +49,8 @@ class MeshDiagnosticsGenerator : public MeshGenerator void checkNonConformalMeshFromAdaptivity(const std::unique_ptr & mesh) const; /// Routine to check whether the Jacobians (elem and side) are not negative void checkLocalJacobians(const std::unique_ptr & mesh) const; + //// Routing to check for non matching edges + void checkNonMatchingEdges(const std::unique_ptr & mesh) const; /** * Utility routine to output the final diagnostics level in the desired mode @@ -82,6 +84,8 @@ class MeshDiagnosticsGenerator : public MeshGenerator const MooseEnum _check_non_conformal_mesh; /// tolerance for detecting when meshes are not conformal const Real _non_conformality_tol; + //// whether to check for intersecting edges + const MooseEnum _check_non_matching_edges; /// whether to check for the adaptivity of non-conformal meshes const MooseEnum _check_adaptivity_non_conformality; /// whether to check for negative jacobians in the domain diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index 2328e5626240..bbfa818a85f4 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -24,4 +24,8 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, const unsigned int num_outputs, const Real conformality_tol, unsigned int & num_nonconformal_nodes); + +bool checkEdgeOverlap(const std::unique_ptr & edge1, + const std::unique_ptr & edge2, + double tol); } diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 016763919660..f52270dd4bf0 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -64,6 +64,9 @@ MeshDiagnosticsGenerator::validParams() params.addParam("examine_non_conformality", chk_option, "whether to examine the conformality of elements in the mesh"); + params.addParam("examine_non_matching_edges", + chk_option, + "Whether to check if there are any intersecting edges"); params.addParam("nonconformal_tol", TOLERANCE, "tolerance for element non-conformality"); params.addParam( "search_for_adaptivity_nonconformality", @@ -93,6 +96,7 @@ MeshDiagnosticsGenerator::MeshDiagnosticsGenerator(const InputParameters & param _check_non_planar_sides(getParam("examine_nonplanar_sides")), _check_non_conformal_mesh(getParam("examine_non_conformality")), _non_conformality_tol(getParam("nonconformal_tol")), + _check_non_matching_edges(getParam("examine_non_matching_edges")), _check_adaptivity_non_conformality( getParam("search_for_adaptivity_nonconformality")), _check_local_jacobian(getParam("check_local_jacobian")), @@ -108,10 +112,18 @@ MeshDiagnosticsGenerator::MeshDiagnosticsGenerator(const InputParameters & param paramError("examine_non_conformality", "You must set this parameter to true to trigger mesh conformality check"); if (_check_sidesets_orientation == "NO_CHECK" && _check_watertight_sidesets == "NO_CHECK" && +<<<<<<< HEAD _check_watertight_nodesets == "NO_CHECK" && _check_element_volumes == "NO_CHECK" && _check_element_types == "NO_CHECK" && _check_element_overlap == "NO_CHECK" && _check_non_planar_sides == "NO_CHECK" && _check_non_conformal_mesh == "NO_CHECK" && _check_adaptivity_non_conformality == "NO_CHECK" && _check_local_jacobian == "NO_CHECK") +======= + _check_watertight_nodesets == "NO_CHECK" && _check_element_volumes == "NO_CHECK" && + _check_element_types == "NO_CHECK" && _check_element_overlap == "NO_CHECK" && + _check_non_planar_sides == "NO_CHECK" && _check_non_conformal_mesh == "NO_CHECK" && + _check_adaptivity_non_conformality == "NO_CHECK" && _check_local_jacobian == "NO_CHECK" && + _check_non_matching_edges == "NO_CHECK") +>>>>>>> a759b2e1ae (Basic outline for checkNonMatchingEdges) mooseError("You need to turn on at least one diagnostic. Did you misspell a parameter?"); } @@ -158,6 +170,9 @@ MeshDiagnosticsGenerator::generate() if (_check_local_jacobian != "NO_CHECK") checkLocalJacobians(mesh); + if (_check_non_matching_edges != "NO_CHECK") + checkNonMatchingEdges(mesh); + return dynamic_pointer_cast(mesh); } @@ -1388,6 +1403,76 @@ MeshDiagnosticsGenerator::checkLocalJacobians(const std::unique_ptr & num_negative_side_qp_jacobians); } +void +MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr & mesh) const +{ + /*Algorithm Overview + 1)Prechecks + a)This algorithn only works for 3D so check for that first + b)Optional->In theory this should only be needed for tetrahedral meshes so checking to ensure the mesh has tet cells could be usefull + 2)Loop + a)Loop through every node + b)For each node get the edges associated with it + c)For each edge check overlap with any edges nearby + d)Have check to make sure the same pair of edges are not being tested twice for overlap + e)Have check to see if the edges are close to eachother before testing for overlap + 3)Overlap check + a)Use algorithm from Paul Bourke + b)Finds shortest line that connects two line segments + c)If length of line is below some threshold, print out info about which two edges are intersecting + */ + unsigned int num_intersecting_edges = 0; + //unsigned int num_elem_check = 0; + for (auto & elem : mesh->active_element_ptr_range()) + { + std::vector> elem_edges(elem->n_edges()); + //num_elem_check++; + for (auto i : elem->edge_index_range()) + elem_edges[i] = elem->build_edge_ptr(i); + for (auto & other_elem : mesh->active_element_ptr_range()) + { + std::vector> other_edges(other_elem->n_edges()); + //num_elem_check++; + for (auto j : other_elem->edge_index_range()) + other_edges[j] = other_elem->build_edge_ptr(j); + for (auto & edge : elem_edges) + { + for (auto & other_edge : other_edges) + { + // Now compare edge with other_edge + double tol = 0.000001; + bool overlap = MeshBaseDiagnosticsUtils::checkEdgeOverlap(edge, other_edge, tol); + if (overlap) + { + num_intersecting_edges++; + const auto & n = edge->get_nodes()[0]; + //const Point * const p = n; + const Real x = n->operator()(0); + const Real y = n->operator()(1); + const Real z = n->operator()(2); + + std::string x_coord = std::to_string(x); + std::string y_coord = std::to_string(y); + std::string z_coord = std::to_string(z); + + std::string message = "Non-matching edges found near (" + x_coord + ", " + y_coord + ", " + z_coord + ")"; + _console << message << std::endl; + } + } + /*for (const auto other_i : other_elem->edge_index_range()) + { + const auto other_edge = other_elem->build_edge_ptr(); + // you now have edge, and other_edge + } + */ + } + } + } + std::string string_message = "Number of intersecting edges: " + std::to_string(num_intersecting_edges); + _console << string_message << std::endl; + diagnosticsLog("Number of intesecting edges: ", _check_non_matching_edges, num_intersecting_edges); +} + void MeshDiagnosticsGenerator::diagnosticsLog(std::string msg, const MooseEnum & log_level, diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 9025e472db32..23229ae113ad 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -67,4 +67,99 @@ checkNonConformalMesh(const std::unique_ptr & mesh, } pl->unset_close_to_point_tol(); } + +bool +checkEdgeOverlap(const std::unique_ptr & edge1, + const std::unique_ptr & edge2, + double tol) +{ + //get node array from two edges + const auto & node_list1 = edge1->get_nodes(); + const auto & node_list2 = edge2->get_nodes(); + + //make sure two edges are not the same and don't share any nodes edge before checking for overlap + auto & n1 = (*node_list1)[0]; + auto & n2 = (*node_list1)[-1]; + auto & n3 = (*node_list2)[0]; + auto & n4 = (*node_list2)[-1]; + //auto num_nodes = node_list1->size(); + //const Point *p1, *p2, *p3, *p4; + //const Point * const p1 = n1; + //const Point * const p2 = n2; + //const Point * const p3 = n3; + //const Point * const p4 = n4; + + auto n1x = n1.operator()(0); + auto n1y = n1.operator()(1); + auto n1z = n1.operator()(2); + auto n2x = n2.operator()(0); + auto n2y = n2.operator()(1); + auto n2z = n2.operator()(2); + auto n3x = n3.operator()(0); + auto n3y = n3.operator()(1); + auto n3z = n3.operator()(2); + auto n4x = n4.operator()(0); + auto n4y = n4.operator()(1); + auto n4z = n4.operator()(2); + + if(std::abs(n1x-n3x) Date: Thu, 19 Sep 2024 16:51:22 -0600 Subject: [PATCH 02/28] Basic working algorithm for detecting intersecting mesh edges -Given a mesh it iterates through each element and compares it to every other element -For each pair of elements it compares their edges to see if any intersect -Checks that edges don't share any nodes before checking for overlap -Prints out where the intersections are and prints out how many intersecting edges were found -Can be improved in future versions so that it doesn't have to loop through every element twice -Tests still need to be added refs #26596 --- .../meshgenerators/MeshDiagnosticsGenerator.C | 44 +++++--- .../src/utils/MeshBaseDiagnosticsUtils.C | 103 ++++++++++++------ 2 files changed, 100 insertions(+), 47 deletions(-) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index f52270dd4bf0..5cdd97bfd689 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1422,14 +1422,15 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr c)If length of line is below some threshold, print out info about which two edges are intersecting */ unsigned int num_intersecting_edges = 0; + std::vector checked_edges; //unsigned int num_elem_check = 0; - for (auto & elem : mesh->active_element_ptr_range()) + for (auto elem : mesh->active_element_ptr_range()) { std::vector> elem_edges(elem->n_edges()); //num_elem_check++; for (auto i : elem->edge_index_range()) elem_edges[i] = elem->build_edge_ptr(i); - for (auto & other_elem : mesh->active_element_ptr_range()) + for (auto other_elem : mesh->active_element_ptr_range()) { std::vector> other_edges(other_elem->n_edges()); //num_elem_check++; @@ -1444,16 +1445,35 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr bool overlap = MeshBaseDiagnosticsUtils::checkEdgeOverlap(edge, other_edge, tol); if (overlap) { - num_intersecting_edges++; - const auto & n = edge->get_nodes()[0]; + const auto node_list1 = edge->get_nodes(); + const auto node_list2 = other_edge->get_nodes(); + auto n1 = *node_list1[0]; + auto n2 = *node_list1[1]; + auto n3 = *node_list2[0]; + auto n4 = *node_list2[1]; + if(std::find(checked_edges.begin(), checked_edges.end(), n1) !=checked_edges.end() && + std::find(checked_edges.begin(), checked_edges.end(), n2) !=checked_edges.end() && + std::find(checked_edges.begin(), checked_edges.end(), n3) !=checked_edges.end() && + std::find(checked_edges.begin(), checked_edges.end(), n4) !=checked_edges.end()) + { + continue; + } + checked_edges.push_back(n1); + checked_edges.push_back(n2); + checked_edges.push_back(n3); + checked_edges.push_back(n4); + num_intersecting_edges+=2; //const Point * const p = n; - const Real x = n->operator()(0); - const Real y = n->operator()(1); - const Real z = n->operator()(2); + auto x1 = n1.operator()(0); + auto y1 = n1.operator()(1); + auto z1 = n1.operator()(2); + auto x2 = n2.operator()(0); + auto y2 = n2.operator()(1); + auto z2 = n2.operator()(2); - std::string x_coord = std::to_string(x); - std::string y_coord = std::to_string(y); - std::string z_coord = std::to_string(z); + std::string x_coord = std::to_string((x1+x2)/2); + std::string y_coord = std::to_string((y1+y2)/2); + std::string z_coord = std::to_string((z1+z2)/2); std::string message = "Non-matching edges found near (" + x_coord + ", " + y_coord + ", " + z_coord + ")"; _console << message << std::endl; @@ -1468,9 +1488,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr } } } - std::string string_message = "Number of intersecting edges: " + std::to_string(num_intersecting_edges); - _console << string_message << std::endl; - diagnosticsLog("Number of intesecting edges: ", _check_non_matching_edges, num_intersecting_edges); + diagnosticsLog("Number of intesecting edges: " + Moose::stringify(num_intersecting_edges), _check_non_matching_edges, num_intersecting_edges); } void diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 23229ae113ad..65fb32066cc7 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -74,14 +74,19 @@ checkEdgeOverlap(const std::unique_ptr & edge1, double tol) { //get node array from two edges - const auto & node_list1 = edge1->get_nodes(); - const auto & node_list2 = edge2->get_nodes(); + const auto node_list1 = edge1->get_nodes(); + const auto node_list2 = edge2->get_nodes(); //make sure two edges are not the same and don't share any nodes edge before checking for overlap - auto & n1 = (*node_list1)[0]; - auto & n2 = (*node_list1)[-1]; - auto & n3 = (*node_list2)[0]; - auto & n4 = (*node_list2)[-1]; + auto n1 = *node_list1[0]; + auto n2 = *node_list1[1]; + auto n3 = *node_list2[0]; + auto n4 = *node_list2[1]; + //auto size1 = sizeof(*node_list1); + //auto size2 = std::size(*node_list1); + //if(size1 < 3){ + //return false; + //} //auto num_nodes = node_list1->size(); //const Point *p1, *p2, *p3, *p4; //const Point * const p1 = n1; @@ -89,27 +94,48 @@ checkEdgeOverlap(const std::unique_ptr & edge1, //const Point * const p3 = n3; //const Point * const p4 = n4; - auto n1x = n1.operator()(0); - auto n1y = n1.operator()(1); - auto n1z = n1.operator()(2); - auto n2x = n2.operator()(0); - auto n2y = n2.operator()(1); - auto n2z = n2.operator()(2); - auto n3x = n3.operator()(0); - auto n3y = n3.operator()(1); - auto n3z = n3.operator()(2); - auto n4x = n4.operator()(0); - auto n4y = n4.operator()(1); - auto n4z = n4.operator()(2); - - if(std::abs(n1x-n3x) & edge1, */ //There's a chance that they overlap. Find shortest line that connects two edges and if its length is close enough to 0 return true - double n13x, n13y, n13z, n21x, n21y, n21z, n43x, n43y, n43z; double d1343, d4321, d1321, d4343, d2121, numerator, denominator, mua, mub; - n13x = n1x - n3x; - n13y = n1y - n3y; - n13z = n1z - n3z; - n21x = n2x - n1x; - n21y = n2y - n1y; - n21z = n2z - n1z; - n43x = n4x - n3x; - n43y = n4y - n3y; - n43z = n4z - n3z; - d1343 = n13x * n43x + n13y * n43y + n13z * n43z; d4321 = n43x * n21x + n43y * n21y + n43z * n21z; d1321 = n13x * n21x + n13y * n21y + n13z * n21z; @@ -142,6 +157,11 @@ checkEdgeOverlap(const std::unique_ptr & edge1, denominator = d2121 * d4343 - d4321 * d4321; numerator = d1343 * d4321 - d1321 * d4343; + if(std::fabs(denominator) < tol) + { + //This indicates that the intersecting line is vertical so they don't intersect + return false; + } mua = numerator/denominator; mub = (d1343 + (mua * d4321)) / d4343; @@ -155,6 +175,21 @@ checkEdgeOverlap(const std::unique_ptr & edge1, nby = n3y + mub * n43y; nbz = n3z + mub * n43z; + //This method assume the two lines are infinite. This check to make sure na and nb are part of their respective line segments + if((nax < std::min(n1x, n2x)) || (nax > std::max(n1x, n2x)) || + (nay < std::min(n1y, n2y)) || (nay > std::max(n1y, n2y)) || + (naz < std::min(n1z, n2z)) || (naz > std::max(n1z, n2z))) + { + return false; + } + + if((nbx < std::min(n3x, n4x)) || (nax > std::max(n3x, n4x)) || + (nby < std::min(n3y, n4y)) || (nay > std::max(n3y, n4y)) || + (nbz < std::min(n3z, n4z)) || (naz > std::max(n3z, n4z))) + { + return false; + } + //Calculate distance between these two nodes double distance = std::sqrt(std::pow(nax - nbx, 2) + std::pow(nay - nby, 2) + std::pow(naz - nbz, 2)); if (distance < tol) From 95d9550bcc7997e9778ae440592ca0981f6f1188 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Mon, 23 Sep 2024 14:14:56 -0600 Subject: [PATCH 03/28] Added more checks to intersecing edges algorithm -Now it checks if the bounding boxes of the two elements being checked intersect -If they do not intersect then they can't have edges that intersect so they are skipped -Moved up edge check so now before checking for overlap, it checks to see if the two edges have already been checked -Removed commented out code, and changed certain variable types from double to auto -console output is now done in MehBaseDiagnosticsUtils.C and prints the coordinates where the edges intersect refs #26596 --- .../include/utils/MeshBaseDiagnosticsUtils.h | 2 +- .../meshgenerators/MeshDiagnosticsGenerator.C | 66 ++++++------- .../src/utils/MeshBaseDiagnosticsUtils.C | 99 ++++++++----------- 3 files changed, 68 insertions(+), 99 deletions(-) diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index bbfa818a85f4..79078dd9dd04 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -27,5 +27,5 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, bool checkEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, - double tol); + const ConsoleStream & console); } diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 5cdd97bfd689..86ac067893f3 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1423,68 +1423,58 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr */ unsigned int num_intersecting_edges = 0; std::vector checked_edges; - //unsigned int num_elem_check = 0; for (auto elem : mesh->active_element_ptr_range()) { std::vector> elem_edges(elem->n_edges()); - //num_elem_check++; for (auto i : elem->edge_index_range()) elem_edges[i] = elem->build_edge_ptr(i); for (auto other_elem : mesh->active_element_ptr_range()) { + //If they're the same element then there's no need to check for overlap + if (elem == other_elem) + { + continue; + } + //Get bounding boxes for both elements. If the bounding boxes don't intersect then no edges will either + auto boundingBox1 = elem->loose_bounding_box(); + auto boundingBox2 = other_elem->loose_bounding_box(); + if(!(boundingBox1.intersects(boundingBox2))) + { + continue; + } std::vector> other_edges(other_elem->n_edges()); - //num_elem_check++; for (auto j : other_elem->edge_index_range()) other_edges[j] = other_elem->build_edge_ptr(j); for (auto & edge : elem_edges) { for (auto & other_edge : other_edges) { + const auto node_list1 = edge->get_nodes(); + const auto node_list2 = other_edge->get_nodes(); + auto n1 = *node_list1[0]; + auto n2 = *node_list1[1]; + auto n3 = *node_list2[0]; + auto n4 = *node_list2[1]; + //Check if the edges have already been added to our check_edges list + if (std::find(checked_edges.begin(), checked_edges.end(), n1) !=checked_edges.end() && + std::find(checked_edges.begin(), checked_edges.end(), n2) !=checked_edges.end() && + std::find(checked_edges.begin(), checked_edges.end(), n3) !=checked_edges.end() && + std::find(checked_edges.begin(), checked_edges.end(), n4) !=checked_edges.end()) + { + continue; + } // Now compare edge with other_edge - double tol = 0.000001; - bool overlap = MeshBaseDiagnosticsUtils::checkEdgeOverlap(edge, other_edge, tol); + bool overlap = MeshBaseDiagnosticsUtils::checkEdgeOverlap(edge, other_edge, _console); if (overlap) { - const auto node_list1 = edge->get_nodes(); - const auto node_list2 = other_edge->get_nodes(); - auto n1 = *node_list1[0]; - auto n2 = *node_list1[1]; - auto n3 = *node_list2[0]; - auto n4 = *node_list2[1]; - if(std::find(checked_edges.begin(), checked_edges.end(), n1) !=checked_edges.end() && - std::find(checked_edges.begin(), checked_edges.end(), n2) !=checked_edges.end() && - std::find(checked_edges.begin(), checked_edges.end(), n3) !=checked_edges.end() && - std::find(checked_edges.begin(), checked_edges.end(), n4) !=checked_edges.end()) - { - continue; - } + // Add the nodes that make up the 2 edges to the vector checked_edges checked_edges.push_back(n1); checked_edges.push_back(n2); checked_edges.push_back(n3); checked_edges.push_back(n4); num_intersecting_edges+=2; - //const Point * const p = n; - auto x1 = n1.operator()(0); - auto y1 = n1.operator()(1); - auto z1 = n1.operator()(2); - auto x2 = n2.operator()(0); - auto y2 = n2.operator()(1); - auto z2 = n2.operator()(2); - - std::string x_coord = std::to_string((x1+x2)/2); - std::string y_coord = std::to_string((y1+y2)/2); - std::string z_coord = std::to_string((z1+z2)/2); - - std::string message = "Non-matching edges found near (" + x_coord + ", " + y_coord + ", " + z_coord + ")"; - _console << message << std::endl; } } - /*for (const auto other_i : other_elem->edge_index_range()) - { - const auto other_edge = other_elem->build_edge_ptr(); - // you now have edge, and other_edge - } - */ } } } diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 65fb32066cc7..33bcb4452c58 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -71,7 +71,7 @@ checkNonConformalMesh(const std::unique_ptr & mesh, bool checkEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, - double tol) + const ConsoleStream & console) { //get node array from two edges const auto node_list1 = edge1->get_nodes(); @@ -82,68 +82,40 @@ checkEdgeOverlap(const std::unique_ptr & edge1, auto n2 = *node_list1[1]; auto n3 = *node_list2[0]; auto n4 = *node_list2[1]; - //auto size1 = sizeof(*node_list1); - //auto size2 = std::size(*node_list1); - //if(size1 < 3){ - //return false; - //} - //auto num_nodes = node_list1->size(); - //const Point *p1, *p2, *p3, *p4; - //const Point * const p1 = n1; - //const Point * const p2 = n2; - //const Point * const p3 = n3; - //const Point * const p4 = n4; - - double n1x = 1.0*(n1.operator()(0)); - double n1y = 1.0*(n1.operator()(1)); - double n1z = 1.0*(n1.operator()(2)); - double n2x = 1.0*(n2.operator()(0)); - double n2y = 1.0*(n2.operator()(1)); - double n2z = 1.0*(n2.operator()(2)); - double n3x = 1.0*(n3.operator()(0)); - double n3y = 1.0*(n3.operator()(1)); - double n3z = 1.0*(n3.operator()(2)); - double n4x = 1.0*(n4.operator()(0)); - double n4y = 1.0*(n4.operator()(1)); - double n4z = 1.0*(n4.operator()(2)); + auto n1x = 1.0*(n1.operator()(0)); + auto n1y = 1.0*(n1.operator()(1)); + auto n1z = 1.0*(n1.operator()(2)); + auto n2x = 1.0*(n2.operator()(0)); + auto n2y = 1.0*(n2.operator()(1)); + auto n2z = 1.0*(n2.operator()(2)); + auto n3x = 1.0*(n3.operator()(0)); + auto n3y = 1.0*(n3.operator()(1)); + auto n3z = 1.0*(n3.operator()(2)); + auto n4x = 1.0*(n4.operator()(0)); + auto n4y = 1.0*(n4.operator()(1)); + auto n4z = 1.0*(n4.operator()(2)); + + //Tolerance should be based on edge length. Here it is arbitrarily set to 1/1000th of the average of the two edge lengths + auto edge1Length = std::sqrt(std::pow(n1x - n2x, 2) + std::pow(n1y - n2y, 2) + std::pow(n1z - n2z, 2)); + auto edge2Length = std::sqrt(std::pow(n3x - n4x, 2) + std::pow(n3y - n4y, 2) + std::pow(n3z - n4z, 2)); + auto tol = ((edge1Length + edge2Length)/2)/1000; //double n13x, n13y, n13z, n21x, n21y, n21z, n43x, n43y, n43z; - double n13x = n1x - n3x; - double n13y = n1y - n3y; - double n13z = n1z - n3z; - double n21x = n2x - n1x; - double n21y = n2y - n1y; - double n21z = n2z - n1z; - double n43x = n4x - n3x; - double n43y = n4y - n3y; - double n43z = n4z - n3z; - - //double n13xfabs = std::fabs(n13x); - //double n13xabs = std::abs(n13x); - if((std::fabs(n1x - n3x) & edge1, //Calculate distance between these two nodes double distance = std::sqrt(std::pow(nax - nbx, 2) + std::pow(nay - nby, 2) + std::pow(naz - nbz, 2)); if (distance < tol) + { + std::string x_coord = std::to_string(nax); + std::string y_coord = std::to_string(nay); + std::string z_coord = std::to_string(naz); + std::string message = "Non-matching edges found near (" + x_coord + ", " + y_coord + ", " + z_coord + ")"; + console << message << std::endl; return true; + } else return false; } From 22663ae19730ea39384191ff4cad03b8b1cf51b7 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Tue, 24 Sep 2024 16:54:35 -0600 Subject: [PATCH 04/28] Added intersecting edge test and added tolerance parameter -New test that generates a mesh with intersecting edges and tests it -Made tolerance parameter adjustable. User can change it by adding intersection_tol= to diagnostic portion of input file -Changed some varible type to auto and made comments more specific closes #26596 --- .../meshgenerators/MeshDiagnosticsGenerator.h | 2 + .../include/utils/MeshBaseDiagnosticsUtils.h | 3 +- .../meshgenerators/MeshDiagnosticsGenerator.C | 23 +++-- .../src/utils/MeshBaseDiagnosticsUtils.C | 96 ++++++++----------- .../mesh_diagnostics_generator/all_at_once.i | 3 +- .../intersecting_edge_test.i | 66 +++++++++++++ .../mesh_diagnostics_generator/tests | 10 +- 7 files changed, 136 insertions(+), 67 deletions(-) create mode 100644 test/tests/meshgenerators/mesh_diagnostics_generator/intersecting_edge_test.i diff --git a/framework/include/meshgenerators/MeshDiagnosticsGenerator.h b/framework/include/meshgenerators/MeshDiagnosticsGenerator.h index 9bc214768562..ca65b6a79143 100644 --- a/framework/include/meshgenerators/MeshDiagnosticsGenerator.h +++ b/framework/include/meshgenerators/MeshDiagnosticsGenerator.h @@ -86,6 +86,8 @@ class MeshDiagnosticsGenerator : public MeshGenerator const Real _non_conformality_tol; //// whether to check for intersecting edges const MooseEnum _check_non_matching_edges; + //// tolerance for detecting when edges intersect + const Real _non_matching_edge_tol; /// whether to check for the adaptivity of non-conformal meshes const MooseEnum _check_adaptivity_non_conformality; /// whether to check for negative jacobians in the domain diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index 79078dd9dd04..b94c4bd69cca 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -27,5 +27,6 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, bool checkEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, - const ConsoleStream & console); + const ConsoleStream & console, + const Real insersection_tol); } diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 86ac067893f3..875942e75ee5 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -67,6 +67,7 @@ MeshDiagnosticsGenerator::validParams() params.addParam("examine_non_matching_edges", chk_option, "Whether to check if there are any intersecting edges"); + params.addParam("intersection_tol", TOLERANCE, "tolerence for intersecting edges"); params.addParam("nonconformal_tol", TOLERANCE, "tolerance for element non-conformality"); params.addParam( "search_for_adaptivity_nonconformality", @@ -97,6 +98,7 @@ MeshDiagnosticsGenerator::MeshDiagnosticsGenerator(const InputParameters & param _check_non_conformal_mesh(getParam("examine_non_conformality")), _non_conformality_tol(getParam("nonconformal_tol")), _check_non_matching_edges(getParam("examine_non_matching_edges")), + _non_matching_edge_tol(getParam("intersection_tol")), _check_adaptivity_non_conformality( getParam("search_for_adaptivity_nonconformality")), _check_local_jacobian(getParam("check_local_jacobian")), @@ -1409,18 +1411,19 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr /*Algorithm Overview 1)Prechecks a)This algorithn only works for 3D so check for that first - b)Optional->In theory this should only be needed for tetrahedral meshes so checking to ensure the mesh has tet cells could be usefull 2)Loop - a)Loop through every node - b)For each node get the edges associated with it - c)For each edge check overlap with any edges nearby + a)Loop through every element + b)For each element get the edges associated with it + c)For each edge check overlap with any edges of nearby elements d)Have check to make sure the same pair of edges are not being tested twice for overlap - e)Have check to see if the edges are close to eachother before testing for overlap 3)Overlap check - a)Use algorithm from Paul Bourke - b)Finds shortest line that connects two line segments - c)If length of line is below some threshold, print out info about which two edges are intersecting + a)Shortest line that connects both lines is perpendicular to both lines + b)A good overview of the math for finding intersecting lines can be found here->paulbourke.net/geometry/pointlineplane/ */ + if(mesh->mesh_dimension() != 3) + mooseError("The edge intersection algorithm only works with 3D meshes"); + if (!mesh->is_serial()) + mooseError("Only serialized/replicated meshes are supported"); unsigned int num_intersecting_edges = 0; std::vector checked_edges; for (auto elem : mesh->active_element_ptr_range()) @@ -1464,7 +1467,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr continue; } // Now compare edge with other_edge - bool overlap = MeshBaseDiagnosticsUtils::checkEdgeOverlap(edge, other_edge, _console); + bool overlap = MeshBaseDiagnosticsUtils::checkEdgeOverlap(edge, other_edge, _console, _non_matching_edge_tol); if (overlap) { // Add the nodes that make up the 2 edges to the vector checked_edges @@ -1478,7 +1481,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr } } } - diagnosticsLog("Number of intesecting edges: " + Moose::stringify(num_intersecting_edges), _check_non_matching_edges, num_intersecting_edges); + diagnosticsLog("Number of intersecting edges: " + Moose::stringify(num_intersecting_edges), _check_non_matching_edges, num_intersecting_edges); } void diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 33bcb4452c58..023d48234b66 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -71,7 +71,8 @@ checkNonConformalMesh(const std::unique_ptr & mesh, bool checkEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, - const ConsoleStream & console) + const ConsoleStream & console, + const Real intersection_tol) { //get node array from two edges const auto node_list1 = edge1->get_nodes(); @@ -83,69 +84,56 @@ checkEdgeOverlap(const std::unique_ptr & edge1, auto n3 = *node_list2[0]; auto n4 = *node_list2[1]; - auto n1x = 1.0*(n1.operator()(0)); - auto n1y = 1.0*(n1.operator()(1)); - auto n1z = 1.0*(n1.operator()(2)); - auto n2x = 1.0*(n2.operator()(0)); - auto n2y = 1.0*(n2.operator()(1)); - auto n2z = 1.0*(n2.operator()(2)); - auto n3x = 1.0*(n3.operator()(0)); - auto n3y = 1.0*(n3.operator()(1)); - auto n3z = 1.0*(n3.operator()(2)); - auto n4x = 1.0*(n4.operator()(0)); - auto n4y = 1.0*(n4.operator()(1)); - auto n4z = 1.0*(n4.operator()(2)); - - //Tolerance should be based on edge length. Here it is arbitrarily set to 1/1000th of the average of the two edge lengths - auto edge1Length = std::sqrt(std::pow(n1x - n2x, 2) + std::pow(n1y - n2y, 2) + std::pow(n1z - n2z, 2)); - auto edge2Length = std::sqrt(std::pow(n3x - n4x, 2) + std::pow(n3y - n4y, 2) + std::pow(n3z - n4z, 2)); - auto tol = ((edge1Length + edge2Length)/2)/1000; - //double n13x, n13y, n13z, n21x, n21y, n21z, n43x, n43y, n43z; - auto n13x = n1x - n3x; - auto n13y = n1y - n3y; - auto n13z = n1z - n3z; - auto n21x = n2x - n1x; - auto n21y = n2y - n1y; - auto n21z = n2z - n1z; - auto n43x = n4x - n3x; - auto n43y = n4y - n3y; - auto n43z = n4z - n3z; - - if(((std::fabs(n1x - n3x) std::max(n1x, n2x)) || @@ -164,7 +152,7 @@ checkEdgeOverlap(const std::unique_ptr & edge1, //Calculate distance between these two nodes double distance = std::sqrt(std::pow(nax - nbx, 2) + std::pow(nay - nby, 2) + std::pow(naz - nbz, 2)); - if (distance < tol) + if (distance < intersection_tol) { std::string x_coord = std::to_string(nax); std::string y_coord = std::to_string(nay); diff --git a/test/tests/meshgenerators/mesh_diagnostics_generator/all_at_once.i b/test/tests/meshgenerators/mesh_diagnostics_generator/all_at_once.i index fe649af20b3a..236a2d75c325 100644 --- a/test/tests/meshgenerators/mesh_diagnostics_generator/all_at_once.i +++ b/test/tests/meshgenerators/mesh_diagnostics_generator/all_at_once.i @@ -15,7 +15,8 @@ examine_nonplanar_sides = INFO examine_sidesets_orientation = WARNING check_for_watertight_sidesets = WARNING - check_for_watertight_nodesets = WARNING + check_for_watertight_nodesets = WARNING + examine_non_matching_edges = WARNING search_for_adaptivity_nonconformality = WARNING check_local_jacobian = WARNING [] diff --git a/test/tests/meshgenerators/mesh_diagnostics_generator/intersecting_edge_test.i b/test/tests/meshgenerators/mesh_diagnostics_generator/intersecting_edge_test.i new file mode 100644 index 000000000000..4f5fdccd1cbf --- /dev/null +++ b/test/tests/meshgenerators/mesh_diagnostics_generator/intersecting_edge_test.i @@ -0,0 +1,66 @@ +[Mesh] + [gmg1] + type = GeneratedMeshGenerator + dim = 3 + nx = 2 + ny = 2 + nz = 2 + xmin = -0.5 + xmax = 0.5 + ymin = -0.5 + ymax = 0.5 + zmin = -1.5 + zmax = -0.5 + [] + [gmg2] + type = GeneratedMeshGenerator + dim = 3 + nx = 2 + ny = 2 + nz = 2 + xmin = -0.5 + xmax = 0.5 + ymin = -0.5 + ymax = 0.5 + zmin = -0.5 + zmax = 0.5 + [] + [block_1] + type = ParsedSubdomainMeshGenerator + input = gmg1 + combinatorial_geometry = 'z < -0.5' + block_id = 1 + [] + [block_2] + type = ParsedSubdomainMeshGenerator + input = gmg2 + combinatorial_geometry = 'z > 0.5' + block_id = 2 + [] + [convert1] + type = ElementsToTetrahedronsConverter + input = 'block_1' + [] + [convert2] + type = ElementsToTetrahedronsConverter + input = 'block_2' + [] + [rotate] + type = TransformGenerator + input = convert2 + transform = 'rotate' + vector_value = '180 0 0' + [] + [cmbn] + type = CombinerGenerator + inputs = 'convert1 rotate' + [] + [diag] + type = MeshDiagnosticsGenerator + input = cmbn + examine_non_matching_edges = INFO + [] +[] +[Outputs] + exodus = true +[] diff --git a/test/tests/meshgenerators/mesh_diagnostics_generator/tests b/test/tests/meshgenerators/mesh_diagnostics_generator/tests index 655af16e3e4e..588aed3a88f9 100644 --- a/test/tests/meshgenerators/mesh_diagnostics_generator/tests +++ b/test/tests/meshgenerators/mesh_diagnostics_generator/tests @@ -13,6 +13,13 @@ expect_out = 'Number of non-conformal nodes: 3' detail = 'element overlapping,' [] + [intersecting_edges] + type = RunApp + input = 'intersecting_edge_test.i' + cli_args = '--mesh-only' + mesh_mode = replicated + expect_out = 'Number of intersecting edges: 4' + [] [inconsistent_sidesets] type = RunApp input = 'consistent_domains.i' @@ -269,7 +276,8 @@ Mesh/diag/examine_nonplanar_sides=NO_CHECK Mesh/diag/examine_sidesets_orientation=NO_CHECK Mesh/diag/check_for_watertight_sidesets=NO_CHECK - Mesh/diag/check_for_watertight_nodesets=NO_CHECK + Mesh/diag/check_for_watertight_nodesets=NO_CHECK + Mesh/diag/examine_non_matching_edges=NO_CHECK Mesh/diag/search_for_adaptivity_nonconformality=NO_CHECK Mesh/diag/check_local_jacobian=NO_CHECK --mesh-only" detail = 'a diagnostics object is created but no diagnostics are requested.' From 833cada78995aabb1fa2b37da922e26aff7103e4 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Tue, 24 Sep 2024 17:14:48 -0600 Subject: [PATCH 05/28] Implemented style changes reccomended by Civet --- .../include/utils/MeshBaseDiagnosticsUtils.h | 2 +- .../meshgenerators/MeshDiagnosticsGenerator.C | 35 +++++++----- .../src/utils/MeshBaseDiagnosticsUtils.C | 56 ++++++++++--------- 3 files changed, 52 insertions(+), 41 deletions(-) diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index b94c4bd69cca..3e291c69dceb 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -25,7 +25,7 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, const Real conformality_tol, unsigned int & num_nonconformal_nodes); -bool checkEdgeOverlap(const std::unique_ptr & edge1, +bool checkEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, const ConsoleStream & console, const Real insersection_tol); diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 875942e75ee5..68385b048518 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -65,8 +65,8 @@ MeshDiagnosticsGenerator::validParams() chk_option, "whether to examine the conformality of elements in the mesh"); params.addParam("examine_non_matching_edges", - chk_option, - "Whether to check if there are any intersecting edges"); + chk_option, + "Whether to check if there are any intersecting edges"); params.addParam("intersection_tol", TOLERANCE, "tolerence for intersecting edges"); params.addParam("nonconformal_tol", TOLERANCE, "tolerance for element non-conformality"); params.addParam( @@ -1418,9 +1418,10 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr d)Have check to make sure the same pair of edges are not being tested twice for overlap 3)Overlap check a)Shortest line that connects both lines is perpendicular to both lines - b)A good overview of the math for finding intersecting lines can be found here->paulbourke.net/geometry/pointlineplane/ + b)A good overview of the math for finding intersecting lines can be found + here->paulbourke.net/geometry/pointlineplane/ */ - if(mesh->mesh_dimension() != 3) + if (mesh->mesh_dimension() != 3) mooseError("The edge intersection algorithm only works with 3D meshes"); if (!mesh->is_serial()) mooseError("Only serialized/replicated meshes are supported"); @@ -1433,15 +1434,16 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr elem_edges[i] = elem->build_edge_ptr(i); for (auto other_elem : mesh->active_element_ptr_range()) { - //If they're the same element then there's no need to check for overlap + // If they're the same element then there's no need to check for overlap if (elem == other_elem) { continue; } - //Get bounding boxes for both elements. If the bounding boxes don't intersect then no edges will either + // Get bounding boxes for both elements. If the bounding boxes don't intersect then no edges + // will either auto boundingBox1 = elem->loose_bounding_box(); auto boundingBox2 = other_elem->loose_bounding_box(); - if(!(boundingBox1.intersects(boundingBox2))) + if (!(boundingBox1.intersects(boundingBox2))) { continue; } @@ -1458,16 +1460,17 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr auto n2 = *node_list1[1]; auto n3 = *node_list2[0]; auto n4 = *node_list2[1]; - //Check if the edges have already been added to our check_edges list - if (std::find(checked_edges.begin(), checked_edges.end(), n1) !=checked_edges.end() && - std::find(checked_edges.begin(), checked_edges.end(), n2) !=checked_edges.end() && - std::find(checked_edges.begin(), checked_edges.end(), n3) !=checked_edges.end() && - std::find(checked_edges.begin(), checked_edges.end(), n4) !=checked_edges.end()) + // Check if the edges have already been added to our check_edges list + if (std::find(checked_edges.begin(), checked_edges.end(), n1) != checked_edges.end() && + std::find(checked_edges.begin(), checked_edges.end(), n2) != checked_edges.end() && + std::find(checked_edges.begin(), checked_edges.end(), n3) != checked_edges.end() && + std::find(checked_edges.begin(), checked_edges.end(), n4) != checked_edges.end()) { continue; } // Now compare edge with other_edge - bool overlap = MeshBaseDiagnosticsUtils::checkEdgeOverlap(edge, other_edge, _console, _non_matching_edge_tol); + bool overlap = MeshBaseDiagnosticsUtils::checkEdgeOverlap( + edge, other_edge, _console, _non_matching_edge_tol); if (overlap) { // Add the nodes that make up the 2 edges to the vector checked_edges @@ -1475,13 +1478,15 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr checked_edges.push_back(n2); checked_edges.push_back(n3); checked_edges.push_back(n4); - num_intersecting_edges+=2; + num_intersecting_edges += 2; } } } } } - diagnosticsLog("Number of intersecting edges: " + Moose::stringify(num_intersecting_edges), _check_non_matching_edges, num_intersecting_edges); + diagnosticsLog("Number of intersecting edges: " + Moose::stringify(num_intersecting_edges), + _check_non_matching_edges, + num_intersecting_edges); } void diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 023d48234b66..c0ac8bd8d444 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -69,22 +69,22 @@ checkNonConformalMesh(const std::unique_ptr & mesh, } bool -checkEdgeOverlap(const std::unique_ptr & edge1, +checkEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, const ConsoleStream & console, const Real intersection_tol) { - //get node array from two edges + // get node array from two edges const auto node_list1 = edge1->get_nodes(); const auto node_list2 = edge2->get_nodes(); - //make sure two edges are not the same and don't share any nodes edge before checking for overlap + // make sure two edges are not the same and don't share any nodes edge before checking for overlap auto n1 = *node_list1[0]; auto n2 = *node_list1[1]; auto n3 = *node_list2[0]; auto n4 = *node_list2[1]; - //get x,y,z coordinates for each point + // get x,y,z coordinates for each point auto n1x = n1.operator()(0); auto n1y = n1.operator()(1); auto n1z = n1.operator()(2); @@ -98,16 +98,21 @@ checkEdgeOverlap(const std::unique_ptr & edge1, auto n4y = n4.operator()(1); auto n4z = n4.operator()(2); - //Check that none of these points are the same - if(((std::fabs(n1x - n3x) & edge1, auto denominator = d2121 * d4343 - d4321 * d4321; auto numerator = d1343 * d4321 - d1321 * d4343; - if(std::fabs(denominator) < intersection_tol) + if (std::fabs(denominator) < intersection_tol) { - //This indicates that the intersecting line is vertical so they don't intersect + // This indicates that the intersecting line is vertical so they don't intersect return false; } - auto mua = numerator/denominator; + auto mua = numerator / denominator; auto mub = (d1343 + (mua * d4321)) / d4343; - //Use these values to solve for the two poits that define the shortest line segment + // Use these values to solve for the two poits that define the shortest line segment auto nax = n1x + mua * (n2x - n1x); auto nay = n1y + mua * (n2y - n1y); auto naz = n1z + mua * (n2z - n1z); @@ -135,29 +140,30 @@ checkEdgeOverlap(const std::unique_ptr & edge1, auto nby = n3y + mub * (n4y - n3y); auto nbz = n3z + mub * (n4z - n3z); - //This method assume the two lines are infinite. This check to make sure na and nb are part of their respective line segments - if((nax < std::min(n1x, n2x)) || (nax > std::max(n1x, n2x)) || - (nay < std::min(n1y, n2y)) || (nay > std::max(n1y, n2y)) || - (naz < std::min(n1z, n2z)) || (naz > std::max(n1z, n2z))) + // This method assume the two lines are infinite. This check to make sure na and nb are part of + // their respective line segments + if ((nax < std::min(n1x, n2x)) || (nax > std::max(n1x, n2x)) || (nay < std::min(n1y, n2y)) || + (nay > std::max(n1y, n2y)) || (naz < std::min(n1z, n2z)) || (naz > std::max(n1z, n2z))) { return false; } - if((nbx < std::min(n3x, n4x)) || (nax > std::max(n3x, n4x)) || - (nby < std::min(n3y, n4y)) || (nay > std::max(n3y, n4y)) || - (nbz < std::min(n3z, n4z)) || (naz > std::max(n3z, n4z))) + if ((nbx < std::min(n3x, n4x)) || (nax > std::max(n3x, n4x)) || (nby < std::min(n3y, n4y)) || + (nay > std::max(n3y, n4y)) || (nbz < std::min(n3z, n4z)) || (naz > std::max(n3z, n4z))) { return false; } - //Calculate distance between these two nodes - double distance = std::sqrt(std::pow(nax - nbx, 2) + std::pow(nay - nby, 2) + std::pow(naz - nbz, 2)); + // Calculate distance between these two nodes + double distance = + std::sqrt(std::pow(nax - nbx, 2) + std::pow(nay - nby, 2) + std::pow(naz - nbz, 2)); if (distance < intersection_tol) { std::string x_coord = std::to_string(nax); std::string y_coord = std::to_string(nay); std::string z_coord = std::to_string(naz); - std::string message = "Non-matching edges found near (" + x_coord + ", " + y_coord + ", " + z_coord + ")"; + std::string message = + "Non-matching edges found near (" + x_coord + ", " + y_coord + ", " + z_coord + ")"; console << message << std::endl; return true; } From eb3bed8d3059be08197662a26e5cd0656f635ed9 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Wed, 25 Sep 2024 17:27:12 -0600 Subject: [PATCH 06/28] Made changes reccomended in PR -To check if two edges should be checked for overlap, node ids are compared to see if they share any nodes -checked_edges_nodes is now in a set -_console error message is now in MeshDiagnosticsGenerator.C and includes element ids -In test the two blocks are stitched so any touching nodes should be merged -Other minor improvements to make code more readable --- .../meshgenerators/MeshDiagnosticsGenerator.h | 2 +- .../include/utils/MeshBaseDiagnosticsUtils.h | 2 +- .../meshgenerators/MeshDiagnosticsGenerator.C | 53 +++++++----- .../src/utils/MeshBaseDiagnosticsUtils.C | 83 +++++++------------ .../intersecting_edge_test.i | 7 +- 5 files changed, 68 insertions(+), 79 deletions(-) diff --git a/framework/include/meshgenerators/MeshDiagnosticsGenerator.h b/framework/include/meshgenerators/MeshDiagnosticsGenerator.h index ca65b6a79143..b61eb5cb4d0d 100644 --- a/framework/include/meshgenerators/MeshDiagnosticsGenerator.h +++ b/framework/include/meshgenerators/MeshDiagnosticsGenerator.h @@ -49,7 +49,7 @@ class MeshDiagnosticsGenerator : public MeshGenerator void checkNonConformalMeshFromAdaptivity(const std::unique_ptr & mesh) const; /// Routine to check whether the Jacobians (elem and side) are not negative void checkLocalJacobians(const std::unique_ptr & mesh) const; - //// Routing to check for non matching edges + //// Routine to check for non matching edges void checkNonMatchingEdges(const std::unique_ptr & mesh) const; /** diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index 3e291c69dceb..f12add10aff7 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -27,6 +27,6 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, bool checkEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, - const ConsoleStream & console, + std::vector &intersection_coords, const Real insersection_tol); } diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 68385b048518..25c5ee440085 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1426,8 +1426,8 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr if (!mesh->is_serial()) mooseError("Only serialized/replicated meshes are supported"); unsigned int num_intersecting_edges = 0; - std::vector checked_edges; - for (auto elem : mesh->active_element_ptr_range()) + std::set checked_edges_nodes; + for (const auto elem : mesh->active_element_ptr_range()) { std::vector> elem_edges(elem->n_edges()); for (auto i : elem->edge_index_range()) @@ -1436,17 +1436,15 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr { // If they're the same element then there's no need to check for overlap if (elem == other_elem) - { continue; - } + // Get bounding boxes for both elements. If the bounding boxes don't intersect then no edges // will either auto boundingBox1 = elem->loose_bounding_box(); auto boundingBox2 = other_elem->loose_bounding_box(); if (!(boundingBox1.intersects(boundingBox2))) - { continue; - } + std::vector> other_edges(other_elem->n_edges()); for (auto j : other_elem->edge_index_range()) other_edges[j] = other_elem->build_edge_ptr(j); @@ -1454,31 +1452,42 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr { for (auto & other_edge : other_edges) { - const auto node_list1 = edge->get_nodes(); - const auto node_list2 = other_edge->get_nodes(); - auto n1 = *node_list1[0]; - auto n2 = *node_list1[1]; - auto n3 = *node_list2[0]; - auto n4 = *node_list2[1]; + // Get nodes from edges + const Node * n1 = edge->get_nodes()[0]; + const Node * n2 = edge->get_nodes()[1]; + const Node * n3 = other_edge->get_nodes()[0]; + const Node * n4 = other_edge->get_nodes()[1]; + // Check if the edges have already been added to our check_edges list - if (std::find(checked_edges.begin(), checked_edges.end(), n1) != checked_edges.end() && - std::find(checked_edges.begin(), checked_edges.end(), n2) != checked_edges.end() && - std::find(checked_edges.begin(), checked_edges.end(), n3) != checked_edges.end() && - std::find(checked_edges.begin(), checked_edges.end(), n4) != checked_edges.end()) + if (checked_edges_nodes.count(*n1) && + checked_edges_nodes.count(*n2) && + checked_edges_nodes.count(*n3) && + checked_edges_nodes.count(*n4)) { continue; } + // Now compare edge with other_edge + std::vector intersection_coords(3); bool overlap = MeshBaseDiagnosticsUtils::checkEdgeOverlap( - edge, other_edge, _console, _non_matching_edge_tol); + edge, other_edge, intersection_coords, _non_matching_edge_tol); if (overlap) { - // Add the nodes that make up the 2 edges to the vector checked_edges - checked_edges.push_back(n1); - checked_edges.push_back(n2); - checked_edges.push_back(n3); - checked_edges.push_back(n4); + // Add the nodes that make up the 2 edges to the vector checked_edges_nodes + checked_edges_nodes.insert(*n1); + checked_edges_nodes.insert(*n2); + checked_edges_nodes.insert(*n3); + checked_edges_nodes.insert(*n4); num_intersecting_edges += 2; + + // Print error message + std::string elem_id = std::to_string(elem->id()); + std::string other_elem_id = std::to_string(other_elem->id()); + std::string x_coord = std::to_string(intersection_coords[0]); + std::string y_coord = std::to_string(intersection_coords[1]); + std::string z_coord = std::to_string(intersection_coords[2]); + std::string message = "Intersecting edges found between elements " + elem_id + " and " + other_elem_id + " near point (" + x_coord + ", " + y_coord + ", " + z_coord + ")"; + _console << message << std::endl; } } } diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index c0ac8bd8d444..c341bc07325e 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -71,42 +71,32 @@ checkNonConformalMesh(const std::unique_ptr & mesh, bool checkEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, - const ConsoleStream & console, + std::vector &intersection_coords, const Real intersection_tol) { - // get node array from two edges - const auto node_list1 = edge1->get_nodes(); - const auto node_list2 = edge2->get_nodes(); - - // make sure two edges are not the same and don't share any nodes edge before checking for overlap - auto n1 = *node_list1[0]; - auto n2 = *node_list1[1]; - auto n3 = *node_list2[0]; - auto n4 = *node_list2[1]; - - // get x,y,z coordinates for each point - auto n1x = n1.operator()(0); - auto n1y = n1.operator()(1); - auto n1z = n1.operator()(2); - auto n2x = n2.operator()(0); - auto n2y = n2.operator()(1); - auto n2z = n2.operator()(2); - auto n3x = n3.operator()(0); - auto n3y = n3.operator()(1); - auto n3z = n3.operator()(2); - auto n4x = n4.operator()(0); - auto n4y = n4.operator()(1); - auto n4z = n4.operator()(2); + // Get nodes from the two edges + const Node * n1 = edge1->get_nodes()[0]; + const Node * n2 = edge1->get_nodes()[1]; + const Node * n3 = edge2->get_nodes()[0]; + const Node * n4 = edge2->get_nodes()[1]; + + // get x,y,z coordinates for each node + auto n1x = n1->operator()(0); + auto n1y = n1->operator()(1); + auto n1z = n1->operator()(2); + auto n2x = n2->operator()(0); + auto n2y = n2->operator()(1); + auto n2z = n2->operator()(2); + auto n3x = n3->operator()(0); + auto n3y = n3->operator()(1); + auto n3z = n3->operator()(2); + auto n4x = n4->operator()(0); + auto n4y = n4->operator()(1); + auto n4z = n4->operator()(2); // Check that none of these points are the same - if (((std::fabs(n1x - n3x) < intersection_tol) && (std::fabs(n1y - n3y) < intersection_tol) && - (std::fabs(n1z - n3z) < intersection_tol)) || - ((std::fabs(n2x - n4x) < intersection_tol) && (std::fabs(n2y - n4y) < intersection_tol) && - (std::fabs(n2z - n4z) < intersection_tol)) || - ((std::fabs(n1x - n4x) < intersection_tol) && (std::fabs(n1y - n4y) < intersection_tol) && - (std::fabs(n1z - n4z) < intersection_tol)) || - ((std::fabs(n2x - n3x) < intersection_tol) && (std::fabs(n2y - n3y) < intersection_tol) && - (std::fabs(n2z - n3z) < intersection_tol))) + if (n1->id() == n3->id() || n1->id() == n4->id() || + n2->id() == n3->id() || n2->id() == n4->id()) return false; /* @@ -124,14 +114,13 @@ checkEdgeOverlap(const std::unique_ptr & edge1, auto numerator = d1343 * d4321 - d1321 * d4343; if (std::fabs(denominator) < intersection_tol) - { - // This indicates that the intersecting line is vertical so they don't intersect + // This indicates that the two lines are parallel so they don't intersect return false; - } + auto mua = numerator / denominator; auto mub = (d1343 + (mua * d4321)) / d4343; - // Use these values to solve for the two poits that define the shortest line segment + // Use these values to solve for the two points that define the shortest line segment auto nax = n1x + mua * (n2x - n1x); auto nay = n1y + mua * (n2y - n1y); auto naz = n1z + mua * (n2z - n1z); @@ -142,29 +131,19 @@ checkEdgeOverlap(const std::unique_ptr & edge1, // This method assume the two lines are infinite. This check to make sure na and nb are part of // their respective line segments - if ((nax < std::min(n1x, n2x)) || (nax > std::max(n1x, n2x)) || (nay < std::min(n1y, n2y)) || - (nay > std::max(n1y, n2y)) || (naz < std::min(n1z, n2z)) || (naz > std::max(n1z, n2z))) - { + if (mua < 0 || mua > 1) return false; - } - - if ((nbx < std::min(n3x, n4x)) || (nax > std::max(n3x, n4x)) || (nby < std::min(n3y, n4y)) || - (nay > std::max(n3y, n4y)) || (nbz < std::min(n3z, n4z)) || (naz > std::max(n3z, n4z))) - { + if (mub < 0 || mub > 1) return false; - } - + // Calculate distance between these two nodes double distance = std::sqrt(std::pow(nax - nbx, 2) + std::pow(nay - nby, 2) + std::pow(naz - nbz, 2)); if (distance < intersection_tol) { - std::string x_coord = std::to_string(nax); - std::string y_coord = std::to_string(nay); - std::string z_coord = std::to_string(naz); - std::string message = - "Non-matching edges found near (" + x_coord + ", " + y_coord + ", " + z_coord + ")"; - console << message << std::endl; + intersection_coords[0] = nax; + intersection_coords[1] = nay; + intersection_coords[2] = naz; return true; } else diff --git a/test/tests/meshgenerators/mesh_diagnostics_generator/intersecting_edge_test.i b/test/tests/meshgenerators/mesh_diagnostics_generator/intersecting_edge_test.i index 4f5fdccd1cbf..73865b00ed57 100644 --- a/test/tests/meshgenerators/mesh_diagnostics_generator/intersecting_edge_test.i +++ b/test/tests/meshgenerators/mesh_diagnostics_generator/intersecting_edge_test.i @@ -51,13 +51,14 @@ transform = 'rotate' vector_value = '180 0 0' [] - [cmbn] - type = CombinerGenerator + [stitch] + type = StitchedMeshGenerator inputs = 'convert1 rotate' + stitch_boundaries_pairs = 'front back' [] [diag] type = MeshDiagnosticsGenerator - input = cmbn + input = stitch examine_non_matching_edges = INFO [] [] From 7252293bd7c29c6b0a6a4a14bd97d29af4b63c09 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Thu, 26 Sep 2024 07:57:47 -0600 Subject: [PATCH 07/28] Style changes to pass Civet tests --- framework/include/utils/MeshBaseDiagnosticsUtils.h | 2 +- .../src/meshgenerators/MeshDiagnosticsGenerator.C | 12 ++++++------ framework/src/utils/MeshBaseDiagnosticsUtils.C | 7 +++---- 3 files changed, 10 insertions(+), 11 deletions(-) diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index f12add10aff7..6396737ba1df 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -27,6 +27,6 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, bool checkEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, - std::vector &intersection_coords, + std::vector & intersection_coords, const Real insersection_tol); } diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 25c5ee440085..eebf7242814f 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1437,7 +1437,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr // If they're the same element then there's no need to check for overlap if (elem == other_elem) continue; - + // Get bounding boxes for both elements. If the bounding boxes don't intersect then no edges // will either auto boundingBox1 = elem->loose_bounding_box(); @@ -1459,10 +1459,8 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr const Node * n4 = other_edge->get_nodes()[1]; // Check if the edges have already been added to our check_edges list - if (checked_edges_nodes.count(*n1) && - checked_edges_nodes.count(*n2) && - checked_edges_nodes.count(*n3) && - checked_edges_nodes.count(*n4)) + if (checked_edges_nodes.count(*n1) && checked_edges_nodes.count(*n2) && + checked_edges_nodes.count(*n3) && checked_edges_nodes.count(*n4)) { continue; } @@ -1486,7 +1484,9 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr std::string x_coord = std::to_string(intersection_coords[0]); std::string y_coord = std::to_string(intersection_coords[1]); std::string z_coord = std::to_string(intersection_coords[2]); - std::string message = "Intersecting edges found between elements " + elem_id + " and " + other_elem_id + " near point (" + x_coord + ", " + y_coord + ", " + z_coord + ")"; + std::string message = "Intersecting edges found between elements " + elem_id + " and " + + other_elem_id + " near point (" + x_coord + ", " + y_coord + + ", " + z_coord + ")"; _console << message << std::endl; } } diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index c341bc07325e..68573c3edc09 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -71,7 +71,7 @@ checkNonConformalMesh(const std::unique_ptr & mesh, bool checkEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, - std::vector &intersection_coords, + std::vector & intersection_coords, const Real intersection_tol) { // Get nodes from the two edges @@ -95,8 +95,7 @@ checkEdgeOverlap(const std::unique_ptr & edge1, auto n4z = n4->operator()(2); // Check that none of these points are the same - if (n1->id() == n3->id() || n1->id() == n4->id() || - n2->id() == n3->id() || n2->id() == n4->id()) + if (n1->id() == n3->id() || n1->id() == n4->id() || n2->id() == n3->id() || n2->id() == n4->id()) return false; /* @@ -135,7 +134,7 @@ checkEdgeOverlap(const std::unique_ptr & edge1, return false; if (mub < 0 || mub > 1) return false; - + // Calculate distance between these two nodes double distance = std::sqrt(std::pow(nax - nbx, 2) + std::pow(nay - nby, 2) + std::pow(naz - nbz, 2)); From b5ede426d04a29620a4267959810143ff0ead24f Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Thu, 26 Sep 2024 11:36:28 -0600 Subject: [PATCH 08/28] Added detail tests file --- test/tests/meshgenerators/mesh_diagnostics_generator/tests | 1 + 1 file changed, 1 insertion(+) diff --git a/test/tests/meshgenerators/mesh_diagnostics_generator/tests b/test/tests/meshgenerators/mesh_diagnostics_generator/tests index 588aed3a88f9..18b8e89ab817 100644 --- a/test/tests/meshgenerators/mesh_diagnostics_generator/tests +++ b/test/tests/meshgenerators/mesh_diagnostics_generator/tests @@ -19,6 +19,7 @@ cli_args = '--mesh-only' mesh_mode = replicated expect_out = 'Number of intersecting edges: 4' + detail = 'edges intersecting' [] [inconsistent_sidesets] type = RunApp From fbec7dd05cc6887656f20a560542683c809628f6 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Thu, 26 Sep 2024 16:38:25 -0600 Subject: [PATCH 09/28] Made more changes for PR -Changed most variable types to const auto -Changed name of overlap check to checkFirstOrderEdgeOverlap -Added mooseAssert to check that elements being checked are of type EDGE2 -Updated comments and console output message --- .../include/utils/MeshBaseDiagnosticsUtils.h | 2 +- .../meshgenerators/MeshDiagnosticsGenerator.C | 10 +-- .../src/utils/MeshBaseDiagnosticsUtils.C | 77 ++++++++++--------- .../mesh_diagnostics_generator/tests | 2 +- 4 files changed, 48 insertions(+), 43 deletions(-) diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index 6396737ba1df..d2d6adbab1bc 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -25,7 +25,7 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, const Real conformality_tol, unsigned int & num_nonconformal_nodes); -bool checkEdgeOverlap(const std::unique_ptr & edge1, +bool checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, std::vector & intersection_coords, const Real insersection_tol); diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index eebf7242814f..2466050796a1 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1432,7 +1432,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr std::vector> elem_edges(elem->n_edges()); for (auto i : elem->edge_index_range()) elem_edges[i] = elem->build_edge_ptr(i); - for (auto other_elem : mesh->active_element_ptr_range()) + for (const auto other_elem : mesh->active_element_ptr_range()) { // If they're the same element then there's no need to check for overlap if (elem == other_elem) @@ -1440,8 +1440,8 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr // Get bounding boxes for both elements. If the bounding boxes don't intersect then no edges // will either - auto boundingBox1 = elem->loose_bounding_box(); - auto boundingBox2 = other_elem->loose_bounding_box(); + const auto boundingBox1 = elem->loose_bounding_box(); + const auto boundingBox2 = other_elem->loose_bounding_box(); if (!(boundingBox1.intersects(boundingBox2))) continue; @@ -1467,7 +1467,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr // Now compare edge with other_edge std::vector intersection_coords(3); - bool overlap = MeshBaseDiagnosticsUtils::checkEdgeOverlap( + bool overlap = MeshBaseDiagnosticsUtils::checkFirstOrderEdgeOverlap( edge, other_edge, intersection_coords, _non_matching_edge_tol); if (overlap) { @@ -1493,7 +1493,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr } } } - diagnosticsLog("Number of intersecting edges: " + Moose::stringify(num_intersecting_edges), + diagnosticsLog("Number of intersecting element edges: " + Moose::stringify(num_intersecting_edges), _check_non_matching_edges, num_intersecting_edges); } diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 68573c3edc09..64c04d7a1f8f 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -69,32 +69,37 @@ checkNonConformalMesh(const std::unique_ptr & mesh, } bool -checkEdgeOverlap(const std::unique_ptr & edge1, +checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, std::vector & intersection_coords, const Real intersection_tol) { + //check that the two elements are of type EDGE2 + mooseAssert(edge1->side_type(*(edge1->side_index_range()).end())==0, + "Elements must be of type EDGE2"); + mooseAssert(edge2->side_type(*(edge2->side_index_range()).end())==0, + "Elements must be of type EDGE2"); // Get nodes from the two edges - const Node * n1 = edge1->get_nodes()[0]; - const Node * n2 = edge1->get_nodes()[1]; - const Node * n3 = edge2->get_nodes()[0]; - const Node * n4 = edge2->get_nodes()[1]; + const Node * const n1 = edge1->get_nodes()[0]; + const Node * const n2 = edge1->get_nodes()[1]; + const Node * const n3 = edge2->get_nodes()[0]; + const Node * const n4 = edge2->get_nodes()[1]; // get x,y,z coordinates for each node - auto n1x = n1->operator()(0); - auto n1y = n1->operator()(1); - auto n1z = n1->operator()(2); - auto n2x = n2->operator()(0); - auto n2y = n2->operator()(1); - auto n2z = n2->operator()(2); - auto n3x = n3->operator()(0); - auto n3y = n3->operator()(1); - auto n3z = n3->operator()(2); - auto n4x = n4->operator()(0); - auto n4y = n4->operator()(1); - auto n4z = n4->operator()(2); - - // Check that none of these points are the same + const auto n1x = n1->operator()(0); + const auto n1y = n1->operator()(1); + const auto n1z = n1->operator()(2); + const auto n2x = n2->operator()(0); + const auto n2y = n2->operator()(1); + const auto n2z = n2->operator()(2); + const auto n3x = n3->operator()(0); + const auto n3y = n3->operator()(1); + const auto n3z = n3->operator()(2); + const auto n4x = n4->operator()(0); + const auto n4y = n4->operator()(1); + const auto n4z = n4->operator()(2); + + // Check that the two edges are not sharing a node if (n1->id() == n3->id() || n1->id() == n4->id() || n2->id() == n3->id() || n2->id() == n4->id()) return false; @@ -103,30 +108,30 @@ checkEdgeOverlap(const std::unique_ptr & edge1, is close enough to 0 return true The shortest line between the two edges will be perpendicular to both. */ - auto d1343 = (n1x - n3x) * (n4x - n3x) + (n1y - n3y) * (n4y - n3y) + (n1z - n3z) * (n4z - n3z); - auto d4321 = (n4x - n3x) * (n2x - n1x) + (n4y - n3y) * (n2y - n1y) + (n4z - n3z) * (n2z - n1z); - auto d1321 = (n1x - n3x) * (n2x - n1x) + (n1y - n3y) * (n2y - n1y) + (n1z - n3z) * (n2z - n1z); - auto d4343 = (n4x - n3x) * (n4x - n3x) + (n4y - n3y) * (n4y - n3y) + (n4z - n3z) * (n4z - n3z); - auto d2121 = (n2x - n1x) * (n2x - n1x) + (n2y - n1y) * (n2y - n1y) + (n2z - n1z) * (n2z - n1z); + const auto d1343 = (n1x - n3x) * (n4x - n3x) + (n1y - n3y) * (n4y - n3y) + (n1z - n3z) * (n4z - n3z); + const auto d4321 = (n4x - n3x) * (n2x - n1x) + (n4y - n3y) * (n2y - n1y) + (n4z - n3z) * (n2z - n1z); + const auto d1321 = (n1x - n3x) * (n2x - n1x) + (n1y - n3y) * (n2y - n1y) + (n1z - n3z) * (n2z - n1z); + const auto d4343 = (n4x - n3x) * (n4x - n3x) + (n4y - n3y) * (n4y - n3y) + (n4z - n3z) * (n4z - n3z); + const auto d2121 = (n2x - n1x) * (n2x - n1x) + (n2y - n1y) * (n2y - n1y) + (n2z - n1z) * (n2z - n1z); - auto denominator = d2121 * d4343 - d4321 * d4321; - auto numerator = d1343 * d4321 - d1321 * d4343; + const auto denominator = d2121 * d4343 - d4321 * d4321; + const auto numerator = d1343 * d4321 - d1321 * d4343; if (std::fabs(denominator) < intersection_tol) // This indicates that the two lines are parallel so they don't intersect return false; - auto mua = numerator / denominator; - auto mub = (d1343 + (mua * d4321)) / d4343; + const auto mua = numerator / denominator; + const auto mub = (d1343 + (mua * d4321)) / d4343; // Use these values to solve for the two points that define the shortest line segment - auto nax = n1x + mua * (n2x - n1x); - auto nay = n1y + mua * (n2y - n1y); - auto naz = n1z + mua * (n2z - n1z); + const auto nax = n1x + mua * (n2x - n1x); + const auto nay = n1y + mua * (n2y - n1y); + const auto naz = n1z + mua * (n2z - n1z); - auto nbx = n3x + mub * (n4x - n3x); - auto nby = n3y + mub * (n4y - n3y); - auto nbz = n3z + mub * (n4z - n3z); + const auto nbx = n3x + mub * (n4x - n3x); + const auto nby = n3y + mub * (n4y - n3y); + const auto nbz = n3z + mub * (n4z - n3z); // This method assume the two lines are infinite. This check to make sure na and nb are part of // their respective line segments @@ -136,8 +141,8 @@ checkEdgeOverlap(const std::unique_ptr & edge1, return false; // Calculate distance between these two nodes - double distance = - std::sqrt(std::pow(nax - nbx, 2) + std::pow(nay - nby, 2) + std::pow(naz - nbz, 2)); + const auto distance = + std::sqrt(Utility::pow<2>(nax - nbx) + Utility::pow<2>(nay - nby) + Utility::pow<2>(naz - nbz)); if (distance < intersection_tol) { intersection_coords[0] = nax; diff --git a/test/tests/meshgenerators/mesh_diagnostics_generator/tests b/test/tests/meshgenerators/mesh_diagnostics_generator/tests index 18b8e89ab817..ace59e1c956c 100644 --- a/test/tests/meshgenerators/mesh_diagnostics_generator/tests +++ b/test/tests/meshgenerators/mesh_diagnostics_generator/tests @@ -18,7 +18,7 @@ input = 'intersecting_edge_test.i' cli_args = '--mesh-only' mesh_mode = replicated - expect_out = 'Number of intersecting edges: 4' + expect_out = 'Number of intersecting element edges: 4' detail = 'edges intersecting' [] [inconsistent_sidesets] From dab7b62665c1123458088d9cc997ed9b460b2516 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Thu, 26 Sep 2024 16:49:18 -0600 Subject: [PATCH 10/28] More style changes needed to pass Civet tests --- .../include/utils/MeshBaseDiagnosticsUtils.h | 6 ++-- .../meshgenerators/MeshDiagnosticsGenerator.C | 3 +- .../src/utils/MeshBaseDiagnosticsUtils.C | 31 +++++++++++-------- 3 files changed, 23 insertions(+), 17 deletions(-) diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index d2d6adbab1bc..42e70a5922d2 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -26,7 +26,7 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, unsigned int & num_nonconformal_nodes); bool checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, - const std::unique_ptr & edge2, - std::vector & intersection_coords, - const Real insersection_tol); + const std::unique_ptr & edge2, + std::vector & intersection_coords, + const Real insersection_tol); } diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 2466050796a1..6f82f5cb3f54 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1493,7 +1493,8 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr } } } - diagnosticsLog("Number of intersecting element edges: " + Moose::stringify(num_intersecting_edges), + diagnosticsLog("Number of intersecting element edges: " + + Moose::stringify(num_intersecting_edges), _check_non_matching_edges, num_intersecting_edges); } diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 64c04d7a1f8f..a0f3715307dc 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -70,14 +70,14 @@ checkNonConformalMesh(const std::unique_ptr & mesh, bool checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, - const std::unique_ptr & edge2, - std::vector & intersection_coords, - const Real intersection_tol) + const std::unique_ptr & edge2, + std::vector & intersection_coords, + const Real intersection_tol) { - //check that the two elements are of type EDGE2 - mooseAssert(edge1->side_type(*(edge1->side_index_range()).end())==0, + // check that the two elements are of type EDGE2 + mooseAssert(edge1->side_type(*(edge1->side_index_range()).end()) == 0, "Elements must be of type EDGE2"); - mooseAssert(edge2->side_type(*(edge2->side_index_range()).end())==0, + mooseAssert(edge2->side_type(*(edge2->side_index_range()).end()) == 0, "Elements must be of type EDGE2"); // Get nodes from the two edges const Node * const n1 = edge1->get_nodes()[0]; @@ -108,11 +108,16 @@ checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, is close enough to 0 return true The shortest line between the two edges will be perpendicular to both. */ - const auto d1343 = (n1x - n3x) * (n4x - n3x) + (n1y - n3y) * (n4y - n3y) + (n1z - n3z) * (n4z - n3z); - const auto d4321 = (n4x - n3x) * (n2x - n1x) + (n4y - n3y) * (n2y - n1y) + (n4z - n3z) * (n2z - n1z); - const auto d1321 = (n1x - n3x) * (n2x - n1x) + (n1y - n3y) * (n2y - n1y) + (n1z - n3z) * (n2z - n1z); - const auto d4343 = (n4x - n3x) * (n4x - n3x) + (n4y - n3y) * (n4y - n3y) + (n4z - n3z) * (n4z - n3z); - const auto d2121 = (n2x - n1x) * (n2x - n1x) + (n2y - n1y) * (n2y - n1y) + (n2z - n1z) * (n2z - n1z); + const auto d1343 = + (n1x - n3x) * (n4x - n3x) + (n1y - n3y) * (n4y - n3y) + (n1z - n3z) * (n4z - n3z); + const auto d4321 = + (n4x - n3x) * (n2x - n1x) + (n4y - n3y) * (n2y - n1y) + (n4z - n3z) * (n2z - n1z); + const auto d1321 = + (n1x - n3x) * (n2x - n1x) + (n1y - n3y) * (n2y - n1y) + (n1z - n3z) * (n2z - n1z); + const auto d4343 = + (n4x - n3x) * (n4x - n3x) + (n4y - n3y) * (n4y - n3y) + (n4z - n3z) * (n4z - n3z); + const auto d2121 = + (n2x - n1x) * (n2x - n1x) + (n2y - n1y) * (n2y - n1y) + (n2z - n1z) * (n2z - n1z); const auto denominator = d2121 * d4343 - d4321 * d4321; const auto numerator = d1343 * d4321 - d1321 * d4343; @@ -141,8 +146,8 @@ checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, return false; // Calculate distance between these two nodes - const auto distance = - std::sqrt(Utility::pow<2>(nax - nbx) + Utility::pow<2>(nay - nby) + Utility::pow<2>(naz - nbz)); + const auto distance = std::sqrt(Utility::pow<2>(nax - nbx) + Utility::pow<2>(nay - nby) + + Utility::pow<2>(naz - nbz)); if (distance < intersection_tol) { intersection_coords[0] = nax; From ab7d1c31fcd61fc2c78faaf8485602882c8da674 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Thu, 26 Sep 2024 17:15:09 -0600 Subject: [PATCH 11/28] Fixed edge assertion to call edge->type() --- framework/src/utils/MeshBaseDiagnosticsUtils.C | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index a0f3715307dc..ca60ef0b3f32 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -75,9 +75,9 @@ checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, const Real intersection_tol) { // check that the two elements are of type EDGE2 - mooseAssert(edge1->side_type(*(edge1->side_index_range()).end()) == 0, + mooseAssert(edge1->type() == 0, "Elements must be of type EDGE2"); - mooseAssert(edge2->side_type(*(edge2->side_index_range()).end()) == 0, + mooseAssert(edge2->type() == 0, "Elements must be of type EDGE2"); // Get nodes from the two edges const Node * const n1 = edge1->get_nodes()[0]; From 72b3258dcc45be59d32f77d8f926523b1b7b98db Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Fri, 27 Sep 2024 09:56:57 -0600 Subject: [PATCH 12/28] Added edge check in MeshDiagnosticsGenerator -If edges are not of type EDGE2 a warning is printed out with edge and cell type -The edge intersection check is then skipped for that element --- .../meshgenerators/MeshDiagnosticsGenerator.C | 18 ++++++++++++++++++ framework/src/utils/MeshBaseDiagnosticsUtils.C | 6 ++---- 2 files changed, 20 insertions(+), 4 deletions(-) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 6f82f5cb3f54..78f3b5ec3f29 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1465,6 +1465,24 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr continue; } + //Check element/edge type + if (edge->type() != 0) + { + std::string element_message = "Edge of type " + Utility::enum_to_string(edge->type()) + " was found in cell " + std::to_string(elem->id()) + + " which is of type " + Utility::enum_to_string(elem->type()) + '\n' + + "The edge intersection check only works for EDGE2 elements"; + _console << element_message << std::endl; + continue; + } + if (other_edge->type() != 0) + { + std::string element_message = "Edge of type " + Utility::enum_to_string(other_edge->type()) + " was found in cell " + std::to_string(other_elem->id()) + + " which is of type " + Utility::enum_to_string(other_elem->type()) + '\n' + + "The edge intersection check only works for EDGE2 elements"; + _console << element_message << std::endl; + continue; + } + // Now compare edge with other_edge std::vector intersection_coords(3); bool overlap = MeshBaseDiagnosticsUtils::checkFirstOrderEdgeOverlap( diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index ca60ef0b3f32..535919ee7a59 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -75,10 +75,8 @@ checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, const Real intersection_tol) { // check that the two elements are of type EDGE2 - mooseAssert(edge1->type() == 0, - "Elements must be of type EDGE2"); - mooseAssert(edge2->type() == 0, - "Elements must be of type EDGE2"); + mooseAssert(edge1->type() == 0, "Elements must be of type EDGE2"); + mooseAssert(edge2->type() == 0, "Elements must be of type EDGE2"); // Get nodes from the two edges const Node * const n1 = edge1->get_nodes()[0]; const Node * const n2 = edge1->get_nodes()[1]; From 08747679f17593388f17ab45e89dca6ff532f751 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Fri, 27 Sep 2024 10:15:00 -0600 Subject: [PATCH 13/28] Changed edge check from type == 0 to type == EDGE2 --- framework/src/meshgenerators/MeshDiagnosticsGenerator.C | 4 ++-- framework/src/utils/MeshBaseDiagnosticsUtils.C | 4 ++-- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 78f3b5ec3f29..e351382f40ff 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1466,7 +1466,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr } //Check element/edge type - if (edge->type() != 0) + if (edge->type() != EDGE2) { std::string element_message = "Edge of type " + Utility::enum_to_string(edge->type()) + " was found in cell " + std::to_string(elem->id()) + " which is of type " + Utility::enum_to_string(elem->type()) + '\n' + @@ -1474,7 +1474,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr _console << element_message << std::endl; continue; } - if (other_edge->type() != 0) + if (other_edge->type() != EDGE2) { std::string element_message = "Edge of type " + Utility::enum_to_string(other_edge->type()) + " was found in cell " + std::to_string(other_elem->id()) + " which is of type " + Utility::enum_to_string(other_elem->type()) + '\n' + diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 535919ee7a59..2737add9c226 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -75,8 +75,8 @@ checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, const Real intersection_tol) { // check that the two elements are of type EDGE2 - mooseAssert(edge1->type() == 0, "Elements must be of type EDGE2"); - mooseAssert(edge2->type() == 0, "Elements must be of type EDGE2"); + mooseAssert(edge1->type() == EDGE2, "Elements must be of type EDGE2"); + mooseAssert(edge2->type() == EDGE2, "Elements must be of type EDGE2"); // Get nodes from the two edges const Node * const n1 = edge1->get_nodes()[0]; const Node * const n2 = edge1->get_nodes()[1]; From 33f1f5b182a179550d7459d7905bfdbd36b493ac Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Fri, 27 Sep 2024 10:23:40 -0600 Subject: [PATCH 14/28] More style changes for Civet tests --- .../meshgenerators/MeshDiagnosticsGenerator.C | 18 +++++++++++------- 1 file changed, 11 insertions(+), 7 deletions(-) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index e351382f40ff..97339abdf6fc 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1465,20 +1465,24 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr continue; } - //Check element/edge type + // Check element/edge type if (edge->type() != EDGE2) { - std::string element_message = "Edge of type " + Utility::enum_to_string(edge->type()) + " was found in cell " + std::to_string(elem->id()) + - " which is of type " + Utility::enum_to_string(elem->type()) + '\n' + - "The edge intersection check only works for EDGE2 elements"; + std::string element_message = + "Edge of type " + Utility::enum_to_string(edge->type()) + " was found in cell " + + std::to_string(elem->id()) + " which is of type " + + Utility::enum_to_string(elem->type()) + '\n' + + "The edge intersection check only works for EDGE2 elements"; _console << element_message << std::endl; continue; } if (other_edge->type() != EDGE2) { - std::string element_message = "Edge of type " + Utility::enum_to_string(other_edge->type()) + " was found in cell " + std::to_string(other_elem->id()) + - " which is of type " + Utility::enum_to_string(other_elem->type()) + '\n' + - "The edge intersection check only works for EDGE2 elements"; + std::string element_message = + "Edge of type " + Utility::enum_to_string(other_edge->type()) + + " was found in cell " + std::to_string(other_elem->id()) + " which is of type " + + Utility::enum_to_string(other_elem->type()) + '\n' + + "The edge intersection check only works for EDGE2 elements"; _console << element_message << std::endl; continue; } From 4d9f6e06c3410d69a64e80fcabb693203c2b9237 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Fri, 27 Sep 2024 15:54:46 -0600 Subject: [PATCH 15/28] Modified wrong edge check so that it only prints message once --- .../meshgenerators/MeshDiagnosticsGenerator.C | 21 +++++++------------ 1 file changed, 7 insertions(+), 14 deletions(-) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 97339abdf6fc..5ab866883010 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1468,24 +1468,17 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr // Check element/edge type if (edge->type() != EDGE2) { - std::string element_message = - "Edge of type " + Utility::enum_to_string(edge->type()) + " was found in cell " + - std::to_string(elem->id()) + " which is of type " + - Utility::enum_to_string(elem->type()) + '\n' + - "The edge intersection check only works for EDGE2 elements"; - _console << element_message << std::endl; + std::string element_message = "Edge of type " + Utility::enum_to_string(edge->type()) + + " was found in cell " + std::to_string(elem->id()) + + " which is of type " + + Utility::enum_to_string(elem->type()) + '\n' + + "The edge intersection check only works for EDGE2 " + "elements.\nThis message will not be output again"; + mooseDoOnce(_console << element_message << std::endl); continue; } if (other_edge->type() != EDGE2) - { - std::string element_message = - "Edge of type " + Utility::enum_to_string(other_edge->type()) + - " was found in cell " + std::to_string(other_elem->id()) + " which is of type " + - Utility::enum_to_string(other_elem->type()) + '\n' + - "The edge intersection check only works for EDGE2 elements"; - _console << element_message << std::endl; continue; - } // Now compare edge with other_edge std::vector intersection_coords(3); From a245b89e264ad9c7add2c571a6316b8bd6f62772 Mon Sep 17 00:00:00 2001 From: DanielYankura Date: Mon, 30 Sep 2024 08:02:44 -0600 Subject: [PATCH 16/28] changed checked_edges_nodes to overlapping_edges_nodes --- .../meshgenerators/MeshDiagnosticsGenerator.C | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 5ab866883010..c0e3275efe23 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1426,7 +1426,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr if (!mesh->is_serial()) mooseError("Only serialized/replicated meshes are supported"); unsigned int num_intersecting_edges = 0; - std::set checked_edges_nodes; + std::set overlapping_edges_nodes; for (const auto elem : mesh->active_element_ptr_range()) { std::vector> elem_edges(elem->n_edges()); @@ -1459,8 +1459,8 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr const Node * n4 = other_edge->get_nodes()[1]; // Check if the edges have already been added to our check_edges list - if (checked_edges_nodes.count(*n1) && checked_edges_nodes.count(*n2) && - checked_edges_nodes.count(*n3) && checked_edges_nodes.count(*n4)) + if (overlapping_edges_nodes.count(*n1) && overlapping_edges_nodes.count(*n2) && + overlapping_edges_nodes.count(*n3) && overlapping_edges_nodes.count(*n4)) { continue; } @@ -1486,11 +1486,11 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr edge, other_edge, intersection_coords, _non_matching_edge_tol); if (overlap) { - // Add the nodes that make up the 2 edges to the vector checked_edges_nodes - checked_edges_nodes.insert(*n1); - checked_edges_nodes.insert(*n2); - checked_edges_nodes.insert(*n3); - checked_edges_nodes.insert(*n4); + // Add the nodes that make up the 2 edges to the vector overlapping_edges_nodes + overlapping_edges_nodes.insert(*n1); + overlapping_edges_nodes.insert(*n2); + overlapping_edges_nodes.insert(*n3); + overlapping_edges_nodes.insert(*n4); num_intersecting_edges += 2; // Print error message From c1704883f9cc3bf55fc16737a23971cdb1c0e3d5 Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Tue, 29 Oct 2024 10:40:12 -0600 Subject: [PATCH 17/28] added string_to_enum.h reference for Civet test --- framework/src/meshgenerators/MeshDiagnosticsGenerator.C | 1 + 1 file changed, 1 insertion(+) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index c0e3275efe23..1b40c7d0aac8 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -20,6 +20,7 @@ #include "libmesh/cell_tet4.h" #include "libmesh/face_quad4.h" #include "libmesh/cell_hex8.h" +#include "libmesh/string_to_enum.h" registerMooseObject("MooseApp", MeshDiagnosticsGenerator); From c9a78e947586c1062a18f420706670be21534847 Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Tue, 29 Oct 2024 14:01:49 -0600 Subject: [PATCH 18/28] Included utility.h for civet test --- framework/src/utils/MeshBaseDiagnosticsUtils.C | 1 + 1 file changed, 1 insertion(+) diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 2737add9c226..8da5a0aaf095 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -16,6 +16,7 @@ #include "libmesh/node.h" #include "libmesh/mesh_base.h" #include "libmesh/point_locator_base.h" +#include "libmesh/utility.h" namespace MeshBaseDiagnosticsUtils { From 9ed4ee9cdd1f38d699cb503b5c3f27720c87783d Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Thu, 14 Nov 2024 16:47:40 -0700 Subject: [PATCH 19/28] changed overlapping nodes set to take an id instead of a node to avoid libmesh issues --- .../src/meshgenerators/MeshDiagnosticsGenerator.C | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 1b40c7d0aac8..9bab900cf385 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1427,7 +1427,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr if (!mesh->is_serial()) mooseError("Only serialized/replicated meshes are supported"); unsigned int num_intersecting_edges = 0; - std::set overlapping_edges_nodes; + std::set overlapping_edges_nodes; for (const auto elem : mesh->active_element_ptr_range()) { std::vector> elem_edges(elem->n_edges()); @@ -1460,8 +1460,8 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr const Node * n4 = other_edge->get_nodes()[1]; // Check if the edges have already been added to our check_edges list - if (overlapping_edges_nodes.count(*n1) && overlapping_edges_nodes.count(*n2) && - overlapping_edges_nodes.count(*n3) && overlapping_edges_nodes.count(*n4)) + if (overlapping_edges_nodes.count(n1->id()) && overlapping_edges_nodes.count(n2->id()) && + overlapping_edges_nodes.count(n3->id()) && overlapping_edges_nodes.count(n4->id())) { continue; } @@ -1488,10 +1488,10 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr if (overlap) { // Add the nodes that make up the 2 edges to the vector overlapping_edges_nodes - overlapping_edges_nodes.insert(*n1); - overlapping_edges_nodes.insert(*n2); - overlapping_edges_nodes.insert(*n3); - overlapping_edges_nodes.insert(*n4); + overlapping_edges_nodes.insert(n1->id()); + overlapping_edges_nodes.insert(n2->id()); + overlapping_edges_nodes.insert(n3->id()); + overlapping_edges_nodes.insert(n4->id()); num_intersecting_edges += 2; // Print error message From fa7d61c127589c411c6c467028e311133aebb1b4 Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Tue, 19 Nov 2024 16:36:49 -0700 Subject: [PATCH 20/28] Made changes to improve efficiency and eliminate risk of false positives/negatives -Bounding boxes are now calculated in their own loop and added to a map -This map gets called later when they are needed -Intersecting edges are stored as a set of arrays with 4 nodes each -This prevents false positives/negatives when nodes may be part of multiple intersecting edges -Made changes to loop so it doesn't calculate distance between the same edges twice -Replaced Node calculations with Point instead to improve readability --- .../include/utils/MeshBaseDiagnosticsUtils.h | 4 +- .../meshgenerators/MeshDiagnosticsGenerator.C | 39 ++++++++----- .../src/utils/MeshBaseDiagnosticsUtils.C | 57 +++++-------------- 3 files changed, 41 insertions(+), 59 deletions(-) diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index 42e70a5922d2..14eae1d71587 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -27,6 +27,6 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, bool checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, - std::vector & intersection_coords, - const Real insersection_tol); + Point & intersection_point, + const Real intersection_tol); } diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 9bab900cf385..c5d0c0f05989 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1411,7 +1411,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr { /*Algorithm Overview 1)Prechecks - a)This algorithn only works for 3D so check for that first + a)This algorithm only works for 3D so check for that first 2)Loop a)Loop through every element b)For each element get the edges associated with it @@ -1427,22 +1427,31 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr if (!mesh->is_serial()) mooseError("Only serialized/replicated meshes are supported"); unsigned int num_intersecting_edges = 0; - std::set overlapping_edges_nodes; + + // Create map of element to bounding box to avoing reinitializing the same bounding box multiple times + std::unordered_map bounding_box_map; + for (const auto elem : mesh->active_element_ptr_range()) + { + const auto boundingBox = elem->loose_bounding_box(); + bounding_box_map.insert({elem, boundingBox}); + } + + std::set> overlapping_edges_nodes; for (const auto elem : mesh->active_element_ptr_range()) { std::vector> elem_edges(elem->n_edges()); for (auto i : elem->edge_index_range()) elem_edges[i] = elem->build_edge_ptr(i); + const auto boundingBox1 = bounding_box_map[elem]; for (const auto other_elem : mesh->active_element_ptr_range()) { // If they're the same element then there's no need to check for overlap - if (elem == other_elem) + if (elem->id() >= other_elem->id()) continue; // Get bounding boxes for both elements. If the bounding boxes don't intersect then no edges // will either - const auto boundingBox1 = elem->loose_bounding_box(); - const auto boundingBox2 = other_elem->loose_bounding_box(); + const auto boundingBox2 = bounding_box_map[other_elem]; if (!(boundingBox1.intersects(boundingBox2))) continue; @@ -1459,9 +1468,12 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr const Node * n3 = other_edge->get_nodes()[0]; const Node * n4 = other_edge->get_nodes()[1]; + // Create array to check against set + std::array node_id_array = {n1->id(), n2->id(), n3->id(), n4->id()}; + std::sort(node_id_array.begin(), node_id_array.end()); + // Check if the edges have already been added to our check_edges list - if (overlapping_edges_nodes.count(n1->id()) && overlapping_edges_nodes.count(n2->id()) && - overlapping_edges_nodes.count(n3->id()) && overlapping_edges_nodes.count(n4->id())) + if (overlapping_edges_nodes.count(node_id_array)) { continue; } @@ -1482,24 +1494,21 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr continue; // Now compare edge with other_edge - std::vector intersection_coords(3); + Point intersection_coords; bool overlap = MeshBaseDiagnosticsUtils::checkFirstOrderEdgeOverlap( edge, other_edge, intersection_coords, _non_matching_edge_tol); if (overlap) { // Add the nodes that make up the 2 edges to the vector overlapping_edges_nodes - overlapping_edges_nodes.insert(n1->id()); - overlapping_edges_nodes.insert(n2->id()); - overlapping_edges_nodes.insert(n3->id()); - overlapping_edges_nodes.insert(n4->id()); + overlapping_edges_nodes.insert(node_id_array); num_intersecting_edges += 2; // Print error message std::string elem_id = std::to_string(elem->id()); std::string other_elem_id = std::to_string(other_elem->id()); - std::string x_coord = std::to_string(intersection_coords[0]); - std::string y_coord = std::to_string(intersection_coords[1]); - std::string z_coord = std::to_string(intersection_coords[2]); + std::string x_coord = std::to_string(intersection_coords(0)); + std::string y_coord = std::to_string(intersection_coords(1)); + std::string z_coord = std::to_string(intersection_coords(2)); std::string message = "Intersecting edges found between elements " + elem_id + " and " + other_elem_id + " near point (" + x_coord + ", " + y_coord + ", " + z_coord + ")"; diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 8da5a0aaf095..950f96971dd7 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -72,34 +72,20 @@ checkNonConformalMesh(const std::unique_ptr & mesh, bool checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, - std::vector & intersection_coords, + Point & intersection_point, const Real intersection_tol) { // check that the two elements are of type EDGE2 mooseAssert(edge1->type() == EDGE2, "Elements must be of type EDGE2"); mooseAssert(edge2->type() == EDGE2, "Elements must be of type EDGE2"); // Get nodes from the two edges - const Node * const n1 = edge1->get_nodes()[0]; - const Node * const n2 = edge1->get_nodes()[1]; - const Node * const n3 = edge2->get_nodes()[0]; - const Node * const n4 = edge2->get_nodes()[1]; - - // get x,y,z coordinates for each node - const auto n1x = n1->operator()(0); - const auto n1y = n1->operator()(1); - const auto n1z = n1->operator()(2); - const auto n2x = n2->operator()(0); - const auto n2y = n2->operator()(1); - const auto n2z = n2->operator()(2); - const auto n3x = n3->operator()(0); - const auto n3y = n3->operator()(1); - const auto n3z = n3->operator()(2); - const auto n4x = n4->operator()(0); - const auto n4y = n4->operator()(1); - const auto n4z = n4->operator()(2); + const Point & p1 = edge1->point(0); + const Point & p2 = edge1->point(1); + const Point & p3 = edge2->point(0); + const Point & p4 = edge2->point(1); // Check that the two edges are not sharing a node - if (n1->id() == n3->id() || n1->id() == n4->id() || n2->id() == n3->id() || n2->id() == n4->id()) + if (p1 == p3 || p1 == p4 || p2 == p3 || p2 == p4) return false; /* @@ -107,16 +93,11 @@ checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, is close enough to 0 return true The shortest line between the two edges will be perpendicular to both. */ - const auto d1343 = - (n1x - n3x) * (n4x - n3x) + (n1y - n3y) * (n4y - n3y) + (n1z - n3z) * (n4z - n3z); - const auto d4321 = - (n4x - n3x) * (n2x - n1x) + (n4y - n3y) * (n2y - n1y) + (n4z - n3z) * (n2z - n1z); - const auto d1321 = - (n1x - n3x) * (n2x - n1x) + (n1y - n3y) * (n2y - n1y) + (n1z - n3z) * (n2z - n1z); - const auto d4343 = - (n4x - n3x) * (n4x - n3x) + (n4y - n3y) * (n4y - n3y) + (n4z - n3z) * (n4z - n3z); - const auto d2121 = - (n2x - n1x) * (n2x - n1x) + (n2y - n1y) * (n2y - n1y) + (n2z - n1z) * (n2z - n1z); + const auto d1343 = (p1-p3)*(p4-p3); + const auto d4321 = (p4-p3)*(p2-p1); + const auto d1321 = (p1-p3)*(p2-p1); + const auto d4343 = (p4-p3)*(p4-p3); + const auto d2121 = (p2-p1)*(p2-p1); const auto denominator = d2121 * d4343 - d4321 * d4321; const auto numerator = d1343 * d4321 - d1321 * d4343; @@ -129,13 +110,8 @@ checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, const auto mub = (d1343 + (mua * d4321)) / d4343; // Use these values to solve for the two points that define the shortest line segment - const auto nax = n1x + mua * (n2x - n1x); - const auto nay = n1y + mua * (n2y - n1y); - const auto naz = n1z + mua * (n2z - n1z); - - const auto nbx = n3x + mub * (n4x - n3x); - const auto nby = n3y + mub * (n4y - n3y); - const auto nbz = n3z + mub * (n4z - n3z); + const auto pa = p1 + mua * (p2 - p1); + const auto pb = p3 + mub * (p4 - p3); // This method assume the two lines are infinite. This check to make sure na and nb are part of // their respective line segments @@ -145,13 +121,10 @@ checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, return false; // Calculate distance between these two nodes - const auto distance = std::sqrt(Utility::pow<2>(nax - nbx) + Utility::pow<2>(nay - nby) + - Utility::pow<2>(naz - nbz)); + const auto distance = (pa - pb).norm(); if (distance < intersection_tol) { - intersection_coords[0] = nax; - intersection_coords[1] = nay; - intersection_coords[2] = naz; + intersection_point = pa; return true; } else From cee9baee132288636e92a49ffc3fcf80943d35ca Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Tue, 19 Nov 2024 16:49:56 -0700 Subject: [PATCH 21/28] style changes for Civet --- .../src/meshgenerators/MeshDiagnosticsGenerator.C | 5 +++-- framework/src/utils/MeshBaseDiagnosticsUtils.C | 10 +++++----- 2 files changed, 8 insertions(+), 7 deletions(-) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index c5d0c0f05989..4882a00f8586 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1428,13 +1428,14 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr mooseError("Only serialized/replicated meshes are supported"); unsigned int num_intersecting_edges = 0; - // Create map of element to bounding box to avoing reinitializing the same bounding box multiple times + // Create map of element to bounding box to avoing reinitializing the same bounding box multiple + // times std::unordered_map bounding_box_map; for (const auto elem : mesh->active_element_ptr_range()) { const auto boundingBox = elem->loose_bounding_box(); bounding_box_map.insert({elem, boundingBox}); - } + } std::set> overlapping_edges_nodes; for (const auto elem : mesh->active_element_ptr_range()) diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 950f96971dd7..6ce683aab8b6 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -93,11 +93,11 @@ checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, is close enough to 0 return true The shortest line between the two edges will be perpendicular to both. */ - const auto d1343 = (p1-p3)*(p4-p3); - const auto d4321 = (p4-p3)*(p2-p1); - const auto d1321 = (p1-p3)*(p2-p1); - const auto d4343 = (p4-p3)*(p4-p3); - const auto d2121 = (p2-p1)*(p2-p1); + const auto d1343 = (p1 - p3) * (p4 - p3); + const auto d4321 = (p4 - p3) * (p2 - p1); + const auto d1321 = (p1 - p3) * (p2 - p1); + const auto d4343 = (p4 - p3) * (p4 - p3); + const auto d2121 = (p2 - p1) * (p2 - p1); const auto denominator = d2121 * d4343 - d4321 * d4321; const auto numerator = d1343 * d4321 - d1321 * d4343; From b8adcbdb49573482aca5846adcaa788fda3cfa7d Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Wed, 20 Nov 2024 08:39:43 -0700 Subject: [PATCH 22/28] Removed utility.h because it is no longer needed for calculating distance --- framework/src/utils/MeshBaseDiagnosticsUtils.C | 1 - 1 file changed, 1 deletion(-) diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 6ce683aab8b6..b10d787ca733 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -16,7 +16,6 @@ #include "libmesh/node.h" #include "libmesh/mesh_base.h" #include "libmesh/point_locator_base.h" -#include "libmesh/utility.h" namespace MeshBaseDiagnosticsUtils { From ced523374bc6a6bbca69732e9dd12aae985b7d93 Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Wed, 27 Nov 2024 16:05:09 -0700 Subject: [PATCH 23/28] Odd issue with merge conflicts. This should fix it --- framework/src/meshgenerators/MeshDiagnosticsGenerator.C | 9 +-------- 1 file changed, 1 insertion(+), 8 deletions(-) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 4882a00f8586..b6f401ec329f 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -115,18 +115,11 @@ MeshDiagnosticsGenerator::MeshDiagnosticsGenerator(const InputParameters & param paramError("examine_non_conformality", "You must set this parameter to true to trigger mesh conformality check"); if (_check_sidesets_orientation == "NO_CHECK" && _check_watertight_sidesets == "NO_CHECK" && -<<<<<<< HEAD _check_watertight_nodesets == "NO_CHECK" && _check_element_volumes == "NO_CHECK" && _check_element_types == "NO_CHECK" && _check_element_overlap == "NO_CHECK" && _check_non_planar_sides == "NO_CHECK" && _check_non_conformal_mesh == "NO_CHECK" && - _check_adaptivity_non_conformality == "NO_CHECK" && _check_local_jacobian == "NO_CHECK") -======= - _check_watertight_nodesets == "NO_CHECK" && _check_element_volumes == "NO_CHECK" && - _check_element_types == "NO_CHECK" && _check_element_overlap == "NO_CHECK" && - _check_non_planar_sides == "NO_CHECK" && _check_non_conformal_mesh == "NO_CHECK" && - _check_adaptivity_non_conformality == "NO_CHECK" && _check_local_jacobian == "NO_CHECK" && + _check_adaptivity_non_conformality == "NO_CHECK" && _check_local_jacobian == "NO_CHECK" && _check_non_matching_edges == "NO_CHECK") ->>>>>>> a759b2e1ae (Basic outline for checkNonMatchingEdges) mooseError("You need to turn on at least one diagnostic. Did you misspell a parameter?"); } From b4b02cefddfa7f9d3bcd2afcd8a29a6e7dd7f48d Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Wed, 27 Nov 2024 16:11:07 -0700 Subject: [PATCH 24/28] Whitespace issue fixed for civet --- .../meshgenerators/mesh_diagnostics_generator/all_at_once.i | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/tests/meshgenerators/mesh_diagnostics_generator/all_at_once.i b/test/tests/meshgenerators/mesh_diagnostics_generator/all_at_once.i index 236a2d75c325..803cd0e2367d 100644 --- a/test/tests/meshgenerators/mesh_diagnostics_generator/all_at_once.i +++ b/test/tests/meshgenerators/mesh_diagnostics_generator/all_at_once.i @@ -15,7 +15,7 @@ examine_nonplanar_sides = INFO examine_sidesets_orientation = WARNING check_for_watertight_sidesets = WARNING - check_for_watertight_nodesets = WARNING + check_for_watertight_nodesets = WARNING examine_non_matching_edges = WARNING search_for_adaptivity_nonconformality = WARNING check_local_jacobian = WARNING From 77358b5709fa5a536b8ea9c375b0594074824f77 Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Mon, 2 Dec 2024 10:48:53 -0700 Subject: [PATCH 25/28] Improved efficiency of nonmatching edge check -Should run in O(n log n) -Nested for loop looks for nearby elements and loops through only those -No longer had to check each element against every other element --- .../include/utils/MeshBaseDiagnosticsUtils.h | 2 +- .../meshgenerators/MeshDiagnosticsGenerator.C | 45 ++++++++++--------- .../src/utils/MeshBaseDiagnosticsUtils.C | 2 +- 3 files changed, 27 insertions(+), 22 deletions(-) diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index 14eae1d71587..d427428636d5 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -26,7 +26,7 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, unsigned int & num_nonconformal_nodes); bool checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, - const std::unique_ptr & edge2, + const std::unique_ptr & edge2, Point & intersection_point, const Real intersection_tol); } diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index b6f401ec329f..4b01accb388b 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -21,6 +21,7 @@ #include "libmesh/face_quad4.h" #include "libmesh/cell_hex8.h" #include "libmesh/string_to_enum.h" +#include "libmesh/enum_point_locator_type.h" registerMooseObject("MooseApp", MeshDiagnosticsGenerator); @@ -1430,26 +1431,28 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr bounding_box_map.insert({elem, boundingBox}); } + std::unique_ptr point_locator = PointLocatorBase::build(TREE_ELEMENTS, *mesh); std::set> overlapping_edges_nodes; for (const auto elem : mesh->active_element_ptr_range()) { + // loop through elem's nodes and find nearby elements with it + std::set candidate_elems; + std::set nearby_elems; + for (unsigned int i = 0; i < elem->n_nodes(); i++) + { + (*point_locator)(elem->point(i), candidate_elems); + nearby_elems.insert(candidate_elems.begin(), candidate_elems.end()); + } std::vector> elem_edges(elem->n_edges()); for (auto i : elem->edge_index_range()) elem_edges[i] = elem->build_edge_ptr(i); - const auto boundingBox1 = bounding_box_map[elem]; - for (const auto other_elem : mesh->active_element_ptr_range()) + for (const auto other_elem : nearby_elems) { // If they're the same element then there's no need to check for overlap if (elem->id() >= other_elem->id()) continue; - // Get bounding boxes for both elements. If the bounding boxes don't intersect then no edges - // will either - const auto boundingBox2 = bounding_box_map[other_elem]; - if (!(boundingBox1.intersects(boundingBox2))) - continue; - - std::vector> other_edges(other_elem->n_edges()); + std::vector> other_edges(other_elem->n_edges()); for (auto j : other_elem->edge_index_range()) other_edges[j] = other_elem->build_edge_ptr(j); for (auto & edge : elem_edges) @@ -1496,17 +1499,19 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr // Add the nodes that make up the 2 edges to the vector overlapping_edges_nodes overlapping_edges_nodes.insert(node_id_array); num_intersecting_edges += 2; - - // Print error message - std::string elem_id = std::to_string(elem->id()); - std::string other_elem_id = std::to_string(other_elem->id()); - std::string x_coord = std::to_string(intersection_coords(0)); - std::string y_coord = std::to_string(intersection_coords(1)); - std::string z_coord = std::to_string(intersection_coords(2)); - std::string message = "Intersecting edges found between elements " + elem_id + " and " + - other_elem_id + " near point (" + x_coord + ", " + y_coord + - ", " + z_coord + ")"; - _console << message << std::endl; + if (num_intersecting_edges < _num_outputs) + { + // Print error message + std::string elem_id = std::to_string(elem->id()); + std::string other_elem_id = std::to_string(other_elem->id()); + std::string x_coord = std::to_string(intersection_coords(0)); + std::string y_coord = std::to_string(intersection_coords(1)); + std::string z_coord = std::to_string(intersection_coords(2)); + std::string message = "Intersecting edges found between elements " + elem_id + " and " + + other_elem_id + " near point (" + x_coord + ", " + y_coord + + ", " + z_coord + ")"; + _console << message << std::endl; + } } } } diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index b10d787ca733..8dd88e672b6b 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -70,7 +70,7 @@ checkNonConformalMesh(const std::unique_ptr & mesh, bool checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, - const std::unique_ptr & edge2, + const std::unique_ptr & edge2, Point & intersection_point, const Real intersection_tol) { From b1453218de42a98007be359940f9f8f5b1079edb Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Mon, 2 Dec 2024 10:55:16 -0700 Subject: [PATCH 26/28] Style changes for Civet --- framework/src/meshgenerators/MeshDiagnosticsGenerator.C | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 4b01accb388b..818dcf0e35ce 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1507,9 +1507,9 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr std::string x_coord = std::to_string(intersection_coords(0)); std::string y_coord = std::to_string(intersection_coords(1)); std::string z_coord = std::to_string(intersection_coords(2)); - std::string message = "Intersecting edges found between elements " + elem_id + " and " + - other_elem_id + " near point (" + x_coord + ", " + y_coord + - ", " + z_coord + ")"; + std::string message = "Intersecting edges found between elements " + elem_id + + " and " + other_elem_id + " near point (" + x_coord + ", " + + y_coord + ", " + z_coord + ")"; _console << message << std::endl; } } From 815a28f52a45a6466430105810c4e121594ebf40 Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Mon, 27 Jan 2025 17:51:25 -0700 Subject: [PATCH 27/28] Minor changes for PR - point_locator is now initialized with sub_point_locator() - edge1 is now const to match edge2 --- framework/include/utils/MeshBaseDiagnosticsUtils.h | 2 +- framework/src/meshgenerators/MeshDiagnosticsGenerator.C | 5 +++-- framework/src/utils/MeshBaseDiagnosticsUtils.C | 2 +- 3 files changed, 5 insertions(+), 4 deletions(-) diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index d427428636d5..17494fe0278e 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -25,7 +25,7 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, const Real conformality_tol, unsigned int & num_nonconformal_nodes); -bool checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, +bool checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, Point & intersection_point, const Real intersection_tol); diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index 818dcf0e35ce..a8196f8eed52 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1431,7 +1431,8 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr bounding_box_map.insert({elem, boundingBox}); } - std::unique_ptr point_locator = PointLocatorBase::build(TREE_ELEMENTS, *mesh); + std::unique_ptr point_locator = + mesh->sub_point_locator(); // PointLocatorBase::build(TREE_ELEMENTS, *mesh); std::set> overlapping_edges_nodes; for (const auto elem : mesh->active_element_ptr_range()) { @@ -1443,7 +1444,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr (*point_locator)(elem->point(i), candidate_elems); nearby_elems.insert(candidate_elems.begin(), candidate_elems.end()); } - std::vector> elem_edges(elem->n_edges()); + std::vector> elem_edges(elem->n_edges()); for (auto i : elem->edge_index_range()) elem_edges[i] = elem->build_edge_ptr(i); for (const auto other_elem : nearby_elems) diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 8dd88e672b6b..3960797cdac3 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -69,7 +69,7 @@ checkNonConformalMesh(const std::unique_ptr & mesh, } bool -checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, +checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, const std::unique_ptr & edge2, Point & intersection_point, const Real intersection_tol) From 4d278c6bced09292dfb8c4cd2049906373a79c6b Mon Sep 17 00:00:00 2001 From: Daniel Yankura Date: Tue, 28 Jan 2025 10:20:55 -0700 Subject: [PATCH 28/28] Change for PR -Removed commented out code -checkFirstOrderEdgeOverlap now takes a const Elem instead of a uniqu_ptr --- .../include/utils/MeshBaseDiagnosticsUtils.h | 4 ++-- .../meshgenerators/MeshDiagnosticsGenerator.C | 5 ++--- framework/src/utils/MeshBaseDiagnosticsUtils.C | 16 ++++++++-------- 3 files changed, 12 insertions(+), 13 deletions(-) diff --git a/framework/include/utils/MeshBaseDiagnosticsUtils.h b/framework/include/utils/MeshBaseDiagnosticsUtils.h index 17494fe0278e..830f46c36c7f 100644 --- a/framework/include/utils/MeshBaseDiagnosticsUtils.h +++ b/framework/include/utils/MeshBaseDiagnosticsUtils.h @@ -25,8 +25,8 @@ void checkNonConformalMesh(const std::unique_ptr & mesh, const Real conformality_tol, unsigned int & num_nonconformal_nodes); -bool checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, - const std::unique_ptr & edge2, +bool checkFirstOrderEdgeOverlap(const Elem & edge1, + const Elem & edge2, Point & intersection_point, const Real intersection_tol); } diff --git a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C index a8196f8eed52..3d626c150fe5 100644 --- a/framework/src/meshgenerators/MeshDiagnosticsGenerator.C +++ b/framework/src/meshgenerators/MeshDiagnosticsGenerator.C @@ -1431,8 +1431,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr bounding_box_map.insert({elem, boundingBox}); } - std::unique_ptr point_locator = - mesh->sub_point_locator(); // PointLocatorBase::build(TREE_ELEMENTS, *mesh); + std::unique_ptr point_locator = mesh->sub_point_locator(); std::set> overlapping_edges_nodes; for (const auto elem : mesh->active_element_ptr_range()) { @@ -1494,7 +1493,7 @@ MeshDiagnosticsGenerator::checkNonMatchingEdges(const std::unique_ptr // Now compare edge with other_edge Point intersection_coords; bool overlap = MeshBaseDiagnosticsUtils::checkFirstOrderEdgeOverlap( - edge, other_edge, intersection_coords, _non_matching_edge_tol); + *edge, *other_edge, intersection_coords, _non_matching_edge_tol); if (overlap) { // Add the nodes that make up the 2 edges to the vector overlapping_edges_nodes diff --git a/framework/src/utils/MeshBaseDiagnosticsUtils.C b/framework/src/utils/MeshBaseDiagnosticsUtils.C index 3960797cdac3..9c7fff5d5cdf 100644 --- a/framework/src/utils/MeshBaseDiagnosticsUtils.C +++ b/framework/src/utils/MeshBaseDiagnosticsUtils.C @@ -69,19 +69,19 @@ checkNonConformalMesh(const std::unique_ptr & mesh, } bool -checkFirstOrderEdgeOverlap(const std::unique_ptr & edge1, - const std::unique_ptr & edge2, +checkFirstOrderEdgeOverlap(const Elem & edge1, + const Elem & edge2, Point & intersection_point, const Real intersection_tol) { // check that the two elements are of type EDGE2 - mooseAssert(edge1->type() == EDGE2, "Elements must be of type EDGE2"); - mooseAssert(edge2->type() == EDGE2, "Elements must be of type EDGE2"); + mooseAssert(edge1.type() == EDGE2, "Elements must be of type EDGE2"); + mooseAssert(edge2.type() == EDGE2, "Elements must be of type EDGE2"); // Get nodes from the two edges - const Point & p1 = edge1->point(0); - const Point & p2 = edge1->point(1); - const Point & p3 = edge2->point(0); - const Point & p4 = edge2->point(1); + const Point & p1 = edge1.point(0); + const Point & p2 = edge1.point(1); + const Point & p3 = edge2.point(0); + const Point & p4 = edge2.point(1); // Check that the two edges are not sharing a node if (p1 == p3 || p1 == p4 || p2 == p3 || p2 == p4)