Solve 1D Poiseuille flow defining the total flow rate

I want to solve the 1D Poiseuille flow in cylindrical coordinates defining the total flow rate Q and not the pressure gradient \frac{dp}{dz} . The problem is defined as:

\tau + r\frac{d\tau}{dr} = -r\frac{dp}{dz} (Eq.1) with the constraint \int_{Ω} u r \,dr = 1/2 (Eq.2)

One can impose the axial pressure gradient \frac{dp}{dz} which is the driven force of the problem, but I want to impose the flow rate Q (Eq.2). Now the unknowns of the problem are the axial velocity u(r) and the pressure gradient \frac{dp}{dz} . I know that I can express the constraint of Eq.2 via a Lagrange multiplier in accordance with the example of demo but I have a different integral to impose. How can I do it ?

Thanks

Hi, @giannoko. I solved this problem for the 1D Poiseulle flow in cartesian coordinates, and I think you can easily extend to cylindrical coordinates. Instead of the flow rate Q, I have prescribed the mean velocity U. The problem is given by

\mu \frac{d^2u}{dy^2} = \frac{dp}{dx}, with u(0) = u(H) = 0 and \int u dx/\int dx = U.

The result is in good agreement with the analytical equation \frac{\Delta p}{L} = \frac{12U\mu}{H^2},

from fenics import *
import matplotlib.pyplot as plt

U = 2.0
mu = 1.0
H = 1.0

mesh = IntervalMesh(100, 0.0, H)
V = FiniteElement('CG', mesh.ufl_cell(), 1)
R = FiniteElement('R', mesh.ufl_cell(), 0)
W = FunctionSpace(mesh, MixedElement([V, R]))

bc = DirichletBC(W.sub(0), Constant(0.0), "on_boundary")

u, c = TrialFunctions(W)
v, d = TestFunctions(W)

F1 = dot(mu*grad(u), grad(v))*dx + v*c*dx
F2 = d*(u - U)*dx
F = F1 + F2
a, L = lhs(F), rhs(F)
w = Function(W)
solve(a == L, w, bc)
u, c = w.split()

print("Pressure gradient - Analytical: %.4f Numerical: %.4f" % (-12.0*U*mu/(H**2),c(0.5)))

plot(u)
plt.ylabel("u")
plt.xlabel("y")
plt.show()

Thanks you so much @fforteneto . It works !