Re: [firedrake] Unexpected asymmetry of solution of Poisson equation
Miklos is exactly right. What you see in the very corners is the following: if you have a quad with 3 zeroes and a nonzero (let's say a 1), if the quad is split along a 0-0 diagonal the you get a 0-0-0 triangle and a 0-1-0 triangle, but if the quad is split a 0-1 diagonal then you get two 0-0-1 triangles. It's even worse if the quad is +1 -1 -1 +1 in which case you get a 'valley' or a 'ridge', even though the actual function is symmetric about 0. To answer your last question, the triangular mesh isn't symmetric in that way. On 1 August 2016 at 11:08, George Ovchinnikov <lives9@gmail.com> wrote:
Dear all,
I'm solving Poisson equation and I found some irregularity. Here the code:
from firedrake import *
mesh = UnitSquareMesh(8, 8, quadrilateral=True)
V = FunctionSpace(mesh, "CG", 1) W = FunctionSpace(mesh, "DG", 0) w = Function(W) w.interpolate(Expression(1)) #Initial material distribution
bc0 = DirichletBC(V, Expression(("0")), 1) bc1 = DirichletBC(V, Expression(("0")), 2) bc2 = DirichletBC(V, Expression(("0")), 3) bc3 = DirichletBC(V, Expression(("0")), 4) bcs = [bc0, bc1, bc2, bc3]
u = TrialFunction(V) v = TestFunction(V)
f = Function(V) f.interpolate(Expression("1")) L = f * v * dx rhs = assemble(L)
a = (w * dot(grad(u), grad(v))) * dx frm = assemble(a, bcs = bcs) sol = Function(V) solve(frm, sol, rhs, solver_parameters = {'ksp_type': 'preonly', 'pc_type': 'lu'}) File("u.pvd") << sol
The asymmetry of solution (see the figure attached) is clearly visible: for example, upper left element is not the same as upper right.
I'm wondering what causes this and can this be remedied?
Also, I'm getting exactly same result (with the same asymmetry) on triangular, not quadrilateral mesh. Why solutions of triangular and quadrilateral meshes are the same?
Sincerely, George
Andrew, Miklos, thank you very much! On 01.08.2016 13:35, Andrew McRae wrote:
Miklos is exactly right. What you see in the very corners is the following: if you have a quad with 3 zeroes and a nonzero (let's say a 1), if the quad is split along a 0-0 diagonal the you get a 0-0-0 triangle and a 0-1-0 triangle, but if the quad is split a 0-1 diagonal then you get two 0-0-1 triangles.
It's even worse if the quad is +1 -1 -1 +1 in which case you get a 'valley' or a 'ridge', even though the actual function is symmetric about 0.
To answer your last question, the triangular mesh isn't symmetric in that way.
On 1 August 2016 at 11:08, George Ovchinnikov <lives9@gmail.com <mailto:lives9@gmail.com>> wrote:
Dear all,
I'm solving Poisson equation and I found some irregularity. Here the code:
from firedrake import *
mesh = UnitSquareMesh(8, 8, quadrilateral=True)
V = FunctionSpace(mesh, "CG", 1) W = FunctionSpace(mesh, "DG", 0) w = Function(W) w.interpolate(Expression(1)) #Initial material distribution
bc0 = DirichletBC(V, Expression(("0")), 1) bc1 = DirichletBC(V, Expression(("0")), 2) bc2 = DirichletBC(V, Expression(("0")), 3) bc3 = DirichletBC(V, Expression(("0")), 4) bcs = [bc0, bc1, bc2, bc3]
u = TrialFunction(V) v = TestFunction(V)
f = Function(V) f.interpolate(Expression("1")) L = f * v * dx rhs = assemble(L)
a = (w * dot(grad(u), grad(v))) * dx frm = assemble(a, bcs = bcs) sol = Function(V) solve(frm, sol, rhs, solver_parameters = {'ksp_type': 'preonly', 'pc_type': 'lu'}) File("u.pvd") << sol
The asymmetry of solution (see the figure attached) is clearly visible: for example, upper left element is not the same as upper right.
I'm wondering what causes this and can this be remedied?
Also, I'm getting exactly same result (with the same asymmetry) on triangular, not quadrilateral mesh. Why solutions of triangular and quadrilateral meshes are the same?
Sincerely, George
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
participants (2)
- 
                
                Andrew McRae
- 
                
                George Ovchinnikov