matlab 计算方法的总结

此文章是我的大学课程《计算方法》的总结,所选用的代码是matlab的形式,因为内容都是个人的总结,大部分都是只是一个函数+事例搞定。。所以要问我基础的东西其实我也不是很懂。。

 

所以下面这文章不讲理论,只讲函数怎么用。。能力有限哈。。

 

里面涉及到的定义我是百科直接搬运过来的~~

 

 

目录
一元线性线性方程的求解
   什么是一元线性方程,什么是一元非线性方程?
    二分法
    牛顿法
    弦截法

线性方程的计算方法
       相关命令的基础知识
       高斯消去法
      LT分解
        

插值法 
     线性插值
     拉格朗日插值
     牛顿插值
     样条插值
     其他命令

数值积分
   牛顿-柯西公式
   数值积分命令

常微分方程的数值解法
   (缺课了)


最优方法
    没理论,只有命令


二分法 function [x,y,Result]=erfenfa(a,b) %erfenfa是函数名(a,b)是参数名,x,y,result是返回值
%二分法,抄老师的。已加上注释
%输入a,b--区间端点,要求该区间必须有一个跟
%输出x---近似跟
% y--fun(x)
%


y1=fun(a); %fun--函数名称,输入一个x得到他的函数值
y2 = fun(b);

if y1*y2 >=0
  x=[]; y=[] ; Result=[];
  return ;
end

delta = b-a; %得到他的第一次距离差是否符合精度

Result=zeros(1,2) ; %Result是迭代的试验结果,这句话的意思是初始化Result矩阵
k=1; %k=1也是初始化值。k是表示迭代多少次

while abs(delta) >1e-8 %当迭代到的初值未我想追求的精度的时候
  xm=(a+b)/2; %继续缩小a+b /2
  ym=fun(xm); %得到函数值
  Result(k,1) = xm; %Result是保存每一次的试验结果的矩阵
  Result(k,2) = ym;
  k = k+1; %k记录迭代多少次
  if ym==0 %如果函数值是0,就说明我找到了那个点,函数就退出
    x=xm;
    return;
  elseif %否则,如果 ym*y1>0 ,所以ym*y2<0 ,说明函数点在区间右侧 所以a=xm;
    ym*y1>0
a=xm;
  else
    b=xm; % 说明 ym*y2>0 所以ym*y1<0 ,说明函数点在区间左侧 所以b=xm;
  end
  delta = b-a; %得到下一次距离差,看是否得到精度的差
end

x = (a+b)/2 ; %这里得到的是我们最终点的值
y = fun(x) ; %得到的最终结果的函数值
Result(k,1) = x; %保存最终结果
Result(k,2) = y;

function y = fun(x) %测试一个简单的函数
y = x*(x*x-1) -1;
end

 

一元线性线性方程的求解
 什么是线性方程,什么是非线性方程?
首先我们要理解线性方程也称一次方程。指未知数都是一次的方程。其一般的形式是ax+by+...+cz+d=0线性方程也称为一次方程。
因为在笛卡尔坐标系上任何一个一次方程的表示都是一条直线。组成一次方程的每个项必须是常数或者是一个常数和一个变量的乘积。
且方程中必须包含一个变量,因为如果没有变量只有常数的式子是算数式而非方程式。一元也就是一个未知量了。

非线性方程所谓的非线性方程,那就是未知的次方数不一定是一次的方程咯。非线性方程,就是因变量与自变量之间的关系不是线性的关系,
这类方程很多,例如平方关系、对数关系、指数关系、三角函数关系等等。求解此类方程往往很难得到精确解,经常需要求近似解问题。
相应的求近似解的方法也逐渐得到大家的重视。

一元是只有一个未知量,也就是一个x

 

牛顿法,切线法

基本思想:(日后再补充)

 

 

%割线法

弦截法
弦截法是一种求方程根的基该方法,在计算机编程中常用。

 

他的思路是这样的:任取两个数x1、x2,求得对应的函数值f(x1)、f(x2)。
如果两函数值同 号,则重新取数,直到这两个函数值异号为止。
连接(x1,f(x1))与(x2,f(x2))这两点形成的直线与x轴相交于一点x,求得对应的f(x),
 判断其与f(x1)、f(x2)中的哪个值同号。如f(x)与f(x1)同号,则f(x)为新的f(x1)。
将新的f(x1)与f(x2)连接,如此循 环。体现的是极限的思想。




function[x,y]=gexianfa(x0,x1)

%割线法程序
%输入 x0,x1 ---初始迭代点
%输出  x---近似跟
%  y---fun(x)


%{

%}

