• 设为首页
  • 点击收藏
  • 手机版
    手机扫一扫访问
    迪恩网络手机版
  • 关注官方公众号
    微信扫一扫关注
    迪恩网络公众号

Matlab数值计算示例: 牛顿插值法、LU分解法、拉格朗日插值法、牛顿插值法 ...

原作者: [db:作者] 来自: [db:来源] 收藏 邀请

本文源于一次课题作业,部分自己写的,部分借用了网上的demo

  • 牛顿迭代法(1)
x=1:0.01:2;
y=x.^3-x.^2+sin(x)-1;
plot(x,y,\'linewidth\',2);grid on;%由图像可知 根在1.05到1.15之间

syms x
s0=diff(x^3-x^2+sin(x)-1,x,1);
% 得到s0= cos(x) - 2*x + 3*x^2
% 迭代方程为 y=x-(x.^3-x.^2+sin(x)-1)/(cos(x) - 2.*x + 3*x.^2)
clear;
x=1.15;
for i=1:30
    x=x-(x.^3-x.^2+sin(x)-1)/(cos(x) - 2.*x + 3*x.^2)%根据牛顿迭代法公式。一直迭代计算30次。
end
%可得  x取值为1.0935;
  • 牛顿迭代法(2)
%%   绘制图形。判断跟的大概位置。
x=1:0.01:2;
f=x.^3-x.^2+sin(x)-1;
plot(x,f,\'linewidth\',2);grid on;%由图像可知 根在1.05到1.15之间
%%
clc,clear
syms x 
f=x.^3-x.^2+sin(x)-1;%所求函数
df=diff(f,x); %求取一阶导数
eps=1e-4; %误差判断
x0=1.15; %迭代初始值。
cnt=0; 
MAXCNT=20; %最大循环次数 
while cnt<MAXCNT %防止无限循环 
       x1=x0-subs(f,x,x0)/subs(df,x,x0); %去掉这个分号,可以看到迭代过程.
      if (abs(x1-x0)<eps) 
      break; 
      end 
        x0=x1; 
        cnt=cnt+1; 
end 
if cnt==MAXCNT 
    disp \'不收敛\' 
else 
    vpa(x1,8) 
end
  • LU分解法

被调函数:


function [L,U]=lufj(A)
%  利用紧凑格式法原理 编写的LU 分解
[n,m]=size(A); % 获取A矩阵的行和列
if m~=n    %判断行列相等与否
    error(\'Not a squared matrix1\');
else
    A(2:n)=A(2:n)/A(1,1);
    for k=2:n-1
        A(k,k:n)=A(k,k:n)-A(k,1:k-1)*A(1:k-1,k:n);
        A(k+1:n,k)=(A(k+1:n,k)-A(k+1:n,1:k-1)*A(1:k-1,k))/A(k,k);
    end    %都是根据定义进行循环计算。
end
L=A;U=A;
for i=1:n
    L(i,i)=1;
    L(i,i+1:n)=0;
    U(i,1:i-1)=0;
end

主函数:

%%  需要调用lufj函数;
A=[-2 -2 3 5
    1 2 1 -2
    2 5 3 -2
    1 3 2 3]
b=[-1
    4
    7
    0]
% x=A\b  %左除法求解
[L,U]=lufj(A);
x0=L\b;
x=U\x0%求出的x即为解   

  • 拉格朗日插值法
    被调函数:
function y=lagrange(x0,y0,x);
% 根据拉格朗日插值定义编写
n=length(x0);m=length(x);
for i=1:m
  z=x(i);
  s=0.0;%给s的初值
  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

主函数:

x=[0,1,2,4];
y=[1,9,23,3];
y0=lagrange(x,y,1.5)
  • 牛顿插值
    被调函数:
function yi=New_Int(x,y,xi)
%Newton基本插值公式
%x为向量,全部的插值节点
%y为向量,差值节点处的函数值
%xi为标量,是自变量
%yi为xi出的函数估计值
n=length(x);
m=length(y);
if n~=m
    error(\'The lengths of X ang Y must be equal!\');
    return;
end
%计算均差表Y
Y=zeros(n);
Y(:,1)=y\';
for k=1:n-1
    for i=1:n-k
        if abs(x(i+k)-x(i))<eps
            error(\'the DATA is error!\');
            return;
        end
        Y(i,k+1)=(Y(i+1,k)-Y(i,k))/(x(i+k)-x(i));
    end
end
%计算牛顿插值公式
yi=0;
for i=1:n
    z=1;
    for k=1:i-1
        z=z*(xi-x(k));
    end
    yi=yi+Y(1,i)*z;
end

主函数:

clear all 
clc
x0=[0.4 0.55 0.65 0.80 0.90 1.05];
y0=[0.41075 0.57815 0.69675 0.88811 1.0265 1.25382]; 
x1=0.596;  % 待插值点。
y1=New_Int(x0,y0,x1)% y1即为待插值点的函数值。

TIP:主函数和被调函数要放在一个文件夹内。否则会引起调用错误

NOTE:本文对基本方法做了总结,你可以结合理论知识再来看代码,希望对你有所帮助


鲜花

握手

雷人

路过

鸡蛋
该文章已有0人参与评论

请发表评论

全部评论

专题导读
上一篇:
Rust枚举和模式匹配发布时间:2022-07-18
下一篇:
牛顿插值法,matlab编程计算 - a往南向北发布时间:2022-07-18
热门推荐
阅读排行榜

扫描微信二维码

查看手机版网站

随时了解更新最新资讯

139-2527-9053

在线客服(服务时间 9:00~18:00)

在线QQ客服
地址:深圳市南山区西丽大学城创智工业园
电邮:jeky_zhao#qq.com
移动电话:139-2527-9053

Powered by 互联科技 X3.4© 2001-2213 极客世界.|Sitemap