zoutendijk_可行方向法的matlab实现

2012-2013（1）专业课程实践论文

可行方向法

夏文春，0818180226，R数学08-2班

姜皓缤，0618180213，R数学08-2班

肖云龙，0818180214，R数学08-2班

算法理论

(一) 、基本思想是：

给定一个可行点)(k x 之后，用某种方法确定一个改进的可行方向k d ，然后沿方向

k d ，求解一个有约束的线搜索问题，得极小点k λ，令k

k k k d x

x

λ+=+)

()

1(，如果)

1(+k x 不是最优解，则重复上述步骤。可行方向法就是利用线性规划方法来确定k d 的。

1) 、线性约束问题：

设x 是问题

???

??∈=≤n R

f x e Ex b,Ax x s.t.)

( min 的一个可行解，假定11b x A =，22b x A <，其中

??????=21A A A ，??

?

???=21b b b 则一个非零向量d 是在点x 点的一个可行方向，当且仅当0d A 1≤，0=Ed ；如果0(T

(二) 、算法

1) 、线性不等式约束的Zoutendijk 方法的计算步骤：

1.求一初始可行解0x 。，令k ＝1，转2。

2.对于可行点k x ，设11k A x b =，22k A x b <，12(,)T T T A A A =， 12(,)T T T b b b =求解

问题()11m in . 0 0

11,1,,T j z f x d s t A d P Ed d j n

?=??

≤??=??-≤≤=?

，得最优解k d ，如果()T

k f x d ?=0，计算结束，

k x 是K —T 点；否则转3。

3.求解线搜索问题

max m in ()

. 0k k f x d s t λλλ+??

≤≤?

(a) 其中

{}____m ax

_22_2m in |0,0

,,0

k k b d d d b b A x d A d d λ?>≤≠?==-=?+∞≤?? 设k λ为(a)式最优解，令1,1k k k k x x d k k λ+=+=+，返回2。

二、算法框图

五 流程图

Y

Y

N

Y

N

begin

误差21,εε初始点

,0=

k X Φ

=)(k

X J ？

1

2

)

