New to Fenics


277
views
0
6 months ago by
Hi,
I have a rather simple question. For the code that I want to write I need to have a vector-matrix multiplication. I am basically trying to change the hyperelasticity code on the website. Here is what I have

# Create mesh and define function space
mesh = UnitCubeMesh(24, 16, 16)
V = VectorFunctionSpace(mesh, "Lagrange", 1)

# Mark boundary subdomians
left =  CompiledSubDomain("near(x[0], side) && on_boundary", side = 0.0)
right = CompiledSubDomain("near(x[0], side) && on_boundary", side = 1.0)

# Define Dirichlet boundary (x = 0 or x = 1)
c = Expression(("0.0", "0.0", "0.0"))
r = Expression(("scale*0.0",
                "scale*(y0 + (x[1] - y0)*cos(theta) - (x[2] - z0)*sin(theta) - x[1])",
                "scale*(z0 + (x[1] - y0)*sin(theta) + (x[2] - z0)*cos(theta) - x[2])"),
                scale = 0.5, y0 = 0.5, z0 = 0.5, theta = pi/3)

bcl = DirichletBC(V, c, left)
bcr = DirichletBC(V, r, right)
bcs = [bcl, bcr]

# Define functions
du = TrialFunction(V)            # Incremental displacement
v  = TestFunction(V)             # Test function
u  = Function(V)                 # Displacement from previous iteration
B  = Constant((0.0, -0.5, 0.0))  # Body force per unit volume
T  = Constant((0.1,  0.0, 0.0))  # Traction force on the boundary

# Kinematics
d = u.geometric_dimension()
I = Identity(d)             # Identity tensor
F = I + grad(u)             # Deformation gradient
C = F.T*F                   # Right Cauchy-Green tensor
​

Now suppose I have a constant vector n=(1,1,1) and I want to have the invariant I4=n.C.n^T (for ^T meaning the transpose and the multiplication being the usual row by column matrix multiplication). I tried the code:
n = as_vector([1.0, 1.0, 1.0])
I4 = dot(n, dot(n, C))​

but it is always zero. Am I doing the multiplication wrong?
Community: FEniCS Project

2 Answers


2
6 months ago by
n = Constant((1.0, 1.0, 1.0))
I4 = dot(n, dot(n, C.T))


I have not tried it but this is the formula you describe.

Shouldnt this be..
n = Constant((1.0, 1.0, 1.0))
I4 = dot(n, dot(C, n.T))
written 6 months ago by Ovais  
Thank you so much. It worked. You are right about C being transposed since I had conjugated the dot product. Also I defined n using the constant command. Now my code does not converge which is a whole other problem. But thank you.
written 6 months ago by Navid Mirzaei  
0
6 months ago by
Ovais  
May be this support material on ufl can help you : https://aspace.repository.cam.ac.uk/handle/1810/243981?show=full.
There is demo on products and indexing which may help you. Also consider  consulting  Unified Form Language (UFL) Documentation (Release 2017.2.0.dev0)
Thank you so much. But I had checked them before hand and they made me confused.
written 6 months ago by Navid Mirzaei  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »