GI_GAUSS_EXTENDED was originally added by me, and the poupouse was to be used by the solid-shell. My fault, because after that we have been using it as catch-all (cajón de sastre). My suggestion is to define proper names and keep GI_GAUSS_EXTENDED with current behaviour. I think the problem is that this overloads the size of the geometry class make it more inefficient, so maybe we should rethink completely the integration strategies .

Since it's labelled a general discussion, I'll go ahead and aim a suggestion at @KratosMultiphysics/nurbs-breps and @KratosMultiphysics/mpm.

In my opinion, both the IGA and MPM apps struggle to align their abstractions with Kratos' layout, partially due to the relatively rigid quadrature system.

I'm not exactly sure why, but both applications treat integration points as elements even though they'd have more conventional options. I'm not too familiar with either field, but I'd imagine that:

• using the background elements with a dynamic quadrature rule would fit MPM better. A quadrature rule could be constructed at each solution step, defaulting to the current implementation, but allowing the construction of a reduced integration in case of too many particles in a single background element.
• using knotspans as elements would be much more intuitive for IGA in my mind. Of course its integration rule has to depend on its polynomial order and trimming.

In any case, I think it would be neat if each element/geometry instance had the possibility to own its own integration object (kinda similar to Properties - in most cases all elements share the same Properties but if necessary, they could share subsets or just have their own).

GI_GAUSS_EXTENDED was originally added by me, and the poupouse was to be used by the solid-shell. My fault, because after that we have been using it as cath-all (cajón de sastre). My suggestion is to define proper names and keep GI_GAUSS_EXTENDED with current behaviour.

After a brief discussion yesterday, @RiccardoRossi @AlejandroCornejo and myself reached consensus on removing them if are of course not being used. If they are being used, we'd implement them consistently for all the geometries, not necessarily for order 1 to 5.

I think the problem is that this overloads the size of the geometry class make it more inefficient, so maybe we should rethink completely the integration strategies .

This is what we also spoke about. The point is that, as you already know, everything must happen at compile time, something that turns into an extremely tight design constrain.

* using the background elements with a dynamic quadrature rule would fit MPM better. A quadrature rule could be constructed at each solution step, defaulting to the current implementation, but allowing the construction of a reduced integration in case of too many particles in a single background element.

This I don't think is possible. Also I don't think it is a good idea. Note that integration quadratures cannot be dynamic for the sake of efficiency. Right now, everything happens at compile time, so calling the integration points has no cost.

* using knotspans as elements would be much more intuitive for IGA in my mind. Of course its integration rule has to depend on its polynomial order and trimming.

Ideally yes, but note that our IGA is not "classical" IGA but IBRA, hence, we will always have the trimming. Having said this, we will always have to deal with the trimming, which necessarily needs to happen at runtime and hence goes against the need of setting up the integration points at compile time.

Small comment regarding IGA/IBRA: We plan to have the possibility of creating one layer abstraction by adding KnotSpanGeometry, which is going to use knotspans as elements. The idea is to collect the list of integration points and assign it as one geometry so that we have one element/condition per set of integration points. Note that this will not depreciate the idea of QuadraturePointGeometry as this still might be used in some cases, such as point on a surface/curve or the MPM application. In a couple of weeks, we will prepare the PR.

CC: @manuelmessmer

Having said all this, let me try to put in short the bullet points that we need to address (some of them are of course wishes):

• We need the quadratures "to happen" at compile time
• We need to fix the catch-all (cajón desastre) issue with the GI_EXTENDED quadratures
• When adding a quadrature for one geometry, we would like to avoid the need to define it for all the other geometries.
• If the one above is not possible, would be nice to have a mechanism to throw a runtime no-implementation error (rather than putting anything just for the sake of having it defined as we did so far).

Let me note that some of these complexities and inconveniences come from the fact that we are relying on enum(s) wrap the integration methods. This is the unique thing that was available in ancient days to accomplish with the compile time requirement, but, maybe, with modern c++ we can do something else (pinging @roigcarlo for an imaginative solution).

Mentioning @pablobecker as expert in dynamic allocated quadratures

We need the quadratures "to happen" at compile time

It's unlikely that I'll personally touch kratos quadrature in the near future so I'm not requesting any changes here, but I'll just share my unsolicited thoughts on it:

I agree that storing shape function values at predetermined locations instead of evaluating them on the fly works very well for low order elements/geometries. That said, making this a requirement for every geometry and forbidding other approaches locks us out of whole families of methods (FCM comes to mind first, but I really mean any approach that relies on tweaking its shape functions or quadrature dynamically). This is quite a shame because kratos is primarily a research code and most other parts of our infrastructure would be able to support such methods.

Note that integration quadratures cannot be dynamic for the sake of efficiency. Right now, everything happens at compile time, so calling the integration points has no cost.

The only cost we're sparing is the calculation of shape functions at each integration point. We still have to fetch cached shape function values from memory, which is absolutely not a negligible cost.

We can think alternative approaches, but I think I will discuss with @pooyan-dadvand first, because I think my idea woul have a performance penalty.

