How to Make a Custom Tensor Expression


197
views
1
5 months ago by
Saro  
I am trying to make 2x2 tensor expression, but it is not working. How can I make a tensor expression?
from fenics import *
class MarkerTest(Expression):
    def eval(self, values, x):
        values[0,0] = 1.0
        values[1,0] = 2.0
        values[0,1] = 3.0
        values[1,1] = 4.0
    def value_shape(self):
        return (2,2)

marker = MarkerTest(degree=1)
print(marker(0))​

Here is the type of error message encountered:

---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
<ipython-input> in <module>()
     10 
     11 marker = MarkerTest(degree=1)
---> 12 print(marker(0))

/usr/lib/python3/dist-packages/dolfin/functions/expression.py in __call__(self, *args, **kwargs)
    863 
    864         # The actual evaluation
--> 865         self.eval(values, x)
    866 
    867         # If scalar return statement, return scalar value.

<ipython-input> in eval(self, values, x)
      2 class MarkerTest(Expression):
      3     def eval(self, values, x):
----> 4         values[0,0] = 1.0
      5         values[1,0] = 2.0
      6         values[0,1] = 3.0

IndexError: too many indices for array
Community: FEniCS Project

1 Answer


3
5 months ago by
As the error message says, there ar too many indices. For an array, there should be only one:
class MarkerTest(Expression):
    def eval(self, values, x):
        values[0] = 1.0
        values[1] = 2.0
        values[2] = 3.0
        values[3] = 4.0
    def value_shape(self):
        return (2,2)​

Thank you for your reply.

So would the tensor be as follows?

\begin{bmatrix}
values[0] & values[2]  \\
values[1] & values[3] \\
\end{bmatrix}
written 5 months ago by Saro  
It is the inverse, i.e.,
\[ \begin{bmatrix} \mathrm{values}[0] & \mathrm{values}[1] \\ \mathrm{values}[2] & \mathrm{values}[3] \end{bmatrix} \]
You can check it with something like
mesh = UnitSquareMesh(1, 1)
n1 = Constant(("1.0", "0.0"))
n2 = Constant(("0.0", "1.0"))

print(assemble(dot(n1, marker*n1)*dx(domain=mesh)))
print(assemble(dot(n1, marker*n2)*dx(domain=mesh)))
print(assemble(dot(n2, marker*n1)*dx(domain=mesh)))
print(assemble(dot(n2, marker*n2)*dx(domain=mesh)))​
written 5 months ago by Adam Janecka  
Thank you. That solves my problem.
written 5 months ago by Saro  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »