FEniCS热弹性力学模拟遇ArityMismatch报错的解决咨询

我来帮你解决这个问题——你遇到的ArityMismatch错误其实是变分形式的结构问题，结合代码里的小拼写错误导致的，下面一步步说明修复方案：

错误原因分析 #

你看到的Adding expressions with non-matching form arguments () vs ('v_1',)错误，核心原因是你把热弹性的常数应力项错误地放到了左侧的双线性形式中：
FEniCS要求变分问题的左侧a必须是双线性形式（同时依赖试函数u和测试函数v），而右侧L必须是线性形式（只依赖测试函数v）。

当你在sigma(u)中加入-Kt*Identity(d)时，inner(sigma(u), epsilon(v))*dx会被拆成两部分：

• 双线性项：来自和u相关的弹性部分，依赖u和v；
• 线性项：来自常数Kt的部分，只依赖v。

把这两类形式混合到左侧的a中，就触发了参数不匹配的错误。另外你代码里还有个小拼写错误：Kt定义里的alph应该是alpha。

修复后的完整代码 #

from dolfin import *
from fenics import *
from ufl import nabla_div
import matplotlib.pyplot as plt
import numpy as np

# 材料参数
E = Constant(100*10**9)
nu = Constant(0.3)
rho = Constant(2200)
lambda_ = (E*nu)/((1+nu)*(1-2*nu))
mu = E/(2+2*nu)
alpha = 1.2*(10**(-8))
deltaT = Constant(50)
# 修正拼写错误：alph -> alpha
Kt = E*alpha*deltaT/(1-2*nu)
g = 9.81
tol = 1E-14

# 创建网格和函数空间
mesh = RectangleMesh(Point(-2., 0.),Point(2., 10.), 80, 200)
V = VectorFunctionSpace(mesh, "P", 1)
d = V.geometric_dimension()  # 提前定义几何维度，避免作用域问题

# 边界条件定义
def clamped_boundary(x, on_boundary):
    return on_boundary and x[1] < tol

class UpFace(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary and (x[1] > 10 - tol)

ueN = UpFace()
boundaries = MeshFunction("size_t", mesh, mesh.topology().dim()-1, 0)
ueN.mark(boundaries, 1)
ds = Measure("ds", domain=mesh, subdomain_data=boundaries)

bc = DirichletBC(V, Constant((0, 0)), clamped_boundary)

# 应变张量定义
def epsilon(u):
    return 0.5*(nabla_grad(u) + nabla_grad(u).T)

# 修改为纯弹性应力张量（移除热项）
def sigma(u):
    return lambda_*nabla_div(u)*Identity(d) + 2*mu*epsilon(u)

# 变分问题定义
u = TrialFunction(V)
v = TestFunction(V)
f = Constant((0,-rho*g))
T = Constant((0, 0))
Pr = Constant((0, -2*10**9))

# 左侧保持纯双线性形式
a = inner(sigma(u), epsilon(v))*dx
# 右侧添加热弹性等效载荷项
L = dot(f, v)*dx + dot(T, v)*ds + dot(Pr,v)*ds(1) + inner(Kt*Identity(d), epsilon(v))*dx

# 求解并可视化
u = Function(V)
solve(a == L, u, bc)

plot(u, mode="displacement", color= "red")
plt.colorbar(plot(u))
plt.show()

关键修复点 #

1. 修正拼写错误：把Kt定义中的alph改为alpha，避免变量未定义错误。
2. 拆分热弹性项：
   • 让sigma(u)只保留纯弹性的双线性部分，确保左侧a是标准的双线性形式；
   • 将热应力的常数项转化为等效载荷，添加到右侧的L中，保持L为纯线性形式。
3. 提前定义几何维度：把d = V.geometric_dimension()移到函数空间定义后，避免函数调用时的作用域问题。

内容的提问来源于stack exchange，提问作者MEGalex