we had a first round of discussion on this one, taking also into account the feedback of @matekelemen

as of today IT IS possible to store in a geometry a "flexible" GeometryData (chancing the pointer at construction time). This needs to be better documented but should take into account the use case of @matekelemen

per regards the EXTENSION_GAUSS our impression is that they should be eventually removed to prevent usage as catch all

this is very preliminary and we will iterate it further, nevertheless our current impression is that the integration rules should have:

• GAUSS, order
• LOBATTO, order
• COLLOCATION, order

anyhowthis would be a breaking change so it is just a working idea.

please note however that EXTENSION_:GAUSS with different meaning should NOT be allowed

we had a first round of discussion on this one, taking also into account the feedback of @matekelemen

as of today IT IS possible to store in a geometry a "flexible" GeometryData (chancing the pointer at construction time). This needs to be better documented but should take into account the use case of @matekelemen

per regards the EXTENSION_GAUSS our impression is that they should be eventually removed to prevent usage as catch all

this is very preliminary and we will iterate it further, nevertheless our current impression is that the integration rules should have:

* GAUSS, order

* LOBATTO, order

* COLLOCATION, order

anyhowthis would be a breaking change so it is just a working idea.

please note however that EXTENSION_:GAUSS with different meaning should NOT be allowed

COLLOCATION is not known at builiding time, am I right?

well, you could have a collocation known at building time, for example dividing a triangle in 4, 16 and so on

So, in practical terms, did we reach any consensus regarding this issue?@KratosMultiphysics/technical-committee

Even though I note @loumalouomega 's initial intentions, I'd start by removing the extended Gauss quadratures as these are only correctly implemented for a few geometries and, as far as I know, these are not being used.

We have been further discussing this in @KratosMultiphysics/technical-committee

as a short term fix we shall change this list:

    enum class IntegrationMethod {
        GI_GAUSS_1,
        GI_GAUSS_2,
        GI_GAUSS_3,
        GI_GAUSS_4,
        GI_GAUSS_5,
        GI_EXTENDED_GAUSS_1,
        GI_EXTENDED_GAUSS_2,
        GI_EXTENDED_GAUSS_3,
        GI_EXTENDED_GAUSS_4,
        GI_EXTENDED_GAUSS_5,
        GI_LOBATTO_1,
        NumberOfIntegrationMethods // Note that this entry needs to be always the last to be used as integration methods counter
    };

to

    enum class IntegrationMethod {
        GI_GAUSS_1,
        GI_GAUSS_2,
        GI_GAUSS_3,
        GI_GAUSS_4,
        GI_GAUSS_5,
        GI_LOBATTO_1,
        GI_EXTENDED_GAUSS_1,
        GI_EXTENDED_GAUSS_2,
        GI_EXTENDED_GAUSS_3,
        GI_EXTENDED_GAUSS_4,
        GI_EXTENDED_GAUSS_5,
        NumberOfIntegrationMethods // Note that this entry needs to be always the last to be used as integration methods counter
    };

This implies the requirement of implementing the first order lobatto quadrature for all available geometries.

In the mid term, the EXTENDED* quadratures should be fixed to be consisten across geometries, which would imply ideally naming them in a more self descriptive way.

@rubenzorrilla will suprvise this syncronizing with @loumalouomega and @AlejandroCornejo

We have been further discussing this in @KratosMultiphysics/technical-committee

as a short term fix we shall change this list:

enum class IntegrationMethod {
    GI_GAUSS_1,
    GI_GAUSS_2,
    GI_GAUSS_3,
    GI_GAUSS_4,
    GI_GAUSS_5,
    GI_EXTENDED_GAUSS_1,
    GI_EXTENDED_GAUSS_2,
    GI_EXTENDED_GAUSS_3,
    GI_EXTENDED_GAUSS_4,
    GI_EXTENDED_GAUSS_5,
    GI_LOBATTO_1,
    NumberOfIntegrationMethods // Note that this entry needs to be always the last to be used as integration methods counter
};

to

enum class IntegrationMethod {
    GI_GAUSS_1,
    GI_GAUSS_2,
    GI_GAUSS_3,
    GI_GAUSS_4,
    GI_GAUSS_5,
    GI_LOBATTO_1,
    GI_EXTENDED_GAUSS_1,
    GI_EXTENDED_GAUSS_2,
    GI_EXTENDED_GAUSS_3,
    GI_EXTENDED_GAUSS_4,
    GI_EXTENDED_GAUSS_5,
    NumberOfIntegrationMethods // Note that this entry needs to be always the last to be used as integration methods counter
};

This implies the requirement of implementing the first order lobatto quadrature for all available geometries.

In the mid term, the EXTENDED* quadratures should be fixed to be consisten across geometries, which would imply ideally naming them in a more self descriptive way.

@rubenzorrilla will suprvise this syncronizing with @loumalouomega and @AlejandroCornejo

I am OK, the place where I use my extended gauss points is quite hardcorded, so it can be a little bit more hardcoded.