数值计算方法 Matlab实题训练(内附程序,模型)

《数值计算方法训练》

实习报告

题目：6－A组

院系：上海电力学院数理学院

专业年级：信息与计算科学专业2009级

学生姓名：XX远学号：20092426

2011年7月8日

第1题：含炭量与时间的关系

在某冶炼过程中，钢的含炭量y 与时间t 的统计数据如下 t 0 5 10 15 20 25 30 35 40 45 50 55 y 0 1.27 2.16 2.86 3.44 3.87 4.15 4.37 4.51 4.5

8

4.02

4.64

（1）画出原始数据分布趋势图；

（2）用最小二乘法求钢的含炭量y 与时间t 的拟合曲线32ct bt at y ++=； （3）打印出拟合曲线；

（4）另外选用b at y =进行拟合，比较二种拟合的效果。

解：分析：使用到曲线拟合的最小二乘法，对于拟合函数，尽量转化为可以方便提炼出基函数的方程。在明确基函数的基础上，通过计算，得到各个系数，得到法方程组 （1），程序：function yuan(y)

t=[0:5:55]； plot(t,y,'*')

legend('原始数据分布趋势图')

运行结果：yuan([0 1.27 2.16 2.86 3.44 3.87 4.15 4.37 4.51 4.58 4.02 4.64])

图1 原始数据分布趋势

（2），使用最小二乘法，就必须先取基函数，对于该题流程如下：

①：取基函数为：t =0? 21t =? 32t =? ②：由基函数和y 求法方程组的系数：

??????

???

?

???

?=?=?=?=?=?=?=?=?=∑∑∑∑∑∑∑∑∑=========)(),()(),()(),()

()(),()()(),()()(),()

