- Thank you received: 0
Forum Message
We have moved the forum to https://forum.ngsolve.org . This is an archived version of the topics until 05/05/23. All the topics were moved to the new forum and conversations can be continued there. This forum is just kept as legacy to not invalidate old links. If you want to continue a conversation just look for the topic in the new forum.
Notice
The forum is in read only mode.
Stokes Equtions with pure Dirichlet boundary
5 years 6 months ago - 5 years 6 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: 5 years 6 months ago by noname.
5 years 6 months ago - 5 years 6 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: 5 years 6 months ago by plederer.
Time to create page: 0.125 seconds
