I encountered the error message below. I am wondering is there a problem with the dimension of
ds not being the same. I would greatly appreciate any suggestions or guidance on how to resolve this issue.
Can't add expressions with different shapes.
ERROR:UFL:Can't add expressions with different shapes.
Traceback (most recent call last):
File "/Users/victoriachan/Documents/CityU Doc/FYP/python_files/ProgramV2.py", line 61, in <module>
+ ufl.inner(w_h - ufl.grad(u_h), ufl.grad(kappa * ufl.nabla_div(p))) * ds \
File "/Users/victoriachan/opt/anaconda3/envs/fenicsdolxconda/lib/python3.10/site-packages/ufl/exproperators.py", line 223, in _sub
return Sum(self, -o)
File "/Users/victoriachan/opt/anaconda3/envs/fenicsdolxconda/lib/python3.10/site-packages/ufl/algebra.py", line 41, in __new__
error("Can't add expressions with different shapes.")
File "/Users/victoriachan/opt/anaconda3/envs/fenicsdolxconda/lib/python3.10/site-packages/ufl/log.py", line 135, in error
ufl.log.UFLException: Can't add expressions with different shapes.
The source code (without the import packages):
If you want to receive help, you need to supply a minimal code (in text form) that reproduces the error.
The code should be wrapped as
#add rest of code here
General comment to your code. If you want to solve your problem for both V and W, you need to use a MixedElement. See:
or any of the approaches further up if you want a more technical (But better perfoming) code.
Thanks. I am reading this example. Could I ask how can I change the boundary to a unit square boundary for the biharmonic equation? As I don’t need the velocity term now and I don’t know how to revise it. Thanks for your kind help.
The original boundary code from the link:
# Create mesh
msh = create_rectangle(MPI.COMM_WORLD,
[np.array([0, 0]), np.array([1, 1])],
# Function to mark x = 0, x = 1 and y = 0
return np.logical_or(np.logical_or(np.isclose(x, 0.0),
# Function to mark the lid (y = 1)
return np.isclose(x, 1.0)
# Lid velocity
return np.stack((np.ones(x.shape), np.zeros(x.shape)))
# No-slip boundary condition for velocity field (`V`) on boundaries
# where x = 0, x = 1, and y = 0
noslip = np.zeros(msh.geometry.dim, dtype=PETSc.ScalarType)
facets = locate_entities_boundary(msh, 1, noslip_boundary)
bc0 = dirichletbc(noslip, locate_dofs_topological(V, 1, facets), V)
# Driving velocity condition u = (1, 0) on top boundary (y = 1)
lid_velocity = Function(V)
facets = locate_entities_boundary(msh, 1, lid)
bc1 = dirichletbc(lid_velocity, locate_dofs_topological(V, 1, facets))
# Collect Dirichlet boundary conditions
bcs = [bc0, bc1]
I am not sure what you are asking for here.