[ARCHIVED] Experiment lhs rhs calculators

[rfaasse]

This pull request archives an experimental prototype branch to be able to start discussions. There is no intention of merging this to master

[RiccardoRossi]

hi @rfaasse could you describe a bit what is the objective of this?

[rfaasse]

Hi @RiccardoRossi,
In our team we had some discussions, since we saw that currently, we have to create new elements to take into account transient vs steady-state behavior. For example, the water pressures in our calculations are governed by:
$$C \dot{p} + Hp = q$$,
where $C$ is the compressibility, $\dot{p}$ the time derivative of $p$ (water pressure), $H$ is the permeability and $q$ the RHS. However, if we run a steady state computation (i.e. $\dot{p} = 0$) we don't want to do a redundant calculation of the compressibility matrix. For 'displacement-related' contributions (like the mass and damping matrix) this of course can be handled by the schemes, but the interfaces don't support the more domain specific contributions (like permeability, compressibility, conductivity etc for water pressure/heat calculations).

Currently that's solved in geomechanics by creating a new element where we override the calculation of the compressibility to be a zero matrix, but this way of working can lead to a 'wildgrowth' of elements, which is hard to maintain.

Our idea is to provide the element with a list of abstract calculators which can calculate the contributions. By extracting the calculations in this way, we can very easily decide what contributions are taken into account (e.g. compressibility, permeability, coupling to displacement, fluid body flow etc).

Concretely, this means that we create a transient Pw element like this:

    const auto list_of_contributions = {CalculationContribution::Permeability, CalculationContribution::Compressibility}
    const PwElement<2, 2> mTransientPwElement2D2N{
        0, Kratos::make_shared<Line2D2<NodeType>>(Element::GeometryType::PointsArrayType(2)),
        list_of_calculators};

And the only difference in creating a steady state Pw element would be to not add the compressibility contribution to the list_of_contributions (instead of creating a new element).

The CalculateLocalSystem of the element then would look something like this:

        rRightHandSideVector = CalculateFluidBodyVector();
        rLeftHandSideMatrix  = ZeroMatrix{TNumNodes, TNumNodes};

        for (const auto& rContribution : mContributions) {
            const auto calculator = CreateCalculator(rContribution, rCurrentProcessInfo);
            const auto [LHSContribution, RHSContribution] = calculator->CalculateLeftAndRightHandSide();
            rLeftHandSideMatrix += LHSContribution;
            rRightHandSideVector += RHSContribution;
        }

We have done something conceptually similar for the stress states, where we inject a 3D, PlaneStrain or Axisymmetric stress state policy into an element, such that all the 'dimensional-specific' things are handled there and we don't need to create special 'axisymmetric' elements, we just create a 'normal' element and inject an axisymmetric stress state into it.

We plan on elaborating on this way of writing our elements/components at the Kratos User Meeting this November as well and if you want to stay updated on this specific idea (see #12779 for the first implementation), or if you have any further questions, please let me or anyone from the @KratosMultiphysics/geomechanics team know!

[RiccardoRossi]

hi @rfaasse i think I understand wht you want to do. A few years back @pooyan-dadvand has in mind something (which never came to life) and that was tentatively called "Formulation".

My impression however is that we are going a long way because of the limitations imposed by the builder and solver.

imagine you were able to do something like

auto element_function = 
for(auto& geom: mp.Geometries()){
    eq_ids = get_eq_ids_function(geom)
    lhs,rhs = element_function(geom, ... )
    rA.assemble(lhs, eq_ids)
    rb.assemble(b,eq_ids)
}

your code would look probably much more cleaner.

Our idea (Still in the works) is to leverage the Kratos csr matrix to allow this.

see for example link

for an idea of what i would like to achieve in a new version of the scheme

[rfaasse]

Hi @RiccardoRossi,
Indeed 'Formulations' look very similar to the calculator workflow we are working on! If I understand it correctly the main difference is that we are implementing it on an element-matrix level. If this can be done at a higher level at some point, that might be a nice improvement!

I'll try to find some time to look into the concept of the CSR matrix, it sounds interesting!