受约束的 Newton-Raphson 根发现

计算科学 牛顿法
2021-12-03 11:04:41

原始问题

我有一组非线性方程,我需要找到我的解向量的子集被限制为大于或等于 0 的根。我已经实现了 Newton-Raphson 算法,但我发现其中一些我需要保持积极的数量正在变为负数。我熟悉使用拉格朗日乘数在优化问题中强制执行约束,但我不确定如何为求根做同样的事情。

我突然想到,我可以将必须​​为正的变量映射到对数空间(即求解log(xi)而不是xi) 但由于零是一个重要的情况,我担心这可能会导致数值问题。

先前的答案看来,内点法可能是合适的,但它会涉及修改我的求解器,而不仅仅是我想避免的残差方程。

进一步说明

这些方程是由包含不可恢复变形的材料模型产生的一系列非线性方程(Ep)。通常实现这一点的方式是通过内部状态变量(ξ) 由一些演化方程控制,这些方程由一些演化速率 (γ˙)。

演化方程采用以下形式:

Ep˙=E˙(Ep,ξ,γ˙)

ξ˙=ξ˙(ξ,γ˙)

受限于一些额外的起始条件(屈服函数):

F=F(EEp)

在哪里E是总变形,它是一个 <= 0 的函数。该问题符合 Kuhn-Tucker 条件,可以概括为:

γ˙F=0

暗示进化率为零,如果F<0γ˙>0如果F=0.

我的未知向量是这样的:

x=[Ep,ξ,γ]

我的残差是这样的:

R=[Ep,expectedEp,ξexpectedξ,γ˙F+F]

标记的值在哪里()expected是进化方程的结果和是麦考利括号,定义为:

x=12(x+abs(x))

我包括了麦考利项,否则,F > 0 的情况和γ˙=0将被接受这是不正确的。我尝试在残差计算中使用 if 语句,但这会导致其他问题,而且这似乎表现得更好。

1个回答

我想我现在会回到这个问题,因为我有一个有效的答案。这个问题是多方面的,在最终确定解决方案之前,我尝试了几种不同程度的解决方案。

我尝试过但无法自行解决的事情:

  1. 试图实现障碍的同伦求解器

    g(x,s)=(11a(s))R(x)+1a(s)b(x,s)

    在哪里

    b(x,s)=esa(bx)1

    a(s)=eAs

    在哪里R是原始残差,s, 是伪时间,并且A是一个足够大的值以确保屏障成立。我发现如果 A 至少为 10 左右,那么从障碍函数到真实函数的变化非常平滑。这适用于玩具问题,但不适用于我的实际问题。

  2. 牛顿同伦求解器:

    g(x,s)=R(x)+(1s)R(x0)

    我喜欢这个同伦并最终将它用于我的最终非线性方程求解。在解决方案中,我首先尝试 s = 1,然后根据需要进行缩减。

  3. 执行多个嵌套的 Newton-Raphson 求解。而不是试图解决Ep,ξ, 和γ˙同时我解决了ξ分析和Ep数值上对于给定的值γ˙作为每次迭代求解γ˙. 这极大地帮助了我解决许多问题,但并没有解决更困难的情况。

  4. 通过我的代码工作,发现了一个可能导致意外行为的错误。同样,这有帮助,但没有解决我的问题。

  5. 尝试将 Newton Raphson 中的线搜索限制为始终位于正确的域中。这有帮助,但有时不允许收敛。

现在实际工作:

我无意中设置了我的残差方程,如果误差太大,它会进入一个新的吸引力盆地,将它拉向错误的方向。我需要修改我的残差方程以防止这种行为。请注意,在原始残差中γ˙

Rγ˙=F+γ˙F

问题在于,如果F变得非常大,最终雅可比行列式可能会因为第二项将尝试驱动而变为正数γ˙消极的。所以最后我需要:

  1. 使用牛顿同伦求解器(现在不一定需要,但在问题变得僵硬时会有所帮助)。
  2. 将整体解决方案分解为子问题。
  3. 修改我的残差以排除有问题的行为。

在原始方程的上下文中,我认为 3 可能是关键点。