Fix Multigrid Solver in Forest #1035
Conversation
|
hm, I managed to break something in the last couple commits. Will get back to it on Monday. Likely an issue related to |
now fixed, was due to a bug introduced by a merge. |
|
Is it safe to change the target branch to develop now given that #1028 has been merged? |
|
@pgrete: The target branch is now |
There was a problem hiding this comment.
I looked at the code and followed as much as I could (which, admittedly is not too much given that I haven't really internalized all the logic and associated magic).
Anyway, negative lines are always welcome and I just stumbled across two numbers that I weren't able to immediate connect.
| // lx3() >> nlevel); for most cases | ||
| const int norig = 1LL << std::max(level(), 0); | ||
| const int nparent = 1LL << std::max(level() - nlevel, 0); | ||
| constexpr int nmax_tree_offset = 5; |
There was a problem hiding this comment.
Is this a magic number/should be be obvious why it's 5?
There was a problem hiding this comment.
I think it could be anything that doesn't overflow when multiplied by norig.
| } | ||
| } | ||
| } | ||
| } | ||
| } | ||
|
|
||
| // Build in negative levels | ||
| for (int l = -20; l < 0; ++l) { |
There was a problem hiding this comment.
Is -20 just a large enough number that we should never reach it?
There was a problem hiding this comment.
yeah, for that to not be enough would imply we had a meshblock size of (2^20)^D.
| block_bcs[GetOuterBoundaryFace(dir)] = mesh_bcs[GetOuterBoundaryFace(dir)]; | ||
| } | ||
| } | ||
| block_size = forest.GetBlockDomain(loc); |
There was a problem hiding this comment.
No, block_size is passed by reference. The input to this routine is loc and the output is block_size and block_bcs. I agree that this is a weird design, but it is what was originally in the code.
PR Summary
This PR updates the geometric multigrid infrastructure to work with the new forest of trees based mesh. Generally, I think this has resulted in a much cleaner implementation.
The forest neighbor finding algorithm is extended to find neighbors in different types of meshes (i.e. the standard leaf mesh and two-level composite meshes at different levels). Additionally,
gids are assigned to internal nodes in the forest whenForest::GetMeshBlocksAndResolveGidsis called. The way the level of gmg grids is defined is also simplified so that it just corresponds to the logical level of finer blocks that are members of the grid (which means we store thegmg_block_listsin a map instead of a vector since the index can be negative). All of these changes should be internal and not be visible to downstream codes. Because we no longer have to worry about periodic root grids (with periodicity being taken care of by the forest connectivity) and because we use the forest information explicitly for gmg neighbor finding, a lot of the complex neighbor finding code inLogicalLocationis no longer necessary. This PR also removes a lot of other code that died in switching to a forest.PR Checklist
LogicalLocationcode #1028 is merged into develop