[Tutor] newton's method for system of nonlinear equations
kgotlelelok at galmail.co.za
kgotlelelok at galmail.co.za
Tue Jun 30 12:05:59 CEST 2009
Hi can someone help me with my program demostrating Newton's method for
system of nonlinear equations below. Line 45 gives an error saying return
is outside the function and I don't know how to rectify that. I neewd this
program to print differnt values of x until the differnce between previous
x and present x is less than tol.
Please help
Thank you
from scipy import*
tol=0.01
def swapRows(v,i,j):
if len(v.getshape()) == 1: v[i],v[j] = v[j],v[i]
else:
temp = v[i].copy()
v[i] = v[j]
v[j] = temp
def swapCols(v,i,j):
temp = v[:,j].copy()
v[:,j] = v[:,i]
v[:,i] = temp
def gaussPivot(a,b,tol=1.0e-12):
n = len(b)
# Set up scale factors
s = zeros((n),float)
for i in range(n):
s[i] = max(abs(a[i,:]))
for k in range(0,n-1):
# Row interchange, if needed
p = int(argmax(abs(a[k:n,k])/s[k:n])) + k
if abs(a[p,k]) < tol: error.err('Matrix is singular')
if p != k:
swap.swapRows(b,k,p)
swap.swapRows(s,k,p)
swap.swapRows(a,k,p)
# Elimination
for i in range(k+1,n):
if a[i,k] != 0.0:
lam = a[i,k]/a[k,k]
a[i,k+1:n] = a [i,k+1:n] - lam*a[k,k+1:n]
b[i] = b[i] - lam*b[k]
if abs(a[n-1,n-1]) < tol: error.err('Matrix is singular')
for k in range(n-1,-1,-1):
b[k] = (b[k] - dot(a[k,k+1:n],b[k+1:n]))/a[k,k]
return b
def newtonRaphson2(f,x,tol=1.0e-9):
def jacobian(f,x):
h = 1.0e-4
n = len(x)
jac = zeros((n,n),float)
f0 = f(x)
for i in range(n):
temp = x[i]
x[i] = temp + h
f1 = f(x)
x[i] = temp
jac[:,i] = (f1 - f0)/h
return jac,f0
for i in range(30):
jac,f0 = jacobian(f,x)
if sqrt(dot(f0,f0)/len(x)) < tol: return x
dx = gaussPivot(jac,-f0)
x = x + dx
if sqrt(dot(dx,dx)) < tol*max(max(abs(x)),1.0): return x
print 'Too many iterations'
def f(x):
f = zeros((len(x)),float)
f[0] = sin(x[0]) + x[1]**2 + log(x[2]) - 7.0
f[1] = 3.0*x[0] + 2.0**x[1] - x[2]**3 + 1.0
f[2] = x[0] + x[1] + x[2] - 5.0
return f
x = array([1.0, 1.0, 1.0])
print newtonRaphson2(f,x)
raw_input ("\nPress return to exit")
-----------------------------------------
This email was sent using SquirrelMail.
"Webmail for nuts!"
http://squirrelmail.org/
More information about the Tutor
mailing list