Hello, trying to solve linear elasticity equation, but it looks there is a problem with the sparsity pattern (from what could find, but i don’t understand why). Can some one explain my error?
MWE:
import numpy as np
from dolfinx import fem, cpp, io, mesh
from mpi4py import MPI
from petsc4py import PETSc
import ufl
YOUNG=1.
NU=0.4
G=9.81
MASS=3.5
WEIGHT=G*MASS
L, W, H=0.303, 0.19, 0.064
EMP=0.06
AXE=0.025
NX, NY, NZ = 30, 19, 6
points=[np.array([0, 0, 0]), np.array([L, W, H])]
domain=mesh.create_box(
MPI.COMM_WORLD,
points,
[NX,NY,NZ],
cell_type=cpp.mesh.CellType.hexahedron,
ghost_mode=mesh.GhostMode.shared_facet
)
dim=domain.topology.dim
def axes(x):
return np.logical_and(
np.logical_or(np.sqrt( (x[0]-EMP)**2 + (x[2]-H/2)**2 ) <=AXE,
np.sqrt( (x[0]-L+EMP)**2 + (x[2]-H/2)**2 ) <=AXE),
np.logical_or(np.isclose(x[1], 0, atol=1e-3),
np.isclose(x[1], W, atol=1e-3))
)
def bottom(x):
return np.isclose(x[2], 0, atol=1e-3)
bottom_facets=mesh.locate_entities(domain, dim-1, bottom)
facet_indices=[bottom_facets]
facet_indices=np.array(np.hstack(facet_indices), dtype=np.int32)
facet_markers=[np.full(len(bottom_facets), 1)]
facet_markers=np.array(np.hstack(facet_markers), dtype=np.int32)
sorted_facets=np.argsort(facet_indices)
facet_tag=mesh.meshtags(
domain,
dim-1,
facet_indices[sorted_facets],
facet_markers[sorted_facets]
)
U=fem.VectorFunctionSpace(domain, ("CG", 1))
u=fem.Function(U, name="Déplacement")
LAMBDA=YOUNG*NU/(1+NU)/(1-2*NU)
MU=YOUNG/(2*1+NU)
bcdofs=fem.locate_dofs_geometrical(U, bottom)
u_0=np.array((0,)*dim, dtype=PETSc.ScalarType)
bc=fem.dirichletbc(u_0, bcdofs, U)
ds= ufl.Measure("ds", domain=domain, subdomain_data=facet_tag)
BOTTOM_SURF=fem.assemble_scalar(fem.form(fem.Constant(domain, 1.)*ds(1)))
def eps(v):
return ufl.sym(ufl.grad(v))
def sigma(v):
return (LAMBDA*ufl.nabla_div(v)*ufl.Identity(dim)+2*MU*eps(v))
def dE(u, v):
return ufl.inner(sigma(u), eps(v))
u_=ufl.TestFunction(U)
du=ufl.TrialFunction(U)
a=fem.form(ufl.inner(sigma(u_), eps(du))*ufl.dx)
b=fem.form(ufl.dot((WEIGHT/BOTTOM_SURF)*ufl.FacetNormal(domain), u_)*ds(1))
A=fem.assemble_matrix(a, bcs=[bc])
A.assemble()
The error message:
Traceback (most recent call last):
File "/root/opti copy.py", line 81, in <module>
A=fem.assemble_matrix(a, bcs=[bc])
File "/usr/lib/python3.10/functools.py", line 889, in wrapper
return dispatch(args[0].__class__)(*args, **kw)
File "/usr/local/dolfinx-real/lib/python3.10/dist-packages/dolfinx/fem/assemble.py", line 273, in _assemble_matrix_form
A: la.MatrixCSRMetaClass = create_matrix(a)
File "/usr/local/dolfinx-real/lib/python3.10/dist-packages/dolfinx/fem/assemble.py", line 97, in create_matrix
return la.matrix_csr(sp, dtype=a.dtype)
File "/usr/local/dolfinx-real/lib/python3.10/dist-packages/dolfinx/la.py", line 81, in matrix_csr
return matrixcls(sp)
File "/usr/local/dolfinx-real/lib/python3.10/dist-packages/dolfinx/la.py", line 54, in __init__
super().__init__(sp)
RuntimeError: Block size not yet supported