← Back to team overview

dolfin team mailing list archive

Re: [Bug 531207] [NEW] Problem with DOLFIN_EPS

 

Defining boundary conditions based on coordinates is prone to
round-off errors. That's why DOLFIN_EPS is used in many of the demos.
But it's not a fool-proof solution.

Try scaling the margin of error by your stretch parameter and see if
that helps:

  def boundary1(x):
      return x[0] > stretch - stretch*DOLFIN_EPS

--
Anders


On Wed, Mar 03, 2010 at 09:25:39AM -0000, Chris Richardson wrote:
> Public bug reported:
>
> I am having problems defining boundary conditions using DOLFIN_EPS.
> See example below, which works when stretch=1.0, but fails if stretch> 4.5
> I suppose this could be something to do with machine precision, so would
> be interesting to see if it is reproducable on different systems...
>
>
> from dolfin import *
>
> stretch = 4.5
>
> mesh = UnitSquare(40,20)
> for x in mesh.coordinates():
>     x[0] *= stretch
>
> Vs = VectorFunctionSpace(mesh, "CG", 1)
>
> E  = 10.0
> nu = 0.35
>
> mu    = E / (2.0*(1.0 + nu))
> lmbda = E*nu / ((1.0 + nu)*(1.0 - 2.0*nu))
>
> def sigma(v):
>     return 2.0*mu*sym(grad(v)) + lmbda*tr(sym(grad(v)))*Identity(v.cell().d)
>
>
> # Define Dirichlet boundaries (x = 0 or x = stretch)
> def boundary0(x):
>     return x[0] < DOLFIN_EPS
>
> def boundary1(x):
>     return x[0] > stretch - DOLFIN_EPS
>
>
> # Define boundary condition
> u0 = Expression(("0.0","0.0"))
> bc0 = DirichletBC(Vs, u0, boundary0)
> u1 = Expression(("0.1","0.0"))
> bc1 = DirichletBC(Vs, u1, boundary1)
> bcs = [bc0,bc1]
>
> # Define variational problem
> v = TestFunction(Vs)
> u = TrialFunction(Vs)
>
> f = Constant((0.0,0.0))
> a = inner(grad(v), sigma(u))*dx
> L = inner(v, f)*dx
>
> # Compute solution
> problem = VariationalProblem(a, L, bcs)
> u = problem.solve()
>
> # Plot solution
> plot(u, mode="displacement", interactive=True)
>
> mesh.move(u)
> plot(mesh,interactive=True)
>
> ** Affects: dolfin
>      Importance: Undecided
>          Status: New
>

-- 
Problem with DOLFIN_EPS
https://bugs.launchpad.net/bugs/531207
You received this bug notification because you are a member of DOLFIN
Team, which is subscribed to DOLFIN.

Status in DOLFIN: New

Bug description:
I am having problems defining boundary conditions using DOLFIN_EPS.
See example below, which works when stretch=1.0, but fails if stretch> 4.5
I suppose this could be something to do with machine precision, so would
be interesting to see if it is reproducable on different systems...


from dolfin import *

stretch = 4.5

mesh = UnitSquare(40,20)
for x in mesh.coordinates():
    x[0] *= stretch

Vs = VectorFunctionSpace(mesh, "CG", 1)

E  = 10.0
nu = 0.35

mu    = E / (2.0*(1.0 + nu))
lmbda = E*nu / ((1.0 + nu)*(1.0 - 2.0*nu))

def sigma(v):
    return 2.0*mu*sym(grad(v)) + lmbda*tr(sym(grad(v)))*Identity(v.cell().d)


# Define Dirichlet boundaries (x = 0 or x = stretch)
def boundary0(x):
    return x[0] < DOLFIN_EPS
  
def boundary1(x):
    return x[0] > stretch - DOLFIN_EPS


# Define boundary condition
u0 = Expression(("0.0","0.0"))
bc0 = DirichletBC(Vs, u0, boundary0)
u1 = Expression(("0.1","0.0"))
bc1 = DirichletBC(Vs, u1, boundary1)
bcs = [bc0,bc1]

# Define variational problem
v = TestFunction(Vs)
u = TrialFunction(Vs)

f = Constant((0.0,0.0))
a = inner(grad(v), sigma(u))*dx
L = inner(v, f)*dx

# Compute solution
problem = VariationalProblem(a, L, bcs)
u = problem.solve()

# Plot solution
plot(u, mode="displacement", interactive=True)

mesh.move(u)
plot(mesh,interactive=True)





References