1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134
|
#!/usr/bin/env python
#coding=utf8
from __future__ import print_function
import sys
import math
try:
import numpy
except ImportError:
print("WARNING: This test requires NumPy: http://http://www.numpy.org/")
sys.exit(0)
import vtk
# Create a simple hexahedron
points = vtk.vtkPoints()
points.SetNumberOfPoints(8)
points.InsertPoint(0, 0, 0, 0)
points.InsertPoint(1, 1, 0, 0)
points.InsertPoint(2, 1, 1, 0)
points.InsertPoint(3, 0, 1, 0)
points.InsertPoint(4, 0, 0, 1)
points.InsertPoint(5, 1, 0, 1)
points.InsertPoint(6, 1, 1, 1)
points.InsertPoint(7, 0, 1, 1)
#for k_point in xrange(8): print points.GetPoint(k_point)
hexahedron = vtk.vtkHexahedron()
for k_point in xrange(8):
hexahedron.GetPointIds().SetId(k_point, k_point)
cell_array = vtk.vtkCellArray()
cell_array.InsertNextCell(hexahedron)
farray_disp = vtk.vtkFloatArray()
farray_disp.SetName("displacement")
farray_disp.SetNumberOfComponents(3)
farray_disp.SetNumberOfTuples(8)
ugrid_hex = vtk.vtkUnstructuredGrid()
ugrid_hex.SetPoints(points)
ugrid_hex.SetCells(hexahedron.GetCellType(), cell_array)
ugrid_hex.GetPointData().AddArray(farray_disp)
ugrid_hex.GetPointData().SetActiveVectors("displacement")
# Test a few deformation gradients
F_lst = []
# Simple extentions
F_lst += [numpy.array([[2,0,0],
[0,1,0],
[0,0,1]])]
F_lst += [numpy.array([[1,0,0],
[0,2,0],
[0,0,1]])]
F_lst += [numpy.array([[1,0,0],
[0,1,0],
[0,0,2]])]
# Simple shears
F_lst += [numpy.array([[1,1,0],
[0,1,0],
[0,0,1]])]
F_lst += [numpy.array([[1,0,1],
[0,1,0],
[0,0,1]])]
F_lst += [numpy.array([[1,0,0],
[0,1,1],
[0,0,1]])]
F_lst += [numpy.array([[1,0,0],
[1,1,0],
[0,0,1]])]
F_lst += [numpy.array([[1,0,0],
[0,1,0],
[1,0,1]])]
F_lst += [numpy.array([[1,0,0],
[0,1,0],
[0,1,1]])]
# Rotation
alpha = math.pi/2
F_lst += [numpy.array([[+math.cos(alpha),-math.sin(alpha),0],
[+math.sin(alpha),+math.cos(alpha),0],
[ 0, 0,1]])]
for F in F_lst:
print("F = " + str(F))
# Assuming an homogeneous deformation gradient, compute the displacement and displacement gradient: u = x(X)-X = F.X-X = (F-1).X = GU.X with GU = F-1
GU = F - numpy.eye(3)
print("GU = " + str(GU))
for k_point in xrange(8):
farray_disp.SetTuple(k_point, numpy.dot(GU, points.GetPoint(k_point)))
#print farray_disp.GetTuple(k_point)
# Compute the small (linearized) strain tensor: e = (GU)_sym = (GU+GU^T)/2
e = (GU + numpy.transpose(GU))/2
print("e = " + str(e))
# Compute the Green-Lagrange strain tensor: E = (F.F^T-1)/2 = (GU+GU^T+GU^T.GU)/2 = e + GU^T.GU/2
E = e + numpy.dot(numpy.transpose(GU), GU)/2
print("E = " + str(E))
# Compute displacement gradient using vtkCellDerivatives
cell_derivatives = vtk.vtkCellDerivatives()
cell_derivatives.SetVectorModeToPassVectors()
cell_derivatives.SetTensorModeToComputeGradient()
cell_derivatives.SetInputData(ugrid_hex)
cell_derivatives.Update()
vector_gradient = numpy.reshape(cell_derivatives.GetOutput().GetCellData().GetArray("VectorGradient").GetTuple(0), (3,3))
print("VectorGradient = " + str(vector_gradient))
assert numpy.allclose(vector_gradient, GU)
# Compute small strain tensor using vtkCellDerivatives
cell_derivatives = vtk.vtkCellDerivatives()
cell_derivatives.SetVectorModeToPassVectors()
cell_derivatives.SetTensorModeToComputeStrain()
cell_derivatives.SetInputData(ugrid_hex)
cell_derivatives.Update()
strain = numpy.reshape(cell_derivatives.GetOutput().GetCellData().GetArray("Strain").GetTuple(0), (3,3))
print("Strain = " + str(strain))
assert numpy.allclose(strain, e)
# Compute Green-Lagrange strain tensor using vtkCellDerivatives
cell_derivatives = vtk.vtkCellDerivatives()
cell_derivatives.SetVectorModeToPassVectors()
cell_derivatives.SetTensorModeToComputeGreenLagrangeStrain()
cell_derivatives.SetInputData(ugrid_hex)
cell_derivatives.Update()
green_lagrange_strain = numpy.reshape(cell_derivatives.GetOutput().GetCellData().GetArray("GreenLagrangeStrain").GetTuple(0), (3,3))
print("GreenLagrangeStrain = " + str(green_lagrange_strain))
assert numpy.allclose(green_lagrange_strain, E)
|