()(),()()(),()()(),(212121121101210212

1222212

11211121111212

1

02011210100121000i i i i i i i i i i

i i i i i i i i i i i i i i i i i x y f x y f x y f x x x x x x x x x x x x ?????????????????????????????? ③：由这些系数，确定法方程组：

B

X A =???

??

?

?????

????

???????=????

?

?????=),(),(),(),(),(),(),(),(),(),(),(),(210222120121110020100?????????????????????f f f B A ④：解这个法方程组：

??

??

?

?????==?c b a X B X A ，得到拟合函数：3

2ct bt at y ++= 程序：function [a,b,c]=xian(y0)

t0=[0:5:55]; k1=t0;

k2=t0.*t0; k3=t0.*t0.*t0;

A=[sum(k1.*k1) sum(k2.*k1) sum(k3.*k1);sum(k1.*k2) sum(k2.*k2)

sum(k3.*k2);sum(k1.*k3) sum(k2.*k3) sum(k3.*k3)]; B=[sum(k1.*y0);sum(k2.*y0);sum(k3.*y0)]; x=pinv(A)*B; a=x(1,1); b=x(2,1); c=x(3,1); t=0:55;

y=a.*t+b.*t.^2+c.*t.^3; plot(t,y,'--') hold on

plot(t0,y0,'*')

legend('y=a*t+b*t^2+c*t^3拟合效果','真实值')

运行结果：[a,b,c]=xian([0 1.27 2.16 2.86 3.44 3.87 4.15 4.37 4.51 4.58 4.02 4.64]) a =

0.2657 b =

-0.0053 c =

3.5168e-005

（3）拟合的图形，即上一题显示的图像

图2 拟合函数3

2ct bt at y ++=效果

（4），用于这种非线性模型的拟合

①：把其化作线性：b t y ?=α→两边同时取以e 为底的对数→ t b y ln ln )ln(?+=α

②：重复上面第二题的步骤进行，其中需要强调的是（0,0）的点需要另外输入，因为)0ln(不存在，为了在同图出现，故对第二条拟合函数，取b n a m ==

程序：function [m,n,a,b,c]=fei(y2,y0)%y2=y0除了0以外的数

y1=log(y2);

t1=[5:5:55];

n=length(t1);

k1=ones(1,n);

k2=log(t1);

A=[sum(k1.*k1) sum(k2.*k1);sum(k1.*k2) sum(k2.*k2)];

B=[sum(k1.*y1);sum(k2.*y1)];

x=pinv(A)*B;

m=exp(x(1,1));

n=x(2,1);

t=0:55;

y=m*t.^n;

plot(t,y,'-')

hold on

[a,b,c]=xian(y0)

plot(t,y,'--')

hold on

plot(t1,y2,'*',0,0,'*')

legend('y=m*t.^n拟合效果','y=a*t+b*t^2+c*t^3拟合效果','真实值')

得到的拟合图像：

图3 两种拟合函数拟合效果对比

结论：在实际生活当中，不免需要对一组数据进行拟合，通过采用最佳的拟

合，找到一个近似的函数来研究数据的共性。通过这一道题目，发现不同的函数，拟合效果差别也是蛮大的。

第2题：特征值与特征向量

用幂法求下列矩阵的主特征值与相应的特征向量

(1) ??????????=3616415936421A (2) ??

??

?

?????--=1333643432A 解：利用幂法求矩阵A 的主特征值与相应的特征向量，首先要给一个初始向量：

①：定义一个和A 行数一致的1列全一矩阵，即1

01111??????

???????????=n M ν

②：010νν?=A ,为了方便计算，减少计算量，需要求出1ν中按模最大的那个分量的值1P ,同时得到向量1

1

1u P v =

,由此可知00v u = ③：重复第二步，得到22v P

④：计算12P P -，看其是否大于给定的误差，如果大于，则令21v v =，21P P =并从第三步开始重复；如果小于，则2P 为所求的按模最大的特征值，即主特征值，2v 为其对应的特征向量

流程图：

表1 幂法求主特征值的流程图

程序：function [v2,p2]=maxtr(A,err) n=size(A,2); v0=ones(n,1); v1=A*v0; p1=max(v1); u1=v1/p1; v2=A*u1; p2=max(v2); i=1;

while i>0

if abs(p2-p1)<=err break end v1=v2; p1=p2; u1=v1/p1; v2=A*u1; p2=max(v2); i=i+1; end

fprintf('求得矩阵A 的按模最大的特征值p2=:%10.8f\n 对应的特征向量为v2=\n',p2) disp(v2)

对问题（1），??

??

?

?????=3616415936421A ，误差=1err 0.0001 运行结果：[v2,p2]=maxtr([2 4 6;3 9 15;4 16 36],0.0001);

求得矩阵A 的按模最大的特征值p2=:43.87998819 对应的特征向量为v2= 8.1559 19.5719 43.8800

对问题（2），??

??

?

?????--=1333643432A ，误差=2err 0.0001 运行结果：[v2,p2]=maxtr([3 -4 3;-4 6 3;3 3 1],0.0001);

求得矩阵A 的按模最大的特征值p2=:8.86979621 对应的特征向量为v2= -5.3602 8.8698

1.3380

结论：对一个向量进行单位化，可以为下面的计算节省很多的计算步骤。我用到的程序，就是循环的赋予1ν，1P 的值，并通过1ν，1P 求出22v P ，最后达到误差要求，就求出了矩阵的主特征值与相应的特征向量。

第3题：微分方程求解

用四阶龙格-库塔法求解初值问题

??

?≤≤=-=')50(2

)0(2x y xy y

取步长h=0.25，并画出 50≤≤x 数值解的曲线。

解：对于四阶龙格-库塔法，我们在求第k+1个点的值时，必须知道第k 个点的值，在通过中间的一些变量，使结果更加精确。 （1）：由题意可以确定初始值，即20

00==y x

（2）：分别计算出：???????++?=++?=++?=?=)

,(),()

,()

,(34223222121k y h x f h k y x f h k y x f h k y x f h k n n k n h

n k n h n n n （3）:h x x k k k k y y n n n n +=++++=++14321611)

