Conditional Expression


157
views
0
3 months ago by
Birkan  

Dear all,

I need to assign the distance from the center on a circle with radius 0.25 within a 1x1 square mesh. I want to assign 0.25 everywhere out of circle. I tried the following but does not work.

from fenics import *
from mshr import *
import ufl

parameters["form_compiler"]["representation"] = 'uflacs'

# Generate mesh
mesh = RectangleMesh(Point(0.0,0.0), Point(1.0,1.0), 30, 30)

# Function space
P=FiniteElement("Lagrange", mesh.ufl_cell(), 2)
V = FunctionSpace(mesh,P)

# Plot file
cfile = File("cfile.pvd")

# Expression
c=Expression("sqrt((x[0]-0.5)*(x[0]-0.5) + (x[1]-0.5)*(x[1]-0.5))", degree=2)
c=conditional(ge(sqrt((x[0]-0.5)*(x[0]-0.5) + (x[1]-0.5)*(x[1]-0.5)),0.25),Constant(0.25),c)

# Save results
cplot=interpolate(c,V)
cplot.rename("c", "c"); cfile << cplot


The error is the following:

c=conditional(ge(sqrt((x[0]-0.5)*(x[0]-0.5) + (x[1]-0.5)*(x[1]-0.5)),0.25),Constant(0.25),c)
NameError: name 'x' is not defined
Aborted (core dumped)

Other than conditional statement, the code is working. What I need to get is the following (plotted solution obtained w/o conditional statement by scaling result between 0-0.25):



Do you know how to define a conditional expression for this?

Thank you
Birkan

Community: FEniCS Project

1 Answer


0
3 months ago by
Just define x as x = SpatialCoordinate(mesh)
Thank you Adam,

If I do this, I get the following error:

cplot=interpolate(c,V)
File "/usr/lib/python2.7/dist-packages/dolfin/fem/interpolation.py", line 76, in interpolate
Pv.interpolate(v)
TypeError: in method 'Function_interpolate', argument 2 of type 'dolfin::GenericFunction const &'
Aborted (core dumped)

Do you know how to fix this?

Thanks
written 3 months ago by Birkan  
I fixed by cplot=project(c,V)
written 3 months ago by Birkan  
Yeah, using project instead of interpolate resolves it. I'm not sure about the degree though since conditional is an UFL operator, not an Expression.
written 3 months ago by Adam Janecka  
Hi Adam,

I wanted to use this conditional expression (will be defined as exact solution of the problem inside circle) to compute error norm. I realized errornorm function of fenics does not accept conditional statements. Thus, do you know how I can define this as an expression without using conditional UFL operator? What I need to have an exact solution inside circle and a constant outside.

Thank you
written 3 months ago by Birkan  
Most probably, the problem is that errornorm is implemented using interpolate, see https://github.com/FEniCS/dolfin/blob/master/site-packages/dolfin/fem/norms.py. As before, I'd try to project the conditional to express the norm, i.e., something like this (from the source code):
# Interpolate functions into finite element space
pi_u = interpolate(u, V)
pi_uh = project(uh, V)

# Compute the difference
e = Function(V)
e.assign(pi_u)
e.vector().axpy(-1.0, pi_uh.vector())

# Compute norm
norm(e, norm_type=norm_type, mesh=mesh)​
where uh is your conditional statement and V is an appropriate function space.
written 3 months ago by Adam Janecka  
It worked, thank you
written 12 weeks ago by Birkan  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »