[WIP] Mortar-based periodic boundary conditions#62
Draft
rcarson3 wants to merge 29 commits into
Draft
Conversation
Largely had Claude drive this port as it the original SNLS trust region dogleg solver was relatively straightforward. Most of it would map over to MFEM's framework and wasn't overly complicated. The big boon though from having Claude do this is I also had it add the necessary mult transpose operators for the nonlinearform integrators and in particular this was done for the PA forms as the EA was trivial and the full matrix version already implemented it. I was also surprised that it was able to do the BBar PA Grad implementation of everything as well. The math behind it was actually more straightforward than I was expecting, and it was cool to see how it could derive all of it using different relationship it had discovered from looking at the full integration form.
… mfem::Array<bool> Was getting some odd GPU failures here at one point and just moved over to std::array<bool, 3> to get rid of them as I didn't want to deal with the odd memory issues I was hitting. The MFEM stuff was still fine. I just think I was hitting an odd bug in some GPU kernel, but moving to the std::array is ultimately the better approach as it honestly makes more sense for this stuff...
Had codex help work out what changes were needed to get ExaConstit updated to work with MFEM v4.9+ as I wanted to make sure we could use newer versions of Hypre as apparently that works better with ROCm v7+
Found a couple bugs in the post-processing: First one is that volume average values were getting outputted every time step designated for the viz files and the viz files were getting outputted every time step. So caught and fixed that issue. Next found in certain cases there could be segfaults due to some models not having a variable defined and thus a vdim equal to 0 causing segfaults which was a fun bug to run down... Finally this was one I had codex help dissect and chase down, we had some fun MPI stalls occurring and it was due to how data was being solved requiring the global communicator rather than a region defined communicator... Didn't notice this in earlier testing as our regions were usually on every rank... I'd still blame this one though on how MFEM is handling the communicators for meshes / par finite element spaces and in turn how this relates to the data collections as well...
Was running into a nasty GPU-bug where on newer versions of MFEM if we used > 1 GPUs we got different answers than our CPU runs as certain terms were just 0.0. I could not for the life of me figure it out other than it was likely due to some thing in the velocity field near the time the boundary conditions were being applied was not getting set. I threw codex at the problem and it was able to over a couple iterations of debugging work out the error and find a suitable new MFEM API that we could use that fixed the GPU error and kept our answers on the CPU the same.
Probably should have captured more of the intermediate updates through-out this journey of "vibe-coding"... However, I used Claude Opus 4.7 and Claude's project / github integration with its web chat to drive an initial implementation of a mortar based periodic boundary conditions. I initially had it work things out in pyMFEM and create a plan to build things up from 2D serial to parallel to creating a ton of tests to then moving onto 3D with more tests. Some of these tests involved multi-material elastic simulations to break up in homogeneity and to actually test the methods. Outside of the tests here, I did validate what it was doing with visualizations making sure the solutions at each stage of the way made sense and if not would iterate with Claude until we got something decent. After the python implementation looked more or less good enough, we moved over to C++. Here, I had it work on the initial port and verify that things matched up with the Python version. After we did that, I had it move on to a more scaleable set-up and assembly phase for the constraint matrix as it was initially doing everything on a single rank or broadcasting everything to all ranks. Once that seemed like in a good state, we moved onto an element assembly formulation as that seemed like a good potential path to get us on the GPU and that was verified against the existing sparse matrix / hypre implementation. An initial stab at GPU support also happened. I know this will need to be expanded on. Finally, I had it add support for 3D non-conforming mesh faces. It needed Axom's BVH and similar spatial search algorithms for this which is why we know support Axom. The solutions there seem to line up with the expected results and when the meshes and as the mesh gets closer to conformal the distortional displacement required in the solution approached 0 and matched the conformal variation which was a good indication that the method was doing what we would expect. Also, this work added a ton of new tests which is good but yeah there are a ton of new tests to work through...
For nonconformal tests added the checkerboard and stripe test. Outside that relaxed the F_tol test condition so the tests would pass. Then the mortar saddle point API was updated to make it more useable for when ExaConstit would havee to use it.
Claude created a manager that controls all of the fun-ness that we'll need to manage for the Mortar PBCs in both the SystemDriver and elsewhere.
…PBCs and add test
…t matrix Had Claude add some functionality to create the boundary sub-mesh, add a simple function to compute the current deformation gradient, and then finally some other things to get out the Fbar and Fdot_bar terms.
… fact The mortar method should now have several different capabilities to check the fluctuation field or the hill-mandel condition and a number of other things kinda related to this stuff.
…thod into SystemDriver and nonlinear solvers This was over a couple different iterations with Claude but the above has most of the framework in-place to where we should be able to in the very near term test the new mortar-based periodic boundary conditions on real problems.
After lots of debugging finally got a working example with periodic boundary conditions when using ExaConstit directly. Initial tests are pretty much just an isotropic elastic examples using ExaCMech but it's a start. The good things is we do see converged results which is super awesome.
…ic boundary conditions Had Claude add a number of useful validation / diagnostic type fields as part of the periodic boundary conditions. Also had it add all of the affine and fluctuation fields for the velocity field so we can see how the PBCs have driven things to be different.
Pushing Claude to move towards "semi-periodic" BCSs so we don't end up with overly constrained systems and we can relax the constraint matrix to allow for more natural relaxation of faces for like in monotonic deformation modes
So, I found in testing that the default method of clamping all of the corner DOFs to be way to restrictive and wouldn't really allow me to do like monotonic tension tests easily. Therefore, I had Claude drive things so that we could have what is like semi-periodic BCs by allowing the constraint matrix to be set up based on which DOFs you actually want. Overall, the whole thing works and we now have something that should be useable by the broader community.
…ewton/Krylov wiring
Add the Phase 5.11 saddle-residual scaling stack for mortar PBC solves,
including the residual scaler, scaled saddle/Jacobian/preconditioner
wrappers, trust-region dogleg support for scaled coordinates, and
per-Newton-iteration diagnostic logging.
This change also wires the scaling path through SystemDriver and the
mortar manager, adds sub-block partition support for lambda rows, and
extends the options path to configure saddle scaling from TOML.
While debugging the new path, fix several solver-state and lifecycle
issues that could make first-step behavior nondeterministic or unsafe:
- zero Newton/Krylov correction buffers before inner solves so the
linear solver never consumes stale data as an initial guess
- propagate iterative_mode through wrapped solvers so zero-initial-guess
semantics survive the wrapper stack
- remove a duplicate physical residual evaluation in the scaling
pre-attempt path; the residual callback is stateful and must not be
probed twice before Newton starts
- refresh scaled wrapper state, TRDOG offsets, and diagnostic wiring
after periodic-BC spec changes that resize the lambda block
- remove the temporary inspecting iterative solver used during
root-cause analysis once the underlying initialization bug was fixed
At a high level this commit introduces the new saddle scaling
infrastructure, preserves diagnostic visibility for the scaled solve,
and fixes the wrapper/solver initialization bugs that were causing
inconsistent cycle-1 residual behavior and occasional MINRES NaNs.
…aling checks and fix a cmake issue
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
So this has been completely driven Claude Opus 4.7 writing all of the code and tests for this feature as I'm in no way an expert on mortar methods, and I've been very much interested in supporting periodic boundary conditions in ExaConstit for a while. The code was largely inspired by https://doi.org/10.1016/j.cma.2021.113930 which was a really interesting paper I came across on periodic boundary conditions and a way to enforce periodic boundary conditions through mortar methods which allows you to have non-conforming mesh faces or heterogeneous materials across various faces.
So, I largely "vibe" coded in the sense that I did not write anything but directed where things should go, what assumptions it was making were bad, when it should be doing a deeper research of the field to ensure we were doing things correctly as well as many, many different things. The initial prototype of this code started as a python miniapp using pyMFEM where I could quickly iterate with it and verify what it was doing worked starting from 2D and moving up to 3D with a linear elastic problem or a bi-material elastic material with a soft and hard phase. After doing that, I felt comfortable moving onto the C++ port and bringing over all of the python tests to C++ to verify the port was working as expected. The C++ side of things took a lot of work as we needed scalable algorithms (largely for the constraint systems), GPU-capable algorithms (still getting there), usable preconditioners, support for non-conforming mesh faces, proper translation from the total lagrangian framework these methods were written in to the updated lagrangian framework, and then like a ton of tests to make sure things actually work and help to track bugs down later.
This is still a work in progress as I need to get it working with higher-order elements as we're currently limited to just linear elements. Outside of that, I'm also working on a number of extensions to the nonlinear / linear solvers so that we see faster convergences in some polycrystal simulations I've been running with only a handful of grains and linear tets, which are definitely a harder test for these models. Once I get a lot of this done, I will work on having Claude create some better / extensive standalone documentation to go with this new feature as this is not a small feature addition and there's a lot of things that need to be documented on how it works and what Claude was pulling from while working on this.