(ε≤k

X f

over k X

m in

()(),()11,1,2,...k T k T k j i f X g X j J X d i n

ηηη??≤??-≤∈??-≤≤=?

2k ηε≤

over

{}0:min (max (0,1,2,...k k

k k

k

j f X D g X D j n λλ

λλλλλ≤≤+=+≥=）

）

1

:1

k k

k

X

X

D

k k λ+=+=+

算法代码

1)、主程序

function [X0,f_val]=zoutendijk(A,b,x0,Aeq,beq,label)

%自定义函数diff_val(x0)作用是求所给函数在x0出的偏导数%自定义函数fval(x0)作用是求所给函数在x0出的函数值format long;

eps=1.0e-6;

x0=transpose(x0);%刚开始给的x0为行向量

func

sz=length(x0);

if label==1

[m,n]=size(A);

%把A分解为A1,A2，其中A1为起作用约束

for k=1:1:100

A1=A;

A2=A;

b1=b;

b2=b;

for i=m:-1:1

if abs(A2(i,:)*x0-b2(i,:)) < 0.1

A2(i,:)=[];

b2(i,:)=[];

end

end

for i=m:-1:1

if abs(A1(i,:)*x0-b1(i,:))>=0.1

A1(i,:)=[];

b1(i,:)=[];

end

end

A1;

A2;

b1;

b2;

i2=rank(A2);

AE=[A1;Aeq];

[i1,j1]=size(AE);

r=rank(AE);

if r

'行不满秩'

return

end

if i2==0

'无效'

return

end

%求解线性规划问题得到可行下降方向d0

s=diff_val(x0);

c=double(s);

lb=-1*ones(sz,1);

ub=ones(sz,1);

k1=length(b1);

k2=length(beq);

p=zeros(k1,1);

q=zeros(k2,1);

[d0,mn,m1,m2,m3]=linprog(c,A1,p,Aeq,q,lb,ub);

d0;mn;

df=abs(s*d0);

if df<0.1

'最优解为：'

'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@

@@@@@@@@@@@@@@@@@@@@@'

x0

f_val=fval(x0)

k

'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ @@@@@@@@@@@@@@@@@@@@@'

return

else

%进行一维搜索，求f(x(k+1))的最小值

b_=b2-A2*x0;

d_=A2*d0;

[dh,dl]=size(d_);

ul=1;

for i=1:1:dh

if d_(i,:)>=0

u=1;

else

u=0;

end

ul=ul*u;

end

ul;b_;d_;

vmax=inf;

if ul==0

vmax=inf;

else

for i=1:1:dh

if d_(i,:)>0

v=b_(i,:)/d_(i,:);

if v

vmax=v;

end

end

end

end

end

vmax;

h=fmin(x0,d0,vmax);

a=x0+h*d0;

f_val=fval(a);

x0=x0+h*d0;

'****************'

X0=x0

f_val=fval(x0)

k

'****************'

end

end

if label==0

for k=1:1:100

%确定指标集

'f(x)梯度：'

sf=diff_val(x0)

sf=eval(sf)

'f(x)梯度长度：'

norm_s=norm(sf)

GL=length(G);

G_copy=G;

for i=1:1:GL

g=subs(G(i,:),x,x0);

G(i,:)=g;

end

G_zero=eval(G);

for i=GL:-1:1

if abs(G_zero(i,:))>0.1

G_zero(i,:)=[];

G_copy(i,:)=[];

I=length(G_zero);

end

end

'x0时为零的g(x)：'

G_copy

'指标集I(x)：'

I

add=-ones(I,1);

%根据指标集确定不同情况下的搜索方向

if I==0

if norm_s<=10

'最优解为：'

'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ @@@@@@@@@@@@@@@@@@@@@'

x0

f_val=fval(x0)

k

'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ @@@@@@@@@@@@@@@@@@@@@'

return

else

d0=-sf;

%线搜索问题

vmax=100;

d0=transpose(d0);

h=fmin(x0,d0,vmax);

x0=x0+h*d0;

end

else

%线性规划问题

grad=jacobian(G_copy,x);

G_zero=subs(grad,x,x0);

G_zero=[G_zero,add];

sf=[sf,-1];

'线性规划问题A矩阵：'

Ac=[sf;G_zero]

lb=-1*ones(sz,1);

ub=ones(sz,1);

p=zeros(I+1,1);

c=zeros(1,sz);

c=[c,1];

[dz,mn,m1,m2,m3]=linprog(c,Ac,p,[],[],lb,ub);

dz;

mn;

sd=length(dz);

d1=dz(1:sd-1,1:1);

z0=dz(sd,1);

z0=abs(z0)

if z0<0.01

'最优解为：'

'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ @@@@@@@@@@@@@@@@@@@@@'

x0

f_val=fval(x0)

k

'@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@ @@@@@@@@@@@@@@@@@@@@@'

return

else

d0=d1;

%线搜索问题

vmax=10000;

h=fmin(x0,d0,vmax);

a=x0+h*d0;

x0=x0+h*d0;

end

end

k

end

end

2)、回调函数

●func 记录函数及其自变量信息，如：

syms x1 x2;

f=2*x1^2+2*x2^2-2*x1*x2-4*x1-6*x2;

x=[x1,x2];

●fval 计算函数在x0处得函数值

function f_val=fval(x0)

x0=transpose(x0);

func;

f_val=subs(f,x,x0);

●diff_val(x0)计算函数在x0处得导数值function s=diff_val(x0)

func

grad=jacobian(f,x);

s=subs(grad,x,x0);

●fmin(x0,d0,vmax)求函数在[0,vmax]上的最小值function h=fmin(x0,d0,vmax)

func

syms h;

a=x0+h*d0;

f_val=inline(subs(f,x,a));

if vmax==inf

min_h=fminbnd(f_val,0,10000);

else

min_h=fminbnd(f_val,0,vmax);

end

h=min_h;

四、算法实现

()22

12121212121

2m in 22246. 2 55

0f x x x x x x x s t x x x x x x ?=+---?

+≤??

+≤??-≤?

?-≤?

解：输入

A=[1 1;1 5;-1 0;0 -1] b=[2 5 0 0]' x0=[0 0] Aeq=[] beq=[] label=1

最后运行程序：zoutendijk(A,b,x0,Aeq,beq,label)

得到结果：

()22

1212121

2m in 4. 1 151012 0

0f x x x s t x x x x x x ?=+?

+≥??

+≥??≥?

?≥?

解：输入

A=[-1 -1;-15 -10;-1 0;0 -1] b=[-1 -2 0 0]' x0=[0 2] Aeq=[] beq=[] label=1;

最后运行程序：zoutendijk(A,b,x0,Aeq,beq,label) 得到结果：

三次样条插值---matlab实现

计算方法实验—三次样条插值 机电学院075094-19 苏建加 20091002764 题目：求压紧三次样条曲线，经过点（-3,2），（-2,0），（1,3），（4,1），而且一阶导 数边界条件S'(-3)=-1;S'（4）=1。 解：首先计算下面的值： 记 1--=j j j x x h ； 1++=j j j j h h h u ；1=+j j u λ ； ?? ????????---+=-++++-j j j j j j j j j j j h y y h y y h h x x x f 1111 111],,[ ；M j =)(''j x s ;],,[611+-=j j j j x x x f d ; h1=-2-(-3)=1;h2=1-(-2)=3;h3=4-1=3; u1=1/4;u2=3/6; d1=6/4*(3/3-(-2)/1)=4.5;d2=6/6*(-2/3-3/3)=-5/3; 由于边界条件S'(-3)=-1;S'（4）=1,得到如下 式子： d0=6/1*(-2/1-(-1))=-6； d3=6/3*(1-(-2)/3)=10/3; 所以得到4个含参数m0~m3 的线性代数方程组为： 2.0000 1.0000 0 0 m0 0.2500 2.0000 0.7500 0 m1 0 0.5000 2.0000 0.5000 m2 0 0 1.0000 2.0000 m3 利用matlab 求解方程得： m = -4.9032 3.8065 -2.5161 2.9247 所以 S1(x)=-0.8172*(-2-x)^3+ 0.6344*(x+3)^3+2.8172*(-2-x)-0.6344*(x+3) x ∈[-3，-2] S2(x)=0.2115*(1-x)^3 -0.1398*(x+2)^3- 1.9032*(1-x)+ 2.2581*(x+2) x ∈[-2，1] S3(x)=-0.1398*(4-x)^3+0.1625(x-1)^3+ 2.2581*(4-x)-1.1290*(x-1) x ∈[1，4] 化简后得：S1(x)=1.4516*x^3 + 10.6128*x^2 + 23.4836*x + 16.1288 x ∈[-3，-2] S2(x)=-0.3513x^3-0.2043x^2+1.8492x+1.7061 x ∈[-2，1] S3(x)=0.3023x^3-2.1651x^2+3.8108x+1.0517 x ∈[1，4] 画图验证：

时域有限差分法的Matlab仿真

时域有限差分法的Matlab仿真 关键词: Matlab 矩形波导时域有限差分法 摘要：介绍了时域有限差分法的基本原理，并利用Matlab仿真，对矩形波导谐振腔中的电磁场作了模拟和分析。 关键词：时域有限差分法；Matlab；矩形波导；谐振腔 目前，电磁场的时域计算方法越来越引人注目。时域有限差分（Finite Difference Time Domain，FDTD）法［1］作为一种主要的电磁场时域计算方法，最早是在1966年由K. S. Yee提出的。这种方法通过将Maxwell旋度方程转化为有限差分式而直接在时域求解，通过建立时间离散的递进序列，在相互交织的网格空间中交替计算电场和磁场。经过三十多年的发展，这种方法已经广泛应用到各种电磁问题的分析之中。 Matlab作为一种工程仿真工具得到了广泛应用［2］。用于时域有限差分法，可以简化编程，使研究者的研究重心放在FDTD法本身上，而不必在编程上花费过多的时间。 下面将采用FDTD法，利用Matlab仿真来分析矩形波导谐振腔的电磁场，说明了将二者结合起来的优越性。 1FDTD法基本原理 时域有限差分法的主要思想是把Maxwell方程在空间、时间上离散化，用差分方程代替一阶偏微分方程，求解差分方程组，从而得出各网格单元的场值。FDTD 空间网格单元上电场和磁场各分量的分布如图1所示。 电场和磁场被交叉放置，电场分量位于网格单元每条棱的中心，磁场分量位于网格单元每个面的中心，每个磁场（电场）分量都有4个电场（磁场）分量环绕。这样不仅保证了介质分界面上切向场分量的连续性条件得到自然满足，而且

还允许旋度方程在空间上进行中心差分运算，同时也满足了法拉第电磁感应定律和安培环路积分定律，也可以很恰当地模拟电磁波的实际传播过程。 1.1Maxwell方程的差分形式 旋度方程为： 将其标量化，并将问题空间沿3个轴向分成若干网格单元，用Δx，Δy和Δz 分别表示每个网格单元沿3个轴向的长度，用Δt表示时间步长。网格单元顶点的坐标（x，y，z）可记为： 其中：i，j，k和n为整数。 同时利用二阶精度的中心有限差分式来表示函数对空间和时间的偏导数，即可得到如下FDTD基本差分式： 由于方程式里出现了半个网格和半个时间步，为了便于编程，将上面的差分式改写成如下形式：

MATLAB三次样条插值之三弯矩法

MATLAB三次样条插值之三弯矩法 首先说这个程序并不完善，为了实现通用（1，2，…，n）格式解题，以及为调用追赶法程序，没有针对节点数在三个以下的情况进行分类讨论。希望能有朋友给出更好的方法。 首先，通过函数 sanwanj得到方程的系数矩阵，即追赶法方程的四个向量参数，接下来调用 追赶法（在intersanwj函数中），得到三次样条分段函数系数因子，然后进行多项式合并得 到分段函数的解析式，程序最后部分通过判断输入值的区间自动选择对应的分段函数并计算改 点的值。附：追赶法程序 chase %%%%%%%%%%%%%% function [newv,w,newu,newd]=sanwj(x,y,x0,y0,y1a,y1b) % 三弯矩样条插值 % 将插值点分两次输入，x0 y0 单独输入 % 边值条件a的二阶导数 y1a 和b的二阶导数 y1b，这里建议将y1a和y1b换成y2a和 y2b，以便于和三转角代码相区别 n=length(x);m=length(y); if m~=n error('x or y 输入有误，再来'); end v=ones(n-1,1);u=ones(n-1,1);d=zeros(n-1,1); w=2*ones(n+1); h0=x(1)-x0; h=zeros(n-1,1); for k=1:n-1 h(k)=x(k+1)-x(k); end v(1)=h0/(h0+h(1)); u(1)=1-v(1); d(1)=6*((y(2)-y(1))/h(1)-(y(1)-y0)/h0)/(h0+h(1)); % for k=2:n-1 v(k)=h(k-1)/(h(k-1)+h(k)); u(k)=1-v(k); d(k)=6*((y(k+1)-y(k))/h(k)-(y(k)-y(k-1))/h(k-1))/(h(k-1)+h(k)); end newv=[v;1]; newu=[1;u]; d0=6*((y(1)-y0)/h0-y1a)/h0;

各种BP学习算法MATLAB仿真

3.3.2 各种BP学习算法MATLAB仿真 根据上面一节对BP神经网络的MATLAB设计，可以得出下面的通用的MATLAB程序段，由于各种BP学习算法采用了不同的学习函数，所以只需要更改学习函数即可。 MATLAB程序段如下： x=-4:0.01:4; y1=sin((1/2)*pi*x)+sin(pi*x); %trainlm函数可以选择替换 net=newff(minmax(x),[1,15,1],{'tansig','tansig','purelin'},'trainlm'); net.trainparam.epochs=2000; net.trainparam.goal=0.00001; net=train(net,x,y1); y2=sim(net,x); err=y2-y1; res=norm(err); %暂停，按任意键继续 Pause %绘图，原图（蓝色光滑线）和仿真效果图（红色+号点线） plot(x,y1); hold on plot(x,y2,'r+'); 注意：由于各种不确定因素，可能对网络训练有不同程度的影响，产生不同的效果。如图3-8。 标准BP算法（traingd）

图3-8 标准BP算法的训练过程以及结果（原图蓝色线，仿真图+号线）增加动量法(traingdm) 如图3-9。 图3-9 增加动量法的训练过程以及结果（原图蓝色线，仿真图+号线）弹性BP算法（trainrp）如图3-10 图3-10 弹性BP算法的训练过程以及结果（原图蓝色线，仿真图+号线）

动量及自适应学习速率法（traingdx）如图3-11。 图3-11 动量及自适应学习速率法的训练过程以及结果（原图蓝色线，仿真图+号线）共轭梯度法(traincgf)如图3-12。

三次样条插值的Matlab实现(自然边界和第一边界条件)

（第一边界条件）源代码：function y=yt1(x0,y0,f_0,f_n,x)_____________(1) %第一类边界条件下三次样条插值； %xi所求点； %yi所求点函数值； %x已知插值点； %y已知插值点函数值； %f_0左端点一次导数值； %f_n右端点一次导数值； n = length(x0); z = length(y0); h = zeros(n-1,1); k=zeros(n-2,1); l=zeros(n-2,1); S=2*eye(n); fori=1:n-1 h(i)= x0(i+1)-x0(i); end fori=1:n-2 k(i)= h(i+1)/(h(i+1)+h(i)); l(i)= 1-k(i);

end %对于第一种边界条件： k = [1;k];_______________________(2) l = [l;1];_______________________(3) %构建系数矩阵S： fori = 1:n-1 S(i,i+1) = k(i); S(i+1,i) = l(i); end %建立均差表： F=zeros(n-1,2); fori = 1:n-1 F(i,1) = (y0(i+1)-y0(i))/(x0(i+1)-x0(i)); end D = zeros(n-2,1); fori = 1:n-2 F(i,2) = (F(i+1,1)-F(i,1))/(x0(i+2)-x0(i)); D(i,1) = 6 * F(i,2); end %构建函数D： d0 = 6*(F(1,2)-f_0)/h(1);___________(4)

三次样条插值的MATLAB实现

MATLAB 程序设计期中考查 在许多问题中，通常根据实验、观测或经验得到的函数表或离散点上的信息，去研究分析函数的有关特性。其中插值法是一种最基本的方法，以下给出最基本的插值问题——三次样条插值的基本提法： 对插值区间[]b a ,进行划分：b x x x a n ≤

内点法matlab仿真doc资料

编程方式实现： 1.惩罚函数 function f=fun(x,r) f=x(1,1)^2+x(2,1)^2-r*log(x(1,1)-1); 2.步长的函数 function f=fh(x0,h,s,r) %h为步长 %s为方向 %r为惩罚因子 x1=x0+h*s; f=fun(x1,r); 3. 步长寻优函数 function h=fsearchh(x0,r,s) %利用进退法确定高低高区间，利用黄金分割法进行求解h1=0;%步长的初始点 st=0.001; %步长的步长 h2=h1+st; f1=fh(x0,h1,s,r); f2=fh(x0,h2,s,r); if f1>f2 h3=h2+st; f3=fh(x0,h3,s,r); while f2>f3 h1=h2; h2=h3; h3=h3+st; f2=f3; f3=fh(x0,h3,s,r); end else st=-st; v=h1; h1=h2; h2=v; v=f1; f1=f2; f2=v; h3=h2+st; f3=fh(x0,h3,s,r); while f2>f3 h1=h2; h2=h3; h3=h3+st; f2=f3;

f3=fh(x0,h3,s,r); end end %得到高低高的区间 a=min(h1,h3); b=max(h1,h3); %利用黄金分割点法进行求解 h1=1+0.382*(b-a); h2=1+0.618*(b-a); f1=fh(x0,h1,s,r); f2=fh(x0,h2,s,r); while abs(a-b)>0.0001 if f1>f2 a=h1; h1=h2; f1=f2; h2=a+0.618*(b-a); f2=fh(x0,h2,s,r); else b=h2; h2=h1; f2=f1; h1=a+0.382*(b-a); f1=fh(x0,h1,s,r); end end h=0.5*(a+b); 4. 迭代点的寻优函数 function f=fsearchx(x0,r,epson) x00=x0; m=length(x0); s=zeros(m,1); for i=1:m s(i)=1; h=fsearchh(x0,r,s); x1=x0+h*s; s(i)=0; x0=x1; end while norm(x1-x00)>epson x00=x1; for i=1:m s(i)=1; h=fsearchh(x0,r,s);

PID控制算法的matlab仿真

PID 控制算法的matlab 仿真 PID 控制算法就是实际工业控制中应用最为广泛的控制算法,它具有控制器设计简单,控制效果好等优点。PID 控制器参数的设置就是否合适对其控制效果具有很大的影响,在本课程设计中一具有较大惯性时间常数与纯滞后的一阶惯性环节作为被控对象的模型对PID 控制算法进行研究。被控对象的传递函数如下: ()1d s f Ke G s T s τ-= + 其中各参数分别为30,630,60f d K T τ===。MATLAB 仿真框图如图1所示。 图1 2 具体内容及实现功能 2、1 PID 参数整定 PID 控制器的控制参数对其控制效果起着决定性的作用,合理设置控制参数就是取得较好的控制效果的先决条件。常用的PID 参数整定方法有理论整定法与实验整定法两类,其中常用的实验整定法由扩充临界比例度法、试凑法等。在此处选用扩充临界比例度法对PID 进行整定,其过程如下: 1) 选择采样周期 由于被控对象中含有纯滞后,且其滞后时间常数为 60d τ=,故可选择采样周期1s T =。 2) 令积分时间常数i T =∞,微分时间常数0d T =,从小到大调节比例系数K , 使得系统发生等幅震荡,记下此时的比例系数k K 与振荡周期k T 。 3) 选择控制度为 1.05Q =,按下面公式计算各参数:

