Skip to content
New issue

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

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

Already on GitHub? Sign in to your account

Fix of loading boundary elements of ParMesh from files #4128

Merged
merged 11 commits into from
Apr 15, 2024

Conversation

najlkin
Copy link
Contributor

@najlkin najlkin commented Feb 9, 2024

This PR fixes the problem recognized in #4083 . When ParMesh is loaded from files/streams through its constructor or Load() method with refine or fix_orientation parameters set true (as by default), boundary elements are generated for every boundary face including the shared ones if the local mesh does not have by chance at least one boundary element already loaded from the file, in which case the global boundary is preserved and no boundary elements are generated. This is a completely inconsistent behavior causing spurious appearance of internal boundaries depending on partitioning of the ParMesh.

The problem is/was in method Mesh::GenerateBoundaryElements(), which does not check anyhow if the boundary is shared or not. Here, it is overridden in ParMesh and no boundary elements are generated on shared faces.

PR Author Editor Reviewers Assignment Approval Merge
#4128 @najlkin @tzanio @sebastiangrimberg + @sohailreddy + @jandrej 2/10/24 4/12/24 4/15/24
PR Checklist
  • Code builds.
  • Code passes make style.
  • Update CHANGELOG:
    • Is this a new feature users need to be aware of? New or updated example or miniapp?
    • Does it make sense to create a new section in the CHANGELOG to group with other related features?
  • Update INSTALL:
    • Had a new optional library been added? If so, what range of versions of this library are required? (Make sure the external library is compatible with our BSD license, e.g. it is not licensed under GPL!)
    • Have the version ranges for any required or optional libraries changed?
    • Does make or cmake have a new target?
    • Did the requirements or the installation process change? (rare)
  • Update continuous integration server configurations if necessary (e.g. with new version requirements for each of MFEM's dependencies)
    • .github
    • .appveyor.yml
  • Update .gitignore:
    • Check if make distclean; git status shows any files that were generated from the source by the project (not an IDE) but we don't want to track in the repository.
    • Add new patterns (just for the new files above) and re-run the above test.
  • New examples:
    • All sample runs at the top of the example source file work.
    • Update examples/makefile:
      • Add the example code to the appropriate SEQ_EXAMPLES and PAR_EXAMPLES variables.
      • Add any files generated by it to the clean target.
      • Add the example binary and any files generated by it to the top-level .gitignore file.
    • Update examples/CMakeLists.txt:
      • Add the example code to the ALL_EXE_SRCS variable.
      • Make sure THIS_TEST_OPTIONS is set correctly for the new example.
    • List the new example in doc/CodeDocumentation.dox.
    • If new examples directory (e.g.examples/pumi), list it in doc/CodeDocumentation.conf.in
    • Companion pull request for documentation in mfem/web repo:
      • Update or add example-specific documentation, see e.g. the src/examples.md.
      • Add the description, labels and screenshots in src/examples.md and src/img.
      • In examples.md, list the example under the appropriate categories, add new categories if necessary.
      • Add a short description of the example in the "Extensive Examples" section of features.md.
  • New miniapps:
    • All sample runs at the top of the miniapp source file work.
    • Update top-level makefile and makefile in corresponding miniapp directory.
    • Add the miniapp binary and any files generated by it to the top-level .gitignore file.
    • Update CMake build system:
      • Update the CMakeLists.txt file in the miniapps directory, if the new miniapp is in a new directory.
      • Add/update the CMakeLists.txt file in the new miniapp directory.
      • Consider adding a new test for the new miniapp.
    • List the new miniapp in doc/CodeDocumentation.dox
    • If new miniapps directory (e.g.miniapps/nurbs), add it to MINIAPP_SUBDIRS in the makefile.
    • If new miniapps directory (e.g.miniapps/nurbs), list it in doc/CodeDocumentation.conf.in
    • Companion pull request for documentation in mfem/web repo:
      • Update or add miniapp-specific documentation, see e.g. the src/meshing.md and src/electromagnetics.md files.
      • Add the description, labels and screenshots in src/examples.md and src/img.
      • The miniapps go at the end of the page, and are usually listed only under a specific "Application (PDE)" category.
      • Add a short description of the miniapp in the "Extensive Examples" section of features.md.
  • New capability:
    • All new public, protected, and private classes, methods, data members, and functions have full Doxygen-style documentation in source comments. Documentation should include descriptions of member data, function arguments and return values, template parameters, and prerequisites for calling new functions.
    • Pointer arguments and return values must specify whether ownership is being transferred or lent with the call.
    • Any new functions should include descriptions of their intended use e.g. for internal use only, user-facing, etc., along with references to example code whenever possible/appropriate.
    • Consider adding new sample runs in existing examples to highlight the new capability.
    • Consider saving cool simulation pictures with the new capability in the Confluence gallery (LLNL only) or submitting them, via pull request, to the gallery section of the mfem/web repo.
    • If this is a major new feature, consider mentioning it in the short summary inside README (rare).
    • List major new classes in doc/CodeDocumentation.dox (rare).
  • Update this checklist, if the new pull request affects it.
  • Run make unittest to make sure all unit tests pass.
  • Run the tests in tests/scripts.
  • (LLNL only) After merging:
    • Update internal tests to include the new features.

mesh/pmesh.cpp Outdated Show resolved Hide resolved
mesh/pmesh.cpp Outdated Show resolved Hide resolved
@tzanio
Copy link
Member

tzanio commented Feb 11, 2024

This PR is now under review (see the table in the PR description). To help with the review process, please do not force push to the branch.

@tzanio tzanio added this to Review Now in Pull Requests via automation Feb 11, 2024
@tzanio tzanio added this to the mfem-4.7 milestone Feb 11, 2024
@sebastiangrimberg
Copy link
Contributor

This is a good bug fix and resolves the issue in #4083. I would suggest also including the following change to mesh.cpp as part of this PR to avoid calling GenerateBoundaryElements in parallel situations when it is unintended:

diff --git a/mesh/mesh.cpp b/mesh/mesh.cpp
index ce8dc7cf5..b830c499b 100644
--- a/mesh/mesh.cpp
+++ b/mesh/mesh.cpp
@@ -3014,7 +3014,7 @@ void Mesh::FinalizeTetMesh(int generate_edges, int refine, bool fix_orientation)
    FinalizeCheck();
    CheckElementOrientation(fix_orientation);

-   if (NumOfBdrElements == 0)
+   if (ReduceInt(NumOfBdrElements) == 0)
    {
       GetElementToFaceTable();
       GenerateFaces();
@@ -3056,7 +3056,7 @@ void Mesh::FinalizeWedgeMesh(int generate_edges, int refine,
    FinalizeCheck();
    CheckElementOrientation(fix_orientation);

-   if (NumOfBdrElements == 0)
+   if (ReduceInt(NumOfBdrElements) == 0)
    {
       GetElementToFaceTable();
       GenerateFaces();
@@ -3093,7 +3093,7 @@ void Mesh::FinalizeHexMesh(int generate_edges, int refine, bool fix_orientation)
    GetElementToFaceTable();
    GenerateFaces();

-   if (NumOfBdrElements == 0)
+   if (ReduceInt(NumOfBdrElements) == 0)
    {
       GenerateBoundaryElements();
    }
@@ -3166,7 +3166,7 @@ void Mesh::FinalizeTopology(bool generate_bdr)
    {
       GetElementToFaceTable();
       GenerateFaces();
-      if (NumOfBdrElements == 0 && generate_bdr)
+      if (ReduceInt(NumOfBdrElements) == 0 && generate_bdr)
       {
          GenerateBoundaryElements();
          GetElementToFaceTable(); // update be_to_face
@@ -3186,7 +3186,7 @@ void Mesh::FinalizeTopology(bool generate_bdr)
       if (Dim == 2)
       {
          GenerateFaces(); // 'Faces' in 2D refers to the edges
-         if (NumOfBdrElements == 0 && generate_bdr)
+         if (ReduceInt(NumOfBdrElements) == 0 && generate_bdr)
          {
             GenerateBoundaryElements();
          }
@@ -3200,7 +3200,7 @@ void Mesh::FinalizeTopology(bool generate_bdr)
    if (Dim == 1)
    {
       GenerateFaces();
-      if (NumOfBdrElements == 0 && generate_bdr)
+      if (ReduceInt(NumOfBdrElements) == 0 && generate_bdr)
       {
          // be_to_face will be set inside GenerateBoundaryElements
          GenerateBoundaryElements();
@@ -3264,7 +3264,7 @@ void Mesh::Finalize(bool refine, bool fix_orientation)
    }
    if (refine)
    {
-      MarkForRefinement();   // may change topology!
+      MarkForRefinement(); // may change topology!
    }

    if (may_change_topology)
@@ -3277,7 +3277,8 @@ void Mesh::Finalize(bool refine, bool fix_orientation)
       }
       else
       {
-         FinalizeTopology(); // Re-computes some data unnecessarily.
+         bool generate_bdr = false;
+         FinalizeTopology(generate_bdr); // re-computes some data unnecessarily
       }

       // TODO: maybe introduce Mesh::NODE_REORDER operation and FESpace::

@najlkin
Copy link
Contributor Author

najlkin commented Feb 12, 2024

I am not sure about bool generate_bdr = false;. It is an intended behavior perhaps 🤔

@sebastiangrimberg
Copy link
Contributor

I am not sure about bool generate_bdr = false;. It is an intended behavior perhaps 🤔

Because Finalize already assumes FinalizeTopology has been called, there are two options:

  • Either FinalizeTopology was called with generate_bdr = true, in which case calling it here with generate_bdr = false has no effect
  • FinalizeTopology was called with generate_bdr = false, in which case calling it here with generate_bdr = true as is currently the behavior on master is producing unintended consequences for the user

So, I think in either of these cases adding the explicit generate_bdr = false argument is not going to introduce unwelcome changes. But the other reviewers should chime in.

@najlkin
Copy link
Contributor Author

najlkin commented Feb 27, 2024

To move things forward a little bit, @jandrej , @sohailreddy , what are your opinions on generate_bdr = false? Note that ParMesh::Finalize() is a public method and Mesh::Finalize() generates boundaries automatically (when refine or fix_orientation are true).

@tzanio tzanio mentioned this pull request Mar 4, 2024
46 tasks
@najlkin
Copy link
Contributor Author

najlkin commented Apr 9, 2024

Yes, that is where we have started 😄 . I edited the header of #4240 to mention it avoids this as well 😉 . Please keep this PR minimal, so it can get finally merged. The overhead is only in the loader, so practically negligible, and it is not caused by this PR.

@sebastiangrimberg
Copy link
Contributor

sebastiangrimberg commented Apr 9, 2024

Yes, that is where we have started 😄 . I edited the header of #4240 to mention it avoids this as well 😉 . Please keep this PR minimal, so it can get finally merged. The overhead is only in the loader, so practically negligible and not caused by this PR.

Just need to reemphasize that my suggestion to remove the boundary element generation in the first part of this comment would be one way of remove that limitation and make the fixes in this PR even more minimal. I don't think there was clear consensus that this was a necessary feature to add though maybe I am missing something.

Another option would be to use the approach suggested in this comment and avoid accessing sedge_ledge or sface_lface inside of ParMesh::GenerateBoundaryElements. You can get the same information from the (already constructed) shared_edges, shared_trias, and shared_quads for your boundary element generation.

@najlkin
Copy link
Contributor Author

najlkin commented Apr 9, 2024

Just need to reemphasize that my suggestion to remove the boundary element generation in the first part of #4128 (comment) would be one way of remove that limitation and make the fixes in this PR even more minimal. I don't think there was clear consensus that this was a necessary feature to add though maybe I am missing something.

This bugfix is not adding a feature, but recovers the intended behavior, which is consistent with the behavior of serial meshes. I believe majority is on the side of keeping things consistent, or at least for a bugfix. We may return to this discussion in another PR if you want to open it for this.

Another option would be to use the approach suggested in #4128 (comment) and avoid accessing sedge_ledge or sface_lface inside of ParMesh::GenerateBoundaryElements. You can get the same information from the (already constructed) shared_edges, shared_trias, and shared_quads for your boundary element generation.

If we do not want to go beyond the bugfix and touch calls of FinalizeParTopo(), there is no need for the change. We would be just risking breaking the implementation. However, it is a good suggestion, just I am not convinced it is worth of it, but we may discuss it within #4240 😉 .

@sebastiangrimberg
Copy link
Contributor

OK, I've voiced my concerns but can only defer to editors/other reviewers at this point.

@v-dobrev
Copy link
Member

v-dobrev commented Apr 10, 2024

@najlkin, I think the decision at the meeting was to NOT include the parallel boundary generation but just make the checks (NumOfBdrElements == 0) global and leave the parallel version of GenerateBoundaryElements empty, at least for now.

This way we postpone the controversial changes for further discussion but get the important bugfix in.

@najlkin
Copy link
Contributor Author

najlkin commented Apr 10, 2024

@v-dobrev , I thought the decision was to split this PR, as I did, moving the potentially controversial change of the initialization of some structures to #4240 . The implementation of Mesh::GenerateBoundaryElements() in ParMesh is not controversial in the sense we all believe it is working correctly. It has been tested on the code of the user and the user also confirmed it is working for him in #4083 . @sebastiangrimberg is proposing two other things:

Both are API changes and should be discussed in #4240 , but not here (!). Please someone just confirm this is working, so it can go to v4.7 asap 🙏 .

@v-dobrev
Copy link
Member

I'll let @tzanio handle this -- he's the editor.

@najlkin
Copy link
Contributor Author

najlkin commented Apr 10, 2024

@jandrej , have you tested this PR is working? That is all we need now, the discussion related to the API changes should go to #4240 .

@tzanio
Copy link
Member

tzanio commented Apr 11, 2024

As a compromise: is it just possible to just ParMesh::GenerateBoundaryElements = delete which will explicitly expose this as a bug to the users?

@sebastiangrimberg, @v-dobrev and @najlkin -- what do you think?

@sebastiangrimberg
Copy link
Contributor

sebastiangrimberg commented Apr 11, 2024

As a compromise: is it just possible to just ParMesh::GenerateBoundaryElements = delete which will explicitly expose this as a bug to the users?

@sebastiangrimberg, @v-dobrev and @najlkin -- what do you think?

I like this idea but am not sure you can do that in C++ if I understand correctly:

./mfem/linalg/../fem/../mesh/pmesh.hpp:402:9: error: deleted function 'GenerateBoundaryElements' cannot override a non-deleted function
  402 |    void GenerateBoundaryElements() override = delete;

Leaving it defined as {} or as { MFEM_ABORT(" ..."); } would achieve similar effect at runtime (either silently or not)?

@najlkin
Copy link
Contributor Author

najlkin commented Apr 11, 2024

@sebastiangrimberg , it must be without virtual (and override)

@sebastiangrimberg
Copy link
Contributor

@sebastiangrimberg , it must be without virtual (and override)

Oh I see, but doesn't that mean you could get this behavior for a valid ParMesh?

ParMesh pmesh(...);
pmesh.GenerateBoundaryElements();  // Compiler error
Mesh *mesh = &pmesh;
mesh->GenerateBoundaryElements();  // Not a compiler error, bug-prone?

@najlkin
Copy link
Contributor Author

najlkin commented Apr 11, 2024

I think that behavior is correct in the sense that when you explicitly use the local sub-mesh, you can do many harmful things to the distributed mesh for sure, but that is responsibility of the user and might be desired for some master hacker 🤓

@sebastiangrimberg
Copy link
Contributor

sebastiangrimberg commented Apr 11, 2024

I think it means that when Mesh::FinalizeTopology eventually gets called from ParMesh::Load, you will get a call to Mesh::GenerateBoundaryElements (not ParMesh::GenerateBoundaryElements, which is correctly deleted) which is what I thought was what we wanted to avoid.

@tzanio
Copy link
Member

tzanio commented Apr 11, 2024

OK, then I vote to implement it as { MFEM_ABORT(" ..."); } and note that users shouldn't be calling that function.

@najlkin
Copy link
Contributor Author

najlkin commented Apr 11, 2024

Hmm, well, it is not critical, because that happens only when you do not have any boundary elements globally, which was not the case for the user and most other users. Anyway, I think a working implementation would be a better option, but to follow this compromise we may make FinalizeTopology() virtual and call the serial implementation with generate_bdr=false from ParMesh.

@najlkin
Copy link
Contributor Author

najlkin commented Apr 11, 2024

...or abort, if generate_bdr==true

@najlkin
Copy link
Contributor Author

najlkin commented Apr 12, 2024

Ok, so followed the pathway of @tzanio 's comment 😉 .

@tzanio
Copy link
Member

tzanio commented Apr 13, 2024

Merged in next for testing...

@tzanio tzanio added the in-next label Apr 13, 2024
@tzanio tzanio merged commit 12ed616 into master Apr 15, 2024
27 checks passed
Pull Requests automation moved this from Review Now to Merged Apr 15, 2024
@tzanio tzanio deleted the najlkin/pmesh-load-fix branch April 15, 2024 12:25
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
Pull Requests
  
Merged
Development

Successfully merging this pull request may close these issues.

None yet

6 participants