matlab拉格朗日插值演算法

① 拉格朗日插值的matlab代碼

1、給出一列數據之後，作圖如下：aa= randn(100,1);plot(aa);。

② 在MATLAB編程實驗中，用拉格朗日插值法跟牛頓插值法運行之後計算的結果為什麼是一樣的

根據插值多項式的唯一性，兩種方法的結果應該是一樣的。條條道路通羅馬，只是方法不同而已，牛頓法要比拉格朗日法優越簡單。

Matlab函數M文件Lagrange程序function yy=lagrange(x,y,xi) m=length(x)上面是拉格朗日插值法，其中xi為要計算的數值比如 x=[0 3 5 9 31];Q

clear all;clc

x0=1:5;

y0=sin(x0);

x=1:0.2:2;

y0=lagrange(x0,y0,x)

命令窗口輸這個就沒有問題。

(2)matlab拉格朗日插值演算法擴展閱讀：

如果這特定函數是多項式，就稱它為插值多項式。利用插值基函數很容易得到拉格朗日插值多項式，公式結構緊湊，在理論分析中甚為方便，但當插值節點增減時全部插值基函數均要隨之變化，整個公式也將發生變化，這在實際計算中是很不方便的，為了克服這一缺點，提出了牛頓插值。

③ 怎麼用matlab利用拉格朗日插值計演算法的原理編寫並計算函數所在節點的近似值。 謝謝。

.m文件
function yy=lagrange(x1,y1,xx)
%本程序為Lagrange1插值，其中x1,y1
%為插值節點和節點上的函數值，輸出為插值點xx的函數值，
%xx可以是向量。
syms x
n=length(x1);
for i=1:n
t=x1;t(i)=[];L(i)=prod((x-t)./(x1(i)-t));% L向量用來存放插值基函數
end
u=sum(L.*y1);
p=simplify(u) % p是簡化後的Lagrange插值函數（字元串）
yy=subs(p,x,xx);
clf
plot(x1,y1,'ro',xx,yy,'*')
====================================
x=[ 0.4 0.5 0.6 0.7 0.8];
y=[-0.916291;-0.693147;-0.510826;-0.356675;-0.223144]';
yy=lagrange(x,y,0.54)

p =

- (14363668061545223*x^4)/6755399441055744 + (229230406283396627*x^3)/33776997205278720 - (6086876668119665137*x^2)/675539944105574400 + (23595121244981107513*x)/3377699720527872000 - 186390055565518223/70368744177664000

yy =

-0.6161

④ MATLAB中如何利用拉格朗日插值法作圖

計算衛星定位么？？ 要幾次插值的 給你一次和二次的吧
1.n個節點分段Lagrange插值多項式；
%2.使用格式y=lagrange(x0,y0,x,k)；
%3.輸入項x0為n維插值節點向量，y0為n維被插函數值向量；
%4.x為m維插值點向量，k為分段插值多項式次數，不超過3，預設為k=1；
%5.輸出y為插值點x處的函數值；
%6.本程序於2002.4.21.編寫?
function y=lagrange(x0,y0,x,k)
if nargin<4,k=1;end
if k>3,error('分段次數過高，容易產生Runge現象，請重新選擇次數k.'),end
n=length(x0);
m=length(x);
nn=1;
for i=1:m
u=x(i);
switch k
%---------------------------針對不同的k判斷插值區間
case 1
if u<=x0(2)
t=1;
elseif u>x0(n-1)
t=n-1;
else
for j=nn+1:n-2
if u>x0(j)&u<=x0(j+1)
t=j;nn=t-1;break
end
end
end
%---------------------------
case 2
if u<=(x0(2)+x0(3))/2
t=1;
elseif u>(x0(n-2)+x0(n-1))/2
t=n-2;
else
for j=nn+1:n-3
if u>(x0(j)+x0(j+1))/2 & u<=(x0(j+1)+x0(j+2))/2
t=j;nn=t-1;break
end
end
end
%---------------------------
case 3
if u<=x0(3)
t=1;
elseif u>x0(n-2)
t=n-3;
else
for j=nn+2:n-3
if u>x0(j) & u<=x0(j+1)
t=j-1;nn=t;break
end
end
end
end
%--------------------------------主程序
s=0;
for j=t:t+k
L=1;
for p=t:t+k
if p~=j
L=L*(u-x0(p))/(x0(j)-x0(p));
end
end
s=s+L*y0(j);
end
y(i)=s;
end

