Смекни!
smekni.com

Решение обратных задач динамики (стр. 5 из 5)

Рис. 8. Графики выходных сигналов системы

Таким образом, можно построить следующий алгоритм решения задачи синтеза входного сигнала нелинейной системы:

1) задается эталонный выходной сигнал;

2) из (1) находится сигнал на выходе нелинейного элемента, который на выходе системы обеспечивает требуемый эталонный процесс;

3) найденный в предыдущем пункте сигнал представляется как сигнал на входе нелинейного элемента и находится реальный сигнал на выходе нелинейного элемента и уточняется эталонный сигнал на выходе системы;

4) поскольку известны сигналы на входе нелинейного элемента и на выходе системы, то, представив искомый входной сигнал в виде (4), строится невязка (6) и функционал (7);

5) минимизируя полученный функционал, находятся числовые значения искомых коэффициентов

,
;

6) проводится анализ полученных результатов.

5. Результаты расчёта

1. Эталонный закон изменения угла teta(t)

Число точек квантования по времени: Nt = 499;

Шаг квантования: h_t = 0.020000 c;

Время поражения цели: T = 9.960000 c;

2. Числовые значения параметров системы самонаведения

Krp = 1.000000;

Trp = 0.330000, с;

Xmax = 0.418879, рад;

Ksn = 0.283000, рад/с;

Tsn = 0.155000, с;

DZsn = 0.052000;

V = 686.700000, м/с;

G = 9.810000, м/с^2;

Kdy = 0.140000;

Kv = 1.200000, c;

mu = 0.115000, с;

Tc = 3.050000, с;

3. Базис - функции Уолша

Число элементов: Nl = 64;

Оператор интегрирования Ai размерностью 5x5

+4.980000e+000 +2.490000e+000 +0.000000e+000 +1.245000e+000 +0.000000e+000

-2.490000e+000 +0.000000e+000 +1.245000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 -1.245000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000

-1.245000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000 +6.225000e-001

+0.000000e+000 +0.000000e+000 +0.000000e+000 -6.225000e-001 +0.000000e+000

Оператор дифференцирования Ad размерностью 5x5

+0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000

4. Матричные операторы системы

Arp размерностью 5x5

+9.668675e-001 +3.313252e-002 -3.310223e-002 +3.310224e-002 -3.171612e-002

-3.313252e-002 +9.006024e-001 +9.930669e-002 +3.310223e-002 -3.171610e-002

-3.310223e-002 -9.930669e-002 +8.343980e-001 +3.307196e-002 -3.168711e-002

-3.310224e-002 +3.310223e-002 -3.307196e-002 +7.682541e-001 +2.220418e-001

-3.171612e-002 +3.171610e-002 -3.168711e-002 -2.220418e-001 +7.048218e-001

Asn размерностью 5x5

+2.824568e-001 -1.904545e-003 -3.561384e-003 +6.907620e-003 -4.945520e-003

+1.904545e-003 +2.862659e-001 +1.403039e-002 +3.561384e-003 -6.087807e-003

-3.561384e-003 -1.403039e-002 +2.791431e-001 +1.571978e-002 -7.488732e-003

-6.907620e-003 +3.561384e-003 -1.571978e-002 +2.477036e-001 +3.501836e-002

-4.945520e-003 +6.087807e-003 -7.488732e-003 -3.501836e-002 +2.378125e-001

Aos1 размерностью 5x5

-6.831527e+001 +0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 -6.831527e+001 +0.000000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 -6.831527e+001 +0.000000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 +0.000000e+000 -6.831527e+001 +0.000000e+000

+0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000 -6.831527e+001

Aos2 размерностью 5x5

+9.800000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 +9.800000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 +9.800000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 +0.000000e+000 +9.800000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000 +9.800000e+000

Apr размерностью 5x5

+2.730338e-001 +9.500096e-003 -1.194782e-002 +1.128111e-002 -7.250387e-003

-9.500096e-003 +2.540337e-001 +3.517674e-002 +1.194782e-002 -8.567413e-003

-1.194782e-002 -3.517674e-002 +2.301380e-001 +1.306212e-002 -5.530241e-003

-1.128111e-002 +1.194782e-002 -1.306212e-002 +2.040138e-001 +5.461586e-002

-7.250387e-003 +8.567413e-003 -5.530241e-003 -5.461586e-002 +1.895130e-001

Aos размерностью 5x5

-5.851527e+001 +0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 -5.851527e+001 +0.000000e+000 +0.000000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 -5.851527e+001 +0.000000e+000 +0.000000e+000

