[GeoMechanicsApplication] Extract biot coefficient calculation #12398

rfaasse · 2024-05-23T10:23:02Z

📝 Description
This PR extracts pre-calculated lists for the biot coefficients and creates a utility function (used by diff order and non-diff order) to calculate the biot modulus inverse for a single integration point (including unit tests). In a next PR, the extraction of the element-wide list for the biot modulus inverse will be done, that is out of scope for this one.

Changes:

Removing redundant hasBiotCoefficient flag
Extracted lists for the biot coefficients
Added utility function for biot modulus inverse with unit tests and used in diff and non-diff order elements
Extracting Bulk Modulus calculation with tests
Extract utility function for calculating the biot coefficients with tests and used in diff and non-diff order elements

…ficients` and renamed it `CalculateBiotModulusInverse`

…tead of ElementVariables

… instead of ElementVariables and input the properties as an argument

… to diff order (i.e. extract lists for biot coefficients and use utility function)

…iot-coefficient-calculation

rfaasse · 2024-05-23T11:44:05Z

applications/GeoMechanicsApplication/tests/cpp_tests/test_transport_equation_utilities.cpp

It might be a point of discussion if we want all transport_equation_utilities tests in a single file

avdg81

A clean and clear extraction of functionality, which makes the code more readable and easier to understand. Well done! I have a few minor suggestions and a question for you.

applications/GeoMechanicsApplication/custom_elements/U_Pw_small_strain_element.cpp

applications/GeoMechanicsApplication/custom_elements/U_Pw_small_strain_element.hpp

applications/GeoMechanicsApplication/custom_elements/small_strain_U_Pw_diff_order_element.cpp

applications/GeoMechanicsApplication/custom_utilities/transport_equation_utilities.hpp

avdg81 · 2024-05-23T12:47:55Z

applications/GeoMechanicsApplication/tests/cpp_tests/test_transport_equation_utilities.cpp

+    const double derivative_of_saturation = 0.2;
+
+    const double expected_value = -0.09999992485;
+    KRATOS_EXPECT_DOUBLE_EQ(GeoTransportEquationUtilities::CalculateBiotModulusInverse(


Wouldn't KRATOS_EXPECT_NEAR be more robust?

Depends on the tolerance of course, under the hood this macro does: KRATOS_EXPECT_NEAR(a,b,std::numeric_limits<double>::epsilon(). Do we prefer a different tolerance? Problem is, with the currnet numbers we need to be quite precise to check that the results are actually okay (no 1e-16, but 1e-9 is needed).

applications/GeoMechanicsApplication/tests/cpp_tests/test_transport_equation_utilities.cpp

avdg81 · 2024-05-23T13:04:00Z

applications/GeoMechanicsApplication/tests/cpp_tests/test_transport_equation_utilities.cpp

+                           biot_coefficient, degree_of_saturation, derivative_of_saturation, properties) > large_number)
+}
+
+KRATOS_TEST_CASE_IN_SUITE(CalculateBiotModulusInverse_DoesNotThrow_ForEmptyProperties, KratosGeoMechanicsFastSuite)


Is this desired behavior, since it would require me to test the returned value before actually using it ...?

This question is exactly the reason I wanted to make the current behavior (this has just moved) explicit in this test. I don't know exactly what the desired behavior is. E.g. do we want to always check if all the properties are set with physical values inside this function, or do we assume it's okay? It is user input though (via the MaterialParameters), so it could indeed result in NaN now if the wrong inputs are given.

I can of course add checks and throw (and UT for that), what's your opinion?

I would suggest to not change the behavior yet, since that would be beyond the scope of this PR. However, I would like you to add an issue to the product backlog about it, so we can revisit it later. Now it has been documented through a unit test, which is a first step in the right direction, so thanks for that.

WPK4FEM

Dear Richard,

The purpose of creating lists for the BiotCoefficient alpha and Biot Modulus Inverse 1/Q is fulfilled here. Good. However, I strongly think that the computation should be done directly on given material input and not on a constitutive matrix.

We can discuss whether that should be in this PR or in a next one.

Regards, Wijtze Pieter

WPK4FEM · 2024-05-23T13:12:18Z

applications/GeoMechanicsApplication/custom_elements/small_strain_U_Pw_diff_order_element.cpp

-            (rVariables.BiotCoefficient - rProp[POROSITY]) / rProp[BULK_MODULUS_SOLID] +
-            rProp[POROSITY] / rProp[BULK_MODULUS_FLUID];
+        // calculate Bulk modulus from stiffness matrix
+        return 1.0 - CalculateBulkModulus(rConstitutiveMatrix) / rProp[BULK_MODULUS_SOLID];


Contents of CalculateBUlkModulus:
return ConstitutiveMatrix(0, 0) - (4.0 / 3.0) * ConstitutiveMatrix(IndexG, IndexG);
The ConstitutiveMatrix is computed (relevant entries only)as:
const Properties& r_material_properties = rValues.GetMaterialProperties();
const auto E = r_material_properties[YOUNG_MODULUS];
const auto NU = r_material_properties[POISSON_RATIO];
const double c1 = (1.0 - NU) * c0;
const double c3 = this->GetConsiderDiagonalEntriesOnlyAndNoShear() ? 0.0 : (0.5 - NU) * c0;
C(INDEX_2D_PLANE_STRAIN_XX, INDEX_2D_PLANE_STRAIN_XX) = c1;
C(INDEX_2D_PLANE_STRAIN_XY, INDEX_2D_PLANE_STRAIN_XY) = c3;

I strongly think that the input of CalculateBulkModulus should be rProp and no a constitutive matrix. That changed the argument lists of all functions calling CalculateBulkModulus. Further CalculateBulkModulus should be in a utility thing and functionality like CalculateBiotCoefficient as well. E.g. both in the transport equation utilities.

You may argue that it is out of scope for this PR, but the above indicates where I would like to end.

I completely agree that here it is a lot of overkill to go through all the calculations (the b matrices, deformation gradients, strains and constitutive matrix), in some cases just to calculate a bulk modulus that we could derive from the material.

I will move the CalculateBulkModulus and CalculateBiotCoefficient(s) to the transport equation utilities (and add some tests), I agree we should make this more generic and not have more than one implementation and it shouldn't be a lot of work.

However, taking a simpler approach for deriving the bulk modulus out of material parameters and also do the testing that this works for all constitutive laws (and not only for this linear elastic example) is in my opinion indeed out of scope for this refactoring PR.

I took the elastic constitutive thing, because I think that it is always the elastic matrix and not some plastic or otherwise nonlinear thing. That is why I'm hammering on taking it from material parameters. There are no objections against doing that in a different PR.

WPK4FEM

Thank you for all efforts, I can live with this implementation

avdg81

This one is good to go.

avdg81 · 2024-05-23T14:42:05Z

applications/GeoMechanicsApplication/tests/cpp_tests/test_transport_equation_utilities.cpp

+                           biot_coefficient, degree_of_saturation, derivative_of_saturation, properties) > large_number)
+}
+
+KRATOS_TEST_CASE_IN_SUITE(CalculateBiotModulusInverse_DoesNotThrow_ForEmptyProperties, KratosGeoMechanicsFastSuite)


I would suggest to not change the behavior yet, since that would be beyond the scope of this PR. However, I would like you to add an issue to the product backlog about it, so we can revisit it later. Now it has been documented through a unit test, which is a first step in the right direction, so thanks for that.

…s yet).

WPK4FEM

Remains good to go.

avdg81

I'm very satisfied with the result of this PR. Let's merge it.

rfaasse added 11 commits May 21, 2024 14:33

Removing redundant hasBiotCoefficient flag

beb415e

Removing redundant hasBiotCoefficient flag from non-diff-order element

7d8713d

Made constitutive matrix direct input instead of via ElementVariables

5d5f060

Extracted lists for the biot coefficients for diff order

48dd387

Moved calculation of the biot coefficients out of `InitializeBiotCoef…

cc94954

…ficients` and renamed it `CalculateBiotModulusInverse`

Changed signature of CalculateBiotModulusInverse to take inputs ins…

f831784

…tead of ElementVariables

Changed signature of CalculateBiotModulusInverse to use return type…

b18af44

… instead of ElementVariables and input the properties as an argument

Added utility function for biot modulus inverse with unit tests

c9f39b8

Minor refactoring of new function and added a test for empty properties

f1a5df8

Used utility in non-diff order element and did comparable refactoring…

955f75e

… to diff order (i.e. extract lists for biot coefficients and use utility function)

Merge remote-tracking branch 'origin/master' into geo/12319-extract-b…

0ab8637

…iot-coefficient-calculation

rfaasse linked an issue May 23, 2024 that may be closed by this pull request

[GeoMechanicsApplication] Split the calculation of the inputs for the utility functions for element matrices and the calling of the functions #12319

Closed

rfaasse added 3 commits May 23, 2024 12:27

Removed unused parameter

0d3b600

Removed some redundant comments

4eef777

Removed some redundant lines

4433ab6

rfaasse marked this pull request as ready for review May 23, 2024 11:42

rfaasse commented May 23, 2024

View reviewed changes

rfaasse requested review from avdg81 and WPK4FEM May 23, 2024 11:49

avdg81 reviewed May 23, 2024

View reviewed changes

rfaasse self-assigned this May 23, 2024

rfaasse added GeoMechanics Issues related to the GeoMechanicsApplication Refactor When code is moved or rewrote keeping the same behavior labels May 23, 2024

WPK4FEM requested changes May 23, 2024

View reviewed changes

Processing review comments Anne

2c95436

WPK4FEM previously approved these changes May 23, 2024

View reviewed changes

avdg81 previously approved these changes May 23, 2024

View reviewed changes

rfaasse added 3 commits May 23, 2024 17:00

Extracting Bulk Modulus calculation with tests

8fc7a57

Extract utility function for calculating the biot coefficients (no UT…

558791b

…s yet).

Added unit tests for the biot coefficient calculations

8e8cc45

rfaasse dismissed stale reviews from avdg81 and WPK4FEM via 8e8cc45 May 23, 2024 15:47

WPK4FEM approved these changes May 24, 2024

View reviewed changes

rfaasse requested a review from avdg81 May 24, 2024 07:10

avdg81 approved these changes May 24, 2024

View reviewed changes

rfaasse merged commit 281c6d3 into master May 24, 2024
11 checks passed

rfaasse deleted the geo/12319-extract-biot-coefficient-calculation branch May 24, 2024 08:42

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

[GeoMechanicsApplication] Extract biot coefficient calculation #12398

[GeoMechanicsApplication] Extract biot coefficient calculation #12398

rfaasse commented May 23, 2024 •

edited

Loading

rfaasse May 23, 2024

avdg81 left a comment

avdg81 May 23, 2024

rfaasse May 23, 2024

avdg81 May 23, 2024

rfaasse May 23, 2024 •

edited

Loading

avdg81 May 23, 2024

WPK4FEM left a comment

WPK4FEM May 23, 2024

rfaasse May 23, 2024

WPK4FEM May 23, 2024

WPK4FEM left a comment

avdg81 left a comment

avdg81 May 23, 2024

WPK4FEM left a comment

avdg81 left a comment

[GeoMechanicsApplication] Extract biot coefficient calculation #12398

[GeoMechanicsApplication] Extract biot coefficient calculation #12398

Conversation

rfaasse commented May 23, 2024 • edited Loading

Choose a reason for hiding this comment

avdg81 left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

rfaasse May 23, 2024 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

WPK4FEM left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

WPK4FEM left a comment

Choose a reason for hiding this comment

avdg81 left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

WPK4FEM left a comment

Choose a reason for hiding this comment

avdg81 left a comment

Choose a reason for hiding this comment

rfaasse commented May 23, 2024 •

edited

Loading

rfaasse May 23, 2024 •

edited

Loading