Stokes alpha-model

Asked by Coraline Faurie

Hello all !

I want to solve the stokes alpha equation, which is basically the stokes equation (momentum equation + incompressibility for the velocity v) with in addition the definition of a filtered velocity u: v = (1-alpha² nabla²)u where alpha in a constant.

To do so, I defined 3 trial and 3 test functions corresponding to the velocity v, the filtered velocity u, and the pressure p:

V = VectorFunctionSpace(mesh, "CG", 2, constrained_domain=pbc)
U = VectorFunctionSpace(mesh, "CG", 2, constrained_domain=pbc)
Q = FunctionSpace(mesh, "CG", 1, constrained_domain=pbc)
W = MixedFunctionSpace([V, U, Q])

1) When I defined W = V * U * Q in stead of with MixedFunctionSpace, I have a problem with the dirichlet boundary condition (no-slip condition for u and v) which I defined with the following lines:

zero = Expression(("0.0","0.0"))
bc1 = DirichletBC(W.sub(0), zero, top_bottom)
bc2 = DirichletBC(W.sub(1), zero, top_bottom)

I get:

*** Error: Unable to create Dirichlet boundary condition.
*** Reason: Illegal value dimension (2), expecting (4).
*** Where: This error was encountered inside DirichletBC.cpp.
*** Process: 0

Which I can solve with:

zero = Expression(("0.0","0.0"))
bc1 = DirichletBC(W.sub(0).sub(1), zero, top_bottom)

But then, what W.sub(0).sub(1) does ? Is the no slip condition applied to v and u ?

2) To solve my problem, i set:

(v, u, p) = TrialFunction(W)
(w, m, q) = TestFunction(W)

and after defining a and L, I assemble my system and solve it for a function S =Function(W). Then, I want to split S to obtain v, u and p. Unfortunately, I have this message error:

    (v, u, p) = TrialFunction(W)
ValueError: too many values to unpack

Is it possible to define 3 Trial/Test Functions ?

I would welcome any answers or suggestions,

Thanks a lot for your help !



Question information

English Edit question
DOLFIN Edit question
No assignee Edit question
Solved by:
Coraline Faurie
Last query:
Last reply:
Revision history for this message
Coraline Faurie (linesintheplace) said :

For those who might be interested I just forgot the 's' to trial and test functions. It works with:
(v, u, p) = TrialFunctions(W)
(w, m, q) = TestFunctions(W)