getfem-users
[Top][All Lists]

## [Getfem-users] Query for contact bricks use example [trying again since

 From: SIMON AMEYE Subject: [Getfem-users] Query for contact bricks use example [trying again since I have a GetFem account] Date: Mon, 22 Jan 2018 09:46:32 +0000

Hi all,

I have a new request about contact bricks.

(At the end of the mail, there is also the answer to a problem I faced using sliding conditions with dirichlet bricks.)

I need your help again ! J

I have a lot of pain using the nonmatching meshes contact brick.

Do someone have any example using it ? I use Matlab, but Scilab or Python are easy to translate.

I have two meshes mesh_1 and mesh_2 that I need to assemble with contact (and no friction).

For now the answer I get is :

“SuperLU solve failed: info=3041”

Here is the code I use with getfem 2010:

% Numerical parameters

E = 21e6; Nu = 0.3;

lambda = E*Nu/((1+Nu)*(1-2*Nu));

mu = E/(2*(1+Nu));

mesh_1= gfMesh('from string', StringMesh);

mesh_2 = gfMesh('from string', StringMesh2);

mim=gfMeshIm(mesh_1);  set(mim, 'integ',gfInteg('IM_TRIANGLE(6)'));

mfu=gfMeshFem(mesh_1,2); set(mfu, 'fem',gfFem('FEM_PK(2,1)'));

mfd=gfMeshFem(mesh_1); set(mfd, 'fem',gfFem('FEM_PK(2,1)'));

mf0=gfMeshFem(mesh_1); set(mf0, 'fem',gfFem('FEM_PK(2,0)'));

mfdu=gfMeshFem(mesh_1); set(mfdu, 'fem',gfFem('FEM_PK_DISCONTINUOUS(2,1)'));

mim2=gfMeshIm(mesh_2);  set(mim2, 'integ',gfInteg('IM_TRIANGLE(6)'));

mfu2=gfMeshFem(mesh_2,2); set(mfu2, 'fem',gfFem('FEM_PK(2,1)'));

mfd2=gfMeshFem(mesh_2); set(mfd2, 'fem',gfFem('FEM_PK(2,1)'));

mf02=gfMeshFem(mesh_2); set(mf02, 'fem',gfFem('FEM_PK(2,0)'));

mfdu2=gfMeshFem(mesh_2); set(mfdu2, 'fem',gfFem('FEM_PK_DISCONTINUOUS(2,1)'));

% Boundary calculation

P=get(mesh_1, 'pts');

pidlow=find(abs(P(2,:))<1e-6);

flow =get(mesh_1,'faces from pid',pidlow);

