Repository navigation
Conversation
…ixed_mimetic_discretization
herve-gross
left a comment
There was a problem hiding this comment.
Excellent, thank you Omar. A Sphinx documentation page would be great. I recommend using your very clear GitHub PR description for this documentation.
|
Will be reviewed by @jafranc, @joshua-white |
| FaceManager & faceManager = mesh.getFaceManager(); | ||
| arrayView1d< integer > const isPresBcFace = faceManager.getField< flow::isBoundaryFace >(); | ||
|
|
||
| fsManager.apply< FaceManager >( 0.0, |
There was a problem hiding this comment.
This initializes the pressure-boundary mask at time 0, but I'm wondering if field specifications can activate and expire later. Then, should we clear and rebuild the mask together with the pressure values at the current assembly time? Might not be an issue, but just wanted to bring up
| } ); | ||
|
|
||
| // evaluate the boundary face pressure values used in the constitutive rows | ||
| applyFacePressureBCValues( time_n + dt, domain ); |
There was a problem hiding this comment.
Looks like boundary values are evaluated during step setup, which is not repeated when the nonlinear retry loop reduces dt. A retry can therefore use the original endpoint's pressure. Could we refresh boundary data before flux assembly using the actual attempt dt, without overwriting beginning-of-step snapshots?
| bool const valid0 = ( m_elemRegionList[kf][0] >= 0 && m_elemSubRegionList[kf][0] >= 0 && m_elemList[kf][0] >= 0 ); | ||
| bool const valid1 = ( m_elemRegionList[kf][1] >= 0 && m_elemSubRegionList[kf][1] >= 0 && m_elemList[kf][1] >= 0 ); | ||
| bool const onBoundary = !( valid0 && valid1 ); | ||
| stack.isNoFlowFace[i] = ( onBoundary && m_isPresBcFace[kf] == 0 ) ? 1 : 0; |
There was a problem hiding this comment.
Looks like this counts both physical neighbors, while the condensed kernel counts only target cells at lines 564–583. A target/non-target interface is therefore interior here but a boundary there; the condensed path then returns assuming an identity closure was assembled. Should both paths use target-relative adjacency, with matching orientation and a partial-region regression?
| vtkArray->GetName(), src, cellIdx ), | ||
| InputError ); | ||
| } | ||
| val = static_cast< DstType >( src ); |
There was a problem hiding this comment.
I'm wondering if we should check that this casting makes sense before actually doing it. For example, a value such as 1e20 can reach an out-of-range floating-to-integer cast
|
|
||
| // the same one-sided conductance and floor as the diagonal entries of TPFAInnerProduct::computeM | ||
| real64 const areaTolerance = m_lengthTolerance * m_lengthTolerance; | ||
| real64 const weightTolerance = 1e-30 * m_lengthTolerance; |
There was a problem hiding this comment.
I'm wondering if 1e-30 should be a variable instead
|
|
||
| FieldIdentifiers fieldsToBeSync; | ||
| fieldsToBeSync.addElementFields( { mixedMimetic::mfdFlag::key(), mixedMimetic::consistencyIndicator::key() }, regionNames ); | ||
| CommunicationTools::getInstance().synchronizeFields( fieldsToBeSync, mesh, neighbors, false ); |
There was a problem hiding this comment.
Looks like the marking kernel writes these fields on the device, but synchronization selects host packing. Should we make memory-space transitions explicit using coherent host views or device packing?
| subRegion.faceList().toViewConst(), | ||
| subRegion.getElementCenter(), | ||
| subRegion.getElementVolume(), | ||
| permeability[er][esr], |
There was a problem hiding this comment.
Looks like host indexing of permeability[er][esr] follows a device projection that captures the whole nested accessor
victorapm
left a comment
There was a problem hiding this comment.
Thanks for the great work @chauj96 and @OmarDuran ! I left a few comments and I'm also working on a follow-up PR with some additional improvements built on top of this branch
This PR adds a mixed mimetic finite difference (MFD) discretization of single-phase flow, solved as a saddle-point problem. It adapts the consistency, each cell uses the TPFA inner product where TPFA is consistent and a consistent
MFD product elsewhere; fluxes between two TPFA cells are condensed. The associated theory can be found here: https://arxiv.org/abs/2607.23568.
It adds the
mixedMimeticmodule, theSinglePhaseMixedMFDsolver, a mixed-form RT inner product that requires no stabilization on simplexes, as it is the exact lower-order RT element, the ConsistencyAdaptation cell classification, a MGR strategy, and unit tests.This is a work in collaboration with @chauj96 , @castelletto1 , @victorapm and myself.
The discretization block will look like this:
innerProductType: the consistent product used in the MFD cells.adaptiveConsistency: 1 selects TPFA or MFD per cell, 0 uses MFD everywhere.consistencyTolerance: a cell becomes MFD when the relative two-point flux error, probed withnominalGradient, exceeds it.degeneracyTolerance: a cell is degenerate, and stays TPFA, when its volume is smaller than this percentage ofthe combined volume of all cells sharing a vertex with it. This controls the conditioning of the system on meshes with cells with extreme aspect ratio.
Alternatively, the choice can be prescribed per cell from the VTK mesh with the integer cell array
prescribedMfdFlag: 0 = TPFA, 1 = MFD, −1 = letconsistencyToleranceanddegeneracyTolerancedecide. Prescribed cells with 0 = TPFA, 1 = MFD are never altered. The result of the classification is written asmfdFlag.The following is an toy example of consistency adaptation on the full SPE10, where the green regions (MFD subregion) represent a simplex region connected with a large portion of hexahedral cells (TPFA subregion).
The following figure shows the MFD sub-region automatically identified employing a
consistencyTolerance="0.01"The following figures show the pressure colormaps for adaptive-MFD (left panel), TPFA on a hexahedral mesh (middle panel), and TPFA on the hybrid mesh (right panel). The TPFA approximation on the hybrid mesh exhibits a visible effect due to the consistency error localized in the simplex cells. The adaptive-MFD approximation is quite similar to the TPFA approximation on the hexahedral mesh, where it is exact. The small differences are due to the different types of meshes.
The following table was constructed using 8 ranks to show the following features: