Dear Yves, this time, it's perfect, the test passed. Thanks a lot!
The only question i have left now is the possibility to compute the divergence of a tensor field ( arising from the gradient of an existing variable in a model), any suggestions? Best regards, David. Le mar. 14 juin 2022 à 13:59, Renard Yves <yves.ren...@insa-lyon.fr> a écrit : > Dear Danan, > > I made an additional correction. Normally the tangent term of matrix_j2 is > ok now. > > Best regards, > > Yves > > > Le 13/06/2022 à 17:37, David Danan a écrit : > > Dear Yves, > > i gave it a try, recompiled the master version available at > https://git.savannah.nongnu.org/cgit/getfem.git/commit/ with your last fix > but it still fails to converge, unfortunately. > > I have tested this with the program enclosed in the third email in this > thread, it seems there is still an issue as far as i can tell. > > Best regards, > David. > > Le lun. 13 juin 2022 à 15:52, Renard Yves <yves.ren...@insa-lyon.fr> a > écrit : > >> Dear Danan and Kostas, >> >> I had a rapid check to the implementation of Mooney-Rivlin hyperelastic >> law and matrix_j2 operators in getfem_nonlinear_elasticity.cc >> >> I do not see any difference between the implementation and the expression >> in the documentation ( >> https://getfem.org/userdoc/model_nonlinear_elasticity.html#some-recalls-on-finite-strain-elasticity) >> for the Mooney-Rivlin law. However, it seems to me that a (det) is missing >> in the matrix_j2 gradient expression. I committed a fix. Could you try >> again if there is still a difference between the use of Mooney-Rivlin law >> and its expression using matrix_j2 ? >> >> Best regards, >> >> Yves >> >> >> Le 05/06/2022 à 20:01, David Danan a écrit : >> >> Dear Kostas, >> >> once again, thank you for pointing out where the problem is! >> As soon as i have some time, it would be a pleasure to have a look and >> trying to fix this. Thanks a lot! >> >> Regarding the question in my second message in this thread about the >> computation of the divergence of a tensor field, arising from the gradient >> of an existing variable in a model, do you have any idea? >> >> BR, >> David. >> >> >> >> Le sam. 28 mai 2022 à 22:10, Konstantinos Poulios < >> logar...@googlemail.com> a écrit : >> >>> I see now that Yves has left a comment "to be verified": >>> >>> // Second derivative >>> void second_derivative(const arg_list &args, size_type, size_type, >>> base_tensor &result) const { // To be verified >>> >>> in the second derivative of matrix_j2_operator in >>> getfem_nonlinear_elasticity.cc, so maybe you can help by checking/fixing >>> this. >>> >>> BR >>> Kostas >>> >>> >>> >>> On Sat, May 28, 2022 at 10:05 PM Konstantinos Poulios < >>> logar...@googlemail.com> wrote: >>> >>>> Dear David, >>>> >>>> It works for me with >>>> model.add_nonlinear_generic_assembly_brick(mim, >>>> "paramsIMR(1)* ( Matrix_j1(F*F') - 3 )\ >>>> + paramsIMR(2)* ( >>>> Matrix_i2(F*F')*pow(Det(F*F'),-2/3) - 3 )") >>>> so you must have hit a bug in Matrix_j2. >>>> >>>> BR >>>> Kostas >>>> >>>> On Sat, May 28, 2022 at 4:42 PM David Danan <daviddanan9...@gmail.com> >>>> wrote: >>>> >>>>> Dear Yves, >>>>> >>>>> first, thanks for your detailed answers! >>>>> >>>>> For the first question, at this point, i guess it's much more >>>>> convenient at this point to share a script containing a minimal (not) >>>>> working test regarding this difference. >>>>> I have kept the incompressibility brick in both cases. >>>>> >>>>> For the second question, it seems clear, thanks. >>>>> For the third question, providing the expression of the stress tensor >>>>> directly instead is fine by me, it's not very difficult using the >>>>> invariants already available in the GWFL combined with some macro. >>>>> >>>>> While i am at it, i have another question. >>>>> It is somehow related to the one in the following thread but for the >>>>> divergence of a tensor field >>>>> >>>>> https://lists.nongnu.org/archive/html/getfem-users/2015-03/msg00009.html >>>>> >>>>> I would like to check the satisfaction of the equilibrium equation >>>>> (i.e. div(P) + f0, where Pi is a stress tensor and f0 is a volumic force). >>>>> To do so, i gave it a try without volumic forces >>>>> value=gf.asm_generic(mim=mim, order=0, expression= >>>>> "Norm_sqr(Div(clambda*Trace(Green_Lagrange_E)*Id(2)+2*cmu*Green_Lagrange_E))", >>>>> model=model,region=-1) >>>>> where Green_Lagrange_E is a macro defined by >>>>> model.add_macro("epsilon", "(Grad_u+(Grad_u)')/2") >>>>> model.add_macro("Green_Lagrange_E", "epsilon+ (Grad_u)'*Grad_u/2") >>>>> but it failed (Error in macro expansion. Only variable name are >>>>> allowed for macro parameter preceded by Grad_ Hess_ Test_ or Test2_ >>>>> prefixes) >>>>> >>>>> Next, i tried to check whether i could at least apply the Div operator >>>>> to Grad_u directly instead, like that >>>>> value=gf.asm_generic(mim=mim, order=0, expression= >>>>> "Norm_sqr(Div(Grad_u))", model=model,region=-1) >>>>> but it also failed, Grad_u was not recognized. >>>>> >>>>> Is there a way to do that? I see that it's possible to have access to >>>>> the hessian of the model variables; is rebuilding div(P) from the hessian >>>>> the only way to compute this quantity? >>>>> Best regards, >>>>> >>>>> David. >>>>> >>>>> >>>>> Le jeu. 26 mai 2022 à 10:17, <yves.ren...@insa-lyon.fr> a écrit : >>>>> >>>>>> Dear David, >>>>>> >>>>>> I do not see any reason why your expression for the Mooney-Rivlin >>>>>> hyper elastic law would give a different result than the brick. At least >>>>>> if >>>>>> of cours you keep the >>>>>> incompressible brick >>>>>> md.add_finite_strain_incompressibility_brick(mim, 'u', 'p') >>>>>> (note that using sqr(expr) instead of pow(expr,2) is slightly more >>>>>> efficient). >>>>>> >>>>>> Concerning the plane strain approximation, it is a priori possible to >>>>>> do this in GWFL, yes. You can define a 3D matrix from your 2D gradient >>>>>> and >>>>>> use the 3D hyperelastic law. >>>>>> >>>>>> For the computation of the Von Mises stress from the potential, it is >>>>>> a little bit more tricky because the automatic differentiation of Getfem >>>>>> will give you the directional derivative (in the direction of the test >>>>>> function). So may be it is possible to extract the components with a >>>>>> specific choice of test functions, but it could be not so easy ... I >>>>>> never >>>>>> done that. Unfortunately it is preferable to give the expression of the >>>>>> law >>>>>> in term of stress tensor directly. >>>>>> >>>>>> Best regards, >>>>>> >>>>>> Yves >>>>>> >>>>>> ------------------------------ >>>>>> *De: *"David Danan" <daviddanan9...@gmail.com> >>>>>> *À: *"getfem-users" <getfem-users@nongnu.org> >>>>>> *Envoyé: *Lundi 23 Mai 2022 13:41:14 >>>>>> *Objet: *[Getfem-users] Questions about user-defined strain energy >>>>>> >>>>>> Dear Getfem-users, >>>>>> >>>>>> to check whether i understand correctly the implementation (and >>>>>> because it's actually much clearer in my code that way), i am trying to >>>>>> replace some predefined bricks for strain energy by the equivalent >>>>>> version >>>>>> using GWFL. >>>>>> For that, i have several questions >>>>>> >>>>>> *First question* >>>>>> For instance >>>>>> lawname = 'SaintVenant Kirchhoff' >>>>>> clambda,cmu = params["clambda"],params["cmu"] >>>>>> model.add_initialized_data('paramsSVK', [clambda, cmu]) >>>>>> idBrick=model.add_finite_strain_elasticity_brick(mim, lawname, 'u', >>>>>> 'paramsSVK') >>>>>> becomes >>>>>> clambda,cmu = params["clambda"],params["cmu"] >>>>>> model.add_initialized_data('paramsSVK', [clambda, cmu]) >>>>>> idBrick=model.add_nonlinear_generic_assembly_brick(mim, >>>>>> "(paramsSVK(1)/2)*pow(Trace(Green_Lagrange_E),2)+paramsSVK(2)*Trace(Green_Lagrange_E'*Green_Lagrange_E)" >>>>>> ) >>>>>> >>>>>> And the results seems to be the same. >>>>>> However, for the incompressible Mooney Rivlin strain energy case, it >>>>>> does not behave as i expected >>>>>> >>>>>> Using the example there as a basis in a 3D case >>>>>> >>>>>> https://github.com/getfem-doc/getfem/blob/master/interface/tests/python/demo_nonlinear_elasticity.py >>>>>> and after some simplifications >>>>>> >>>>>> I have tried to replace the last line in >>>>>> lawname = 'Incompressible Mooney Rivlin' >>>>>> model.add_initialized_data('paramsIMR', [1,1]) >>>>>> model.add_finite_strain_elasticity_brick(mim, lawname, 'u', >>>>>> 'paramsIMR') >>>>>> By this >>>>>> model.add_nonlinear_generic_assembly_brick(mim, >>>>>> "Incompressible_Mooney_Rivlin_potential(Grad_u, >>>>>> [paramsIMR(1);paramsIMR(2)])") >>>>>> Which was exactly the same, of course. Next, i have tried to replace >>>>>> it by these lines >>>>>> model.add_macro("F", "Grad_u+Id(3)") >>>>>> model.add_nonlinear_generic_assembly_brick(mim, "paramsIMR(1)* ( >>>>>> Matrix_j1(Right_Cauchy_Green(F)) - 3 )+ paramsIMR(2)* ( >>>>>> Matrix_j2(Right_Cauchy_Green(F)) - 3 )") >>>>>> Which failed to converge >>>>>> >>>>>> I thought this expression was consistent with the implementation of >>>>>> Mooney_Rivlin_hyperelastic_law::strain_energy and the compute_invariants >>>>>> here >>>>>> >>>>>> http://download-mirror.savannah.gnu.org/releases/getfem/doc/getfem_reference/getfem__nonlinear__elasticity_8cc_source.html >>>>>> and the documentation >>>>>> >>>>>> https://getfem.org/userdoc/model_nonlinear_elasticity.html?highlight=finite%20strain >>>>>> But, clearly, i am missing something. Could you explain what I am >>>>>> doing wrong? >>>>>> >>>>>> Second question >>>>>> In a 2D case, i would like to be able to use either a plane strain >>>>>> approximation based on a given strain energy expression in 3D. >>>>>> In the implementation, it is nicely done in scalar_type >>>>>> plane_strain_hyperelastic_law::strain_energy there >>>>>> >>>>>> http://download-mirror.savannah.gnu.org/releases/getfem/doc/getfem_reference/getfem__nonlinear__elasticity_8cc_source.html >>>>>> >>>>>> Is it possible to do the same using the GWFL? >>>>>> >>>>>> *Third question* >>>>>> I would like to be able to compute the von-mises field for any strain >>>>>> energy function. >>>>>> >>>>>> If It is an existing brick, the method >>>>>> md.compute_finite_strain_elasticity_Von_Mises(lawname, 'u', 'params', >>>>>> mfdu) >>>>>> will do the trick just fine. >>>>>> >>>>>> If it's not the case, i guess i can use something akin to the actual >>>>>> implementation of compute_finite_strain_elasticity_Von_Mises >>>>>> <http://download-mirror.savannah.gnu.org/releases/getfem/doc/getfem_reference/namespacegetfem.html#a9166bc2065a4f2b5633e2e933ec3e693> >>>>>> >>>>>> std::string expr = >>>>>> "sqrt(3/2)*Norm(Deviator(Cauchy_stress_from_PK2(" >>>>>> + adapted_lawname + "_PK2(Grad_" + varname + "," + params + >>>>>> "),Grad_" >>>>>> + varname + ")))"; >>>>>> ga_interpolation_Lagrange_fem(md, expr, mf_vm, VM, rg); >>>>>> } >>>>>> >>>>>> combined with local_projection to get the value at the elements. >>>>>> The question is: is it possible to compute Piola Kirchhoff 2 from the >>>>>> strain energy within the GWFL expression given to local_projection? >>>>>> I have the impression it's the only thing i need to be able to do >>>>>> what i want. >>>>>> Let W be a strain energy function and E be the Green-Lagrange tensor >>>>>> (defined as macros, let's say), is Diff(W, E) doing exactly what i >>>>>> am expecting for this purpose? >>>>>> If it's not the case, is there another way to do it? >>>>>> >>>>>> Thank you in advance, >>>>>> kind regards, >>>>>> David. >>>>>> >>>>>>
