Problems interpolating a UserExpression and plotting it

1

Hi everyone,

I’m using macOS Mojave (10.14.6), miniconda, and FEniCS 2018 with Python. I want to define have an expression/function e: \mathbb{R}^2 \rightarrow \mathbb{R} which is defined as such that f(x, y) = 1 for x \leq 0.2, f(x, y) = 2 for x \geq 0.8, f(x, y) = 3 for y \geq 0.9 and 0 else. So I thought that this

from fenics import *
N = 16
mesh = UnitSquareMesh(N, N)
V = FunctionSpace(mesh, 'P', 1)

class MyExpression(UserExpression):
  def eval(self, value, x):
    if x[0] <= 0.2:
        value[0] = 1
    elif x[0] >= 0.8:
        value[0] = 2
    elif x[1] >= 0.9:
        value[0] = 3
    else:
        value[0] = 0

  def value_shape(self):
    return (1,)

e = MyExpression(degree = 2)
Ie = interpolate(e, V)

should work. But I get the error

RuntimeError: 

*** -------------------------------------------------------------------------
*** DOLFIN encountered an error. If you are not able to resolve this issue
*** using the information listed below, you can ask for help at
***
***     fenics-support@googlegroups.com
***
*** Remember to include the error message listed below and, if possible,
*** include a *minimal* running example to reproduce the error.
***
*** -------------------------------------------------------------------------
*** Error:   Unable to interpolate function into function space.
*** Reason:  Rank of function (1) does not match rank of function space (0).
*** Where:   This error was encountered inside FunctionSpace.cpp.
*** Process: 0
*** 
*** DOLFIN version: 2018.1.0
*** Git changeset:  
*** -------------------------------------------------------------------------

I do not know what the reason wants me to tell. Projecting does also not work. Could there be a problem with if/else test in my expression? Thanks for help.

2

If I omit

  def value_shape(self):
    return (1,)

then I get the warning:
WARNING: user expression has not supplied value_shape method or an element. Assuming scalar element.
Could someone explain this to me?

3

UFL distinguishes between a scalar, and a vector of rank 1 and dimension 1. You can get rid of the warning with

  def value_shape(self):
    #return (1,)
    return ()
4 Likes
4

Thank you very much.