数值分析编程总结

时间:2024.4.20

1. LU分解:

[L U]=lu(A);

2. 追赶法

function x=zhuiganfa(A,b)

[n,n]=size(A);

for i=1:n

    if(i==1)

        l(i)=A(i,i);

        y(i)=b(i)/l(i);

    else

        l(i)=A(i,i)-A(i,i-1)*u(i-1);

        y(i)=(b(i)-y(i-1)*A(i,i-1))/l(i);

    end

    if(i

        u(i)=A(i,i+1)/l(i);

    end

end

x(n)=y(n)

for j=n-1:-1:1

    x(j)=y(j)-u(j)*x(j+1);

end

数值试验:

n=101;

关键是如何定义上述矩阵:

>>n=101;

c1=ones(1,n-1);

a1=diag(c1,-1);   这个-1说明行位置-1

c2=12*ones(1,n);

a2=diag(c2); 

c3=ones(1,n-1);

a3=diag(c3,1); 

a=a1+a2+a3;

3. 拉格朗日插值

function yh=lage(x,y,xh)

n=length(x);

m=length(xh);

yh=zeros(1,m);

c1=ones(n-1,1);

c2=ones(1,m);

for i=1:n

  xp=x([1:i-1 i+1:n]);

  yh=yh+y(i)*prod((c1*xh-xp'*c2)./(x(i)-xp'*c2));

end

end

>> x=[11,12];

>> y=[2,4];

>> xh=[11.75];

>> lage(x,y,xh)

ans =

3.5000

4 最小二乘法

1.       最小二乘的xi和yi为:

要拟合的函数为:y=a+bx-cxy   注意不是多项式

2.       编程函数为:

function z = erchen(x,y)

x1=ones(5,1);         

A=[x1,x,-x.*y];       注意点乘

z=A\y;              注意左除

a=z(1);

b=z(2);

c=z(3);

end        

输入:

?? x=[1953 1964 1982 1990 2000]';

?? y=[5.82 6.95 10.08 11.34 12.66]';

?? erchen(x,y)

ans =

    2.9456         =a

   -0.0014         =b

   -0.0005         =c

1.       最小二乘的xi和yi为:

要拟合的函数为:y=a+bx+cx2   是多项式

2.       编程函数为:

function z = erchen2(x,y)

x1=ones(5,1);         

A=[x1,x,x.^2];

z=A\y;

a=z(1);

b=z(2);

c=z(3);

end

输入:

?? x=[0 0.25 0.5 0.75 1.00]';

?? y=[1.00 1.284 1.6487 2.1170 2.7183]';

?? erchen2(x,y)

ans =

    1.0051

    0.8642

        0.8437

最小二乘多项式拟合的简单函数方法:

?? x=[0 0.25 0.5 0.75 1.00]';

?? y=[1.00 1.284 1.6487 2.1170 2.7183]';

?? P=polyfit(x,y,2)    要拟合成4次,则2改成4就可以了

P =

0.8437    0.8642    1.0051   注意此内置函数输出的结果c,b,a是反的

5 复合辛普森公式求解积分

先定义函数:

function v=f(x)

 v=sin(x);  “若定义有除数要点除,分母有0时要特殊定义”

end

定义程序:

function  I=fsps(f,a,b,n)

h=(b-a)/n;

x=linspace(a,b,2*n+1);

y=feval(f,x);

I=(h/6)*(y(1)+2*sum(y(3:2:2*n-1))+4*sum(y(2:2:2*n))+y(2*n+1));

end

>>  fsps('f ',0,1,4)

ans =

    0.4597

6.不动点迭代思路

不动点迭代常常有好几个迭代的不动点函数,所以要分别定义这些函数是很困难的,如是乎使用SWITCH内置函数进行切换,叫切换函数.

1.先定义函数后进行编程的方法

先需要定义不动点函数:

function v=f(x)

v=x^3-x-1;

end

再定义编程:

function  [it,x]=fixpnt1(f,a,maxit,tol)

it=0;

x=feval(f,a);

while it<=maxit & abs(x-a)>tol,

      it=it+1;

      a=x;

      x=feval(f,a);

end

此函数的调用:

>>  fixpnt1('f',2,100,1e-5)

ans =

     1

3.利用切换函数SWITCH的方法(多个不动点迭代函数)

function [x,it]=fixpnt(np,a,maxit,tol)

switch  np,

    case 1,

        phi=inline('(3*x+10)^(1/5)');

    case 2,

        phi=inline('sin(10*x)+2*cos(x)-3');

    case 3,

        phi=inline('3-atan(x)');

    case 4,

        phi=inline('-2-1/log(x^2+x+1)');

end

it=0;

x=phi(a);

while it<=maxit & abs(x-a)>tol,

      it=it+1;

      a=x;

      x=phi(a);

 end

使用与输入:

>>  fixpnt(2,1,100,1e-5)

ans =

   -4.2696

7. 雅可比迭代

function [x it]=jacobi(A,b,tol)

D=diag(diag(A));

L=D-tril(A);

U=D-triu(A);

x=zeros(size(b));

for it=1:500

    x=D\(b+L*x+U*x);

error=norm(b-A*x)/norm(b);

if (error

    break;

end

end

8. 高斯迭代

function [x it]=gaosi(A,b,tol)

D=diag(diag(A));

L=D-tril(A);

U=D-triu(A);

x=zeros(size(b));

for it=1:500

    x=(D-L)\(b+U*x);

error=norm(b-A*x)/norm(b);

if (error

    break;

end

end

9. SOR迭代

function [x it]=SOR(A,b,w,tol)

D=diag(diag(A));

L=D-tril(A);

U=D-triu(A);

x=zeros(size(b));

for it=1:500

    x=(D-w*L)\(w*b+(1-w)*D*x+w*U*x);

error=norm(b-A*x)/norm(b);

if (error

    break;

end

end

10.二分法:

1. 先要定义所求的函数:

function v=f(x)

v=x^3-x-1;

end

2. 二分法程序如下:

function [x,it]=erfengfa(a,b,f,tol)  

fa=feval(f,a);      

fb=feval(f,b);

it=0;

while abs(b-a)>tol,      

    it=it+1;

    x=a/2+b/2;

    fx=feval(f,x);

if   sign(fx)==sign(fa),    

a=x;fa=fx;

else

    b=x;fb=fx;

end

 end

11. 牛顿法:

1.先定义函数后进行编程的方法

先需要定义不动点函数

需要计算的函数f

function v = f( x )

v=x^5-3*x-10;

end

需要计算的函数的导数g

function v = g( x )

v=5*x^4-3;

end

2. 再定义编程:

function v = newton(a,f,g,maxit,tol)

it=0;

x=a;

while   it<=maxit & abs(feval(f,x))>tol,

        it=it+1;

        x=x-feval(f,x)/feval(g,x);

end

v=[x,it];

end

12. 牛顿下山法:

1.先定义函数后进行编程的方法

先需要定义不动点函数

需要计算的函数f

function v = f(x)

v=x^2+sin(10*x)-1

end

需要计算的函数的导数g

function v =g(x)

v = 2*x+10*cos(10*x)

end

2. 再定义编程1:

function v = newtonxiashang(x0,f,g,maxit,tol)

x=x0;

it=0;

while it<=maxit & abs(feval(f,x))>tol,

    it=it+1;

    d=-feval(f,x)/feval(g,x);

    lambda=1;

    isdone=0;

    while ~isdone,

        xn=x+lambda*d;

        if abs(feval(f,xn))

            isdone=1;

        else lambda=lambda*0.5;

        end

    end

   x=xn;

end

v=[x,it];

end

3. 再定义编程2:

function v = newtonxiashang2(x0,f,g,maxit,tol)

x=x0;

it=0;

while it<=maxit & abs(feval(f,x))>tol,

    it=it+1;

    d=-feval(f,x)/feval(g,x);

    lambda=1;

   

    while   abs(feval(f,x+lambda*d))>=abs(feval(f,x)),

            lambda=0.5*lambda;

    end

   x=x+lambda*d;

end

   v=[x,it];

end

13. 割线法

1.先定义函数后进行编程的方法

先需要定义函数

需要计算的函数f

function v = f( x )

v=x^5-3*x-10;

end

2. 再定义编程:

function v = gexian(a,b,f,maxit,tol)

it=0;

x0=a;

x=b;

while it<=maxit & abs(feval(f,x))>tol,

    it=it+1;

    xt=x-feval(f,x)/(feval(f,x)-feval(f,x0))*(x-x0);

    x0=x;

    x=xt; 

end

v=[x,it];

end

输入:

>>  gexian(1,2,'f',100,1e-5)

ans =

    1.7226    7.0000

14. 乘幂法

function [t,y] = chenmifa( a,xinit,ep )

v0=xinit;

[tv ti]= max(abs(v0));

lam0=v0(ti);

u0=v0/lam0;

flag=0;

while (flag==0)

    v1=a*u0;

    [tv ti]= max(abs(v1));

    lam1=v1(ti);

    u0=v1/lam1;

    if (abs(lam0-lam1)<=ep)

        flag=1;

    end

    lam0=lam1;

end

        t=lam1;

        y=u0;

end

运行及其结果:

>> a=[2,3,2;10,3,4;3,6,1];

>> xinit=[1 1 1]';

>> ep=0.0001;

>> [t y]=chenmifa(a,xinit,ep)

t =

   11.0000

y =

    0.5000

    1.0000

    0.7500

15. 反幂法

function [t y] =fanmifa( a,xinit,ep )

v0=xinit;

[tv ti]= max(abs(v0));

lam0=v0(ti);

u0=v0/lam0;

flag=0;

while (flag==0)

    v1=a^-1*u0;

    [tv ti]= max(abs(v1));

    lam1=v1(ti);

    u0=v1/lam1;

    if (abs(lam0^-1-lam1^-1)<=ep)

        flag=1;

    end

    lam0=lam1;

end

        t=lam1^-1;

        y=u0;

end

运行提示:

>> a=[12 6 -6;6 16 2;-6 2 16];

>>xinit=[1 -0.5 0.5]';

>> [t y] =fanmifa( a,xinit,0.0001 )

t =

    4.4560

y =

    1.0000

   -0.6287

0.6287

16. 结合原点平移反幂法

function [t y] =yuandian( a,p,xinit,ep )

[n n]=size(a);

v0=xinit;

[tv ti]= max(abs(v0));

lam0=v0(ti);

u0=v0/lam0;

flag=0;

while (flag==0)

    v1=(a-p*eye(n))^-1*u0;

    [tv ti]= max(abs(v1));

    lam1=v1(ti);

    u0=v1/lam1;

    if (abs(lam0^-1-lam1^-1)<=ep)

        flag=1;

    end

    lam0=lam1;

end

        t=lam1^-1+p;

        y=u0;

end

运行结果:

>> a=[6 2 1;2 3 1;1 1 1];

>> p=6;

>> ep=0.0001;

>> xinit=[1 1 1]';

>> [t y]=yuandian(a,p,xinit,ep)

t =

    7.2880

y =

    1.0000

    0.5229

    0.2422

17. 欧拉公式初值问题

1. 先定义导数函数f:

function v = f( x,y)

v=-y+x+1;

end

2. 定义欧拉公式编程:

function [x,y] = oula( f,y0,a,b,n )

y(1)=y0;

h=(b-a)/n;

x=a:h:b;

for i=1:n;

    y(i+1)=y(i)+h*feval(f,x(i),y(i));

end

end

18. 改进的欧拉公式初值问题

1. 先定义导数函数f:

function v = f( x,y)

v=-y+x+1;

end

要先运行一下。

2. 定义欧拉公式编程:

function [x,y] =gaijingoula(f,y0,a,b,n)

y(1)=y0;

h=(b-a)/n;

x=a:h:b;

for i=1:n;

    yp=y(i)+h*feval(f,x(i),y(i));

    yc=y(i)+h*feval(f,x(i+1),yp);

    y(i+1)=0.5*(yc+yp);

end

end

19. 梯形公式初值问题

1. 先定义导数函数f:

function v = f( x,y)

v=-y+x+1;

end

要先运行一下。

2. 定义欧拉公式编程:

function [x,y] = tixing( f,y0,a,b,n )

y(1)=y0;

h=(b-a)/n;

x=a:h:b;

for i=1:n;

    y(i+1)=1/(0.5*h+1)*(y(i)+0.5*h*(feval(f,x(i),y(i))+x(i+1)+1));

end

end

20. 标准四阶四段龙格库塔公式初值问题

1. 先定义导数函数f:

function v = f( x,y)

v=y+x;

end

要先运行一下。

2. 定义欧拉公式编程:

function [x,y] = longgekuta(f,y0,a,b,n)

y(1)=y0;

h=(b-a)/n;

x=a:h:b;

for i=1:n,

    k1=h*feval(f,x(i),y(i));

    k2=h*feval(f,x(i)+0.5*h,y(i)+0.5*k1);

    k3=h*feval(f,x(i)+0.5*h,y(i)+0.5*k2);

    k4=h*feval(f,x(i)+h,y(i)+k3);

    y(i+1)=y(i)+(k1+2*k2+2*k3+k4)/6;

end

end

21. 基本问题

1. 特征值:

eig(A)

2. 数值显示位数:

>> format long

>> format short

3. 左除和右除:

AX=B,则X=A\B;  左除\

XB=A,则X=A/B;  右除/   (右除是真除)

4. 特殊符号

In(x)=log(x)

根号用sqrt(x)

5.绘图功能

Plot(x,y,style)

>>x=0:0.5:10;

>>y=x;

>>hold on;

>>plot(x,y,’g+:’)

后面style的意义:

y,g,r,b分别指颜色黄色、绿色、红色、蓝色;

o * + p 代表节点的形状,圆形、*号、+号、五角星;

--、-、:代表线型,分别为虚线、线段、点线

更多相关推荐:
数值分析学习总结感想

数值分析学习感想摘要数值分析主要介绍现代科学计算中常用的数值计算方法及其基本原理研究并解决数值问题的近似解是数学理论与计算机和实际问题的有机结合随着科学技术迅速发展运用数学方法解决工程技术领域中的实际问题已经得...

数值分析总结

数值分析复习总结第二章数值分析基本概念教学内容1误差与有效数字误差误差限相对误差相对误差限和有效数字的定义及相互关系误差的来源和误差的基本特性误差的计算估计的基本方法2算法的适定性问题数值分析中的病态和不稳定性...

数值分析课程总结

课程内容1误差了解误差的来源与分类及误差的基本概念与性质;熟悉绝对误差及绝对误差限、相对误差及相对误差限和有效数字之间的关系;掌握一元和二元函数的误差估计式并会应用;熟悉减小误差的积累和传播应注意的几大原则和通…

数值分析考试复习总结

第一章1误差相对误差和绝对误差得概念例题当用数值计算方法求解一个实际的物理运动过程时一般要经历哪几个阶段在哪些阶段将有哪些误差产生答实际问题数学模型数值方法计算结果在这个过程中存在一下几种误差建立数学模型过程中...

数值分析总结

1章引论2章非线性方程求根3章解线性方程组的直接法4章解线性方程组的迭代法5章插值法6章数值积分7章常微分方程的数值解法第2章非线性方程的迭代法方程求根与二分法迭代法迭代收敛的加速方法牛顿法弦截法第3章解线性代...

数值分析总结

一1数值分析的特点1首先要有可靠的理论分析以确保算法在理论上的收敛性和数值上的稳定性2其次要对计算的结果进行误差估计以确定其是否满足精度3还要考虑算法的运行效率即算法的运算量和存储量2数值分析的误差种类1截断误...

数值分析总结

数值分析复习提要一纲要数值积分与数值微分一章中主要的要点如下数值积分的提法插值型求积公式的导出及其余项估计低阶数值积分公式及其余项的估计数值积分的加速过程Romberg算法与埃特金方法高精度求积公式Gauss求...

数值分析(颜庆津)第7章 学习小结

第7章常微分方程初值问题的数值解法--------学习小结一、本章学习体会本章的主要内容是要掌握如何用数值解代替其精确解,这对于一些特殊的微分方程,特别是一些不好解其通解方程是非常有用的。对于本章我总结如下几点…

西南交通大学研究 数值分析上机实习报告20xx

数值分析上机实习报告要求1应提交一份完整的实习报告具体要求如下1要有封面封面上要标明姓名学号专业和联系电话2要有序言说明所用语言及简要优特点说明选用的考量3要有目录指明题目程序计算结果图标和分析等内容所在位置作...

数值分析考试总结

1数值分析深度总结大兵第一章绪论1数值分析的特点面向计算机算法包括和逻辑运算都是计算机能够直接处理的有可靠的理论分析能任意逼近并达到精读要求对近似算法要保证收敛性和数值稳定性还要对误差进行分析有好的计算复杂性有...

数值分析实验报告

机电工程学院机械工程陈星星6720xx0109数值分析课程设计实验报告实验一函数插值方法一问题提出对于给定的一元函数yfx的n1个节点值yjfxjj01n试用Lagrange公式求其插值多项式或分段二次Lagr...

数值分析实验报告

课程实验报告课程实验报告课程实验报告课程实验报告

数值分析总结(39篇)