Boundary condition in uniaxial compression case

I’m trying to simulate uniaxial compression of a 3D cylindrical specimen. the mesh is imported with markers. Loading is in the -z direction and the deformation is expected in the x and y directions for both the top and bottom surfaces (kind of assuming zero friction). However, I get results showing a bias deformation in y direction. Attached is a part of code where boundary conditions are applied. I believe that I’m making a mistake in bcs. Looking forward to some advise. Thank you in advance.

dofs_fixed = fem.locate_dofs_topological(V.sub(2), mesh.topology.dim - 1, facet_tags.indices[facet_tags.values == 7])
dofs_disp = fem.locate_dofs_topological(V.sub(2), mesh.topology.dim - 1, facet_tags.indices[facet_tags.values == 8])
dofs_sym = fem.locate_dofs_topological(V.sub(0), mesh.topology.dim - 1, facet_tags.indices[facet_tags.values == 9])


u_fixed = fem.Constant(mesh, default_scalar_type(0.0))
u_disp = fem.Constant(mesh, default_scalar_type(0.0))
u_sym = fem.Constant(mesh, default_scalar_type(0.0))

#initial bcs

bsym = fem.dirichletbc(u_sym, dofs_sym, V.sub(0))
bcf = fem.dirichletbc(u_fixed, dofs_fixed, V.sub(2))
bcp = fem.dirichletbc(u_disp, dofs_disp, V.sub(2))

bcs = [bcf, bsym, bcp]

Your lines of code look Ok to me, although it is impossible to say anything definitive without the context in which they are written.

Of course, with these boundary conditions, a rigid-body mode (dispacement in y) remains unconstrained. It could simply be that that mode is somehow triggered in the iterative solve. An easy fix is to constrain V.sub(1) in a single point.

Something like:

def node_y0(self, x: NDArray[np.float64]) -> NDArray[np.bool_]:
   return np.isclose(x[0], 0) & np.isclose(x[1], 0) & np.isclose(x[2],0)
dof_y0 = fem.locate_dofs_topological(V.sub(1), mesh.topology.dim - 2, node_y0)
u_y0 = fem.Constant(mesh, default_scalar_type(0.0))
bcy0 = fem.dirichletbc(u_y0 , dof_y_zero, V.sub(1))

I ran into more errors with the suggested modification of the code. Specifically,

’ ’ ’
TypeError: int() argument must be a string, a bytes-like object or a real number, not 'function
’ ’ ’
So, I opted to model a quarter of a cylinder and constrain the x and y directions. It works as expected. Thanks.

I do some similar simulation as you.

u1_bc = fem.Function(V)

u1 = lambda x: np.zeros_like(x, dtype=PETSc.ScalarType)

u1_bc.interpolate(u1)

left_dofs = fem.locate_dofs_topological((state_space.sub(0), V), facet_tag.dim, facet_tag.find(3))

bc_left = fem.dirichletbc(u1_bc, left_dofs, state_space.sub(0))

bcs = [bc_left]

I did what I found in others’ posts and worked for my uniaxial compression.