22(,由0=n 开始，直至5=x

程序：function [x,y]=rk(x0,y0,a,b,h)%x0,y0为方程的初值 x(1)=x0; y(1)=y0; n=(b-a)/h; for j=1:n

k1=h*f(x(j),y(j));

k2=h*f(x(j)+h/2,y(j)+k1/2); k3=h*f(x(j)+h/2,y(j)+k2/2); k4=h*f(x(j)+h,y(j)+k3);

y(j+1)=y(j)+1/6*(k1+2*k2+2*k3+k4); x(j+1)=x(j)+h; end

plot(x,y,'r:') hold on

xx=[0:0.2:5]; yy=2./(xx.^2+1); plot(xx,yy,'g--')

legend('h=0.25的数值解','解析解')

function z=f(x,y)

z=-x*y^2;

运行结果：[x,y]=rk(0,2,0,5,0.25)

x =

Columns 1 through 6

0 0.2500 0.5000 0.7500 1.0000 1.2500

Columns 7 through 12

1.5000 1.7500

2.0000 2.2500 2.5000 2.7500

Columns 13 through 18

3.0000 3.2500 3.5000 3.7500

4.0000 4.2500

Columns 19 through 21

4.5000 4.7500

5.0000

y =

Columns 1 through 6

2.0000 1.8823 1.5999 1.2799 1.0000 0.7806

Columns 7 through 12

0.6155 0.4924 0.4001 0.3299 0.2759 0.2336

Columns 13 through 18

0.2000 0.1730 0.1510 0.1328 0.1177 0.1049

Columns 19 through 21

0.0941 0.0849 0.0769

结论：四阶龙格-库塔的精度为四阶，相对而言，拟合的效果非常好。但在计算中，计算量比较大，借助计算机才可能实现对数据的处理。

第4题：求方程的根

分别应用二分法和牛顿法求解下面方程的根 1sin =x x

要求满足精度4*102

1

-?<-k x x

解：（1）二分法：

①：取上下限a 、b,使0)(()(

b

a c +=

③：看err a b >-是否成立，成立则继续进行，不成立，则直接进行第四步。计算)(c f ，并判断是否有0)()(

④：输出数值解b

流程图：

表2 二分法流程图

程序：function erfen(err) a=1; b=2; n=0;

while abs(b-a)>=err n=n+1; c=(a+b)/2;

if (a*sin(a)-1)*(c*sin(c)-1)<0 b=c; else a=c; end

fprintf('第%2d 次迭代值为:%10.8f\n',n,a) end

fprintf('第%2d 次求得的近似值x*=:%10.8f\n',n,a)

运行结果：erfen(0.00005) 第 1次迭代值为:1.00000000 第 2次迭代值为:1.00000000 第 3次迭代值为:1.00000000 第 4次迭代值为:1.06250000 第 5次迭代值为:1.09375000 第 6次迭代值为:1.10937500 第 7次迭代值为:1.10937500 第 8次迭代值为:1.11328125 第 9次迭代值为:1.11328125 第10次迭代值为:1.11328125 第11次迭代值为:1.11376953 第12次迭代值为:1.11401367 第13次迭代值为:1.11413574 第14次迭代值为:1.11413574 第15次迭代值为:1.11413574

第15次求得的近似值x*=:1.11413574 （2），牛顿法：

①：

2,1,0)

()(0'1=??

???

-=+k x x f x f x x k k

k k

②：由0=k 开始，不断计算1+k x ，直到err x x k k <-+1，最后的数值解为1+k x

程序：function [x1,n]=mynewton(x0,err) x1=x0-f(x0)/df(x0);

n=1;

fprintf('第%2d次迭代值为:%10.8f\n',n,x1)

while abs(x1-x0)>=err

x0=x1;

x1=x0-f(x0)/df(x0);

n=n+1;

fprintf('第%2d次迭代值为:%10.8f\n',n,x1)

end

fprintf('第%2d次求得的近似值x*=:%10.8f\n',n,x1)

function z=f(x)

z=x*sin(x)-1;

function z=df(x)

z=sin(x)+x*cos(x);

运行结果：[x1,n]=mynewton(1,0.00005);

z =

-0.1585

z =

1.3818

第 1次迭代值为:1.11472867

z =

7.9373e-004

z =

1.3887

第 2次迭代值为:1.11415713

z =

-1.9493e-008

z =

1.3888

第 3次迭代值为:1.11415714

第 3次求得的近似值x*=:1.11415714

结论：对同一题而言，采用不同的迭代方法，在相同的精度要求下，得到的近似解是相差不大的！但，不同的迭代，所需要的计算量是相差很大的，而且在特定的一些场合，一部分迭代式不合适的，所以，在面对同一问题的求解时，需要考虑多种因素。

云模型matlab程序

1.绘制云图 Ex=18 En=2 He=0.2 hold on for i=1:1000 Enn=randn(1)*He+En; x(i)=randn(1)*Enn+Ex; y(i)=exp(-(x(i)-Ex)^2/(2*Enn^2)); plot(x(i),y(i),'*') end Ex=48.7 En=9.1 He=0.39 hold on for i=1:1000 Enn=randn(1)*He+En; x(i)=randn(1)*Enn+Ex; y(i)=exp(-(x(i)-Ex)^2/(2*Enn^2)); plot(x(i),y(i),'*')

end 2.求期望、熵及超熵 X1=[51.93 52.51 54.70 43.14 43.85 44.48 44.61 52.08]; Y1=[0.91169241573 0.921875 0.96032303371 0.75737359551 0.76983848315 0.7808988764 0.78318117978 0.9143258427]; m=8; Ex=mean(X1) En1=zeros(1,m); for i=1:m En1(1,i)=abs(X1(1,i)-Ex)/sqrt(-2*log(Y1(1,i))); end En=mean(En1); He=0; for i=1:m He=He+(En1(1,i)-En)^2; end En=mean(En1) He=sqrt(He/(m-1)) 3.平顶山so2环境: X1=[0.013 0.04 0.054 0.065 0.07 0.067 0.058 0.055 0.045]; Y1=[0.175675676 0.540540541 0.72972973 0.878378378

DEA的Matlab程序(数据包络分析)

模型((P C2R)的MATLAB程序 clear X=[]; %用户输入多指标输入矩阵X Y=[]; %用户输入多指标输出矩阵Y n=size(X',1); m=size(X,1); s=size(Y,1); A=[-X' Y']; b=zeros(n, 1); LB=zeros(m+s,1); UB=[]; for i=1:n; f= [zeros(1,m) -Y(:,i)']; Aeq=[X(:,i)' zeros(1,s)]; beq=1; w(:,i)=LINPROG(f,A,b,Aeq,beq,LB,UB); %解线性规划，得DMU;的最佳权向量w; E(i, i)=Y(:,i)'*w(m+1:m+s,i); %求出DMU i的相对效率值E ii end w %输出最佳权向量 E %输出相对效率值E ii Omega=w(1:m,:) %输出投入权向量。 mu=w(m+1:m+s,:) %输出产出权向量。 模型(D C2R)的MATLAB程序 clear X=[]; %用户输入多指标输入矩阵X Y=[]; %用户输入多指标输出矩阵Y n=size(X',1); m=size(X,1); s=size(Y,1); epsilon=10^-10; %定义非阿基米德无穷小 =10-10 f=[zeros(1,n) -epsilon*ones(1,m+s) 1]; %目标函数的系数矩阵： 的系数为0,s-,s+的系数为- e, 的系数为1； A=zeros(1,n+m+s+1); b=0; %<=约束； LB=zeros(n+m+s+1,1); UB=[]; %变量约束； LB(n+m+s+1)= -Inf; %-Inf表示下限为负无穷大。 for i=1:n; Aeq=[X eye(m) zeros(m,s) -X(:,i) Y zeros(s,m) -eye(s) zeros(s,1)]; beq=[zeros(m, 1 ) Y(:,i)]; w(:,i)=LINPROG (f,A,b,Aeq,beq,LB,UB); %解线性规划，得DMU的最佳权向量w; end w %输出最佳权向量 lambda=w(1:n,:) %输出 s_minus=w(n+1:n+m,:) %输出s- s_plus=w(n+m+1:n+m+s,:) %输出s+ theta=w(n+m+s+1,:) %输出

云模型简介及个人理解matlab程序

云模型简介及个人理解m a t l a b程序 集团档案编码：[YTTR-YTPT28-YTNTL98-UYTYNN08]

随着不确定性研究的深入，越来越多的科学家相信，不确定性是这个世界的魅力所在，只有不确定性本身才是确定的。在众多的不确定性中，和是最基本的。针对和在处理不确定性方面的不足，1995年我国工程院院士教授在概率论和模糊数学的基础上提出了云的概念，并研究了模糊性和随机性及两者之间的关联性。自李德毅院士等人提出云模型至今，云模型已成功的应用到、、、智能控制、等众多领域. 设是一个普通集合。 , 称为论域。关于论域中的模糊集合，是指对于任意元素都存在一个有稳定倾向的随机数，叫做对的隶属度。如果论域中的元素是简单有序的，则可以看作是基础变量，隶属度在上的分布叫做隶属云；如果论域中的元素不是简单有序的，而根据某个法则，可将映射到另一个有序的论域上，中的一个且只有一个和对应，则为基础变量，隶属度在上的分布叫做隶属云[1] 。 数字特征

云模型表示自然语言中的基元——语言值，用云的数字特征——期望Ex，熵En和超熵He表示语言值的数学性质 [3] 。 期望 Ex：云滴在论域空间分布的期望，是最能够代表定性概念的点，是这个概念量化的最典型样本。 熵 En：“熵”这一概念最初是作为描述热力学的一个状态参量，此后又被引入统计物理学、信息论、复杂系统等，用以度量不确定的程度。在云模型中，熵代表定性概念的可度量粒度，熵越大，通常概念越宏观，也是定性概念不确定性的度量，由概念的随机性和模糊性共同决定。一方面, En是定性概念随机性的度量，反映了能够代表这个定性概念的云滴的离散程度；另一方面，又是定性概念亦此亦彼性的度量，反映了在论域空间可被概念接受的云滴的取值范围。用同一个数字特征来反映随机性和模糊性，也必然反映他们之间的关联性。 超熵 He：熵的不确定性度量，即熵的熵，由熵的随机性和模糊性共同决定。反映了每个数值隶属这个语言值程度的凝聚性，即云滴的凝聚程度。超熵越大，云的离散程度越大，隶属度的随机性也随之增大，云的厚度也越大。

云模型简介及个人理解matlab程序文件

随着不确定性研究的深入，越来越多的科学家相信，不确定性是这个世界的魅力所在，只有不确定性本身才是确定的。在众多的不确定性中，随机性和模糊性是最基本的。针对概率论和模糊数学在处理不确定性方面的不足，1995年我国工程院院士李德毅教授在概率论和模糊数学的基础上提出了云的概念，并研究了模糊性和随机性及两者之间的关联性。自李德毅院士等人提出云模型至今，云模型已成功的应用到自然语言处理、数据挖掘、 设是一个普通集合。 , 称为论域。关于论域中的模糊集合，是指对于任意元素都存在一个有稳定倾向的随机数，叫做对的隶属度。如果论域中的元素是简单有序的，则可以看作是基础变量，隶属度在上的分布叫做隶属云；如果论域中的元素不是简单有序的，而根据某个法则，可将映射到另一个有序的论域上，中的一个且只有一个和对应，则为基础变量，隶属度在上的分布叫做隶属云[1] 。 数字特征 云模型表示自然语言中的基元——语言值，用云的数字特征

——期望Ex，熵En和超熵He表示语言值的数学性质[3] 。 期望 Ex：云滴在论域空间分布的期望，是最能够代表定性概念的点，是这个概念量化的最典型样本。 熵 En：“熵”这一概念最初是作为描述热力学的一个状态参量，此后又被引入统计物理学、信息论、复杂系统等，用以度量不确定的程度。在云模型中，熵代表定性概念的可度量粒度，熵越大，通常概念越宏观，也是定性概念不确定性的度量，由概念的随机性和模糊性共同决定。一方面, En是定性概念随机性的度量，反映了能够代表这个定性概念的云滴的离散程度；另一方面，又是定性概念亦此亦彼性的度量，反映了在论域空间可被概念接受的云滴的取值范围。用同一个数字特征来反映随机性和模糊性，也必然反映他们之间的关联性。 超熵 He：熵的不确定性度量，即熵的熵，由熵的随机性和模糊性共同决定。反映了每个数值隶属这个语言值程度的凝聚性，即云滴的凝聚程度。超熵越大，云的离散程度越大，隶属度的随机性也随之增大，云的厚度也越大。 1.绘制云图 Ex=18

基于云模型的粒计算方法研究

第6章从云模型理解模糊集合的争论与发展

第1章基于云模型的粒计算方法应用 云模型是一个定性定量转换的双向认知模型，正向高斯云和逆向高斯云算法实现了一个基本概念与数据集合之间的转换关系；本文基于云模型和高斯变换提出的高斯云变换方法给出了一个通用的认知工具，不仅将数据集合转换为不同粒度的概念，而且可以实现不同粒度概念之间的柔性切换，构建泛概念树，解决了粒计算中的变粒度问题，有着广阔的应用前景。 视觉是人类最重要的感觉，人类所感知的外界信息至少有80%以上都来自于视觉[130]。图像分割[131]是一种最基本的计算机视觉技术，是图像分析与理解的基础，一直以来都受到人们的广泛关注。目前图像的分割算法有很多，包括大大小小的改进算法在内不下千种，但大致可以归纳为两类[132]。第一类是采用自顶向下的方式，从数学模型的选择入手，依靠先验知识假定图像中的部分属性特征符合某一模型，例如马尔科夫随机场、引力场等，利用模型描述图像的邻域相关关系，将图像低层的原始属性转换到高层的模型特征空间，进而建模优化求解所采用模型的参数，通常是一个复杂度非常高的非线性能量优化问题。在特征空间对图像建模，其描述具有结构性、分割结果也一般具有语义特征，但是由于对数据的未知性、缺乏足够先验知识的指导，导致模型的参数选择存在一定的困难。第二类是采用自底向上的方式，从底层原始数据入手，针对图像灰度、颜色等属性采用数据聚类的方法进行图像分割，聚类所采用的理论方法通常包括高斯变换、模糊集、粗糙集等；或者预先假设图像的统计特性符合一定的分类准则，通过优化准则产生分割结果，例如Otsu方法的最大方差准则[133][134]、Kapur方法的最大熵准则[135][136]等。这类方法虽然缺乏语义信息表达，但是直接在数据空间建模，方法更具普适性和鲁棒性。 随着计算机视觉研究的深入，简单的图像分割已经不能满足个性化的需求，有时候人们恰恰兴趣的是图像中亦此亦彼的那些不确定性区域，基于云模型的粒计算方法是一种不确定性计算方法，发现图像中存在的不确定性区域是它的一个重要能力。如何模拟人类自然视觉中的认知能力进行图像分割一直以来都是一个难点问题，而基于高斯云变换的可变粒计算正是用来模拟人类认知中的可变粒计算过程，因此可以利用高斯云变换对自然视觉认知能力中选择性注意能力进行形式化。武汉大学秦昆教授等曾基于云综合、云分解等云运算实现图像分割，正如第5章中的分析结果，基于内涵的概念计算方法随着层次的提升，概念脱离原始数据会增加误分率，甚至失效，而且无法实现自适应地概念数量和粒度优化。

正向云发生器代码(matlab)

正向云发生器matlab代码 %正向云算法：由数字特征到定量数据表示 %直接在程序中固定EX/EN/HE的值 Ex=0; En=1; He=0.2; n=2000; X = zeros(1,n); %产生一个1*n型矩阵，其元素都为0 Y = zeros(1,n); X= normrnd ( En, He, 1, n); %产生一个1*n型正态随机数矩阵，EX为期望，ENN为方差for i=1:n Enn=X(1,i); X(1, i) = normrnd ( Ex, Enn, 1) ; %产生一个正态随机数，EX为期望，ENN为方差(1*1型) Y(1, i) = exp ( - (X(1, i) - Ex) ^2 / (2* Enn^2) ) ; end plot(X(1,:),Y(1,:),'r.'); %画图语句 %倘若X(1,i)是确定的随机数时，本代码是自己输入确定值 %保存为.m文件时，文件名要是字母名，不要中文名 disp('- - - - -云发生器程序开始- - - - -'); Ex = input('输入期望值Ex：'); En = input('输入熵值En：'); He = input('输入超熵值He：'); n = input('输入需重复计算次数：'); X = zeros(1,n); %产生一个1*n型矩阵，其元素都为0 Y = zeros(1,n); X= normrnd ( En, He, 1, n); %产生一个1*n型正态随机数矩阵，EX为期望，He为方差Xi = input('输入随机数X（1，i）:'); %手动输入固定随机数X for i=1:n

云模型实现图形-MATLAB程序

一维云模型 程序： clc clear Ex=170;En=5;He=0.5; n=5000; for i=1:n Enn=randn(1)*He+En; x(i)=randn(1)*Enn+Ex; y(i)=exp(-(x(i)-Ex)^2/(2*Enn^2)); end plot(x,y,'.r') title('5000个男生身高的一维云图') ylabel('确定度'); xlabel('身高值'); axis([150,190,0,1]) grid on 一维： clear vars;clc;close all; Ex1=-8; En1=0.7; He1=0.2; n1=200; Ex2=2.2; En2=2; He2=0.5; n2=800; Ex3=18; En3=4; He3=0.7; n3=1500; En1_t = normrnd(En1,He1,n1,1); data1 = normrnd(Ex1,En1_t,n1,1);

mu1 = exp(-0.5*((data1-Ex1)./En1_t).^2); En2_t = normrnd(En2,He2,n2,1); data2 = normrnd(Ex2,En2_t,n2,1); mu2 = exp(-0.5*((data2-Ex2)./En2_t).^2); En3_t = normrnd(En3,He3,n3,1); data3 = normrnd(Ex3,En3_t,n3,1); mu3 = exp(-0.5*((data3-Ex3)./En3_t).^2); figure(1); plot(data1,mu1,'.b',data2,mu2,'*r',data3,mu3,'+k'); axis equal; 二维云模型 程序： clc clear Ex1=170;En1=5;He1=0.5; Ex2=65;En2=3;He2=0.2; n=5000; for i=1:n

灰色系统预测GM(1,1)模型及其Matlab实现

灰色系统预测GM(1,1)模型及其Matlab 实现 预备知识 (1)灰色系统 白色系统是指系统内部特征是完全已知的；黑色系统是指系统内部信息完全未知的；而灰色系统是介于白色系统和黑色系统之间的一种系统，灰色系统其内部一部分信息已知，另一部分信息未知或不确定。 (2)灰色预测 灰色预测，是指对系统行为特征值的发展变化进行的预测，对既含有已知信息又含有不确定信息的系统进行的预测，也就是对在一定范围内变化的、与时间序列有关的灰过程进行 预测。尽管灰过程中所显示的现象是随机的、杂乱无章的，但毕竟是有序的、有界的，因此得到的数据集合具备潜在的规律。灰色预测是利用这种规律建立灰色模型对灰色系统进行预测。 目前使用最广泛的灰色预测模型就是关于数列预测的一个变量、一阶微分的GM(1,1)模型。它是基于随机的原始时间序列，经按时间累加后所形成的新的时间序列呈现的规律可用一阶线性微分方程的解来逼近。经证明，经一阶线性微分方程的解逼近所揭示的原始时间序列呈指数变化规律。因此，当原始时间序列隐含着指数变化规律时，灰色模型GM(1,1)的预测是非常成功的。 1 灰色系统的模型GM(1,1) 1.1 GM(1,1)的一般形式 设有变量X (0)＝{X (0)(i)，i=1,2，...，n}为某一预测对象的非负单调原始数据列，为建立灰色预测模型：首先对X (0)进行一次累加(1—AGO, Acumulated Generating Operator)生成一次累加序列： X (1)＝{X (1)(k )，k ＝1，2，…，n} 其中 X (1) (k )＝ ∑ =k i 1 X (0)(i) ＝X (1)(k －1)+ X (0)(k ) (1) 对X (1)可建立下述白化形式的微分方程： dt dX )1(十) 1(aX ＝u (2) 即GM(1,1)模型。 上述白化微分方程的解为(离散响应)： ∧ X (1)(k +1)＝(X (0)(1)－ a u )ak e -＋a u (3) 或 ∧ X (1)(k )＝(X (0)(1)－ a u ))1(--k a e ＋a u (4)

-word版本hslogic_利用云模型估计车速

1.在隧道内放置多个地感线圈（间距相同），车辆通过时、对通过的相邻两个线圈（或单线圈）的脉冲信号数据进行实时采集，首先利用云模型算法（正向云与逆向云算法结合）或其他，计算得到车速的估计值，将车速估计结果与行车时间作为车辆行驶位置判定的依据，再采用云推理得到车辆行驶位置的估计值，对所估计的结果验证，实现对车辆位置的实时精确估计。 2.最后还想验证一下估计结果的准确度 3.做一下参数寻优，对比结果 4.补充要求：我想用脉冲频率波形数据哦，因为有个原始波形的图更好的。。 （我当时说给你一篇论文，用那个上面的数据。。不知道能不能用。。。见附带的论文） 二、课题解决思路简介 基于"隧道内放置多个地感线圈"，主要是汽车通过多个线圈，产生不同时刻的脉冲，然后计算每个脉冲之间的时间差，来获得车速的计算。 这里，我们主要需要的数据时每个线圈之间的距离参数以及每个脉冲之间的时间间隔，然后我们通过云模型来算法得到车速的估计值。根据得到的车速，我们可以得到最后的位置。 然后，我们可以根据论文最后一章的分析方法来分析最后结果的准确度。 对于参数优化，主要是针对云模型的初始参数，我们使用随机数，然后通过PSO进行迭代优化，从而获得最佳的参数，并估算得到最佳的值。 最后，将普通算法得到的结果和PSO优化之后的结果进行对比，从而验证优化算法的优势。 最后，你需要的是脉冲频率波形，这个，我们在设计的时候，进行处理，可以保证。 三、课题设计介绍和仿真说明 3.1正向云和逆向云 首先介绍一下基本的云模型，正向云和逆向云，其基本的理论如下所示： 云模型的发生器就是指云的生成算法，发生器的形式可以有很多种，一般都采用软件的形式加以实现。云的发生器大体上可以分为正向云发生器和逆向云发生器。正向发生器是指从定性到定量之间的转换模型，即由云的三个数字特征产生云滴的具体过程。图1为正向云模型发生器示意图。

Matlab编程实例

Matlab非线性方程求解器fsolve总结 fsolve是采用最小二乘法来求解非线性方程。它的一般求解方式为： X=FSOLVE(FUN,X0,OPTIONS) 其中，fun是要求解的非线性方程，X0是变量初值，options由optimset函数产生的结构体，用于对优化参数的设置，可以省略（采用默认值）。 Fsolve可以求解简单的一维非线性方程，如： x = fsolve(@myfun,[0.5 2 4],optimset('Display','iter')); %求解在初值分别为0.5，2和4时方程的解 其中，函数myfun的定义为： function F = myfun(x) F = sin(x); Fsolve还可以求解大型的非线性方程组，如 x0 = [51.6;rand;unifrnd(-1,1);rand]; h=optimset; h.MaxFunEvals=20000; h.MaxIter=5000; h.Display='off'; [p,fval] = fsolve(@f,x0,options); 此时，方程组可以写成矩阵形式： function F=f(x) F=[x(1)+x(2)*(1-exp(-(x(3)*(0)^x(4))))-51.61; x(1)+x(2)*(1-exp(-(x(3)*(9.78)^x(4))))-51.91; x(1)+x(2)*(1-exp(-(x(3)*(30.68)^x(4))))-53.27; x(1)+x(2)*(1-exp(-(x(3)*(59.7)^x(4))))-59.68;]; 自编基于龙格库塔法的Matlab数值积分函数 function varargout = rkkt(varargin) %============================================== % 采用四阶龙格库塔法数值积分 % rkkt(F,x0,t0,tfinal,ts) % F为函数名 % x0为积分变量的初始值 % t0为积分初始时刻 % tfinal为积分终止时刻 % ts为积分步长 % example，如果有一个函数名为foo,则，求解指令为 % rkkt(@foo,x0,t0,tfinal,tspan);