- Thank you received: 0
Stokes Equtions with pure Dirichlet boundary
6 years 2 months ago - 6 years 2 months ago #1225
by noname
Stokes Equtions with pure Dirichlet boundary was created by noname
Hello everyone,
I am trying to solve Stokes equations on a unit square with pure Dirichlet boundary conditions. I know in this case I need to modify the bilinear form, otherwise the solution is unstable. I have taken a look to the tutorials but they all consider homogenous Neumann conditions in on of the boundaries, none of them is pure Dirichlet. I am trying something like this (adding a lagrange multiplier to the bilinear form to incorporate pure Dirichlet bcs)
from ngsolve import *
import netgen.gui
But the results makes no sense.. Anyone knows how to impose this kind of boundary condition?
Thanks in advance.
I am trying to solve Stokes equations on a unit square with pure Dirichlet boundary conditions. I know in this case I need to modify the bilinear form, otherwise the solution is unstable. I have taken a look to the tutorials but they all consider homogenous Neumann conditions in on of the boundaries, none of them is pure Dirichlet. I am trying something like this (adding a lagrange multiplier to the bilinear form to incorporate pure Dirichlet bcs)
from ngsolve import *
import netgen.gui
Code:
from netgen.geom2d import SplineGeometry
geo = SplineGeometry()
geo.AddRectangle( (-1, -1), (1, 1), bcs = ("wall" , "outlet", "wall", "inlet"))
mesh = Mesh( geo.GenerateMesh(maxh=0.25))
mesh.Curve(3)
Draw (mesh)
V = VectorH1(mesh, order=2, dirichlet="wall|inlet|outlet")
Q = H1(mesh, order=1)
X = FESpace([V,Q])
u,p = X.TrialFunction()
v,q = X.TestFunction()
r = Q.TrialFunction()
s = Q.TestFunction()
a = BilinearForm(X)
a += SymbolicBFI(InnerProduct(grad(u),grad(v))+div(u)*q+div(v)*p)
a += SymbolicBFI(p*s+q*r)
a.Assemble()
gfu = GridFunction(X)
res = gfu.vec.CreateVector()
res.data = -a.mat * gfu.vec
inv = a.mat.Inverse(freedofs=X.FreeDofs())
gfu.vec.data += inv * res
But the results makes no sense.. Anyone knows how to impose this kind of boundary condition?
Thanks in advance.
Last edit: 6 years 2 months ago by noname.
6 years 2 months ago - 6 years 2 months ago #1244
by plederer
Replied by plederer on topic Stokes Equtions with pure Dirichlet boundary
You implemented Dirchlet BC on the whole domain correctly. I am not sure what's your intention behind adding the term
.
However, of course you need to add a driving force, otherwise your solution is zero. Your code without your mysterious term and an example force:
Code:
a += SymbolicBFI(p*s+q*r)
However, of course you need to add a driving force, otherwise your solution is zero. Your code without your mysterious term and an example force:
Code:
from netgen.geom2d import SplineGeometry
geo = SplineGeometry()
geo.AddRectangle( (-1, -1), (1, 1), bcs = ("wall" , "outlet", "wall", "inlet"))
mesh = Mesh( geo.GenerateMesh(maxh=0.25))
#mesh.Curve(3)
Draw (mesh)
V = VectorH1(mesh, order=2, dirichlet="wall|inlet|outlet")
Q = H1(mesh, order=1)
X = FESpace([V,Q])
u,p = X.TrialFunction()
v,q = X.TestFunction()
r = Q.TrialFunction()
s = Q.TestFunction()
a = BilinearForm(X)
a += SymbolicBFI(InnerProduct(grad(u),grad(v))+div(u)*q+div(v)*p)
#a += SymbolicBFI(p*s+q*r)
a.Assemble()
f = LinearForm(X)
f += SymbolicLFI(-x*v[1])
f.Assemble()
gfu = GridFunction(X)
res = gfu.vec.CreateVector()
res.data = f.vec -a.mat * gfu.vec
inv = a.mat.Inverse(freedofs=X.FreeDofs())
gfu.vec.data += inv * res
Draw(gfu.components[0], mesh, "vel")
Last edit: 6 years 2 months ago by plederer.
Time to create page: 0.120 seconds