Why aren't planes plane?

Question

Hi,

A few days ago I posted a question on how to simulate steady-state heat transfer on a plate and with someone's help I was able to successfully develop a script. However, the plots didn't show a flat plane, but rather a figure with bent/curved edges. Mesh refinements didn't work. I wonder what is wrong, and what should I do to get a perfect plane. The code is attached.

Thanks in advance for any help.

Fausto

from dolfin import *

if has_linear_algebra_backend("Epetra"):
    parameters["linear_algebra_backend"] = "Epetra"

class East(SubDomain): 
    def inside(self, x , on_boundary): 
        return near(x[0], 1.0) 

class West(SubDomain): 
    def inside(self, x , on_boundary): 
        return near(x[0], 0.0) 

class North(SubDomain): 
    def inside(self, x , on_boundary): 
        return near(x[1], 1.0) 

class South(SubDomain): 
    def inside(self, x , on_boundary): 
        return near(x[1], 0.0)

if __name__ == '__main__':

    mesh = UnitSquareMesh(64, 64)
    File("mesh.pvd") << mesh

    V = FunctionSpace(mesh, "Lagrange", 1)

    bcs = [DirichletBC(V, Constant(150), North()), \
           DirichletBC(V, Constant(250), East()), \
           DirichletBC(V, Constant(100), South()), \
           DirichletBC(V, Constant(50), West())] 

    u = Function(V)
    v = TestFunction(V)
    F = inner(grad(u), grad(v))*dx

    # Compute solution
    solve(F == 0, u, bcs, solver_parameters={"newton_solver":
                                            {"relative_tolerance": 1e-6}})

    # Plot solution and solution gradient
    plot(u, title="Solution")
    interactive()

Comments

Given the boundary conditions you have, why would expect a flat profile?

---

Maybe the confusion is just with the visualisation? The built in plot function warps the mesh with the Function value.
Try output to file and view with e.g. Paraview.

---

That worked, thanks! I didn't know the built-in plot() function would behave that way.

As "The FEniCS Book" says:

"The plot function is intended for quick examination of the solution during program development. More in-depth visual investigations of finite element solutions will normally benefit from using highly professional tools such as ParaView and VisIt."

Fausto

---

Hi MiroK, thanks for your comment. The BCs are temperatures; what's so wrong with them, and why would they bend the edges? Would you give me a little hand fixing the code?

Thanks,

Fausto