Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[CORE] Adding new Lobatto quadratures to linear geometries #12918

Merged
merged 17 commits into from
Jan 6, 2025

Conversation

AlejandroCornejo
Copy link
Member

@AlejandroCornejo AlejandroCornejo commented Dec 10, 2024

📝 Description
Adding a new entry to the enum of thequadraturesfor 1st order Lobatto quadratures (integrate exactly polynomials up to order 1).

Added option to triangles, quads, test and hexas.

@AlejandroCornejo AlejandroCornejo requested a review from a team as a code owner December 10, 2024 14:06
Copy link
Member

@loumalouomega loumalouomega left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

FYI @KratosMultiphysics/technical-committee

@AlejandroCornejo
Copy link
Member Author

Now passes! @rubenzorrilla

Copy link
Member

@rubenzorrilla rubenzorrilla left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think there are still some issues.

@AlejandroCornejo
Copy link
Member Author

I think there are still some issues.

I modified them accordingly

@@ -1450,8 +1450,8 @@ template<class TPointType> class Hexahedra3D8 : public Geometry<TPointType>
Quadrature < HexahedronGaussLegendreIntegrationPoints4, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
Quadrature < HexahedronGaussLegendreIntegrationPoints5, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
Quadrature < HexahedronGaussLobattoIntegrationPoints1, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
Quadrature < HexahedronGaussLobattoIntegrationPoints2, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
Quadrature < HexahedronGaussLobattoIntegrationPoints2, 3, IntegrationPoint<3> >::GenerateIntegrationPoints()
Quadrature < HexahedronGaussLobattoIntegrationPoints1, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
Quadrature < HexahedronGaussLobattoIntegrationPoints1, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),

Why repeated?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It is not "repeated" actually. Previously, in the EXTENDED_GI_GAUSS thy were using the Lobatto so I just maintained the previous code.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Still I don't get it. But aren't we calling thrice the same method? In other words, shouldn't Quadrature < HexahedronGaussLobattoIntegrationPoints1, 3, IntegrationPoint<3> >::GenerateIntegrationPoints() appear once? (right now it appears three times).

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I may be totally wrong but, by analogy with other geometries I was implementing, on the one hand:

    static const ShapeFunctionsValuesContainerType AllShapeFunctionsValues()
    {
        ShapeFunctionsValuesContainerType shape_functions_values =
        {
            {
                Hexahedra3D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_GAUSS_1 ),
                Hexahedra3D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_GAUSS_2 ),
                Hexahedra3D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_GAUSS_3 ),
                Hexahedra3D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_GAUSS_4 ),
                Hexahedra3D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_GAUSS_5 ),
                Hexahedra3D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_1 ),
                Hexahedra3D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_EXTENDED_GAUSS_2 ),
                Hexahedra3D8<TPointType>::CalculateShapeFunctionsIntegrationPointsValues( GeometryData::IntegrationMethod::GI_LOBATTO_1 )
            }
        };
        return shape_functions_values;
    }

and then we need to assing an individual quadrature to each os the previous entries:

    static const IntegrationPointsContainerType AllIntegrationPoints()
    {
        IntegrationPointsContainerType integration_points =
        {
            {
                Quadrature < HexahedronGaussLegendreIntegrationPoints1, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
                Quadrature < HexahedronGaussLegendreIntegrationPoints2, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
                Quadrature < HexahedronGaussLegendreIntegrationPoints3, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
                Quadrature < HexahedronGaussLegendreIntegrationPoints4, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
                Quadrature < HexahedronGaussLegendreIntegrationPoints5, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(),
                Quadrature < HexahedronGaussLobattoIntegrationPoints1, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(), // was there already
                Quadrature < HexahedronGaussLobattoIntegrationPoints1, 3, IntegrationPoint<3> >::GenerateIntegrationPoints(), // was there already
                Quadrature < HexahedronGaussLobattoIntegrationPoints1, 3, IntegrationPoint<3> >::GenerateIntegrationPoints() // My contribution! :)
            }
        };
        return integration_points;
    }

Are we now on the same page?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

OK. Now I understand. Your addition is correct. What is weird to me is that GI_EXTENDED_GAUSS_1 and GI_EXTENDED_GAUSS_2 feature the HexahedronGaussLobattoIntegrationPoints1 integration points (something that was already there).

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@loumalouomega together with @AlejandroCornejo we've discussing about this. We think that it might be something weird with the hexa GI_EXTENDED_GAUSS_* integration as it was using the Lobatto integration points (rather than the collocation ones). I think that you can shred some light on this.

If it is actually wrong, or we put is as a boilerplate for a future implementation, my suggestion would be to remove it.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For now, I've removed the extended ones amd maintained the standard gauss and the new Lobatto one. FYI @loumalouomega

@AlejandroCornejo

This comment was marked as outdated.

@markelov208

This comment was marked as resolved.

@AlejandroCornejo

This comment was marked as resolved.

@ipouplana
Copy link
Member

THX @markelov208 for your tips. I kept doing some research and the problem seems to be this code (for the 2d quadrilateral elements), GeoMechanicsApplication\custom_elements\U_Pw_small_strain_interface_element.hpp:

    virtual void CalculateOnLobattoIntegrationPoints(const Variable<array_1d<double, 3>>& rVariable,
                                                     std::vector<array_1d<double, 3>>&    rOutput,
                                                     const ProcessInfo& rCurrentProcessInfo);

    virtual void CalculateOnLobattoIntegrationPoints(const Variable<Matrix>& rVariable,
                                                     std::vector<Matrix>&    rOutput,
                                                     const ProcessInfo&      rCurrentProcessInfo);

    virtual void CalculateOnLobattoIntegrationPoints(const Variable<Vector>& rVariable,
                                                     std::vector<Vector>&    rOutput,
                                                     const ProcessInfo&      rCurrentProcessInfo);

@ipouplana this code seems to be yours, could you give me a hint of what's being done here? I see that you are computing something using Lobatto quadratures (it is using the mThisIntegrationMethod which by default is GeometryData::IntegrationMethod::GI_GAUSS_2 so it is not Lobatto) and then interpolating. Now the Lobatto quadratures can be used via GeometryData::IntegrationMethod::GI_LOBATTO_1

We already discussed this internally, but I will mention it here too for the record.
The interface elements are defined with a special geometry which uses a Lobatto quadrature file. And thus it actually uses a Lobatto quadrature regardless of the label used to call it, e.g. GI_GAUSS_2.

@AlejandroCornejo
Copy link
Member Author

The problem is SOLVED! :)

Copy link
Member

@rubenzorrilla rubenzorrilla left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

#12918 (comment) is not clear to me yet. Aside of this, looks good to me.

@AlejandroCornejo
Copy link
Member Author

#12918 (comment) is not clear to me yet. Aside of this, looks good to me.

For now I've removed the EXTENDED_GAUSS_GI in the hexa since it is not consistent

Copy link
Member

@rubenzorrilla rubenzorrilla left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

All my concerns have been addressed.

@AlejandroCornejo AlejandroCornejo merged commit 885a57e into master Jan 6, 2025
11 checks passed
@AlejandroCornejo AlejandroCornejo deleted the adding-new-lobatto-quadratures-to-geometries branch January 6, 2025 20:44
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

5 participants