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/
_______________________________________________
Tutor maillist - [email protected]
http://mail.python.org/mailman/listinfo/tutor