0.630.490.140.014p k i k d k s k K K T T T T T T ==== 通过仿真可得在1s T =时,0.567,233k k K T ==,故可得: 0.357,114.17,32.62, 3.262p i d s K T T T ==== 0.0053.57 p s i i p d d s K T K T K T K T === = 按此组控制参数得到的系统阶跃响应曲线如图2所示。 01002003004005006007008009001000 0.20.40.60.811.21.41.6 1.8 图2 由响应曲线可知,此时系统虽然稳定,但就是暂态性能较差,超调量过大,且响应曲线不平滑。根据以下原则对控制器参数进行调整以改善系统的暂态过程: 1) 通过减小采样周期,使响应曲线平滑。 2) 减小采样周期后,通过增大积分时间常数来保证系统稳定。 3) 减小比例系数与微分时间常数,以减小系统的超调。 改变控制器参数后得到系统的阶跃响应曲线如图3所示,系统的暂态性能得到明显改善、

实验一 典型环节的MATLAB仿真汇总

实验一 典型环节的MATLAB 仿真 一、实验目的 1．熟悉MATLAB 桌面和命令窗口，初步了解SIMULINK 功能模块的使用方法。 2．通过观察典型环节在单位阶跃信号作用下的动态特性，加深对各典型环节响应曲线的理解。 3．定性了解各参数变化对典型环节动态特性的影响。 二、SIMULINK 的使用 MATLAB 中SIMULINK 是一个用来对动态系统进行建模、仿真和分析的软件包。利用SIMULINK 功能模块可以快速的建立控制系统的模型，进行仿真和调试。 1．运行MATLAB 软件，在命令窗口栏“>>”提示符下键入simulink 命令，按Enter 键或在工具栏单击按钮，即可进入如图1-1所示的SIMULINK 仿真 环境下。 2．选择File 菜单下New 下的Model 命令，新建一个simulink 仿真环境常规模板。 3．在simulink 仿真环境下，创建所需要的系统 三、实验内容 按下列各典型环节的传递函数，建立相应的SIMULINK 仿真模型，观察并记录其单位阶跃响应波形。 ① 比例环节1)(1=s G 和2)(1=s G 实验处理：1)(1=s G SIMULINK 仿真模型

波形图为： 实验处理：2)(1=s G SIMULINK 仿真模型 波形图为： 实验结果分析：增加比例函数环节以后，系统的输出型号将输入信号成倍数放大. ② 惯性环节11)(1+= s s G 和15.01)(2+=s s G 实验处理：1 1 )(1+=s s G SIMULINK 仿真模型

波形图为： 实验处理：1 5.01 )(2+= s s G SIMULINK 仿真模型 波形图为： 实验结果分析：当1 1 )(1+= s s G 时，系统达到稳定需要时间接近5s,当

可行方向法基本算法

可行方向法基本算法 考虑线性规划问题： Min (),{|()0,1,2,...} {n j f X X R E R X g X j l ∈?=≥= 设()k X 是它的一个可行解，但不是要求的极小点，为了求它的极小点或近似极小点，应在()k X 点的可行下降方向中选取某一方向()k D ，并确定步长k λ，使得 (1)()k k k k X X D R λ+=+∈ (1)()()()k k f X f X +< 若满足精度要求，迭代停止，(1)k X +就是所求的点。否则，从(1)k X +出发继续进行迭代。直到满足要求为止。上述这种方法称为可行方向法。 设()k x 点的起作用约束集非空，为求()k x 点的可行下降方向，可由下述不等式组确定响亮D ： ()()()0()0,k T k T i f X D g X D j J ????? ?∈ 这等价于由下面的不等式组求向量D 和实数η： ()()k T f X D η?≤ ()(),i k T g X D j J η-?≤∈ 0η< 现使()()k T f X D ?和()()i k T g X D -?（对所有j J ∈）的最大值极小化（必须同时 限制向量D 的模），即可将上述选取搜索方向的工作，转换为求解下述线性规划问题： Min η ()()k T f X D η?≤ ()()(),()k i k T g X D j J X η-?≤∈ 11,1,2,3,.. i d i n -≤≤=

式中(1,2,3,...,)i d i n =为向量D 的分量。在上式中加入最后一个限制条件，位的是使该线性规划有有限最优解；由于我们的目的在于寻找搜索方向D ，只需知道D 的各分量的相对大小即可。 将上述线性规划的最优解记为()(,)k k D η，如果求出的0k η=，说明在()k X 点不存在可行下降方向，在()()k j g X ?（此处()()k j J X ∈）线性无关的条件下，()k X 为一K-T 点，若解出0k η<，则得到可行下降方向()k D ，这就是我们所要的所搜方向。上述可行方向法的迭代步骤如下： （1）确定允许误差10ε>和20ε>，选取初始近似点(0)X R ∈，并令:0k =。 （2）确定起作用约束指标集()()(){|()0,1}k k j J X j g X j l ==≤≤ ①若()()21()||()||k k J X f X ε=???≤（为空集），而且，停止迭代，得点()k X 。 ②若 ()()21()|()||k k J X f X ε=??>，但，则选搜索方向()()(),k k D f X =-?然后转向第（5）步。 ③若()()k J X ≠?，转下一步 （3）求解线性规划 M i n η ()()k T f X D η?≤ ()()(),()k i k T g X D j J X η-?≤∈ 11,1,2,3,.. i d i n -≤≤= 设它的最优解为()(,)k k D η。 （4）检验是否满足2||k ηε≤，若满足则停止迭代，得到()k X ；否则，以()k D 为搜索方向，并转向下一步。 （5）解下述一维极值问题()()0:()k k k Min f X D λλ λλ-≤≤+，此处 ()()max{|()0,1,2,...,}k k j g X D j l λλλ-=+≥= （6）令 (1)()( k k k k X X D λ+=+ :1k k =+ 转回第（2）步。

LMMSE算法信道均衡MATLAB仿真

一．信道均衡的概念 实际的基带传输系统不可能完全满足无码间串扰传输条件，因而码间串扰是不可避免的。当串扰严重时，必须对系统的传输函数 进行校正，使其达到或接近无码间串扰要求的特性。理论和实践表明，在基带系统中插入一种可调滤波器就可以补偿整个系统的幅频，和相频特性从而减小码间串扰的影响这个对系统校正的过程称为均衡，实现均衡的滤波器称为均衡器。 均衡分为频域均衡和时域均衡。频域均衡是从频率响应考虑，使包括均衡器在内的整个系统的总传输函数满足无失真传输条件。而时域均衡，则是直接从时间响应考虑，使包括均衡器在内的整个系统的冲激响应满足无码间串扰条件。 频域均衡在信道特性不变，且传输低速率数据时是适用的，而时域均衡可以根据信道特性的变化进行调整，能够有效地减小码间串扰，故在高速数据传输中得以广泛应用。 时域均衡的实现方法有多种，但从实现的原理上看，大致可分为预置式自动均衡和自适应式自动均衡。预置式均衡是在实际传数之前先传输预先规定的测试脉冲（如重复频率很低的周期性的单脉冲波形），然后按“迫零调整原理”自动或手动调整抽头增益；自适应式均衡是在传数过程中连续测出距最佳调整值的误差电压，并据此电压去调整各抽头增益。一般地，自适应均衡不仅可以使调整精度提高，而且当信道特性随时间变化时又能有一定的自适应性，因此很受重视。这种均衡器过去实现起来比较复杂，但随着大规模、超大规模集成电路和微处理机的应用，其发展十分迅速。 二．信道均衡的应用 1.考虑如图所示的基带等效数据传输系统，发送信号k x 经过ISI 失真信道传输，叠加高斯加性噪声。 图1基带等效数据传输模型 设发送信号采用QPSK 调制，即(1)k x j =±±ISI 信道的冲击响应以向量的形式表示为h 2211[,,,]T L L L h h h --+=???。典型的ISI 信道响应向量有三种： h [0.04,0.05,0.07,0.21,0.5,0.72,0.36,0,0.21,0.03,0.07]T A =--- h [0.407,0.815,0.407]T B = h [0.227,0.46,0.6888,0.46,0.227]T C = k ω为实部与虚部独立的复高斯白噪声，其均值为零，方差为2 ωσ。 2.实现目的

三次样条插值的Matlab实现(自然边界和第一边界条件)(精)

(第一边界条件源代码: function y=yt1(x0,y0,f_0,f_n,x _____________(1 %第一类边界条件下三次样条插值; %xi 所求点; %yi所求点函数值; %x 已知插值点; %y 已知插值点函数值; %f_0左端点一次导数值; %f_n右端点一次导数值; n = length(x0; z = length(y0; h = zeros(n-1,1; k=zeros(n-2,1; l=zeros(n-2,1; S=2*eye(n; fori=1:n-1 h(i= x0(i+1-x0(i; end fori=1:n-2

k(i= h(i+1/(h(i+1+h(i; l(i= 1-k(i; end %对于第一种边界条件: k = [1;k]; _______________________(2 l = [l;1]; _______________________(3 %构建系数矩阵 S : fori = 1:n-1 S(i,i+1 = k(i; S(i+1,i = l(i; end %建立均差表: F=zeros(n-1,2; fori = 1:n-1 F(i,1 = (y0(i+1-y0(i/(x0(i+1-x0(i; end D = zeros(n-2,1; fori = 1:n-2 F(i,2 = (F(i+1,1-F(i,1/(x0(i+2-x0(i; D(i,1 = 6 * F(i,2;

end %构建函数 D : d0 = 6*(F(1,2-f_0/h(1; ___________(4 dn = 6*(f_n-F(n-1,2/h(n-1; ___________(5 D = [d0;D;dn]; ______________(6 m= S\D; %寻找 x 所在位置,并求出对应插值: fori = 1:length(x for j = 1:n-1 if (x(i<=x0(j+1&(x(i>=x0(j y(i =( m(j*(x0(j+1-x(i^3/(6*h(j+... (m(j+1*(x(i-x0(j^3/(6*h(j+... (y0(j-(m(j*h(j^2/6*(x0(j+1-x(i/h(j+... (y0(j+1-(m(j+1*h(j^2/6*(x(i-x0(j/h(j ; break; else continue; end end end (2 (自然边界条件源代码: 仅仅需要对上面部分标注的位置做如下修改 :

神经网络学习算法matlab仿真

东南大学自动化学院 智能控制概论 神经网络学习算法研究 学院： 姓名： 学号： 日期：

目录 1 任务要求叙述 ..................................................... 错误！未定义书签。 2 系统分析及设计原理 ......................................... 错误！未定义书签。 3 设计实现.............................................................. 错误！未定义书签。4仿真验证.. (6) 5 讨论与分析.......................................................... 错误！未定义书签。

一．任务要求叙述 （1）任务 (a) 运行算法，观察和分析现有学习算法的性能； clear all;close all; nu=20;pi=3.1415926; for i=1:nu p(i)=2*pi*i/nu; t(i)=0.5*(1+cos(p(i))); end minmax=[min(p(:)) max(p(:))] net = newff([ 0 7],[6 1],{'logsig' 'purelin'},'traingd');% traingd traingdm trainlm net.trainParam.epochs = 10000; net.trainParam.goal = 0.0001; net.trainParam.show=200; net.trainParam.lr=0.1; net.trainParam.mc=0.6; %0.9 default value; available for momentum net = train(net,p,t); y1 = sim(net,p); figure(2); plot(p,t,'*-',p,y1,'r--') %************** test data ****************** nu2=nu*3/2; for i=1:(nu2) p2(i)=2*pi*i/(nu2); t2(i)=0.5*(1+cos(p2(i))); end y2 = sim(net,p2); figure(3); plot(t2,'*-');hold on; plot(y2,'r'); xlabel('times');ylabel('outputs'); figure(4); plot(t2-y2); xlabel('times');ylabel('error'); (b) 为了进一步提高学习逼近效果，可以采取那些措施，调节规律如何？根据所提的每种措施，修改算法程序，给出仿真效果验证、过程以及相应的曲线图，给出适当的评述；(c) 联系、结合前向神经网络的算法样本学习、测试等过程，谈谈本人对神经网络系统的一些认识和看法。 （2）要求 提交完整的报告，包括：封面(题目、个人学号姓名等信息)、目录、任务要求叙述、系

优化设计方法学复习资料(修正)

优化设计方法学复习资料 一、数学规划法部分： 1.一维搜索： 当方向k d 给定时，求最佳步长k α就是求一元函数)()()(1k k k k k d x f x f α?α=+=+的极值问题，称为一维搜索。 一维搜索问题的求解方法有两种：解析解法和数值解法。 解析解法思路：利用一元函数的极值条件0)('*=α?求*α。需要指出的是，在用函数)(α?的导数求*α，所用的函数是步长因子α为变量的一元函数，而不是以设计点x 为变量的多元函数)(x f 。如果直接利用)(x f ，此时需要计算k x x =点处的梯度)(x f ?和海赛矩阵G 。 数值解法：利用计算机通过反复迭代计算，求得最佳步长因子的近似值。先确定*α所在的搜索区间，然后根据区间消去法原理不断缩小此区间，从而获得*α的数值近似解。 一维搜索方法的分类：可以分为两大类。一类是试探法，这类方法是按照某种给定规律来确定区间内插入点的位置的，此点位置的确定仅考虑加快区间缩短速度，而不顾及函数值的分布关系。代表算法：黄金分割法，Fibonacci 法等。另一类是插值法，用插值函数的极小点作为区间的插入点。代表算法：牛顿法（切线法），二次插值法，三次插值法等。 2.无约束优化方法 关键：确定搜索方向k d 的构成问题是无约束优化方法的关键。 分类：（1）利用目标函数的一阶或二阶导数的无约束优化方法【解析解法】，如最速下降法、共轭梯度法、牛顿法以及变尺度法等。（2）只利用目标函数值的无约束优化方法【数值解法】，如坐标轮换法、单形替代法及Powell 法等。 3.约束优化方法（修正） 【分类】 直接法：随机方向法、复合形法； 间接法：惩罚函数法，增广乘子法,可行方向法、广义简约梯度法。 二、智能算法部分： 1.遗传算法 ①提出：遗传算法起源于对生物系统所进行的计算机模拟研究。本世纪40年代，就有学者开始研究如何利用计算机进行生物模拟的技术，他们从生物学的角度进行了生物的进化过程模拟、遗传过程模拟等研究工作。进入60年代后，美国密执安大学的Holland 教授及其学生们受到这种生物模拟技术的启发，创造出了一种基于生物遗传和进化机制的适合于复杂系统优化计算的自适应概率优化技术——遗传算法。 ②基本遗传算法： 步骤一：初始化。设置进化代数计算器t ←0；设置最大进化代数T ，随机生成M 个个体作为初始群体P （0）。 步骤二：个体评价。计算群体P （t ）中各个个体的适应度。

《现代设计方法》模拟试题三

现代设计方法模拟试题（一） 一、单项选择题(在每小题的四个备选答案中，选出一个正确的答案，并将其号码写在题干的括号 内。每小题1分，共30分) 1.试判别矩阵1111???? ? ?，它是（ ） ①单位矩阵 ②正定矩阵 ③负定矩阵 ④不定矩阵 2.对于平面桁架中的杆单元，每个节点在整体坐标系中的位移分量个数为（ ） ①1 ②2 ③3 ④4 3.约束极值点的库恩——塔克条件为： -?=?=∑F X g X i i q i ()()* * λ1 ，当约束函数是g i (X)≤0和λi >0 时，则q 应为（ ） ①等式约束数目 ②不等式约束数目 ③起作用的等式约束数目 ④起作用的不等式约束数目 4.机电产品的平均失效率λ(t)，它表征了该产品工作到t 时刻后（ ） ①单位时刻内发生失效的概率 ②单位时刻内发生失效的产品数 ③的累积失效数与受试产品总数之比 ④的累积失效数与仍正常工作的产品数之比 5.应用四节点等参数单元时，由整体坐标系到自然坐标系单元的映射关系是（ ） ①任意四边形→任意四边形 ②正方形→任意四边形 ③任意四边形→正方形 ④正方形→正方形 6.在图示极小化的约束优化问题中，最优点为（ ） ①A ②B ③C ④D 7.图示弹簧系统的总体刚度矩阵为（ ） 8. 则该产品的存活频率R (200)为（ ） ①0.00125 ②0.8 ③0.001 ④0.2 9.下列优化方法中，不需计算迭代点一阶导数和二阶导数的是（ ） ①可行方向法 ②复合形法 ③DFP 法 ④BFGS 法 10.轴对称问题中，值等于零的应变是（ ） ①γr θ ②γrz ③εθ ④εr 11.在任何一个单元内（ ） ①只有节点符合位移模式 ②只有边界点符合位移模式 ③只有边界点和节点符合位移模式 ④单元内任意点均符合位移模式 12.表示机电设备的一般失效曲线(浴盆曲线)中，偶然失效期的失效密度f(t)服从（ ） ①威布尔分布 ②指数分布 ③正态分布 ④泊松分布

时域有限差分法对平面TE波的MATLAB仿真

时域有限差分法对平面TE波的 MATLAB仿真 摘要 时域有限差分法是由有限差分法发展出来的数值计算方法。自1966年Yee 在其论文中首次提出时域有限差分以来，时域有限差分法在电磁研究领域得到了广泛的应用。主要有分析辐射条线、微波器件和导行波结构的研究、散射和雷达截面计算、分析周期结构、电子封装和电磁兼容的分析、核电磁脉冲的传播和散射以及在地面的反射及对电缆传输线的干扰、微光学元器件中光的传播和衍射特性等等。 由于电磁场是以场的形态存在的物质，具有独特的研究方法，采取重叠的研究方法是其重要的特点，即只有理论分析、测量、计算机模拟的结果相互佐证，才可以认为是获得了正确可信的结论。时域有限差分法就是实现直接对电磁工程问题进行计算机模拟的基本方法。在近年的研究电磁问题中，许多学者对时域脉冲源的传播和响应进行了大量的研究，主要是描述物体在瞬态电磁源作用下的理论。另外，对于物体的电特性，理论上具有几乎所有的频率成分，但实际上，只有有限的频带内的频率成分在区主要作用。 文中主要谈到了关于高斯制下完全匹配层的差分公式的问题，通过MATLAB 程序对TE波进行了仿真，模拟了高斯制下完全匹配层中磁场分量瞬态分布。得到了相应的磁场幅值效果图。 关键词：时域有限差分完全匹配层MATLAB 磁场幅值效果图

目录 摘要 (1) 目录 (3) 第一章绪论 (4) 1.1 课题背景与意义 (4) 1.2 时域有限差分法的发展与应用 (4) 2.1 Maxwell方程和Yee氏算法 (7) 2.2 FDTD的基本差分方程 (9) 2.3 时域有限差分法相关技术 (11) 2.3.1 数值稳定性问题 (11) 2.3.2 数值色散 (12) 2.3.3 离散网格的确定 (13) 2.4 吸收边界条件 (13) 2.4.1 一阶和二阶近似吸收边界条件 (14) 2.4.2 二维棱边及角顶点的处理 (17) 2.4.3 完全匹配层 (19) 2.5 FDTD计算所需时间步的估计 (23) 第三章MATLAB的仿真的程序及模拟 (25) 3.1 MATLAB程序及相应说明 (25) 3.2 出图及结果 (28) 3.2.1程序部分 (28) 3.2.2 所出的效果图 (29) 第四章结论 (31) 参考文献 (32)

matlab 牛顿插值法 三次样条插值法

(){} 21 ()(11),5,10,20: 1252 1()1,(0,1,2,,)()2,(0,1,2,,)() ()2 35,20:1100 (i i i i n n k k k Newton f x x n x f x x i i n f x n x y i n Newton N x S x n x k y f x =-≤≤=+=-+====-+ = 题目：插值多项式和三次样条插值多项式。 已知对作、计算函数在点处的值；、求插值数据点 的插值多项式和三次样条插值多项式；、对计算和相应的函数值),()() (1,2,,99)4:()max ()()max ()n k n k n k n k n k n k k k N x S x k E N y N x E S y S x ==-=- 和； 、计算，； 解释你所得到的结果。 算法组织: 本题在算法上需要解决的问题主要是：求出第二问中的Newton 插值多项式 )(x N n 和三次样条插值多项式()n S x 。如此，则第三、四问则迎刃而解。计算两种插值多项式的算法如下： 一、求Newton 插值多项式)(x N n ，算法组织如下： Newton 插值多项式的表达式如下： )())(()()(110010--???--+???+-+=n n n x x x x x x c x x c c x N 其中每一项的系数c i 的表达式如下： 1102110) ,,,(),,,(),,,(x x x x x f x x x f x x x f c i i i i i -???-???= ???=- 根据i c 以上公式，计算的步骤如下： ?? ??? ?? ?????+??????? ???????????----) ,,,,(1) ,,,(),,,,(),(,),,(2)(,),(),(11101111011010n n n n n n n n x x x x f n x x x f x x x f n x x f x x f x f x f x f 、计算、计算、计算、计算 二、求三次样条插值多项式)(x S n ，算法组织如下：

现代设计方法2018年11月复习地的题目

1. 进行旋转坐标变换的变换矩阵???? ? ??? ??-10 0cos sin 0sin cos θθ θθ中，θ角正负号选择原则是什么？ 2. 三维图形变换矩阵T=??? ?? ???? ???s n m l r j i h q f e d p c b a 中，l 、n 表示的变换为什么意思？ 3. 将点P （0, 10）进行?? ?? ? ?? ?? ?-10 0090cos 90 sin 090sin 90cos o o o o 坐标变换，则变换后P 点坐标为多少？ 4. 已知一组数据 用二次多项式拟合。

5. 对第Ⅰ象限中的一个点P 实施-???? ??? ? ??100010001坐标变换，则变换后P 点位于第几象限？ 6. 在CAD 作业中，为了定义图形方便通常采用不同的坐标系，定义域为连续且无界的坐标系是什么坐标系？ 7. CAD 系统的软件根据其用途可分为哪三类？ 8. 与文件系统相比，数据库系统的主要特征有哪些？

9.图形变换矩阵T= 200 020 001 ? ? ? ? ? ? ? ? ? ?，则变换后的图形是原来的几倍？ 10.CAD系统输入设备、输出设备有哪些？ 11.试写出从视区中一点V(x v,y v)到窗口中一点W(x w,y w)的 变换公式及相应参数的意义。 12.请列举出常用的终止准则（三种）。 13.简述梯度法的基本原理和特点。

14. 简述罚函数法的基本原理。罚函数法分为哪几种？ 15. 判断矩阵2014-???? ??，是否为正定矩阵？ 16. 用梯度法求下列无约束优化问题：MinF(X)=x 12+4x 22， 设初始点取为X (0)={2,2}T ，以梯度模为终止迭代准则，其收敛精度为5。 17. 判断函数212 221260)(x x x x x x F +++-=是否为凸函数。

PID控制算法的matlab仿真.doc

PID 控制算法的 matlab 仿真 PID 控制算法是实际工业控制中应用最为广泛的控制算法，它具有控制器设 计简单，控制效果好等优点。 PID 控制器参数的设置是否合适对其控制效果具有很 大的影响，在本课程设计中一具有较大惯性时间常数和纯滞后的一阶惯性环节作为 被控对象的模型对 PID 控制算法进行研究。被控对象的传递函数如下： G(s) Ke d s 1 T f s 其中各参数分别为 K 30, T f 630, d 60 。MATLAB仿真框图如图1所示。 -K- Kp z 30 -K- (z-1) 630s+1 Step Ki Zero-Order Transport Transfer Fcn Hold Delay Add (z-1) -K- z Kd 图 1 2具体内容及实现功能 2.1 PID 参数整定 PID 控制器的控制参数对其控制效果起着决定性的作用，合理设置控制参 数是取得较好的控制效果的先决条件。常用的 PID 参数整定方法有理论整定法和实 验整定法两类，其中常用的实验整定法由扩充临界比例度法、试凑法等。在此处选用扩充临界比例度法对 PID 进行整定，其过程如下： 1)选择采样周期由于被控对象中含有纯滞后，且其滞后时间常数为 d60 ，故可选择采样周期 T s 1。 2)令积分时间常数T i，微分时间常数T d0 ，从小到大调节比例系数K ， 使得系统发生等幅震荡，记下此时的比例系数K k和振荡周期 T k。 1 Out1 3)选择控制度为Q 1.05 ，按下面公式计算各参数：

K p 0.63K k T i 0.49T k T d 0.14T k T s 0.014T k 通过仿真可得在 T s 1时， K k 0.567, T k 233 ，故可得： K p 0.357, T i 114.17, T d 32.62,T s 3.262 K p T s 0.005 K i T i K p T d 3.57 K d T s 按此组控制参数得到的系统阶跃响应曲线如图 2 所示。 1.8 1.6 1.4 1.2 1 0.8 0.6 0.4 0.2 01002003004005006007008009001000 图 2 由响应曲线可知，此时系统虽然稳定，但是暂态性能较差，超调量过大，且响应曲线不平滑。根据以下原则对控制器参数进行调整以改善系统的暂态过程： 1)通过减小采样周期，使响应曲线平滑。 2)减小采样周期后，通过增大积分时间常数来保证系统稳定。 3)减小比例系数和微分时间常数，以减小系统的超调。 改变控制器参数后得到系统的阶跃响应曲线如图3 所示，系统的暂态性能得到明显改善 .