+0.000000e+000 +0.000000e+000 +0.000000e+000 -5.851527e+001 +0.000000e+000

+0.000000e+000 +0.000000e+000 +0.000000e+000 +0.000000e+000 -5.851527e+001

As размерностью 5x5

+3.194591e-001 +1.707523e-001 +5.751752e-004 +8.508857e-002 +5.722004e-004

-1.707523e-001 -2.204553e-002 +8.393822e-002 -5.751752e-004 +5.722004e-004

+5.751752e-004 -8.393822e-002 -2.089518e-002 -5.751662e-004 +5.721915e-004

-8.508857e-002 -5.751752e-004 +5.751662e-004 -1.974485e-002 +3.882646e-002

+5.722004e-004 -5.722004e-004 +5.721915e-004 -3.882646e-002 -1.860045e-002

5. СХ эталонного выхода

Ctheta размерностью 5x1

+2.948462e-001

-7.002572e-002

-4.945100e-002

-5.104576e-002

-1.450117e-002

Начальные значения искомых коэффициентов

Cu_0 размерностью 5x1

+0.000000e+000

+0.000000e+000

+0.000000e+000

+0.000000e+000

+0.000000e+000

Oshibka_0 = 3.145671e-001

Conditioning of Gradient Poor - Switching To LM method

Optimization terminated: directional derivative along

search direction less than TolFun and infinity-norm of

gradient less than 10*(TolFun+TolX).

Оптимальные значения искомых коэффициентов

Cu_opt размерностью 5x1

+5.349004e-001

+5.156158e-001

+3.167675e-001

+3.345843e-001

+3.459092e-002

Oshibka_0 = 1.444098e-004

-------------------------------------------------------------

Время расчета:

0 часов, 0 минут, 34.703 секунд.

Приложение

1) % Программа синтеза управления системы самонаведения (рассматривается часть % системы) методами обратных задач динамики с использованием метода % матричных операторов (линейная модель)

close all;

clear all;

clc;

my_tic;

global Nl;

global U tgl;

global Krp Trp Ksn Tsn DZsn V G Kdy Kv mu Tc Xmax;

%% 1. Эталонный закон изменения угла teta(t)

% Время наведения

fId = fopen('t_navedenija.dat','r');

t_f = fread(fId,inf,'real*8')';

fclose(fId);

Nt_f = length(t_f);

h_t_f = t_f(2)-t_f(1);

T = t_f(Nt_f);

% угол theta(t)

fId = fopen('theta_navedenija.dat','r');

theta_f = fread(fId,[1 Nt_f],'real*8');

fclose(fId);

% расстояние до цели

fId = fopen('r_navedenija.dat','r');

r_f = fread(fId,[1 Nt_f],'real*8');

fclose(fId);

fprintf('1. Эталонный закон изменения угла teta(t)\n');

fprintf('Число точек квантования по времени: Nt = %i;\n',Nt_f);

fprintf('Шаг квантования: h_t = %f c;\n',h_t_f);

fprintf('Время поражения цели: T = %f c;\n',T);

fprintf('\n');

my_plot2(t_f,theta_f,'t, c','theta(t), рад');

my_plot2(t_f,r_f,'t, c','r(t), м');

% пересчет на больший шаг квантования

Nt = 64;

h_t = T/(Nt-1);

t = 0: h_t: T;

theta = spline(t_f,theta_f,t);

r = spline(t_f,r_f,t);

my_plot2(t,theta,'t, c','theta(t), рад');

my_plot2(t,r,'t, c','r(t), м');

%% 2. Параметры системы

% Числовые значения параметров системы самонаведения

Krp = 1; %

Trp = 0.33; % с

Xmax = 24*pi/180; % рад

Ksn = 0.283; % рад/с

Tsn = 0.155; % с

DZsn = 0.052; %

V = 70*9.81; % м/с

G = 9.81; % м/с^2

Kdy = 0.14; %

Kv = 1.2; % c

mu = 0.115; % с

Tc = 3.05; % с

fprintf('2. Числовые значения параметров системы самонаведения\n');

fprintf('Krp = %f;\n',Krp);

fprintf('Trp = %f, с;\n',Trp);

fprintf('Xmax = %f, рад;\n',Xmax);

fprintf('Ksn = %f, рад/с;\n',Ksn);

fprintf('Tsn = %f, с;\n',Tsn);

fprintf('DZsn = %f;\n',DZsn);

fprintf('V = %f, м/с;\n',V);

fprintf('G = %f, м/с^2;\n',G);