分段線性Lagrange插值
% 命令格式：y=lagrange1(x0,y0,x)
% x0為節點向量，y0為對應的函數值向量，
% x為插值點向量，返回值y為x處的函數近似值向量。
function y=lagrange1(x0,y0,x)
[n1,n]=size(x0);[n1,m]=size(x);
for i=1:m
u=x(i);
if u<=x0(2)
y(i)=y0(1)*(u-x0(2))/(x0(1)-x0(2))+...
y0(2)*(u-x0(1))/(x0(2)-x0(1));
elseif u>=x0(n-1)
y(i)=y0(n-1)*(u-x0(n))/(x0(n-1)-x0(n))+...
y0(n)*(u-x0(n-1))/(x0(n)-x0(n-1));
else
for k=2:n-1
if u>=x0(k) & u<=x0(k+1)
y(i)=y0(k)*(u-x0(k+1))/(x0(k)-x0(k+1))+...
y0(k+1)*(u-x0(k))/(x0(k+1)-x0(k));
end
end
end
end
y;

⑤ 大神 求解 拉格朗日插值 matlab法

拉格朗日

function y=lagrange(x0,y0,x)
n=length(x0);m=length(x);
for i=1:m
z=x(i);
s=0.0;
for k=1:n
p=1.0;
for j=1:n
if j~=k
p=p*(z-x0(j))/(x0(k)-x0(j));
end
end
s=p*y0(k)+s;
end
y(i)=s;
end

SOR迭代法的Matlab程序

