MembraneElement and PrestressMembraneElement

[PhilippeBussetta]

I have problem in my test with the MembraneElement, I have overcame it using PrestressMembraneElement. The MembraneElement is still tested? What are the difference between both elements (excepted pre-stress)? Why PrestressMembraneElement does not derive from MembraneElement?

[philbucher]

the MembraneElement was developed by @RiccardoRossi from what I know a newer version (PrestressMembraneElement) was implemented by @longchentum together with FormFinding and is currently further developed by @annasophie91 towards anisotropic prestress and formfinding-improvements

it was discussed already some time ago and from what I remember the MembraneElement will be removed after the development of the PrestressMembraneElement is finished (next few months)

[PhilippeBussetta]

@philbucher thanks

[philbucher]

closing since I thin this is solved

[RiccardoRossi]

@KratosMultiphysics/structural-mechanics let's make a note to remove MembraneElement for the next release

[PhilippeBussetta]

@philbucher , I have a problem with a simple test of PrestressMembraneElement, the results does not feel good and I have no idea why. It is a simple traction of one membrane element. Could you help me please?

Dynamic_Implicit_Bossak.zip

[philbucher]

@PhilippeBussetta I briefly looked at it, you are prescribing the displacement right? How do you compute the force then?

i don't know the impelemtation myself, but @annasophie91 is currently working on it, so you might try it with her branch

[PhilippeBussetta]

@philbucher Yes, I prescribed the displacement and I used the LinearElasticPlaneStress2DLaw I will try with her branch

[annasophie91]

I think it should work regardless which branch you use. How did you compute the force? Which force is it? I cannot run your example in my branch but I reconstructed a similar example in my branch which seems to work nicely. (See attachment)

tensile_test.gid.zip

[PhilippeBussetta]

Your result feel good. The force is sum of the nodal force X computed thanks to output json file.py. I run in the master.

Le ven. 16 mars 2018 17:46, annasophie91 notifications@github.com a écrit :

I think it should work regardless which branch you use. How did you compute the force? Which force is it? I cannot run your example in my branch but I reconstructed a similar example in my branch which seems to work nicely. (See attachment)

tensile_test.gid.zip https://github.com/KratosMultiphysics/Kratos/files/1820253/tensile_test.gid.zip

[image: grafik] https://user-images.githubusercontent.com/28961195/37536177-2342e08e-294a-11e8-93e7-fe65eac12089.png

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/KratosMultiphysics/Kratos/issues/1686#issuecomment-373792190, or mute the thread https://github.com/notifications/unsubscribe-auth/Aa_izKc8ESJGBQfI2dm9OmJ942Q6SLHVks5te_pqgaJpZM4Sr1eh .

[RiccardoRossi]

@PhilippeBussetta i would suggest making a small test taking as example one of the "patch_test" testcases.

That would allow us to debug the problem, and will remain for the future. (make the test very very small please, so that we can run it in the small test suite)

[PhilippeBussetta]

@annasophie91 I reproduce you results with your branch but when I try with more elements, a convergence problem appear. Do you check the tangential stiffness matrix is correctly computed?

[PhilippeBussetta]

@annasophie91 I tested a larger problem and I cannot solve it without numerical stiffness matrix It is the same problem with the truss element, only the numerical stiffness matrix allows us to solve large problem.

Membrane element

QuasiStaticMembrane.zip

Truss element

QuasiStatic.zip

[RiccardoRossi]

guys did we make any advance on this?

[philbucher]

I am not sure what the status of the work of @annasophie91 is We have to discuss this at the chair @inigolcanalejo

@PhilippeBussetta what do you mean by "numerical stiffness matrix"? to me the curves seem correct in a tensile test. I don't fully understand the problem

Regarding the truss element it was implemented by @KlausBSautter and he and others have tested it a lot meanwhile

[PhilippeBussetta]

@philbucher to me the curves seem correct also. Nevertheless, when the number of element increase, some difficulty appear during the Newton algorithm or the convergence is impossible. Then, we can compute the stiffness matrix thanks to perturbation method -- add small displacement and compute the rightHandSideVector thus the variation of this vector is used to fulfil the LeftHandSideMatrix. Of course, this strategy is more expensive but in some case could be used to verify the computation of the stiffness matrix.