y0 = fun(x0); %fun--函数名,输入x可得到其函数值
y1 = fun(x1);  %初始化迭代点
x2 = x1-y1*(x1-x0)/(y1-y0);  % 根据y1-y0 / x1-x0 来定斜率, 通过x1,y1点
                                           %得出方程 y1-y0 / x1-x0 (x-x1) = y-y1 ;再过 这条直y=0的  点,得到这条方程
delta = x2 -x1;

while abs(delta ) >1e-8   %根据规律,先列出一条数列表达式,再算出x(k+1),再发现迭代的变化量
                                    %然后再去搞
     y2=fun(x2);          %这条是迭代式。。我用数列通式将东西表示出来后,发现变量是下面这样赋值
    x0 = x1;
    y0 = y1;
     x1 = x2;
     y1 = y2;
      x2 = x1-y1*(x1-x0)/(y1-y0);
    delta = x2-x1;
end
x= x2;
y = fun(x);
end

function y = fun(x)
y = x-exp(-x);
end

 

线性方程的计算方法
       相关命令的基础知识
       高斯消去法
      LT分解
        
 线性方程方程组
(这些缺课了,百度了一下发现有点儿晕,最后才补充完)

 

插值法 
离散数据的基础上补插连续函数,使得这条连续曲线通过全部给定的离散数据点。插值是离散函数逼近的重要方法,利用它可通过函数在有限个点处的取值状况,估算出函数在其他点处的近似值。插值:用来填充图像变换时像素之间的空隙。
线性插值

X = linspace(-5,5,11)'; %系统自动生成 -5到5的11个自动均分的数字, -5~-1 (5个) 0 (1个) 1~5(5个)
Y = 1./(1+X.^2); % Runge使用的函数
x = linspace(-5.5,5.5,101);
y0 = 1./(1+x.^2); % 用来画原函数的图形
y1 = interp1(X,Y,x);% 线性插值 , 意思是 以X和Y作参考,生成一个线性函数k(X),然后用k(X)来算k(x)得到y1


     拉格朗日插值
拉格朗日插值使用了一个迭代公式
首先,我们给出一个迭代式的例子
x 0 1 2
y 1 2 3

基本例子
P(x) = 1*(x-1)(x-2) / (0-1)(0-2)
+2*(x-0)(x-2) / (1-0)(1-2)
+3*(x-0)(x-1) / (2-0)(2-1)
= x+1

这个例子我认真一看,太有迷惑性,所以我换了一个例子
x 1 2 3 4
y 0 -5 -6 3

p3(x) =0 *(x-2)(x-3)(x-4) / ( 1-2) (1-3) (1-4)
+(-5) *(x-1)(x-3)(x-4) / (2-1)(2-3)(2-4)
+(-6)*(x-1)(x-2)(x-4) /(3-1)(3-2)(3-4)
+3*(x-1)(x-2)(x-3) /(4-1)(4-2)(4-3)
=x^3-4x^2+3



可得到matlab循环式


伪代码
 denominator = 1;
memeber=1;
for i=1...n
for j=1...n

if ( i!=j)
 denominator = denominator*(i-j);
 member = member * f(i)* (k-j);
 end
end
end


我的例子可能有误,老师没有f(i)项

function y = Linterp(X,Y,x)
% 给定数据(X,Y)进行Lagrange插值,并据此返回x点处的函数值y
% X,Y是同型列向量
% x 可以是向量

% 2015.4.14
n = length(X);
L = zeros(n,1);
m = length(x);
y = zeros(m,1);
for k=1:m
    for i=1:n
        tem1 = 1;  tem2 = 1;
        for j=1:n
            if j~=i
                tem1 = tem1*(X(i)-X(j)); %分子
                tem2 = tem2*(x(k)-X(j)); %分母 ,他没有那个f(i)项?
            end
        end
        L(i) = tem2/tem1; % Lagrange基函数
    end
    y(k) = sum(Y.*L);
end
end


牛顿插值
牛顿插值是什么?等后面补充


function y = Ninterp(X,Y,x)
% 给定数据(X,Y)进行Newton插值,并据此返回x点处的函数值y
% X,Y是同型列向量
% x 可以是向量
% 刘群锋
% 2015.4.14
n = length(X);
m = length(x);% 需要预测的点数
y = zeros(m,1);
for t=1:m % 对第t个需要预测的点x(t)...
    y(t) = Y(1);
    for k=2:n
        temp = 0;
        for i=1:k  % 计算f[x1,...,xk],赋值给temp
            temp1 = 1;
            for j=1:k
                if j~=i
                    temp1 = temp1*(X(i)-X(j));
                end
            end
            temp = temp + Y(i)/temp1;
        end
        temp2 = 1;
        for i=1:(k-1)
            temp2 = temp2*(x(t)-X(i));
        end
        y(t) = y(t) + temp*temp2;
    end
