Dear Yves, Thank you very much!! It is working now.
Best regards, Phuoc On 14/07/2018 15:31, Yves Renard wrote:
Dear Phuoc, I think the correct expression with Contract is simply div_sigma = "(Contract(Grad(first_Piola_Kirchhoff_stress), 2, 3))" Contract makes a summation with respect to the given indices, so that Contract(A, 2, 3) = A_{i,j,j} with summation on j. Can your try again with the fix I commited to the master branch ? Best regards, Yves ----- Original Message ----- From: "Huu Phuoc BUI" <[email protected]> To: "yves renard" <[email protected]> Cc: "franz chouly" <[email protected]> Sent: Friday, 13 July, 2018 12:00:55 Subject: Re: [Getfem-users] Divergence of the first Piola-Kirchhoff stress tensor Dear Yves, Thank you very much for your email. Please find enclosed a minimal code. I am not sure that I use correctly the Contract operator, you will find therefore both the original code and the new one for div_sigma: /div_sigma = "[Contract(Grad(first_Piola_Kirchhoff_stress)(1,:,:),1,2), Contract(Grad(first_Piola_Kirchhoff_stress)(2,:,:),1,2), Contract(Grad(first_Piola_Kirchhoff_stress)(3,:,:),1,2)].[1,1,1]"/ Please let me know if you find the problem. Thank you! Best regards, Phuoc On 13/07/2018 08:59, Yves Renard wrote:Dear Phuoc, No, a priori this should work. To simplify a bit your expression, you can use Contract(Grad(...), 2, 3).[1,1,1] instead of the three traces in div_sigma. If you give me a minimal code, I could have a look to the problem. Best regards, Yves Le 12/07/2018 à 15:36, Huu Phuoc BUI a écrit :Dear GetFEM-users, As Grad(expression) is available in the new version of GetFEM++, I compute the divergence of the first Piola-Kirchhoff stress tensor as (python code): /div_sigma = "( Trace(Grad(first_Piola_Kirchhoff_stress)(1,:,:)) + Trace(Grad(first_Piola_Kirchhoff_stress)(2,:,:)) + Trace(Grad(first_Piola_Kirchhoff_stress)(3,:,:)) )"/ in which first_Piola_Kirchhoff_stress is a macro of a model /md/, which is defined as: /md.add_macro("second_Piola_Kirchhoff_stress", "Compressible_Mooney_Rivlin_sigma(Grad_u, [c1;c2;d1])")// //md.add_macro("first_Piola_Kirchhoff_stress", "(Id(meshdim)+Grad_u)*second_Piola_Kirchhoff_stress")/ Then I use generic assembly to assembly this quantity as: /eta_L_1_tmp = gf.asm_generic(mim,1,"(("+div_sigma+").delta)*Test_t",-1// // ,md// // ,"delta", False, mfz, delta// // ,"t", True,mfrj,np.zeros(mfrj.nbdof()) )/ with /delta/ a given vector on a meshfem /mfz/. When executing the code, GetFEM++ only gives the error as /segmentation fault (core dumped)/ If you see any problem in my code, please let me know. Thank you in advance. Best regards, Phuoc-- Yves Renard ([email protected]) tel : (33) 04.72.43.87.08 Pole de Mathematiques, INSA-Lyon fax : (33) 04.72.43.85.29 20, rue Albert Einstein 69621 Villeurbanne Cedex, FRANCE http://math.univ-lyon1.fr/~renard ---------
