Creating vector field using a stress tensor

Question

I want to couple two different problems. One gives me a stress tensor (sigma[i,j]) (For example out of an hyperelastic problem). And then i want to use this stress tensor to calculate a "vector field" for the second problem. I want to calculate the components of the vector field in the following way:

n = Constant((1.0,0.0)
m = Constant((0.0,1.0)
tau_1=m[i]*sigma[i,j]*n[j]
tau_2=n[i]*sigma[i,j]*m[j]


Now i want to create a "vector field" using tau_1 and tau_2:


V = [tau_1,tau_2]


At the end this V will be used to assemble the following K matrix:


Tr = TrialFunction(FunctionSpace(mesh,'CG',1) 
Ts = TestFunction(FunctionSpace(mesh,'CG',1)

K = assemble(-div(Tr*V)*Ts*dV)


To show u what i need: I was able to assemble K with a simplified version of V using "expression" and "interpolation":


V_expression = Expression(('0.1,0.0),degree=2,domain=mesh)
V = interpolate(V_expression, VectorFunctionSpace(mesh,'CG',1)) 


But this simplified version is not sufficient. This simple V is constant in the whole area/mesh and i need it to be dependent on sigma and by that dependent on the location. So is it possible to implement the version above?

Im not sure, if explanation is adequate. Please say if not. Thanks a lot.

Answer 1

Hi, consider

from dolfin import *

mesh = UnitSquareMesh(4, 4)

T = TensorElement('Lagrange', mesh.ufl_cell(), 1)
T = FunctionSpace(mesh, T)

sigma = interpolate(Constant(((1, 2), (4, 5))), T)

e1 = Constant((1, 0))
e2 = Constant((0, 1))

s12 = e1[i]*sigma[i, j]*e2[j] #inner(e1, dot(sigma, e2))
s21 = e2[i]*sigma[i, j]*e1[j] #inner(e2, dot(sigma, e1))
f = as_vector((s12, s21))

V = FunctionSpace(mesh, 'CG', 1)
u = TrialFunction(V)
v = TestFunction(V)

a = -div(f*u)*v*dx
A = assemble(a)

print '|A|', A.norm('linf')

Comments

I did this. It seems to work. But now, each time i use "A = assemble(a)", the FFC just-in-time compiler is called. This makes it very slow.

I did print as_vector((s12, s21)) and the result was for the first iteration:

[sum_{i_1} f_53[i_1] * (sum_{i_0} f_54[i_0] * (({ A | A_{i_8, i_9} = 0.5 * ({ A | A_{i_0, i_1} = 161538461.538 * 0.5 * ((grad(f_33[i_0]))[i_1] + (grad(f_33[i_1]))[i_0]) + I[i_0, i_1] * 121153846.154 * (sum_{i_2} ({ A | A_{i_0, i_1} = 0.5 * ((grad(f_33[i_0]))[i_1] + (grad(f_33[i_1]))[i_0]) })[i_2, i_2] ) })[i_8, i_9] }) + f_26)[i_0, i_1] ) , sum_{i_1} f_54[i_1] * (sum_{i_0} f_53[i_0] * (({ A | A_{i_8, i_9} = 0.5 * ({ A | A_{i_0, i_1} = 161538461.538 * 0.5 * ((grad(f_33[i_0]))[i_1] + (grad(f_33[i_1]))[i_0]) + I[i_0, i_1] * 121153846.154 * (sum_{i_2} ({ A | A_{i_0, i_1} = 0.5 * ((grad(f_33[i_0]))[i_1] + (grad(f_33[i_1]))[i_0]) })[i_2, i_2] ) })[i_8, i_9] }) + f_26)[i_0, i_1] ) ]

and then in the second iteration, it becomes longer and the assembling gets slower:

[sum_{i_1} f_129[i_1] * (sum_{i_0} f_130[i_0] * (({ A | A_{i_{286}, i_{287}} = 0.5 * ({ A | A_{i_0, i_1} = 161538461.538 * 0.5 * ((grad(f_109[i_0]))[i_1] + (grad(f_109[i_1]))[i_0]) + I[i_0, i_1] * 121153846.154 * (sum_{i_2} ({ A | A_{i_0, i_1} = 0.5 * ((grad(f_109[i_0]))[i_1] + (grad(f_109[i_1]))[i_0]) })[i_2, i_2] ) })[i_{286}, i_{287}] }) + ({ A | A_{i_8, i_9} = 0.5 * ({ A | A_{i_0, i_1} = 161538461.538 * 0.5 * ((grad(f_33[i_0]))[i_1] + (grad(f_33[i_1]))[i_0]) + I[i_0, i_1] * 121153846.154 * (sum_{i_2} ({ A | A_{i_0, i_1} = 0.5 * ((grad(f_33[i_0]))[i_1] + (grad(f_33[i_1]))[i_0]) })[i_2, i_2] ) })[i_8, i_9] }) + f_26)[i_0, i_1] ) , sum_{i_1} f_130[i_1] * (sum_{i_0} f_129[i_0] * (({ A | A_{i_{286}, i_{287}} = 0.5 * ({ A | A_{i_0, i_1} = 161538461.538 * 0.5 * ((grad(f_109[i_0]))[i_1] + (grad(f_109[i_1]))[i_0]) + I[i_0, i_1] * 121153846.154 * (sum_{i_2} ({ A | A_{i_0, i_1} = 0.5 * ((grad(f_109[i_0]))[i_1] + (grad(f_109[i_1]))[i_0]) })[i_2, i_2] ) })[i_{286}, i_{287}] }) + ({ A | A_{i_8, i_9} = 0.5 * ({ A | A_{i_0, i_1} = 161538461.538 * 0.5 * ((grad(f_33[i_0]))[i_1] + (grad(f_33[i_1]))[i_0]) + I[i_0, i_1] * 121153846.154 * (sum_{i_2} ({ A | A_{i_0, i_1} = 0.5 * ((grad(f_33[i_0]))[i_1] + (grad(f_33[i_1]))[i_0]) })[i_2, i_2] ) })[i_8, i_9] }) + f_26)[i_0, i_1] ) ]

And so on ... So it gets longer and slower each iteration. Why does it happen? Is there a way to avoid the symbolic representation? Or is there any other reason for the compiler beeing called in each iteration?

Strangely it does not happen, after i once started the simulation and stopped it. Then starting again without changing the code does not lead to a repeated call of the FFC j-i-t compiler.

I feel quite helpless here. Any advice?