我的函数没有遍历循环,找不到问题



我正试图创建一个函数,为多个变量(dx=-Jinv^-1*F(进行Newton-Raphson寻根,但它不会连续迭代一个循环(即输出不会成为下一次迭代的下一个输入(。有人看到这个问题了吗?

PS。我试图在不使用任何内置矩阵函数的情况下对其进行编码,因此这些解决方案将不适用于此。

# Multivariable Root Finding via Newton Raphson

# import necessary functions
import numpy as np
# define array function for f1 and f2
x = np.arange(-3,3,.01)
def f1(x):
return 2*(x[0]**3) - 6*(x[0]*(x[1]**2)) - 1
def f2(x):
return -2*(x[1]**3) + 6*(x[1]*(x[0]**2)) + 3
# define matrix F(x) 2x1
def F(x_vec): 
return [f1(x_vec),f2(x_vec)]
#print(F([2,3]))

# Create Inverse Jacobian 2x2
# InvJacobian Upper Left
def JUL(x):
return (6*x[0]**2 - 6*x[1]**2)/((6*x[0]**2 - 6*x[1]**2)**2 + 144*(x[0]**2)*(x[1]**2))
# InvJacobian Upper Right
def JUR(x):
return (12*x[0]*x[1])/((6*x[0]**2 - 6*x[1]**2)**2 + 144*(x[0]**2)*(x[1]**2))
# InvJacobian Lower Left
def JLL(x):
return (-12*x[0]*x[1])/((6*x[0]**2 - 6*x[1]**2)**2 + 144*(x[0]**2)*(x[1]**2))
# InvJacobian Lower Left
def JLR(x):
return (6*x[0]**2 - 6*x[1]**2)/((6*x[0]**2 - 6*x[1]**2)**2 + 144*(x[0]**2)*(x[1]**2))
# Combine all Jacobian into Matrix
def Jinv(x_vec):
return [JUL(x_vec),JUR(x_vec),JLL(x_vec),JLR(x_vec)]

# Create Newton Raphson Functon
def rf_newton2d(F_system, Jinv_system, x_vec0, tol, maxiter):
# let ou be the starting array x_vec0
ou = x_vec0
i = 0
err = 10**-4

# F*Jinv should output a 2x1 matrix
while err > tol and i <= maxiter:
F_system = F(ou)
Jinv_system = Jinv(ou)
# w is the [0] position of the 2x1 output matrix
w = F_system[0]*Jinv_system[0] + F_system[1]*Jinv_system[1]

# z is the [1] position of the 2x1 output
z = F_system[0]*Jinv_system[2] + F_system[1]*Jinv_system[3]

# u is the output matrix of F*Jinv
u = [w,z]

# out is the new array involving u + ou
out = [ou[0] + u[0],ou[1] +u[1]]

#define nu as out
nu = out

# let nu be the ou for the next iteration
ou = nu

#run continuous iterations
i += 1

return nu
print(rf_newton2d(F([x_vec0]),Jinv([x_vec0]),[2,3],10**-5, 5))
函数中的

return语句会立即停止该函数,并从调用该函数的位置继续运行主代码。在我看来,while循环中的return语句不应该在while循环中。

在最后一行中,您有函数F([x_vec0]),这意味着函数F中有变量x_vec = [x_vec0]。然后将此变量传递给函数f1,因此在f1中为x = [x_vec0]。在x[1]内可以在f1内找到。问题是x永远不会有第二个条目。如果x_vec0 = [0,1](或一些其他数组对象,如np.ndarray(,则不存在[x_vec0] = [[0,1]],即[x_vec0][0] = [0,1][x_vec0][1]

我从回溯中得到了所有信息(显然我必须为x_vec0创造一个值(:

---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
<ipython-input-4-d57a9d10dd94> in <module>
80         return nu
81 x_vec0 = np.array([1,2])
---> 82 print(rf_newton2d(F([x_vec0]),Jinv([x_vec0]),[2,3],10**-5, 5))
<ipython-input-4-d57a9d10dd94> in F(x_vec)
16 # define matrix F(x) 2x1
17 def F(x_vec):
---> 18     return [f1(x_vec),f2(x_vec)]
19 
20 #print(F([2,3]))
<ipython-input-4-d57a9d10dd94> in f1(x)
9 
10 def f1(x):
---> 11     return 2*(x[0]**3) - 6*(x[0]*(x[1]**2)) - 1
12 
13 def f2(x):
IndexError: list index out of range

如果不知道牛顿-拉斐森根的数学原理,我就无法提出解决方案。

相关内容

  • 没有找到相关文章

最新更新