newton raphson算法
哪位能提供一下newton raphson算法的程序阿，求解非线性方程组的算法，我在网上找不到。 一个matlab实例% This m-file takes care of synthetic division.
% By giving one polynomial and one root this function returns
% the polynomial formed with the other roots of the given polynomial excluding the given root.
% Keerthi Venkateswara Rao-M.Tech-NITC-INDIA
function coeff_second=syn_division(coeff_function,fun_root_new)
order_fun=size((coeff_function),2);
coeff_second=0;
for index=1:size((coeff_function),2)-1
if index==1
coeff_second(index)=coeff_function(index);
else
coeff_second(index)=coeff_function(index)+fun_root_new*coeff_second(index-1);
end
end
% This m-file calculates the derivative of the function, the limitation of
% this function is, it can calculate only the derivatives of power(x,n)....
% Keerthi Venkateswara Rao-M.tech-NITC-INDIA.
function coeff_derivative=derivate(coeff_function)
der_order=size((coeff_function),2)-1;
coeff_derivative=0;
for index=1:size((coeff_function),2)-1
coeff_derivative(index)=der_order*coeff_function(index);
der_order=der_order-1;
end
% this m-file calculates the real roots of the given polynomial using
% newton raphson technique.this m-file calls the functions in the two m-files named as syn_division and
% derivate.
% coeff_function is the 1xn matrix conatining the coeff of the polynomial.
% Keerthi Venkateswara Rao-M.Tech-NITC-INDIA
function = newton(coeff_function,initial_guess,error_tolerance,max_iterations)
iterations=0;
max_no_roots=size(coeff_function,2);
final_roots=0;
functionvalue=0;
for no_roots=1:max_no_roots-1
fun_root_new=initial_guess;
flag=1;
coeff_der_function=derivate(coeff_function);
order_fun=size(coeff_function,2);
order_der_fun=size(coeff_der_function,2);
while flag==1
fun_root_old=fun_root_new;
fx=0;
dfx=0;
nonzero=1;
while nonzero==1
powers=order_fun-1;
for index=1:order_fun
fx=fx+coeff_function(index)*fun_root_old^powers;
powers=powers-1;
end
powers=order_der_fun-1;
for index=1:order_der_fun
dfx=dfx+coeff_der_function(index)*fun_root_old^powers;
powers=powers-1;
end
if dfx==0
fun_root_old=fun_root_old+1;
else
nonzero=0;
end
end
iterations = iterations + 1;
fun_root_new = fun_root_old - fx/dfx;
if iterations >= max_iterations
flag=0;
elseifabs(fun_root_new-fun_root_old)<=error_tolerance
flag=0;
final_roots(no_roots)=fun_root_new;
functionvalue(no_roots)=fx;
end
end
coeff_function=syn_division(coeff_function,fun_root_new);
end c语言实例
#include <math.h>
#include <iomanip.h>
#include <iostream.h>
#include <process.h>
void usrfun(double x,double alpha,double beta)
{
int np;
np = 15;
alpha = -2.0 * x;
alpha = -2.0 * x;
alpha = -2.0 * x;
alpha = 1.0;
alpha = 2.0 * x;
alpha = 2.0 * x;
alpha = 2.0 * x;
alpha = 2.0 * x;
alpha = 1.0;
alpha = -1.0;
alpha = 0.0;
alpha = 0.0;
alpha = 0.0;
alpha = 1.0;
alpha = -1.0;
alpha = 0.0;
beta= x*x+x*x+x*x-x;
beta=-x*x-x*x-x*x-x*x+1.0;
beta=-x+x;
beta=-x+x;
}
void main()
{
//program d10r13
//driver for routine mnewt
int i,j,kk,k,ntrial,np,n;
double tolx,tolf,xx;
ntrial = 5;
tolx = 0.000001;
n = 4;
tolf = 0.000001;
np = 15;
double x,alpha,beta;
for (kk = -1; kk<= 1; kk=kk+2)
{
for (k =1; k<=3; k++)
{
xx = 0.2 * k * kk;
cout<< "Starting vector number "<<k<<endl;
for (i = 1; i<=4; i++)
{
x = xx + 0.2 * i;
cout<<"x("<<i<<")= "<<x<<endl;
}
for (j = 1;j <=ntrial; j++)
{
mnewt(1,x,n,tolx,tolf);
usrfun(x,alpha,beta);
cout<<"i x(i) f"<<endl;
for (i = 1; i<=n; i++)
{
cout<< setw(3)<< i;
cout<< setw(14)<< x;
cout<< setw(16)<< -beta<<endl;
}
}
}
}
}
void mnewt(int ntrial,double x,int n,double tolx,double tolf)
{
double alpha,beta;
int k,i,indx;
double errf,d,errx;
for (k = 1; k<=ntrial; k++)
{
usrfun(x,alpha,beta);
errf = 0.0;
for (i = 1; i<=n; i++)
errf = errf + fabs(beta);
if (errf <= tolf)_c_exit();
ludcmp(alpha,n,indx,d);
lubksb(alpha,n,indx,beta);
errx = 0.0;
for (i = 1; i<=n; i++)
{
errx = errx + fabs(beta);
x = x + beta;
}
if (errx <= tolx) _c_exit();
}
} 好资料!
2F的原始来源连接
http://www.mathworks.com/matlabcentral/fileexchange/4313 谢谢分享，已经有些概念了 这也太复杂了吧，我的求解牛顿迭代的时候出现了复数，是初始值不对还是精度不对啊 的确复杂了，有些还是C++
页:
[1]