如何在 FiPy 中求解 Euler-Bernoulli 梁方程？

Question

要了解 FiPy 的工作原理，我想用固定端点解决 Euler–Bernoulli beam equation：

w''''(x) = q(x,t),    w(0) = w(1) = 0,  w'(0) = w'(1) = 0.

为简单起见，让q(x,t) = sin(x)。

如何在 FiPy 中定义和解决它？如何针对方程中唯一的自变量指定源项sin(x)？

from fipy import CellVariable, Grid1D, DiffusionTerm, ExplicitDiffusionTerm
from fipy.tools import numerix

nx = 50
dx = 1/nx
mesh = Grid1D(nx=nx, dx=dx)

w = CellVariable(name="deformation",mesh=mesh,value=0.0)

valueLeft = 0.0
valueRight = 0.0

w.constrain(valueLeft, mesh.facesLeft)
w.constrain(valueRight, mesh.facesRight)
w.faceGrad.constrain(valueLeft, mesh.facesLeft)
w.faceGrad.constrain(valueRight, mesh.facesRight)

# does not work:
eqX = DiffusionTerm((1.0, 1.0)) == numerix.sin(x)
eqX.solve(var=w)

Answer 1

这似乎是您问题的有效版本

from fipy import CellVariable, Grid1D, DiffusionTerm
from fipy.tools import numerix
from fipy.solvers.pysparse.linearPCGSolver import LinearPCGSolver
from fipy import Viewer
import numpy as np


L = 1.
nx = 500
dx = L / nx
mesh = Grid1D(nx=nx, dx=dx)

w = CellVariable(name="deformation",mesh=mesh,value=0.0)

valueLeft = 0.0
valueRight = 0.0

w.constrain(valueLeft, mesh.facesLeft)
w.constrain(valueRight, mesh.facesRight)
w.faceGrad.constrain(valueLeft, mesh.facesLeft)
w.faceGrad.constrain(valueRight, mesh.facesRight)

x = mesh.x
k_0 = 0
k_1 = -1
k_2 = 2 + np.cos(L) - 3 * np.sin(L)
k_3 = -1 + 2 * np.sin(L) - np.cos(L)
w_analytical = numerix.sin(x) + k_3 * x**3 + k_2 * x**2 + k_1 * x + k_0
w_analytical.name = 'analytical'

# does not work:
eqX = DiffusionTerm((1.0, 1.0)) == numerix.sin(x)
eqX.solve(var=w, solver=LinearPCGSolver(iterations=20000))
Viewer([w_analytical, w]).plot()
raw_input('stopped')

运行后，FiPy 的解似乎与解析结果非常接近。

OP 实现的两个重要变化。

• 使用 mesh.x 这是引用空间变量以用于 FiPy 方程的正确方法。

• 指定求解器和迭代次数。这个问题似乎收敛缓慢，因此需要大量迭代。根据我的经验，四阶空间方程通常需要良好的预处理器才能快速收敛。您可以尝试将 Trilinos 求解器包与 Fipy 一起使用，以使这项工作更好，因为它具有更广泛的可用预处理器。

• 为 L 使用显式浮点数以避免 Python 2.7 中的整数数学运算（从评论编辑）