|
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 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.
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
|
[Prev in Thread] | Current Thread | [Next in Thread] |