for example with this algorithm:

    const double deltaDisp = 0.01;
    for (int i = 0; i < number_of_nodes; ++i) {
        int index = i * 3;
        VectorType rightHandSideVector = ZeroVector(mat_size);

        this->GetGeometry()[i].GetDof(DISPLACEMENT_X).GetSolutionStepValue() += deltaDisp;
        this->CalculateRightHandSide(rightHandSideVector, rCurrentProcessInfo);
        for (int k = 0; k < mat_size; ++k)
        {
            rLeftHandSideMatrix(k, index) = (rRightHandSideVector[k] - rightHandSideVector[k]) / deltaDisp;
        }
        this->GetGeometry()[i].GetDof(DISPLACEMENT_X).GetSolutionStepValue() -= deltaDisp;

        this->GetGeometry()[i].GetDof(DISPLACEMENT_Y).GetSolutionStepValue() += deltaDisp;
        this->CalculateRightHandSide(rightHandSideVector, rCurrentProcessInfo);
        for (int k = 0; k < mat_size; ++k)
        {
            rLeftHandSideMatrix(k, index + 1) = (rRightHandSideVector[k] - rightHandSideVector[k]) / deltaDisp;
        }
        this->GetGeometry()[i].GetDof(DISPLACEMENT_Y).GetSolutionStepValue() -= deltaDisp;

        this->GetGeometry()[i].GetDof(DISPLACEMENT_Z).GetSolutionStepValue() += deltaDisp;
        this->CalculateRightHandSide(rightHandSideVector, rCurrentProcessInfo);
        for (int k = 0; k < mat_size; ++k)
        {
            rLeftHandSideMatrix(k, index + 2) = (rRightHandSideVector[k] - rightHandSideVector[k]) / deltaDisp;
        }
        this->GetGeometry()[i].GetDof(DISPLACEMENT_Z).GetSolutionStepValue() -= deltaDisp;

    }

[RiccardoRossi]

@PhilippeBussetta in principle we should get some stiffness in the membrane provided that it is in stressed state. Without any stress the membrane is not statically defined. I am not sure if the technique you propose would solve such type of problems, however there is a major issue with your implementation since it is NOT threadsafe.

an element CANNOT change the value of its nodal variables, since elements are executed in parallel.

this implies that

  this->GetGeometry()[i].GetDof(DISPLACEMENT_Y).GetSolutionStepValue() -= deltaDisp;

is actually a forbidden operation (it is also very slow, since you use a GetDof, it would be faster to use the FastGetSolutionStepValue(DISPLACEMENT) -= DeltaDisp , however it is still incorrect in parallel)

[PhilippeBussetta]

@annasophie91 I cannot compile anymore your branch, I have a link error:

1>------ Build started: Project: KratosStructuralMechanicsApplication, Configuration: Release x64 ------
1>   Creating library D:/Kratos_membrane_analysis/cmake_build/applications/StructuralMechanicsApplication/Release/KratosStructuralMechanicsApplication.lib and object D:/Kratos_membrane_analysis/cmake_build/applications/StructuralMechanicsApplication/Release/KratosStructuralMechanicsApplication.exp
1>formfinding_IO_utility.obj : error LNK2019: unresolved external symbol "public: void __cdecl Kratos::ModelPartIO::WriteDataBlock<class Kratos::PointerVectorSet<class Kratos::Element,class Kratos::IndexedObject,struct std::less<unsigned __int64>,struct std::equal_to<unsigned __int64>,class boost::shared_ptr<class Kratos::Element>,class std::vector<class boost::shared_ptr<class Kratos::Element>,class std::allocator<class boost::shared_ptr<class Kratos::Element> > > > >(class Kratos::PointerVectorSet<class Kratos::Element,class Kratos::IndexedObject,struct std::less<unsigned __int64>,struct std::equal_to<unsigned __int64>,class boost::shared_ptr<class Kratos::Element>,class std::vector<class boost::shared_ptr<class Kratos::Element>,class std::allocator<class boost::shared_ptr<class Kratos::Element> > > > const &,class std::basic_string<char,struct std::char_traits<char>,class std::allocator<char> > const &)" (??$WriteDataBlock@V?$PointerVectorSet@VElement@Kratos@@VIndexedObject@2@U?$less@_K@std@@U?$equal_to@_K@5@V?$shared_ptr@VElement@Kratos@@@boost@@V?$vector@V?$shared_ptr@VElement@Kratos@@@boost@@V?$allocator@V?$shared_ptr@VElement@Kratos@@@boost@@@std@@@5@@Kratos@@@ModelPartIO@Kratos@@QEAAXAEBV?$PointerVectorSet@VElement@Kratos@@VIndexedObject@2@U?$less@_K@std@@U?$equal_to@_K@5@V?$shared_ptr@VElement@Kratos@@@boost@@V?$vector@V?$shared_ptr@VElement@Kratos@@@boost@@V?$allocator@V?$shared_ptr@VElement@Kratos@@@boost@@@std@@@5@@1@AEBV?$basic_string@DU?$char_traits@D@std@@V?$allocator@D@2@@std@@@Z) referenced in function "public: void __cdecl Kratos::FormfindingIOUtility::PrintPrestressData(void)" (?PrintPrestressData@FormfindingIOUtility@Kratos@@QEAAXXZ)
1>D:\Kratos_membrane_analysis\cmake_build\applications\StructuralMechanicsApplication\Release\KratosStructuralMechanicsApplication.pyd : fatal error LNK1120: 1 unresolved externals
1>Done building project "KratosStructuralMechanicsApplication.vcxproj" -- FAILED.