fprintf('Kdy = %f;\n',Kdy);

fprintf('Kv = %f, c;\n',Kv);

fprintf('mu = %f, с;\n',mu);

fprintf('Tc = %f, с;\n',Tc);

fprintf('\n');

%% 3. Формирование ортонормированного базиса

Nl = Nt;

setsize(Nl);

settime(T);

Ai = mkint; % оператор интегрирования

Ad = inv(Ai); % оператор дифференцирования

Ae = eye(Nl); % единичная матрица

fprintf('3. Базис - функции Уолша\n');

fprintf('Число элементов: Nl = %i;\n',Nl);

pr_matrix(Ai,'Оператор интегрирования Ai')

pr_matrix(Ad,'Оператор дифференцирования Ad')

%% 4. Расчет операторов системы

Arp = inv(Trp*Ae+Ai)*(Krp*Ai);

Asn = inv(Tsn^2*Ae+2*DZsn*Tsn*Ai+Ai*Ai)*(Ksn*Ai*Ai);

Aos1 = Kv*mu*Tc*Ad*Ad+Kv*(mu+Tc)*Ad+Kv*Ae;

Aos2 = (Kdy*V/G)*Ae;

Apr = Asn*Arp;

Aos = Aos1+Aos2;

As = inv(Ae+Aos*Apr)*Apr;

As = Ai*As;

fprintf('4. Матричные операторы системы\n');

pr_matrix(Arp,'Arp');

pr_matrix(Asn,'Asn');

pr_matrix(Aos1,'Aos1');

pr_matrix(Aos2,'Aos2');

pr_matrix(Apr,'Apr');

pr_matrix(Aos,'Aos');

pr_matrix(As,'As');

%% 5. Расчет спектральной характеристики эталонного выхода

Ctheta = fwht(theta');

fprintf('5. СХ эталонного выхода\n');

pr_matrix(Ctheta,'Ctheta');

%% 6. Синтез входного сигнала

Cu_0 = zeros(Nl,1);

fprintf('Начальные значения искомых коэффициентов\n');

pr_matrix(Cu_0,'Cu_0');

oshibka = sqrt((As*Cu_0-Ctheta)'*(As*Cu_0-Ctheta));

fprintf('Oshibka_0 = %e\n',oshibka);

my_function = @(Cu)sqrt((As*Cu-Ctheta)'*(As*Cu-Ctheta));

% optimset('Display','iter','NonlEqnAlgorithm','gn','TolFun',1e-8,...

Cu = fsolve(my_function,Cu_0,...

optimset('NonlEqnAlgorithm','gn','TolFun',1e-8,...

'TolX',1e-8,'MaxFunEvals',50000,'MaxIter',50000));

% Cu = inv(As)*Ctheta;

fprintf('Оптимальные значения искомых коэффициентов\n');

pr_matrix(Cu,'Cu_opt');

oshibka = sqrt((As*Cu-Ctheta)'*(As*Cu-Ctheta));

fprintf('Oshibka_0 = %e\n',oshibka);

U = iwht(Cu)';

tgl = t;

my_plot2(t,U,'t, c','U(t)');

%% 7. Анализ полученных результатов (метод Рунге-Кутта (ode45))

[tt,yy] = ode45(@ode_navedenija1,t,[0 0 0 0]);

theta_rr = yy(:,1)';

my_plot2(t,[theta;theta_rr],'t, c','theta(t), рад','',['эталонный ';'реальный ']);

my_toc;

2)второстепенные программы:

function dy = ode_navedenija1(t,y);

global U tgl;

global Krp Trp Ksn Tsn DZsn V G Kdy Kv mu Tc Xmax;

a32 = -1/(Tsn^2);

a33 = -2*DZsn/Tsn;

a3f = Ksn/(Tsn^2);

a42 = -(Krp/Trp)*(Kv-Kv*mu*Tc/(Tsn^2)+Kdy*V/G);

a43 = -(Krp/Trp)*(Kv*(mu+Tc)-2*Kv*mu*Tc*DZsn/Tsn);

a44 = -1/Trp;

a4f = -(Krp/Trp)*Kv*Ksn*mu*Tc/(Tsn^2);

b4 = Krp/Trp;

u = spline(tgl,U,t);

dy = zeros(4,1);

dy(1) = y(2);

dy(2) = y(3);

y4 = y(4);

dy(3) = a32*y(2)+a33*y(3)+a3f*y4;

dy(4) = b4*u+a42*y(2)+a43*y(3)+a44*y(4)+a4f*y4;