I’m really sorry, through yesterday’s questioning I have already managed to apply fixed displacement conditions on one vertex of a built-in grid:https://fenicsproject.discourse.group/t/in-dolfinx-how-to-apply-a-fixed-displacement-condition-to-a-vertex/14886/2?u=french_fries
The reason why I used a built-in grid in the previous post was to present the problem more simply, but when I applied the solution to an externally imported grid, it didn’t seem to work, apart from the grid I didn’t make any other modifications, here is my code (applying fixed displacement conditions on the center of a circle located at the origin):
import numpy as np
from ufl import (tr, grad, Identity, TrialFunction, TestFunction, lhs, rhs, inner, dx)
from dolfinx import (mesh, io, default_scalar_type, fem)
from dolfinx.fem import functionspace, dirichletbc
from dolfinx.fem.petsc import LinearProblem
from mpi4py import MPI
from pathlib import Path
from dolfinx.io import gmshio
# Create mesh
domain, cell_markers, facet_markers = gmshio.read_from_msh("test.msh", MPI.COMM_WORLD, gdim=2)
def center(x):
return np.isclose(x[0], 0) & np.isclose(x[1], 0)
E = 50e3
nu = 0.2
mu = E/2/(1+nu)
lmbda = E*nu/(1+nu)/(1-2*nu)
alpha = 1e-5
def eps(v):
return 0.5 * (grad(v) + grad(v).T)
def sigma(v, dT):
return (lmbda*tr(eps(v))- alpha*(3*lmbda+2*mu)*dT)*Identity(2) + 2.0*mu*eps(v)
VU = functionspace(domain, ("CG", 1, (domain.geometry.dim, )))
dU = TrialFunction(VU)
U_ = TestFunction(VU)
Wint = inner(sigma(dU, 1000), eps(U_)) * dx
aM = lhs(Wint)
LM = rhs(Wint)
fdim = 0
boundary_vertices = mesh.locate_entities_boundary(domain, fdim, center)
print(boundary_vertices)
u_D = np.array([0, 0], dtype = default_scalar_type)
bc = dirichletbc(u_D, fem.locate_dofs_topological(VU, fdim, boundary_vertices), VU)
# Compute solution
problem = LinearProblem(aM, LM, bcs=[bc], petsc_options={"ksp_type": "preonly", "pc_type": "lu"})
Uf = problem.solve()
This is the geo file (a circular domain divided into four parts, ensuring the center falls on a vertex):
// Gmsh project created on Tue May 21 09:56:22 2024
SetFactory("OpenCASCADE");
r = 7.81;
Point(1) = {0, 0, 0, 1.0};
Point(2) = {0, r, 0, 1.0};
Point(3) = {r, 0, 0, 1.0};
Point(4) = {0, -r, 0, 1.0};
Point(5) = {-r, 0, 0, 1.0};
Circle(1) = {2, 1, 3};
Circle(2) = {3, 1, 4};
Circle(3) = {4, 1, 5};
Circle(4) = {5, 1, 2};
Line(5) = {1, 2};
Line(6) = {1, 3};
Line(7) = {1, 4};
Line(8) = {1, 5};
Curve Loop(1) = {5, 1, -6};
Plane Surface(1) = {1};
Curve Loop(2) = {6, 2, -7};
Plane Surface(2) = {2};
Curve Loop(3) = {7, 3, -8};
Plane Surface(3) = {3};
Curve Loop(4) = {8, 4, -5};
Plane Surface(4) = {4};
Physical Curve("surface", 9) = {4, 1, 2, 3};
Physical Surface("circle", 10) = {4, 1, 2, 3};
This is the output (boundary_vertices is an empty list) and the visualization:
Info : Reading 'test.msh'...
Info : 289 nodes
Info : 576 elements
Info : Done reading 'test.msh'
[]