ftop=get(mesh_1,'faces from pid',(TopEdgeNodes'+1));

fall = [1437 403 87;3 1 3];

magn_border = gf_mesh_get(mesh_2, 'outer faces');

mesh_2.set_region(10, magn_border);

% assign boundary numbers

LOW_BOUND = 1; TOP_BOUND = 2; ALL_BOUND = 3;

mesh_1.set_region(LOW_BOUND, flow);

mesh_1.set_region(TOP_BOUND, ftop);

mesh_1.set_region(ALL_BOUND, fall);

for i = 1: TriNumb

ftri_rand = get(mesh_1,'faces from cvid',(i+1));

mesh_1.set_region(20+i, ftri_rand);

end

%% MODEL

md=gf_model('real');

gf_model_set(md, 'add fem variable', 'u', mfu);

gf_model_set(md, 'add initialized data', 'lambda', lambda);

gf_model_set(md, 'add initialized data', 'mu', mu);

gf_model_set(md, 'add isotropic linearized elasticity brick', mim, 'u', 'lambda', 'mu');

% md=gf_model('real');

gf_model_set(md, 'add fem variable', 'u2', mfu2);

gf_model_set(md, 'add initialized data', 'lambda2', lambda);

gf_model_set(md, 'add initialized data', 'mu2', mu);

gf_model_set(md, 'add isotropic linearized elasticity brick', mim2, 'u2', 'lambda2', 'mu2');

for i = 1:TriNumb

gf_model_set(md, 'add initialized data', ['VolumicData' num2str(i)], [Fx_tri(i),Fy_tri(i)]/3);

gf_model_set(md, 'add source term brick', mim, 'u', ['VolumicData' num2str(i)],20+i);

end

% H matrix building : H is a projector (Hat) mathix : H^2=k*h and must be

% symetrical

x4 = 1;x1 = x4/(pi/2/TopAngle);k = x4+x1;x2 = sqrt(x1)*(-sqrt(k-x1));x3 = x2;

gf_model_set(md, 'add initialized data', 'H_LOW', [0 0;0 1]);

gf_model_set(md, 'add initialized data', 'VECTOR_LOW', [0;0]);

gf_model_set(md, 'add initialized data', 'H_TOP', [x1 x2;x3 x4]);

gf_model_set(md, 'add initialized data', 'VECTOR_TOP', [0;0]);

gf_model_set(md, 'add generalized Dirichlet condition with multipliers', mim, 'u', mfu, TOP_BOUND,'VECTOR_TOP', 'H_TOP');

gf_model_set(md, 'add generalized Dirichlet condition with multipliers', mim, 'u', mfu, LOW_BOUND,'VECTOR_LOW', 'H_LOW');

gf_model_set(md, 'add initialized data', 'friction', 155);

gf_model_set(md,'add nonmatching meshes contact brick',mim, mim2,'u','u2','frict','friction',3,10);

%% Solver

gf_model_get(md, 'solve');

U = gf_model_get(md, 'variable', 'u');

VM = gf_model_get(md, 'compute isotropic linearized Von Mises or Tresca', 'u', 'lambda', 'mu', mfdu);

Thank you for your help and for that great work,

Best regards,

Simon

For my last mail about the Free slipping boundary in inclined direction, it finally make it!

Here is the solution I found. I am sure it can help someone. You need to set a matrix for each dirichlet condition, that allow the displacement in a specific dimension.

This matrix H needs (It does no work for me if I don’t do so) to be symmetrical and : H^2 = k*H.

x4 = 1;x1 = x4/(pi/2/TopAngle);k = x4+x1;x2 = sqrt(x1)*(-sqrt(k-x1));x3 = x2;

gf_model_set(md, 'add initialized data', 'H_LOW', [0 0;0 1]);

gf_model_set(md, 'add initialized data', 'VECTOR_LOW', [0;0]);

gf_model_set(md, 'add initialized data', 'H_TOP', [x1 x2;x3 x4]);

gf_model_set(md, 'add initialized data', 'VECTOR_TOP', [0;0]);

gf_model_set(md, 'add generalized Dirichlet condition with multipliers', mim, 'u', mfu, TOP_BOUND,'VECTOR_TOP', 'H_TOP');

gf_model_set(md, 'add generalized Dirichlet condition with multipliers', mim, 'u', mfu, LOW_BOUND,'VECTOR_LOW', 'H_LOW');

Former request :

Dear all getfem users,

First, let me thank you for this great work.

For my apprenticeship, I am currently using getfem to calculate Von Mises constraint in some 2D parts, as it in the fastest way to do so with Matlab.

I need to calculate the constraint in a part with one distributes load F.

 Figure 1 : Dirichlet condition (anchored) Figure 2 : Sliding condition

For now, I am using the dirichlet condition brick to anchor the boundaries (grey), and it works fine (Figure 1).

%Model

md=gf_model('real');

gf_model_set(md, 'add isotropic linearized elasticity brick', mim, 'u''lambda''mu');

gf_model_set(md, 'add source term brick', mim, 'u''VolumicData');

%Boundary Conditions

gf_model_set(md, 'add Dirichlet condition with multipliers', mim, 'u', mfu, LOW_BOUND);

gf_model_set(md, 'add Dirichlet condition with multipliers', mim, 'u', mfu, TOP_BOUND);

%Solver

gf_model_get(md, 'solve');

U = gf_model_get(md, 'variable''u');

VM = gf_model_get(md, 'compute isotropic linearized Von Mises or Tresca''u''lambda''mu', mfdu);

But now, I would like the boundaries to slide in only one direction that is, in one case, not vertical (Figure 2).

I tried the Neumann source term and normal dirichlet conditions but I couldn’t deal with it.

Am I wrong? Do you know how to do it with getfem?

In the hope that it could help some mechanical getfem users,

Best regards,

Simon

 SIMON AMEYE DQI/DRIA Apprenti IFP School Private mail : address@hidden CENTRE TECHNIQUE VELIZY A /