I encountered the error message below. I am wondering is there a problem with the dimension of dx and 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
raise self._exception_type(self._format_raw(*message))
ufl.log.UFLException: Can't add expressions with different shapes.

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])],
[32, 32],
CellType.triangle, GhostMode.none)
# Function to mark x = 0, x = 1 and y = 0
def noslip_boundary(x):
return np.logical_or(np.logical_or(np.isclose(x[0], 0.0),
np.isclose(x[0], 1.0)),
np.isclose(x[1], 0.0))
# Function to mark the lid (y = 1)
def lid(x):
return np.isclose(x[1], 1.0)
# Lid velocity
def lid_velocity_expression(x):
return np.stack((np.ones(x.shape[1]), np.zeros(x.shape[1])))
# 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)
lid_velocity.interpolate(lid_velocity_expression)
facets = locate_entities_boundary(msh, 1, lid)
bc1 = dirichletbc(lid_velocity, locate_dofs_topological(V, 1, facets))
# Collect Dirichlet boundary conditions
bcs = [bc0, bc1]