writing a tensor

Asked by Sungick Kim

Hi.
I want to write a tensor F_ij, where i,j=0,1,2.
Can you give me a hand?

So not working code is
-------------------------------------------------------------
V = VectorFunctionSpace(sphere, "CG", 1)
v = TestFunction(V)

ID = Identity(v.cell().d) # Identity tensor
F= 0*ID #initializing
F[0][0] = 1*u1 # u1 is a scalar field.
F[1][1] = 2*u1
F[2][2] = 3*u1
-------------------------------------------------------------

Question information

Language:
English Edit question
Status:
Solved
For:
DOLFIN Edit question
Assignee:
No assignee Edit question
Solved by:
Johan Hake
Solved:
Last query:
Last reply:
Revision history for this message
Sungick Kim (sungick) said :
#1

https://answers.launchpad.net/dolfin/+question/102501

in the question 102501, I found similar question and it actually solved my problem.
I can write the tensor as the answers in question 102501. FEniCS doesn't complain about it.
But I think that it's not solving problem.
Here are code and the result

<my code>
f00 = Expression("(10.45-9.96)/10.45 *u1c")
f11 = Expression("(5.58-6.05)/5.58 *u1c")
f22 = Expression("(4.86-4.74)/4.86 *u1c")
z = Constant(0.0)
F2= as_matrix( (f00,z,z),(z,f11,z),(z,z,f22) )

I = Identity(v2.cell().d) # Identity tensor
F1 = I + grad(u2c)
F = F1*F2 # Deformation gradient

<result>

(1) computing Transport problem
  Solving nonlinear variational problem
    Starting Newton solve.
      Applying boundary conditions to linear system.
      Applying boundary conditions to linear system.
      Solving linear system of size 7381 x 7381 (PETSc LU solver, umfpack).
      Applying boundary conditions to linear system.
      Newton iteration 1: r (abs) = 0.000e+00 (tol = 1.000e-10) r (rel) = -nan (tol = 1.000e-14)

Revision history for this message
Best Johan Hake (johan-hake) said :
#2

> in the question 102501, I found similar question and it actually solved my
> problem. I can write the tensor as the answers in question 102501. FEniCS
> doesn't complain about it. But I think that it's not solving problem.
> Here are code and the result

No you are probably not solving you problem. u1c is here interpreted as a
scalar value, which is initiated to 0. You can set the value of that scalar
by:

  f00.u1c = some_value

but I guess that ulc really is a Function or another Expression. Then you need
a bit more elaborated Expression. Something like:

code = """
class MyExpr: public Expression
{
public:
  boost::shared_ptr<dolfin::GenericFunction> ulc;
  MyExpr():Expression()
  {}

  void eval(dolfin::Array<double>& values, const dolfin::Array<double>& x)
const
  {
   assert(ulc);
   ulc->eval(values, x);
   values[0] = (10.45-9.96)/10.45*values[0];
  }
};
"""

f00 = Expression(code)
f00.ulc = ulc
ulc.vector()[:]=1
f00(0,0,0)

=> 0.046889952153109919

But... I now see that your expression can be writte directly in ufl which is
much easier:

  f00 = (10.45-9.96)/10.45 *u1c
  f11 = (5.58-6.05)/5.58 *u1c
  f22 = (4.86-4.74)/4.86 *u1c
  z = 0
  F2= as_matrix(((f00,z,z),(z,f11,z),(z,z,f22)))

Johan

Revision history for this message
Sungick Kim (sungick) said :
#3

Thanks Johan Hake, that solved my question.