KratosMultiphysics / Kratos

Kratos Multiphysics (A.K.A Kratos) is a framework for building parallel multi-disciplinary simulation software. Modularity, extensibility and HPC are the main objectives. Kratos has BSD license and is written in C++ with extensive Python interface.
https://kratosmultiphysics.github.io/Kratos/
Other
1.02k stars 245 forks source link

Update of nonlinear element matrices, function calls. #7159

Closed emiroglu closed 3 years ago

emiroglu commented 4 years ago

Hi all,

I have been working on the implementation of modal derivatives for generating a basis to use in @KratosMultiphysics/hromapplication . For this I compute the derivatives of the element matrices with a local finite differencing. In case of a modal derivative with respect to modal coordinates the forward finite differencing formula is simply:

herein:

The implementations can be found in RomApp_ModalDerivatives branch if necessary.

This approach works well for truss and membrane elements but I have not found a solution for doing this with shell and volume elements, which are essential for my simulations.

My open questions are:

Appreciate any help and tips.

@KratosMultiphysics/structural-mechanics

loumalouomega commented 4 years ago

You can try the solid shell element. You can use for shell simulations but the mesh and constitutive laws are the ones of a solid

El mar., 30 jun. 2020 11:19, emiroglu notifications@github.com escribió:

Hi all,

I have been working on the implementation of modal derivatives for generating a basis to use in @KratosMultiphysics/hromapplication https://github.com/orgs/KratosMultiphysics/teams/hromapplication . For this I compute the derivatives of the element matrices with a local finite differencing. In case of a modal derivative with respect to modal coordinates the forward finite differencing formula is simply:

https://render.githubusercontent.com/render/math?math=%5Cdfrac%7BdK_%7Bel%7D%7D%7Bd%5Ceta_%7Bi%7D%7D=%5Cdfrac%7BK_%7Bel%7D(%5Calpha*%5Cphi_%7Bi%7D)-K%7Bel%7D(0)%7D%7B%5Calpha%7D herein:

https://render.githubusercontent.com/render/math?math=%5CK_%7Bel%7D(u):%5Ctext%7Belement%20stiffness%20matrix%20at%20displacement%20u%7D

https://render.githubusercontent.com/render/math?math=%5Ceta_%7Bi%7D:%5Ctext%7Bmodal%20coordinate%7D

https://render.githubusercontent.com/render/math?math=%5Cphi_%7Bi%7D:%5Ctext%7Bmode%20shape%7D

https://render.githubusercontent.com/render/math?math=%5Calpha:%5Ctext%7Bstep%20size%7D The implementations can be found in RomApp_ModalDerivatives RomApp_ModalDerivatives branch if necessary.

This approach works well for truss and membrane elements but I have not found a solution for doing this with shell and volume elements, which are essential for my simulations.

My open questions are:

  • Are there any shell elements that use nonlinear Green-Lagrange strains for the computation of element matrices? Maybe @peterjwilson https://github.com/peterjwilson could have an idea on this?
  • I have not managed to trigger an update for volume elements (currently using TotalLagrangianElement3D8N). Is there a specific function that I should call at the element level? I noticed that calling Element->InitializeNonLinearIteration() and Element->FinalizeNonLinearIteration() functions update the section properties of the corotational shell elements and at least help including the displacements in the element stiffness matrix due to the change in coordinate transformations, but they were not effective for volume elements.

Appreciate any help and tips.

@KratosMultiphysics/structural-mechanics https://github.com/orgs/KratosMultiphysics/teams/structural-mechanics

— You are receiving this because you are on a team that was mentioned. Reply to this email directly, view it on GitHub https://github.com/KratosMultiphysics/Kratos/issues/7159, or unsubscribe https://github.com/notifications/unsubscribe-auth/AEYQZABHUMXN2HTPFZBTBGDRZGTAXANCNFSM4OMAOYWQ .

peterjwilson commented 4 years ago

@emiroglu I'm pretty sure the four shell elements (thin/thick quad/tri) use linear strains. A co-rotational filter can be applied to extend those shells into the large deformation small strain regime.

