### Integrate derivative of known function

127
views
0
4 months ago by
Hi guys,

To calculate the velocity in an advection problem I first calculate pressure and then

vel = (-w**2/(12*mu))*grad(p)​

where w and mu are constants, If I solve the pure advection equation in the case of a constant velocity vector
$\frac{\partial c}{\partial t}+u_i\frac{\partial c}{\partial x_i}=0$ct +uicxi =0
where c is the concentration to solve and u is the velocity vector. After integration by parts to implement Neumann BC, the following code works
q = TestFunction(V)
c = TrialFunction(V)

F = (1/dt)*w*c*q*dx - (1/dt)*w*c_1*q*dx \
+inner(c,dot(vel,n))*q*ds(1)

a, L = lhs(F),rhs(F)

c = Function(V)

t = 0.0

for n in range(num_steps):
solve(a==L,c,bc)
c_1.assign(c)
print('t = %.4f' % (t))
t += dt​

But when I try to solve the case where the velocity changes with x
$\frac{\partial c}{\partial t}+u_i\frac{\partial c}{\partial x_i}+c\frac{\partial u_i}{\partial x_i}=0$ct +uicxi +cuixi =0
the code below pops up the error "Can only integrate scalar expressions. The integrand is a tensor expression with value shape (2,) and free indices with labels ()." How could I improve the code? thank you very much for your comments

F = (1/dt)*w*c*q*dx - (1/dt)*w*c_1*q*dx \
+c*Dx(vel,0)*q*dx+c*Dx(vel,1)*q*dx\
+inner(c,dot(vel,n))*q*ds(1)

#I also try

F = (1/dt)*w*c*q*dx - (1/dt)*w*c_1*q*dx \
+inner(c,dot(vel,n))*q*ds(1)

Community: FEniCS Project

1
4 months ago by
The second line of your first version
+c*Dx(vel,0)*q*dx+c*Dx(vel,1)*q*dx\​

translates to a vectorial expression. Both terms contain derivatives of the whole vector w.r.t. one spatial variable yielding a vectorial expression, which does not represent the divergence you gave in your equation.

Similarly the second line of your second version contains the term grad(vel) which gives a second order tensor.

The divergence  $\frac{\partial u_i}{\partial x_i}$uixi   can be written as div(u) or as u[i].dx(i).