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 ! :)
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
Former request answer about sliding conditions with dirichlet bricks
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.
[cid:[email protected]]
[cid:[email protected]]
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 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');
gf_model_set(md, 'add initialized data', 'VolumicData', F);
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
[cid:[email protected]]
SIMON AMEYE
DQI/DRIA
Apprenti IFP School
Private mail : [email protected]<mailto:[email protected]>
________________________________
CENTRE TECHNIQUE VELIZY A /
