# 关于ti.sym_eig

``````上下cauchy_stress = ti.Matrix([[0.014860193, 0.0], [0.0, 0.014860193]])值相同
``````

# 直接给出

``````    cauchy_stress = ti.Matrix([[0.014860193, 0.0], [0.0, 0.014860193]])
# 计算特征值和特征矩阵
if all(cauchy_stress.transpose() == cauchy_stress):
sym_eigT1, sym_eigT2 = ti.sym_eig(cauchy_stress, ti.f32)
eigen_values = sym_eigT1
eigen_vector = sym_eigT2
print("cauchy_stress:", cauchy_stress, "eig1:", sym_eigT1, "eig2:", sym_eigT2)
else:
T1,T2 = ti.eig(cauchy_stress, ti.f32)
``````

cauchy_stress矩阵是相同的值，想进行特征分解，给定矩阵值能正确分解。

# 不直接给出，

``````    # cauchy_stress = ti.Matrix([[0.014860193, 0.0], [0.0, 0.014860193]])
# 计算特征值和特征矩阵
if all(cauchy_stress.transpose() == cauchy_stress):
sym_eigT1, sym_eigT2 = ti.sym_eig(cauchy_stress, ti.f32)
eigen_values = sym_eigT1
eigen_vector = sym_eigT2
print("cauchy_stress:", cauchy_stress, "eig1:", sym_eigT1, "eig2:", sym_eigT2)
else:
T1,T2 = ti.eig(cauchy_stress, ti.f32)``````

``````   @ti.kernel
def update_d():
# d映射到网格  d：P2G
map_d_to_grid()
# 计算 d 的拉普拉斯算子
central_laplacian()

for p in range(particle_num[None]):
xp = particles[p].position / dx
base = int(xp - 0.5)
fx = xp - base
w = [0.5 * (1.5 - fx) ** 2, 0.75 - (fx - 1) ** 2, 0.5 * (fx - 0.5) ** 2]
# 在网格上计算d的拉普拉斯算子之后，映射回粒子，得到粒子上的拉普拉斯算子 p_laplacian_d
p_laplacian_d = 0.0
for i, j in ti.static(ti.ndrange(3, 3)):  # loop over 3x3 grid node neighborhood
g_laplacian_d = node[base + ti.Vector([i, j])].laplacian_d
weight = w[i][0] * w[j][1]
p_laplacian_d += weight * g_laplacian_d

# 计算阻抗  compute geometric resistance:Dc
Dc = particles[p].d - l0 ** 2 * p_laplacian_d
# print("Dc:", Dc)

# *********计算未退化(g(d) = 1)的柯西应力  compute un-dearaded Cauchy stress***********
# ++++++++ QR分解和符号规定 +++++++++++
Q, R = QR2(particles[p].F)
# sign convention  符号规定
# 之前直接Q[0, 0]报错，取不到值
r00 = R[0, 0]
if r00 < 0:
Q *= -1
R *= -1

# 计算 J
particles[p].J = 1.0
for r in ti.static(range(2)):
particles[p].J *= R[r, r]

# 计算 stress
kx = 0 * mu_0
pk1 = calculate_pk1(Q, R, particles[p].J, kx, 1.0)
# P = Jσ（F-T）（F的逆的转置）     σ = PF/J
cauchy_stress = ti.Matrix.zero(float, 2, 2)
cauchy_stress = pk1 @ particles[p].F / particles[p].J
# print(cauchy_stress)
# ******************************

# 计算sigma plus
sigma_plus = ti.Matrix.zero(float, 2, 2)
T1 = ti.Matrix.zero(float, 2, 2)
T2 = ti.Matrix.zero(float, 4, 2)
sym_eigT1 = ti.Vector.zero(float, 2)
sym_eigT2 = ti.Matrix.zero(float, 2, 2)

# 直接给出
# cauchy_stress = ti.Matrix([[0.014860193, 0.0], [0.0, 0.014860193]])
# 计算特征值和特征矩阵
if all(cauchy_stress.transpose() == cauchy_stress):
sym_eigT1, sym_eigT2 = ti.sym_eig(cauchy_stress, ti.f32)
eigen_values = sym_eigT1
eigen_vector = sym_eigT2
print("cauchy_stress:", cauchy_stress, "eig1:", sym_eigT1, "eig2:", sym_eigT2)
else:
T1,T2 = ti.eig(cauchy_stress, ti.f32)
eigen_values = T1
eigen_vector = T2
``````

Hi @wjy. 欢迎来到太极论坛。

PS：最好能简化一下代码，到最小的复现你问题的程度。目前来看代码有点多，开发者看起来会很费劲。