在 FiPy 中使用扫描函数时的求解器容差和残差

Question

当我意识到命令 sweep 没有按照我想象的方式工作时，我正试图使用 FiPy 求解一组 PDE。这是一个包含我的部分代码的示例：

from pylab import *
import sys
from fipy import *

viscosity = 5.55555555556e-06 

Pe =5.

pfi=100.
lfi=0.01

Ly=1.
Nx =200
Ny=100
Lx=Ly*Nx/Ny
dL=Ly/Ny
mesh = PeriodicGrid2DTopBottom(nx=Nx, ny=Ny, dx=dL, dy=dL)

x, y = mesh.cellCenters

xVelocity = CellVariable(mesh=mesh, hasOld=True,  name='X velocity')

xVelocity.constrain(Pe, mesh.facesLeft)
xVelocity.constrain(Pe, mesh.facesRight)

rad=0.1

var1 = DistanceVariable(name='distance to center', mesh=mesh, value=numerix.sqrt((x-Nx*dL/2.)**2+(y-Ny*dL/2.)**2))

pi_fi= CellVariable(mesh=mesh, value=0.,name='Fluid-interface energy map')
pi_fi.setValue(pfi*exp(-1.*(var1-rad)/lfi), where=(var1 > rad) )
pi_fi.setValue(pfi, where=(var1 <= rad))

xVelocityEq = DiffusionTerm(coeff=viscosity) - ImplicitSourceTerm(pi_fi)

xres=10.
while (xres > 1.e-6) :
        xVelocity.updateOld()
        mySolver = LinearGMRESSolver(iterations=1000,tolerance=1.e-6)
        xres = xVelocityEq.sweep(var=xVelocity,solver=mySolver)
        print 'Result = ', xres
#Thats it

简而言之，我正在声明一个名为 xVelocityEq 的函数，并使用 sweep 求解它。这是我的输出：

Result =  0.0007856742013190237
Result =  6.414470433257661e-07

如您所见，while 循环在两次迭代后结束。我的第一个问题是：为什么我的第一个残差 (=0.0007856742013190237) 高于求解器的容差？我认为，由于 xVelocityEq 对应于一个线性系统，因此求解器公差和残差意味着相同的事情。

如果我增加编号。 mySolver 从 1000 到 10000 的迭代次数，我得到以下输出：

Result =  0.0007856742013190237
Result =  2.4619110931978988e-09

既然第一个残差保持不变，为什么第二个残差发生变化？

如果我将 mySolver 中的公差从 1.e-6 增加到 7.e-4，我会得到以下输出：

Result =  0.0007856742013190237
Result =  6.414470433257661e-07

请注意，这些残差与第一个输出中的残差相同。现在，如果我尝试将公差进一步增加到 8.e-4，这就是我得到的输出：

Result =  0.0007856742013190237
Result =  0.0007856742013190237
Result =  0.0007856742013190237
Result =  0.0007856742013190237
Result =  0.0007856742013190237
...

此时我完全迷失了。 为什么对于小于 7.e-4 的所有求解器公差，残差具有相同的值？以及为什么这些残差是常数且等于 0.0007856742013190237 对于高于 7.e-4 的求解器公差？

如果我将 mySolver 更改为 LinearLUSolver（迭代次数=1000，公差=1.e-6），这就是我得到的结果：

Result =  0.0007856742013190237
Result =  1.6772757200988522e-18

为什么我的第一个残差和以前一样，即使我改变了求解器？

Answer 1

why is my first residual error (=0.0007856742013190237) higher than the solver's tolerance?

.sweep()计算的残差在计算之前调用求解器计算新的解向量。矩阵L和右侧向量b是根据解向量x[=的初始值计算得到的43=]。

残差衡量当前解向量满足非线性 PDE 的程度。求解器公差限制了求解器满足从 PDE 离散化的线性方程组的努力程度。

即使 PDE 是线性的（例如，扩散系数不是解变量的函数），初始值也可能无法求解 PDE，因此残差很大。调用求解器后，x 应求解 PDE，使其在求解器容差范围内。如果 PDE 是非线性的，那么线性代数的收敛性良好的解可能仍然不是 PDE 的好解；这就是清扫的目的。

I thought that, since xVelocityEq corresponds to a linear system, solver tolerance and residual error would mean the same thing.

跟踪两者没有任何用处。除了求解之前的残差和用于终止求解的求解器公差之外，还可以使用不同的归一化，并且许多求解器文档可能有点粗略。 FiPy 使用 |L x - b|_2 作为残差。求解器可能会根据 b 的大小、L 的对角线或月相进行归一化，所有这些都会使其难以直接计算将残差与公差进行比较。

Why did the second residual change, given that the first remained the same?

通过允许 1000 次而不是 100 次迭代，求解器能够达到更严格的公差，这反过来又导致下一次扫描的残差更小。

Why the residuals have the same values for all solver tolerances smaller than 7.e-4? And why these residuals are constant and equal to 0.0007856742013190237 for solver tolerances higher than 7.e-4?

可能是因为解算器失败，所以没有改变解向量的值。一些求解器不报告这一点。在其他情况下，我们应该更好地向您报告该事实。

Why in the world is my first residual the same as before, even though I have changed the solver?

残差不是求解器的属性。它是近似 PDE 的离散化方程组的属性。然后将这些线性代数方程作为求解器的输入。

在 FiPy 中使用扫描函数时的求解器容差和残差

Solver tolerance and residual error when using sweep function in FiPy

python

scipy

fipy