Смекни!
smekni.com

Численные методы интегрирования и оптимизации сложных систем (стр. 5 из 5)

L1(2)=F2;

F3=4*T.^2-1;

L1(3)=F3;

F4=8*T.^3-4*T;

L1(4)=F4;

F5=16*T.^4-12*T.^2+1;

L1(5)=F5;

F6=32*T.^5-32*T.^3+6*T;

L1(6)=F6;

F7=64*T.^6-80*T.^4+24*T.^2-1;

L1(7)=F7;

F8=128*T.^7-192*T.^5+80*T.^3-8*T;

F9=256*T.^8-448*T.^6+240*T.^4-40*T.^2+1;

L1(9)=F9;

F10=512*T.^9-1024*T.^7+672*T.^5-160*T.^3+10*T;

L1(10)=F10;

F1=1;

L1(1)=F1;

G=L'*L1;

In=Kx*G;

r=int(In,T,0,t);

Cx=int(r,t,0,1.5);

In=Ky.*G;

r=int(In,T,0,t);

Cy=int(r,t,0,1.5);

A=((Cx+eye(10)).^-1)*Cy;

Cy=int(L,t,0,1.5);

Cx=A*Cy'

Postr.m

function H=fun(t)

Cx=[3.7672; 1.3134; 0.5181; 0.2065; 0.0819; 0.0323; 0.0127; 0.0491; 0.0189; 0.0723];

F2=2*t;

L(2)=F2;

F3=4*t.^2-1;

L(3)=F3;

F4=8*t.^3-4*t;

L(4)=F4;

F5=16*t.^4-12*t.^2+1;

L(5)=F5;

F6=32*t.^5-32*t.^3+6*t;

L(6)=F6;

F7=64*t.^6-80*t.^4+24*t.^2-1;

L(7)=F7;

F8=128*t.^7-192*t.^5+80*t.^3-8*t;

L(8)=F8;

F9=256*t.^8-448*t.^6+240*t.^4-40*t.^2+1;

L(9)=F9;

F10=512*t.^9-1024*t.^7+672*t.^5-160*t.^3+10*t;

L(10)=F10;

F1=1;

L(1)=F1;

H=(Cx'*L');

t=[0:0.01:5];

plot(t,H)


10. Приложение 5 (Листинг скриптов для оптимизации)

jivs.m

clear

clc

a=0;

b=5;

h=0.1;

Kp1=2; Kd1=1; Ki1=0;

J1=int2(Kp1, Kd1, Ki1);

while h>0.0000001

Kp2=Kp1+h;

J2=int2(Kp2, Kd1, Ki1);

if J2>J1

Kp2=Kp1-h;

J2=int2(Kp2,Kd1,Ki1);

if J2>J1

Kp2=Kp1;

end

end

Kd2=Kd1+h;

J2=int2(Kp2, Kd2, Ki1);

if J2>J1

Kd2=Kd1-h;

J2=int2(Kp2,Kd2,Ki1);

if J2>J1

Kd2=Kd1;

end

end

Ki2=Ki1+h;

J2=int2(Kp2, Kd2, Ki2,h);

if J2>J1

Ki2=Ki1-h;

J2=int2(Kp2,Kd2,Ki2,h);

if J2>J1

Ki2=Ki1;

end

end

h=fibon(a,b,h);

while J2<J1

Kp=Kp1+2*(Kp2-Kp1); Kd=Kd1+2*(Kd2-Kd1); Ki=Ki1+2*(Ki2-Ki1);

J1=J2;

J2=int2(Kp,Kd,Ki,h);

Kp1=Kp2;Kp2=Kp; Kd1=Kd2;Kd2=Kd; Ki1=Ki2;Ki2=Ki;

end

end

disp(Kp)

disp(Kd)

disp(Ki)

int2.m

function J=int2(Kp,Kd,Ki,h)

clc

T=4;

A=[0 1 0 0; 0 0 1 0; 0 0 0 1; -595.23809523*Ki -43.6507936507-595.23809523*Kp -56.34920634920635-595.23809523*Kd -5.59523809523809];

B=[0; 595.23809*Kd; 595.23809*Kp-3330.498866*Kd; 595.23809*Ki-33540.615-354308.277*(Kd)^2-3330.498*Kp-18634.934*Kd];

k=0;

t=0;

while (t < (T-h))

if (t <= 3*h)

K1=A*(A_X(k+1,:))';

K2=A*(A_X(k+1,:))'+1/3*K1;

K3=A*(A_X(k+1,:))'+1/6*K1+1/6*K2;

K4=A*(A_X(k+1,:))'+1/8*K1+3/8*K2;

K5=A*(A_X(k+1,:))'+1/2*K1-3/2*K3+2*K4;

A_X(k+2,:)=(A_X(k+1,:))+h/6*(K1'+4*K4'+K5');

else

h1=h;

t=t+h1;

H=(eye(length(A_X(1,:)))-(9*h1/24)*A);

G=(eye(length(A_X(1,:)))+19*h1/24*A)*(A_X(k+1,:))'+h1/24*A*(-5*(A_X(k,:))'+(A_X(k-1,:))')

+h1/24*B*(9*1+19*1-5*1);

A_X(k+2,:)=(inv(H)*G)';

end

Otr(k+1)=t;

k=k+1;

end

grid on

fibon.m

function h=fibon(a,b,h)

F(1)=1; F(2)=1;n=100;

for i=[1:0.1:n-2]

F(i+2)=F(i+1)+F(i);

end

j=0;

x1=a; x3=b;

L1=x3-x1;

L2=(F(n-1)/F(n))*L1+((-1)^n)/F(n)*eps;

x2=x3-L2;

x4=x1+x3-x2;

while (abs(x3-x1) > eps)

F2=x2;

F4=x4;

if ((x2 < x4)&&(norm(F2) < norm(F4)))

x1=x1; x3=x4;

x4=x1+x3-x2;

elseif ((x2 > x4)&&(norm(F2) < norm(F4)))

x1=x4; x3=x3;

x4=x1+x3-x2;

elseif ((x2 < x4)&&(norm(F2) > norm(F4)))

x1=x2; x3=x3;

x2=x1+x3-x4;

elseif ((x2 > x4)&&(norm(F2) > norm(F4)))

x1=x1; x3=x2;

x2=x1+x3-x4;

end

j=j+1;

la=x1+(x3-x1)/2;

end

l=la;