### TypeError: unsupported operand type(s) for +: 'Sum' and 'Form' Solvi Navier Stokes equation withixed finite elements

82

views

0

Hi,

I am trying to solve a Navier Stokes PDE problem using FEniCS and mixed elements, but I am facing an error that I do not understand:

TypeError: unsupported operand type(s) for +: 'Sum' and 'Form'

Below is presented my complete code. I hope that someone will be able to help me.

Thank you in advance.

I am trying to solve a Navier Stokes PDE problem using FEniCS and mixed elements, but I am facing an error that I do not understand:

TypeError: unsupported operand type(s) for +: 'Sum' and 'Form'

Below is presented my complete code. I hope that someone will be able to help me.

Thank you in advance.

```
from __future__ import print_function
from fenics import *
from dolfin import *
# Define constants and known data
f = Constant((0.0, 0.0))
mu = Constant(0.005)
rho = Constant(1.0)
# Time discretization
T = 1.0 # final time
num_steps = 1000 # number of time steps
dt = T / num_steps # time step size
# Create mesh and define function space
N = 100 # Characteristic number of elements
mesh = RectangleMesh(Point(0,0),Point(1,1),N,N)
# Function space
V = VectorElement('P', mesh.ufl_cell(), 2) # velocity
Q = FiniteElement('P', mesh.ufl_cell(), 1) # pressure
element = MixedElement([V, Q])
W = FunctionSpace(mesh, element)
# Define boundaries
lid = 'near(x[1], 1.0)'
noslip = 'near(x[1], 0.0) || near(x[0], 0.0) || near(x[0], 1.0)'
pressure = 'near(x[1], 0.0)'
# Boundary conditions
bc_noslip = DirichletBC(W.sub(0), Constant((0.0, 0.0)), noslip)
bc_lid = DirichletBC(W.sub(0), Constant((1.0, 0.0)), lid)
bc_pressure = DirichletBC(W.sub(1), Constant(0.0), pressure)
bc = [bc_noslip, bc_lid, bc_pressure]
# Variational formulation
w = TestFunction(W)
(v, q) = split(w)
# Define the initial conditions and previous time step values
# The previous step's solution
w_n = Function(W)
# Split the mixed function
(u_n, p_n) = split(w_n)
# Create and assign the initial conditions
w_init = Constant((0.0,0.0,0.0)) # null velocity field and pressure
w_n.interpolate(w_init)
# Create the current solution.
w_ = Function(W)
(u_, p_) = split(w_)
# Define symmetric gradient
def epsilon(u):
return sym(nabla_grad(u))
# Weak formulation F == 0
Eq1 = dot(u_,v)/dt - dot(u_n,v)/dt \
+ rho*dot(dot(u_, nabla_grad(u_)), v)*dx \
+ 2*mu*inner(epsilon(u_), grad(v))*dx \
- div(v)*p*dx \
- dot(f,v)*dx
Eq2 = q*div(u_)*dx
F = Eq1 + Eq2
# Create the progress bar.
progress = Progress('Time-stepping')
set_log_level(PROGRESS)
# Provide the Jacobian form.
J = derivative(F,w_)
# Create the Newton solver.
problem = NonlinearVariationalProblem(F, w_, bc, J)
solver = NonlinearVariationalSolver(problem)
prm = solver.parameters.newton_solver
prm.error_on_nonconvergence = False
# Save solution to file in VTK format
vtkfile1 = File('Navier Stokes Lid Driven Cavity/velocity.pvd')
vtkfile2 = File('Navier Stokes Lid Driven Cavity/pressure.pvd')
# Form and solve the linear system
for n in range(num_steps):
# Update current time
t += dt
# Solve the variational problem
solver.solve() # for use with the NonlinearVariationalProblem structure
(u_, p_) = split(w_) # Split solution
# Save data
vtkfile1 << (u_, t)
vtkfile << (p_, t)
# Update the previous solution
u_n.assign(u_)
# Update the progress bar
progress.update(t/T)
```

Community: FEniCS Project

### 1 Answer

1

Change the definition of

`Eq1`

to```
Eq1 = dot(u_,v)/dt*dx - dot(u_n,v)/dt*dx \
+ rho*dot(dot(u_, nabla_grad(u_)), v)*dx \
+ 2*mu*inner(epsilon(u_), grad(v))*dx \
- div(v)*p_*dx \
- dot(f,v)*dx
```

Thank you very much for solving my problem !

written
3 months ago by
Matthieu Diaz

Please login to add an answer/comment or follow this question.