feat: Introduce full tensor permeability #2822
Conversation
| void ConstantPermeability::postProcessInput() | ||
| { | ||
| GEOS_ERROR_IF( m_diagonalPermeabilityTensor[0] < 0.0 && m_symmetricFullPermeabilityTensor[0] < 0.0, | ||
| "Either a diagonal permeability tensor or a full tensor must be provided."); |
There was a problem hiding this comment.
| "Either a diagonal permeability tensor or a full tensor must be provided."); | |
| "Either a diagonal permeability tensor or a symmetric full tensor must be provided."); |
|
I confirm it works. I ran an incompressible single-phase flow quarter five spot problem without gravity on two domains. The first case, the domain is aligned with the Cartesian axes and the permeability tensor is diagonal. In the second case, the domain is obtained by rotating the aligned domain and consistently updating the permeability tensor, which becomes full tensor. In both cases the grid is K-orthogonal and TPFA is correct. The computed Jacobian matrices are the same. Below two snapshots with details for the same cell (Id 52) in both domains.
I temporarily added xml and vtk files to reproduce the results in @karimifard, @jafranc, @francoishamon: have a look and share if you think this is headed in the correct direction. |
| m_permeability[ei][q][1] = m_permeabilityComponents[1]; | ||
| m_permeability[ei][q][2] = m_permeabilityComponents[2]; | ||
| for ( int c=0; c < 6; c++) | ||
| m_permeability[ei][q][c] = m_symmetricFullPermeabilityTensor[c]; |
There was a problem hiding this comment.
Shouldn't we enforce positiveness of diagonal component here ?
There was a problem hiding this comment.
We should really throw an error earlier if someone inputs a negative value but, for now, I am using the default negative value to distinguish whether a diagonal or a full tensor was provided so I need to improve that logic fist.
There was a problem hiding this comment.
One option would be to simply provide 2 R1Tensors one for the main diagonal and one for the off diagonal terms and have the off diagonal components default to 0.0 if not provided.
There was a problem hiding this comment.
Just noticed that the order would be imposed by LvArray::Ri_eq_symAijBj implicitly :/
There was a problem hiding this comment.
yes, and it's different from the one used by vtk.
| weight = LvArray::tensorOps::AiBi< 3 >( cellToFace, faceConormal ); | ||
| // Matteo: this one also has to change! | ||
| LvArray::tensorOps::hadamardProduct< 3 >( dWeight_dCoef, cellToFace, faceNormal ); |
There was a problem hiding this comment.
I am not sure about how this one should change.
I let @karimifard and @francoishamon, who are more expert than I am, tell if it is what they had in mind. |
…GEOS into feature/cus-cas/fullTensorPerm
# Conflicts: # src/coreComponents/finiteVolume/CellElementStencilTPFA.hpp
CusiniM
changed the title
rrsettgast
requested review from
rrsettgast,
castelletto1,
paveltomin,
ryar9534,
dkachuma and
jhuang2601
as code owners
@castelletto1 should be usable but I have not tested it yet.