`Expression` with slope discontinuity

1

I am trying to apply a boundary condition for which I have created the following Expression:

expr = Expression(
    f"(x[0] < {ramp_start_x})? {disp_mag} : {disp_mag} * (10 - x[0]) / {10 - ramp_start_x}", degree=1
)

which is psuedo-code is as follows,

if x_coord < x_ramp_start:
    x_disp = fixed_disp
else:
    x_disp = linearly interpolated from fixed_disp to 0

I put degree = 1 since the maximum degree of the expression is 1.

However when I try to compute values from the expression, it’s not correct.

# Expression output for ramped drop in BCs at the top.
2.0, 2.0, 2.0, 2.0, 2.0, 2.0, **2.0, 1.344,** 1.2409, 1.13753, 1.0341, 0.93071, ... , 0.0

You can see that at the values in bold there is sudden drop (this is where the condition changes from true to false.), and then it ramps down to zero.

I have tried playing with the degree parameter for Expressions but that didn’t help. Honestly, I am not sure what the purpose of degree is exactly.

I converted the same thing to UserExpression object and it worked.

# Correct output (from UserExpression)
2.0, 2.0, 2.0, 2.0, 2.0, 2.0, **2.0, 1.9402,** 1.79104, 1.6417, 1.4925, 1.3432, ... , 0.0

I still want to know what’s wrong with Expression since we use it at many other places.

TIA

2

user1,

Can you provide a MWE where you demonstrate this issue? I can’t reproduce this error on fenics 2019.2.0.5 with a built in mesh.

from fenics import *
mesh = IntervalMesh(20, 0, 10)
ramp_start_x = 5
disp_mag = 2
expr = Expression(f"(x[0] < {ramp_start_x})? {disp_mag} : {disp_mag} * (10 - x[0]) / {10 - ramp_start_x}", degree=1)
result = ""
for coord in mesh.coordinates():
    result += "{0}, ".format(expr(coord))

print(result)

Fs = FunctionSpace(mesh, "P", 1)
interp = project(expr, Fs)
result = ""
for coord in mesh.coordinates():
    result += "{0:.3f}, ".format(interp(coord))

print(result)

Which for me results in:

2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 1.8, 1.6, 1.4, 1.2, 1.0, 0.8, 0.6, 0.4, 0.2, 0.0,
2.000, 2.000, 2.000, 2.000, 2.000, 2.000, 2.000, 2.000, 2.000, 2.000, 2.000, 1.800, 1.600, 1.400, 1.200, 1.000, 0.800, 0.600, 0.400, 0.200, 0.000,
4

Thank you for the response. Here’s the MWE:

# User installed module imports
import fenics
import numpy as np

# Displacement settings
top_disp_mag = 2
ramp_start_x = 3.300001   # Computed by different function.

# Define expression for BCs on the top side
expr_top_ramp_to_zero = fenics.Expression(f"(x[0] < {ramp_start_x})? {top_disp_mag} : {top_disp_mag} * (10 - x[0])/({10-ramp_start_x})", degree=1)

# Evaluate expression for points on the top
print("Expression output for ramped drop in BCs at the top..")
for x in np.linspace(0, 10, 21):
    print(expr_top_ramp_to_zero(x, 10), end=', ')

Actually, this printed out the correct value. I find this very strange because in the production code the variable expr_top_ramp_to_zero is not used in any way before the values are printed. I do however define new expressions in between. Can this cause an issue?

5

Turned out there was a silly mistake in the expression. Upon fixing it the issue was resolved. Embarrassing.

Thanks for the help!