Assembly of derivative of Hessian


78
views
0
6 weeks ago by
Hello,

For my problem, I need the derivative of the Hessian (third order tensor), but only one diagonal of it (matrix). As a MWE, if I want to calculate  $A_{ij}=\int_{\Omega}\phi_i\phi_i\phi_jd\Omega$Aij=ΩϕiϕiϕjdΩ (no summation on i), I try the following code
from dolfin import *

mesh = UnitIntervalMesh(5)
V = FunctionSpace(mesh, 'P', 1)

u  = Function(V)
du = TrialFunction(V)
v  = TestFunction(V)
v2 = Argument(V,2)

A = assemble(v*v*du*dx)
A = assemble(v2*v*du*dx)
​

However, I get an error for both of the last two statements. assemble(v*v*du*dx) raises an error that argument number for v is overlapping, whereas assemble(v2*v*du*dx) gives an error because, I believe, FEniCS does not assemble third order tensors. I am able to calculate A only by looping through the degrees of freedom and putting a hat function for v2 one by one, however that is extremely inefficient. Any ideas, how I could achieve this?

Thanks!

Community: FEniCS Project
Where are you going to use this third order tensor? You could assemble its action instead. In the same way matrix free method assemble A*x (for a given x) instead of assembling the whole matrix A.
written 6 weeks ago by Miguel  
Thanks for your response. That works if I assemble the action of $B_{ijk}=\int_{\Omega}\phi_i\phi_j\phi_kd\Omega$Bijk=ΩϕiϕjϕkdΩ on a vector (Function or Coefficient object) as follows
from dolfin import *

mesh = UnitIntervalMesh(5)
V = FunctionSpace(mesh, 'P', 1)

u  = Function(V)
du = TrialFunction(V)
v  = TestFunction(V)
v2 = Argument(V,2)

x = Function(V)
#fill up x with some values
A = assemble(x*v*du*dx)
​

However, for getting  $A_{ij}$Aij I cannot think of a single vector  $x$x  such that  $B\cdot x=A$B·x=A . Since I need the matrix A in my solver, I cannot assemble its action on another vector. I have been trying various ways to link the coefficients of  $x$x with that of $v$v , however cannot figure it out.
written 6 weeks ago by Ankush Aggarwal  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »