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. >>>> >>>>