Dear Kosta, I will try the projection method you mentioned in your previous mail. I am using 'add_integral_large_sliding_contact_brick_raytracing' and then I am adding slave and master brick with boundary conditions for my contact problem. Actually I need to calculate pressure everywhere and not just the contact region. I just tried to calculate in on the contact surface. I will try the projection method you mentioned and will let you know how it goes.
Thanks for your help. Yours sincerely Samyak On Thu, Jan 26, 2017 at 7:21 PM, Konstantinos Poulios < [email protected]> wrote: > Dear Samyak > > The projection method of my previous mail works also for large strains if > you simply replace the expression for the stress tensor with the > corresponding one. Of course, as you suggested you can also do an > interpolation instead of projection, but then you depend on the averaging > that the interpolation function of GetFEM++ will do behind the scene for > the discontinuous field that you are working with. I have never tried it on > a surface, so I cannot say how it works. > > But if you care only about the contact interface, why don't you simply use > the Lagrange multiplier which already describes the traction on this > surface? Which contact algorithm do you use? > > Best regards > Kosta > > > > On Thu, Jan 26, 2017 at 10:47 AM, samyak jain <[email protected]> > wrote: > >> Dear Kostas, >> >> Thanks for the reply. I am using large strain law and using python >> language. >> >> I have a 3D-mesh and mfu1 is the displacement fem for mesh1. >> CONTACT_BOUND is the region where I am trying to calculate the stress >> tensor. >> I am trying something of this sort below. Could you please check and tell >> me if what I am doing is correct or not. I reshape the stress tensor into >> 9*(no. of dof /3) as 3-dimension tensor will have 9 components and md is >> model object. >> >> stress = np.reshape(md.interpolation("lambda*Trace(Grad_u1 + Grad_u1')*Id(3) >> + mu*(Grad_u1 + Grad_u1')", mfu1, CONTACT_BOUND),(9,mfu1.nbdof()/3),'F') >> >> >> Thanking You >> >> Yours sincerely >> Samyak Jain >> >> On Thu, Jan 26, 2017 at 5:07 PM, Konstantinos Poulios < >> [email protected]> wrote: >> >>> Dear Samyak >>> >>> It would be helpful to know which material law you use (small or large >>> strains?) and which programming language (python, matlab, scilab, C++?). >>> >>> E.g for small strains elasticity in python, something like this should >>> work: >>> >>> mass_mat_t = gf.asm_mass_matrix(mim, mft) >>> for i,j in [[1,1],[2,2],[1,2]]: >>> sigmaij = gf.asm_generic(mim, 1, "({Lambda}*Div_u*Id(qdim(u)) + >>> {Mu}*(Grad_u+Grad_u'))({i},{j})*Test_t".format(Lambda=Lambda, Mu=Mu, >>> i=i, j=j), -1, >>> "u", False, mfu, md.variable("u"), >>> "t", True, mft, np.zeros(mft.nbdof())) >>> sigmaij = np.transpose(gf.linsolve_mumps(mass_mat_t, sigma)) >>> >>> Then sigmaij has the ij component of the stress tensor projected on the >>> finite element space defined in mft. If mft is a Lagrangian FEM the values >>> of sigmaij will correspond to the nodes of mft. The projection from gauss >>> points to the finite element space is as Yves mentioned some kind of >>> averaging. >>> >>> I hope this helps. >>> >>> Best regards >>> Kostas >>> >>> >>> >>> On Thu, Jan 26, 2017 at 2:48 AM, samyak jain <[email protected] >>> > wrote: >>> >>>> Dear Yves, >>>> >>>> Thanks a lot for your reply. I am not sure what you mean there. I have >>>> gone through the example and they calculate Von Mises Stress but there is >>>> already a function for that. Could you please tell me which physical >>>> quantity or expression is being interpolated in that example so that I can >>>> use it as a basis for my pressure calculation. It would be really helpful >>>> if you would elaborate a bit on the usage and the example for calculating >>>> any expression as I am quite new to fem and getfem. Thanks a lot. >>>> >>>> Yours sincerely >>>> Samyak >>>> >>>> On Wed, Jan 25, 2017 at 6:44 PM, Yves Renard <[email protected]> >>>> wrote: >>>> >>>>> >>>>> Dear Samyak, >>>>> >>>>> There is no specific function for that, but you can use the >>>>> interpolation facilities of the model object. You can see an example of >>>>> use >>>>> for instance in the demo_thermo_elasticity_electrical_coupling >>>>> example which has a python, Matlab, Scilab and C++ version. You can >>>>> interpolate any expression. Be aware that the stress being discontinuous >>>>> across the element edges/faces, the value will be averaged. >>>>> >>>>> Regards, >>>>> >>>>> Yves >>>>> >>>>> >>>>> Le 25/01/2017 à 06:08, samyak jain a écrit : >>>>> >>>>> Dear getfem-users, >>>>> >>>>> I have started using getfem recently to work on few contact problems >>>>> between elastic bodies. >>>>> >>>>> I am able to solve few simple cases and I am able to compute >>>>> displacement (u) and Von-Mises/Tresca values on my mesh nodes. >>>>> >>>>> The problem I am facing is that I have to calculate pressure at the >>>>> mesh points and for that I believe I need the stress tensor values at >>>>> these >>>>> discrete mesh points. Could someone please guide me or tell me how to get >>>>> the value of these variables. I went through the documentation and checked >>>>> everywhere but I couldn't find a way to calculate stress tensor/pressure. >>>>> >>>>> Thanks a lot. >>>>> >>>>> Samyak >>>>> >>>>> >>>>> >>>>> _______________________________________________ >>>>> Getfem-users mailing >>>>> [email protected]https://mail.gna.org/listinfo/getfem-users >>>>> >>>>> >>>>> >>>>> -- >>>>> >>>>> 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 >>>>> >>>>> --------- >>>>> >>>>> >>>>> _______________________________________________ >>>>> Getfem-users mailing list >>>>> [email protected] >>>>> https://mail.gna.org/listinfo/getfem-users >>>>> >>>>> >>>> >>>> _______________________________________________ >>>> Getfem-users mailing list >>>> [email protected] >>>> https://mail.gna.org/listinfo/getfem-users >>>> >>>> >>> >> >
_______________________________________________ Getfem-users mailing list [email protected] https://mail.gna.org/listinfo/getfem-users