[KlausBSautter]

@PhilippeBussetta

It is the same problem with the truss element, only the numerical stiffness matrix allows us to solve large problem

What do you mean by this? We are using the truss element for large problems which exhibit large deformations without any problem. Can you specify your problem?

[philbucher]

@PhilippeBussetta thanks for the explanation

Unfortunately @annasophie91 is not working at the chair any more @inigolcanalejo and I are taking over her work and are currently contacting her abt the status of it.

Her branch does not compile any more bcs it has not been ported to pybind yet.

Therefore I would ask you to be a bit patient, I will get back to you asap, when we have the infos ourselves.

[PhilippeBussetta]

@philbucher Thank you

@KlausBSautter by large problem I would mean large number of element

[PhilippeBussetta]

@RiccardoRossi I know that this computation of the stiffness is no good regarding numerical performance as well as parallel computing but is useful to verify the good definition in many cases.

I don't understand what do you mean. Without stress, is not possible to define the stiffness of the membrane?

[philbucher]

@PhilippeBussetta there is no stiffness in a membrane if it is not stressed. Think of a piece of cloth

[PhilippeBussetta]

@philbucher I understand the physic of the model but you mean that the total rigidity matrix is null without stress or just the bending part? Did you allow the element to be compressed in-plane? or only buckling? No rigidity in compression? Then did you observe any numerical instability?

I would like a element with only in-plane rigidity -- shell without bending stiffness. It is the case of your membrane element? I would like to model the in-plane deformations.

I try to model composite material with continuous fibres with 2D and 1D elements -- like into the test Membrane_Q4_Truss

With Truss and membrane elements numerical instability appear -- nan into the stiffness matrix. The problem feel solved with Truss and ShellThin elements.

[philbucher]

@PhilippeBussetta the element has NO bending stiffness => it behaves like a piece of cloth, it does not have any "rigidity" in compression! => therefore you can only use it for in-plane deformations if you have tensile loads, it will fail in compression.

What you can do however is to fix the out-of-plane motion (your case is 2D right?). Then it should be able to represent also compression loads. Imagine a rope: If you compress it it fails. If you support is (imagine a rope in a pipe), then it can also take compressive loads

I hope I am explaining myself not too badly

[PhilippeBussetta]

@philbucher Thanks, I have a better understanding of this element. Is it the same with the truss element, no "rigidity" in compression?

[KlausBSautter]

no, the truss can take both compression and tension. If you use the inheriting element CABLE_ELEMENT you will have this feature. But for the membrane element, I am also confused now. I thought it is a plate in membrane action and thus can take also compression.

[longchentum]

@KlausBSautter. Yes, the PrestressMembraneElement can take compression indeed. We don't consider the wrinkling effect in the prestressed membrane element. Thus, when using PrestressMembraneElement for non-linear analysis, one should be aware that the prestress is big enough such that there is no compression.

[philbucher]

@PhilippeBussetta sorry for my confusion I hope @longchentum made it more clear now

[philbucher]

closing since inactive