function [x]=SOR_iterative(A,b)
% 用SOR迭代求解線性方程組,矩陣A是方陣
x0=zeros(1,length(b)); % 賦初值
tol=10^(-2); % 給定誤差界
N=1000; % 給定最大迭代次數
[n,n]=size(A); % 確定矩陣A的階
w=1; % 給定鬆弛因子
k=1;
% 迭代過程
while k<=N
x(1)=(b(1)-A(1,2:n)*x0(2:n)')/A(1,1);
for i=2:n
x(i)=(1-w)*x0(i)+w*(b(i)-A(i,1:i-1)*x(1:i-1)'-A(i,i+1:n)*x0(i+1:n)')/A(i,i);
end
if max(abs(x-x0))<=tol
fid = fopen('SOR_iter_result.txt', 'wt');
fprintf(fid,'\n********用SOR迭代求解線性方程組的輸出結果********\n\n');
fprintf(fid,'迭代次數: %d次\n\n',k);
fprintf(fid,'x的值\n\n');
fprintf(fid, '%12.8f \n', x);
break;
end
k=k+1;
x0=x;
end
if k==N+1
fid = fopen('SOR_iter_result.txt', 'wt');
fprintf(fid,'\n********用SOR迭代求解線性方程組的輸出結果********\n\n');
fprintf(fid,'迭代次數: %d次\n\n',k);
fprintf(fid,'超過最大迭代次數，求解失敗！');
fclose(fid);
end
Matlab中龍格-庫塔(Runge-Kutta)方法原理及實現龍格-庫塔(Runge-Kutta)方法是一種在工程上應用廣泛的高精度單步演算法。由於此演算法精度高，採取措施對誤差進行抑制，所以其實現原理也較復雜。該演算法是構建在數學支持的基礎之上的。龍格庫塔方法的理論基礎來源於泰勒公式和使用斜率近似表達微分，它在積分區間多預計算出幾個點的斜率，然後進行加權平均，用做下一點的依據，從而構造出了精度更高的數值積分計算方法。如果預先求兩個點的斜率就是二階龍格庫塔法，如果預先取四個點就是四階龍格庫塔法。一階常微分方程可以寫作：y'=f(x,y),使用差分概念。
(Yn+1-Yn)/h= f(Xn,Yn)推出（近似等於，極限為Yn'）
Yn+1=Yn+h*f(Xn,Yn)
另外根據微分中值定理，存在0<t<1,使得
Yn+1=Yn+h*f(Xn+th,Y(Xn+th))
這里K＝f(Xn+th,Y(Xn+th))稱為平均斜率，龍格庫塔方法就是求得K的一種演算法。
利用這樣的原理，經過復雜的數學推導（過於繁瑣省略），可以得出截斷誤差為O(h^5)的四階龍格庫塔公式：
K1＝f(Xn,Yn);
K2=f(Xn+h/2,Yn+(h/2)*K1);
K3=f(Xn+h/2,Yn+(h/2)*K2);
K4=f(Xn+h,Yn+h*K3);
Yn+1=Yn+h*(K1+2K2+2K3+K4)*(1/6);
所以，為了更好更准確地把握時間關系，應自己在理解龍格庫塔原理的基礎上，編寫定步長的龍格庫塔函數，經過學習其原理，已經完成了一維的龍格庫塔函數。
仔細思考之後，發現其實如果是需要解多個微分方程組，可以想像成多個微分方程並行進行求解，時間，步長都是共同的，首先把預定的初始值給每個微分方程的第一步，然後每走一步，對多個微分方程共同求解。想通之後發現，整個過程其實很直觀，只是不停的逼近計算罷了。編寫的定步長的龍格庫塔計算函數：
function [x,y]=runge_kutta1(ufunc,y0,h,a,b)%參數表順序依次是微分方程組的函數名稱，初始值向量，步長，時間起點，時間終點（參數形式參考了ode45函數）
n=floor((b-a)/h);%求步數
x(1)=a;%時間起點
y(:,1)=y0;%賦初值，可以是向量，但是要注意維數
for ii=1:n

x(ii+1)=x(ii)+h;

k1=ufunc(x(ii),y(:,ii));

k2=ufunc(x(ii)+h/2,y(:,ii)+h*k1/2);

k3=ufunc(x(ii)+h/2,y(:,ii)+h*k2/2);

k4=ufunc(x(ii)+h,y(:,ii)+h*k3);

y(:,ii+1)=y(:,ii)+h*(k1+2*k2+2*k3+k4)/6;
%按照龍格庫塔方法進行數值求解
end
調用的子函數以及其調用語句：
function dy=test_fun(x,y)
dy = zeros(3,1);%初始化列向量
dy(1) = y(2) * y(3);
dy(2) = -y(1) + y(3);
dy(3) = -0.51 * y(1) * y(2);
對該微分方程組用ode45和自編的龍格庫塔函數進行比較，調用如下：
[T,F] = ode45(@test_fun,[0 15],[1 1 3]);
subplot(121)
plot(T,F)%Matlab自帶的ode45函數效果
title('ode45函數效果')
[T1,F1]=runge_kutta1(@test_fun,[1 1 3],0.25,0,15);%測試時改變test_fun的函數維數，別忘記改變初始值的維數
subplot(122)
plot(T1,F1)%自編的龍格庫塔函數效果
title('自編的 龍格庫塔函數')

⑥ 拉格朗日插值 MATLAB

額，你這是lagrange插值？
這里有個lagrange插值的例子，你參照一下吧。
function y=lagrange_interp(xdata,ydata,x)
% Lagrange插值
% 輸入參數：
% ---xdata：給定的節點橫坐標
% ---ydata：給定的節點縱坐標
% ---x：需要進行插值的節點橫坐標
% 輸出參數：
% ---y：Lagrange插值函數在x處的函數值
n=length(xdata);m=length(ydata);
if n~=m
error('插值數據長度不等！');
end
ii=1:n;y=zeros(size(x));
for i=ii
ij=find(ii~=i);V=1;
for j=1:length(ij)
if abs(xdata(i)-xdata(ij(j)))<eps
error('輸入的n+1個節點不是互異的。');
end
V=V.*(x-xdata(ij(j)));
end
y=y+V*ydata(i)/prod(xdata(i)-xdata(ij));
end

⑦ 用matlab編寫拉格朗日插值演算法的程序

做了一個測試，希望有所幫助。代碼：% 用matlab編寫拉格朗日插值演算法的程序，並以下面給出的函數表為數據基礎，
% 在整個插值區間上採用拉格朗日插值法計算f（0.6），寫出程序源代碼，輸出計算結果
% x -2.15 -1.00 0.01 1.02 2.03 3.25
% y 17.03 7.24 1.05 2.03 17.06 23.05
function main()
clc;
x = [-2.15 -1.00 0.01 1.02 2.03 3.25];
y = [17.03 7.24 1.05 2.03 17.06 23.05 ];
x0 = 0.6;
f = Language(x,y,x0)function f = Language(x,y,x0)
%求已知數據點的拉格朗日插值多項式
%已知數據點的x坐標向量: x
%已知數據點的y坐標向量: y
%插值點的x坐標: x0
%求得的拉格朗日插值多項式或在x0處的插值: fsyms t l;
if(length(x) == length(y))
n = length(x);
else
disp('x和y的維數不相等！');
return; %檢錯
endh=sym(0);
for (i=1:n)
l=sym(y(i));
for(j=1:i-1)
l=l*(t-x(j))/(x(i)-x(j));
end;
for(j=i+1:n)
l=l*(t-x(j))/(x(i)-x(j));
end;

h=h+l;
end
simplify(h);if(nargin == 3)
f = subs (h,'t',x0); %計算插值點的函數值
else
f=collect(h);
f = vpa(f,6); %將插值多項式的系數化成6位精度的小數
end結果：
f = 0.0201>>

⑧ matlab 拉格朗日插值法

function f = Language(x,y,x0)
%求已知數據點的拉格朗日插值多項式
%已知數據點的x坐標向量: x
%已知數據點的y坐標向量: y
%插值點的x坐標: x0
%求得的拉格朗日插值多項式或在x0處的插值: f

syms t;
if(length(x) == length(y))
n = length(x);
else
disp('x和y的維數不相等！');
return; %檢錯
end
f=0.0; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%改為 f=sym(0);
for(i=1:n)
l=y(i); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%改為 l=sym(y(i));
for(j=1:i-1)
l=l*(t-x(j))/(x(i)-x(j));
end;
for(j=i+1:n)
l=l*(t-x(j))/(x(i)-x(j));
end;

f=f+1;
simplify(f);

if(i==n)
if(nargin == 3)
f = subs (f,'t',x0); %計算插值點的函數值
else
f=collect(f);
f = vpa(f,6); %將插值多項式的系數化成6位精度的小數
end
end
end

-------------------------------------------------------
下面的這個應該可以：

function f = Language(x,y,x0)
%求已知數據點的拉格朗日插值多項式
%已知數據點的x坐標向量: x
%已知數據點的y坐標向量: y
%插值點的x坐標: x0
%求得的拉格朗日插值多項式或在x0處的插值: f

syms t l;
if(length(x) == length(y))
n = length(x);
else
disp('x和y的維數不相等！');
return; %檢錯
end

h=sym(0);
for (i=1:n)
l=sym(y(i));
for(j=1:i-1)
l=l*(t-x(j))/(x(i)-x(j));
end;
for(j=i+1:n)
l=l*(t-x(j))/(x(i)-x(j));
end;

h=h+l;
end
simplify(h);

if(nargin == 3)
f = subs (h,'t',x0); %計算插值點的函數值
else
f=collect(h);
f = vpa(f,6); %將插值多項式的系數化成6位精度的小數
end