Re: [Getfem-users] Hyperelastic law Benchmark

[Yves Renard]

Dear David,

No, no idea. Often, if the load is too large, the Newton fails to find the 
solution. It is quite normal. You can decompose into small load steps. 
Alternatively, you can also use the (beautiful) continuations tools of Getfem 
developed by Tomas Ligursky which automatically adapt the step to try to find a 
solution. If the continuation procedure do not find a solution, this may be of 
course due to the fact that there could be no solution to the problem for a 
large load ... This could happen for very large loads.

Best regards,

Yves.

----- Mail original -----
De: "David Danan" <address@hidden>
À: "logari81" <address@hidden>
Cc: address@hidden
Envoyé: Jeudi 22 Juin 2017 12:27:23
Objet: Re: [Getfem-users] Hyperelastic law Benchmark

Dear getfem users,

i encounter stability issues with the second problem of the same benchmark.
 http://rspa.royalsocietypublishing.org/content/471/2184/20150641#sec-15

By using [39], i have tried to improve the stability but the computation
still failed at the first iteration anyway.

So far i can use either:
-The classical version of the hyperelastic law with an incompressibility
condition (lagrange multiplier)
-The classical version of the hyperelastic law with a
quasi-incompressibility condition (penalization)
-A stabilized version of the law (incompressibility with Lagrange
multiplier and penalization)
-An isochoric version that use the isochoric component of the deformation
gradient with either an incompressibility condition, a penalization or both
of them.

Note that the convergence is reached if i decrease the value of the
pressure.

Please, find enclosed:
-the program Ellipsoid.cc and the parameter file Ellipsoid.param
-Two meshes provided by the benchmark. I have also tried with several other
meshes, by using the methods described above, with the same result.

Do you have any idea/advices regarding these issues?
Thanks a lot.

Yours sincerely,
David.




2017-06-15 20:53 GMT+02:00 David Danan <address@hidden>:

> Dear Kostas,
>
> first, thank you for your very fast answer.
>
> So far, i didn't even think about using the high level generic assembly
> language for this (and it's a shame given how easy it is to use in this
> configuration).
> I followed your suggestion and tried your expression (which seems
> perfectly correct to me) and the deformation was, as far as i could tell,
> identical to the one in the Results_nosym.png.
>
> However, your second remark was spot on because the paper actually said it
> explicitly:
> "Please note that in all problems the direction of the pressure boundary
> condition changes with the deformed surface orientation, and its magnitude
> scales with the deformed area."
>
> Therefore, thanks a lot for providing your notes here, it is very
> helpful...And it seems to actually work just fine! There is no visible
> difference between the deformed configuration given by Getfem and the
> reference solution with this modification.
>
> Once again, thank you!
>
> David.
>
>
> 2017-06-14 17:52 GMT+02:00 Konstantinos Poulios <address@hidden>:
>
>> It doesn't seem to be explicitly stated in the paper but since it is
>> about cardiovascular simulation I guess that applying the pressure as a
>> follower load is the standard thing to do, so I am providing you here with
>> my notes on follower loads:
>>
>> [image: Inline image 1]
>>
>>
>> So if your p is the actual blood pressure you need the upper right case
>> of the table with q=-p. If at some point you also need shear stresses from
>> the fluid you can also use the second row of the table.
>>
>> BR
>> Kostas
>>
>>
>> On Wed, Jun 14, 2017 at 4:57 PM, Konstantinos Poulios <
>> address@hidden> wrote:
>>
>>> Dear David
>>>
>>> Have you tried the high level generic assembly language for this?
>>>
>>> In principle you should be able to provide GetFEM with your energy
>>> density function and let GetFEM do the necessary derivations.
>>>
>>> Instead of
>>>
>>> getfem::add_finite_strain_elasticity_brick(...)
>>>
>>> you have to call
>>>
>>>  getfem::add_nonlinear_generic_assembly_brick
>>> (md, mim, "0.5*C*(exp([[bf,bfs,bfs],[bfs,bt,bt],[bfs,bt,bt]]:(Green_La
>>> grangian(Grad_u+Id(3)).*Green_Lagrangian(Grad_u+Id(3))))-1)");
>>>
>>> with C,bf,bfs and bt scalar parameters defined with
>>> md.add_initialized_scalar_data(...).
>>>
>>> I hope I got the expression from the paper right. Can you give it a try?
>>>
>>> Then the other question is how the applied surface pressure p is
>>> distributed, if it is a follower load you need a more complex expression
>>> than
>>>
>>> "-Pressure*Normal.Test_u"
>>>
>>> Because "Normal" is in the undeformed configuration.
>>>
>>> BR
>>> Kostas
>>>
>>>
>>> On Wed, Jun 14, 2017 at 3:58 PM, David Danan <address@hidden>
>>> wrote:
>>>
>>>> Dear Getfem users,
>>>>
>>>> i am trying to implement a new hyperelastic law and, in order to
>>>> validate my results, i am using the following Benchmark
>>>> http://rspa.royalsocietypublishing.org/content/471/2184/20150641#sec-15
>>>>
>>>> There are 3 problems, for now i am working on the first one that is to
>>>> say the deformation of a 3D rectangular beam clamped on one side and with a
>>>> pressure applied to the bottom face.
>>>>
>>>> While the deformed configuration given by Getfem is relatively close to
>>>> the reference(s) solution(s) provided by the benchmark, a visible
>>>> difference between them still remains and i don't understand where it comes
>>>> from.
>>>>
>>>> The material is governed by a transversely isotropic constitutive law
>>>> with an incompressibility constraint, often used in cardiac modelling,
>>>> where the strain energy function is a function of the components of the
>>>> Green–Lagrange strain tensor *E.*
>>>>
>>>> I tried 2 differents implementations of this law:
>>>> -the first use the symmetry of the Green-Lagrange strain tensor to
>>>> simplify the strain energy function
>>>> -The second does not (ergo it is necessary to write the 9 components of
>>>> the second piola Kirchhoff stress tensor and the 81 components of the
>>>> fourth order tensor)
>>>>
>>>> Please find enclosed
>>>> -the comparison in the first case: Results.png
>>>> -the comparison in the second case: Results_nosym.png (slightly better
>>>> results but 15 times as slow as the first version)
>>>> -the python program used to compute the derivative and second
>>>> derivative of the strain energy function in the first case.
>>>> -the implementation of the laws in getfem_nonlinear_elasticity.cc and
>>>> getfem_nonlinear_elasticity.h
>>>> -The program Guccione.cc and Guccione.param used to produce these very
>>>> pictures
>>>>
>>>> in both pictures, the reference solution is in grey.
>>>> The computation uses Q2/Q1 elements (displacement/lagrange multiplier),
>>>> since there is no restrictions regarding these aspects.
>>>> I have tried with a quasi-incompressibility condition instead of the
>>>> Lagrangian multiplier: same result (which was to be expected).
>>>> I have also tried with other meshes (more or less refined) used by
>>>> other teams but in vain.
>>>>
>>>> Could someone have a look and provide some advices regarding this
>>>> case/tell me what i am doing wrong?
>>>>
>>>> Thanks a lot.
>>>>
>>>> Yours sincerely,
>>>> David.
>>>>
>>>>
>>>
>>
>