end


样条插值
X = linspace(-5,5,11)';
Y = 1./(1+X.^2);    % Runge使用的函数
x = linspace(-5.5,5.5,101);
y0 = 1./(1+x.^2);   % 用来画原函数的图形
y1 = interp1(X,Y,x);% 线性插值
y2 = interp1(X,Y,x,'spline'); % 三次样条插值



综合比较程序
% 演示interp1中的插值方法

%
clear;
clc;
X = linspace(-5,5,11)';
Y = 1./(1+X.^2);    % Runge使用的函数
x = linspace(-5.5,5.5,101);
y0 = 1./(1+x.^2);   % 用来画原函数的图形
y1 = interp1(X,Y,x);% 线性插值
y2 = interp1(X,Y,x,'spline'); % 三次样条插值
y3 = Linterp(X,Y,x);  % Lagrange插值 这Linterp是上面的拉格朗日函数
%y4 = Ninterp(X,Y,x);  % Newton插值, 与Lagrange插值图形一样,可略 这个是牛顿插值函数

figure;
plot(X,Y,'bo',x,y0,'k-'); % 画出样本点和函数原图
hold on;
plot(x,y1,'b--',x,y2,'r:'); % 画出线性插值图形、三次样条插值图形
% legend('样本点','函数原图','线性插值','三次样条插值');
 plot(x,y3,'g--'); % 画出Lagrange插值
 legend('样本点','函数原图','线性插值','三次样条插值','Lagrange插值');
hold off;



 

 

数值积分
   牛顿-柯西公式
matlab 命令
 q = quad(@(x) x.^3 ,1,3)
输出
q =

    20
等价于
 3^4*1/4 - 1^4*1/4

ans =

    20
数值积分命令
常微分方程的数值解法
   (缺课了)
odm23






最优方法

什么时候用,怎么用?

没理论,只有命令
fminbnd(极小值)找到极小值(我们现在想找到最小值,所以先画图确定区域,然后去找)

fminbnd(f,1,3)

x=-10:0.1:10;
y=(x-1).^2
plot(x,y) %画一下图,发现极小值的区域(因为fminbnd只能找到极小值,不一定能找到最小值)
%考虑到最小值的区域
在[-10,10]里,找到了有一个极小值点(这里其实还要去证明一下是否存在其他极小值点)

fminbnd(@(x) (x-1).^2 ,-10,10);

%这里也就找到了极小值的点,x=1的点


val 是x点的值,为什么会那么复杂呢?这个问题先留着。。。0 0我也不懂。。


peaks
meshgrid
fminunc(无约束)
1、无约束(无条件)的最优化 
(1)使用fminunc函数   (un-condition) 
(2)可用于任意函数求最小值(局部最小)
(3)将最大、最小问题统一为求最
小值问题(即只能求最小值)。如求最大值,
则要对函数取相反数而变成求最小值问题,
最后再把函数值取反即为函数的最大值。
(求最大值需两次取反(一前一后



是x,y 是用-->x(1),x(2)来表示,而后面的x0=[1,1]是用来初始化x,y的,不是一个点。





对于fminunc。(注意,这个是没有任何约束的)
1.先确定最大值还是最小值
2.先建立.m的函数文件...
3.再定义一个多少函数的初始值。
4运行


fminsearch






解线性规划问题最优解
求解 3x1+4x2+x3;
有约束
{ x1+2x2+3x3<=100
x1+x2-x3<=10
x1,x2,x3>=0}

f=[3,4,1] ; %方程矩阵
A=[1,2,3 ; 1,1,-1] %不等式左边系数矩阵
b=[100;10]%不等式右边参数矩阵
lb=[0; 0 ;0];

linprog(f,A,b,[],[],lb)%求解线性规划问题得到最优解

解非线性规划的问题
fmincon






problem
globlesearch
ga
simulannealbnd(退火算法)
...

目录
一元线性线性方程的求解

命令相关
fzero


线性方程的计算方法 相关命令的基础知识 高斯消去法 LT分解 插值法 线性插值 拉格朗日插值 牛顿插值 样条插值 其他命令 数值积分 牛顿
-柯西公式 数值积分命令 常微分方程的数值解法 (缺课了) 最优方法 没理论,只有命令

 

 

 

posted @ 2015-05-29 15:13  david_way  阅读(1170)  评论(0编辑  收藏  举报