Going from memory, and a quick check of the structural mechanics app now, InitializeNonLinearIteration calls the assigned co-rotation transformation in custom_utilities to do the magic they do. Notice the shells have their own CR transformations defined here (and volume elements don't seem to). I'm not sure how this interfaces with what you want to do, but if you need a CR transformation for your volume element maybe check out adapting/deriving from EICR.hpp (Element Independent CR from Felippa).

emiroglu commented 4 years ago

You can try the solid shell element. You can use for shell simulations but the mesh and constitutive laws are the ones of a solid

Thanks Vicente. These elements seem to be triangular prisms. The element type does not seem to be in the GiD problem type. How do you generate meshes and corresponding mdpa files for it?

emiroglu commented 4 years ago

@emiroglu I'm pretty sure the four shell elements (thin/thick quad/tri) use linear strains. A co-rotational filter can be applied to extend those shells into the large deformation small strain regime.

Going from memory, and a quick check of the structural mechanics app now, InitializeNonLinearIteration calls the assigned co-rotation transformation in custom_utilities to do the magic they do. Notice the shells have their own CR transformations defined here (and volume elements don't seem to). I'm not sure how this interfaces with what you want to do, but if you need a CR transformation for your volume element maybe check out adapting/deriving from EICR.hpp (Element Independent CR from Felippa).

Thanks for the response. I essentially need elements that can handle nonlinear strains. In case of volume elements (total_lagrangian.cpp), I can see that TotalLagrangianElement3D8N is supposed to do the job but InitializeNonLinearIteration function does not seem to do anything for it. Is there any other function that I should call? Is it related to the material law that I am using (currently using LinearElastic3DLaw)? Maybe @loumalouomega has better clue on this?

loumalouomega commented 4 years ago

Thanks for the response. I essentially need elements that can handle nonlinear strains. In case of volume elements (total_lagrangian.cpp), I can see that TotalLagrangianElement3D8N is supposed to do the job but InitializeNonLinearIteration function does not seem to do anything for it. Is there any other function that I should call? Is it related to the material law that I am using (currently using LinearElastic3DLaw)? Maybe @loumalouomega has better clue on this?

TL does not need to did nothing in the InitializeNonLinearIteration, the integration is done in the initial configuration, UL updates the deformation gradient in the FinalizeNonLinearIteration. LinearElastic3DLaw is linear, you can use other laws for NL behavior, like the HyperElastic3DLaw

emiroglu commented 4 years ago

@loumalouomega thanks. It works with the UL element. What puzzles me is the computation of the deformation gradient for TL formulation. The deformation gradient is defined as:

but in the current implementation it returns an identity matrix. This is actually wrong.

loumalouomega commented 4 years ago

@loumalouomega thanks. It works with the UL element. What puzzles me is the computation of the deformation gradient for TL formulation. The deformation gradient is defined as: ![](https://render.githubusercontent.com/render/math?math=\mathbf{F}=\mathbf{I}%2B\dfrac{\partial \mathbf{u}}{\partial \mathbf{X}},) but in the current implementation it returns an identity matrix. This is actually wrong.

TL implementation was originally done by @RiccardoRossi. If it returns identity it is because your displacements are supposed to be zero, can you check it?

emiroglu commented 4 years ago

In the implementation the deformation gradient is computed as in here:

as far as I understood the J is also computed in the initial configuration (see here), so the result returns an identity matrix.

loumalouomega commented 4 years ago

J0 is in the initial configuration, J is in the deformed configuration. If there are displacements should be different of zero

emiroglu commented 4 years ago

Alright just found the problem. It was in the update of the coordinates. Thanks again!

emiroglu commented 4 years ago

@loumalouomega is there a specific way to create an MDPA file for the solid shell elements? As far as I understand from the code these elements have to have a triangular prism geometry. But GiD refuses to create such elements if not chosen in the Kratos problem type and the problem type does not include these elements.

What is your approach to creating a mesh for using these elements?

loumalouomega commented 4 years ago

@loumalouomega is there a specific way to create an MDPA file for the solid shell elements? As far as I understand from the code these elements have to have a triangular prism geometry. But GiD refuses to create such elements if not chosen in the Kratos problem type and the problem type does not include these elements.

What is your approach to creating a mesh for using these elements?

There is a trick I implemented, outside GiD which gives a lot of issues with the prism geometry, it consists in mesh it as a triangular mesh, and in the sprism_process select the option to extrude the geometry and generate a mpda.

emiroglu commented 3 years ago

I forgot this issue and opened a new one #8375 which has a similar problem. The problem persists from my side. Any ideas on how to work this out? The most current implementation is in IMWind branch. So far this approach works with truss and volume elements. It also seems to be working with CrBeam element after the modifications I made in the element. Still, I am getting very weird results when there are mixed element types and/or also it does not work at all with shell elements. Any help or ideas are very much appreciated...

philbucher commented 3 years ago

Might not be related directly to your issue but could be of importance since it is about the design and the intended use

@emiroglu what modifications did you have to do to the beam? Can you elaborate so that I can understand better

emiroglu commented 3 years ago

The problem with the corotational elements in general is the necessity for using incremental deformations. For example the CrBeam element modifies the previous deformation in InitializeNonLinearIteration (here) and gets the current deformation and stores it, instead of reading it from SolutionStepValue everytime. In order to make this flexible I had to add a GetValuesVector call to increment deformation function which then caused many functions including what you mentioned to lose constness. Of course instead of modifying a member variable one can just read from SolutionStepValue to make sure the constness you want remains.

Vector CrBeamElement3D2N::GetIncrementDeformation()
{
    KRATOS_TRY;
    GetValuesVector(mDeformationCurrentIteration);                 <-------------------------------- HERE 
    return mDeformationCurrentIteration - mDeformationPreviousIteration;
    KRATOS_CATCH("")
}

I don't even know where to start with the shell elements though...

philbucher commented 3 years ago

The problem with the corotational elements in general is the necessity for using incremental deformations. For example the CrBeam element modifies the previous deformation in InitializeNonLinearIteration (here) and gets the current deformation and stores it, instead of reading it from SolutionStepValue everytime. In order to make this flexible I had to add a GetValuesVector call to increment deformation function which then caused many functions including what you mentioned to lose constness. Of course instead of modifying a member variable one can just read from SolutionStepValue to make sure the constness you want remains.

Vector CrBeamElement3D2N::GetIncrementDeformation()
{
    KRATOS_TRY;
    GetValuesVector(mDeformationCurrentIteration);                 <-------------------------------- HERE 
    return mDeformationCurrentIteration - mDeformationPreviousIteration;
    KRATOS_CATCH("")
}

I don't even know where to start with the shell elements though...

@KlausBSautter seems like this comes back to haunt us again :D

this was done in #5109 back then @KlausBSautter and I tried a lot to do what was necessary to make the beam work, but we couldn't make it work back then Maybe we can take another look

for details on how it was done like this see #5109

KlausBSautter commented 3 years ago

can we finally have some geometrically non-linear beams and shells? :palm_tree:

emiroglu commented 3 years ago

Please see #8375 . There I computed derivatives of the element matrices with the order of calls given in this comment. It seems to be working with the following elements: