Моделирование системы динамики аппарата
Этот пример показывает нелинейное серо-коробочное моделирование динамики аппарата. Многие новые функции транспортного средства (такие как электронные программы устойчивости (ESP), косвенные системы мониторинга давления в шинах (TPMS), системы мониторинга трения в шинах и т.д.) полагаются на модели базовой динамики аппарата. Так называемая модель велосипедного транспортного средства является довольно простой структурой модели, которая часто используется в динамику аппарата литературе. В этом примере мы начнем с этой структуры модели и попытаемся оценить продольную и боковую жесткость шины. Фактическая работа по моделированию была первоначально выполнена Эриком Нарби в его MSc работе в NIRA Dynamics AB, Швеция.
Моделирование динамики аппарата
Следующий рисунок иллюстрирует ситуацию моделирования транспортного средства, которая должна быть рассмотрена.
Фигура 1: Схематический вид системы динамики аппарата.
При помощи закона движения Ньютона и некоторых основных геометрических отношений продольная скорость v_x (t), поперечная скорость v_y (t) и скорость рыскания r (t), измеренная вокруг центра тяжести ( COG) транспортного средства, могут быть описаны следующими тремя дифференциальными уравнениями:
d
-- v_x(t) = v_y(t)*r(t) + 1/m*( (F_x,FL(t)+F_x,FR(t))*cos(delta(t))
dt - (F_y,FL(t)+F_y,FR(t))*sin(delta(t))
+ F_x,RL(t)+F_x,RR(t)
- C_A*v_x(t)^2)
d
-- v_y(t) = -v_x(t)*r(t) + 1/m*( (F_x,FL(t)+F_x,FR(t))*sin(delta(t))
dt + (F_y,FL(t)+F_y,FR(t))*cos(delta(t))
+ F_y,RL(t)+F_y,RR(t))
d
-- r(t) = 1/J*( a*( (F_x,FL(t)+F_x,FR(t))*sin(delta(t))
dt + (F_y,FL(t)+F_y,FR(t))*cos(delta(t)))
- b*(F_y,RL(t)+F_y,RR(t)))где нижний индекс x используется для обозначения того, что сила F действует в продольном направлении, и y, что она действует в боковом направлении. Сокращения FL, FR, RL и RR маркируют шины: Front Left, Front Right, Rear Left и Rear Right соответственно. Первое уравнение, описывающее продольное ускорение, также содержит член сопротивления воздуха, который принимается квадратичной функцией продольной скорости транспортного средства v_x (t). В сложение дельта (t) (вход) является углом поворота, J - моментом инерции и a и b - расстояниями от центра тяжести до передней и задней осей, соответственно.
Предположим, что силы шин могут быть смоделированы посредством следующих линейных приближений:
F_x,i(t) = C_x*s_i(t)
F_y,i(t) = C_y*alpha_i(t) for i = {FL, FR, RL, RR}где C_x и C_y являются жесткостью продольной и боковой шин, соответственно. Здесь мы приняли, что эти параметры жесткости одинаковы для всех 4 шин. s_i (t) является так называемым (продольным) проскальзыванием шины i и alpha_i (t) углом проскальзывания шины. Для транспортного средства, управляемого передним колесом (как рассматривается здесь), проскальзывания s_FL (t) и s_FR (t) получают из отдельных скоростей колеса (измеряются), принимая, что задние колеса не показывают скольжения (т.е. s_RL (t) = s_RR (t) = 0). Следовательно, скольжения являются входами в нашу структуру модели. Для передних колес углы скольжения шин alpha_Fj (t) могут быть аппроксимированы (когда v_x (t) > 0)
alpha_Fj(t) = delta(t) - arctan((v_y(t) + a*r(t))/v_x(t))
~ delta(t) - (v_y(t) + a*r(t))/v_x(t) for j = {L, R}
Для задних колес углы скольжения alpha_Rj (t) шины аналогично получают и вычисляют как
alpha_Rj(t) = - arctan((v_y(t) - b*r(t))/v_x(t))
~ - (v_y(t) - b*r(t))/v_x(t) for j = {L, R}
С J = 1/( 0,5 * (a + b)) ^ 2 * m) мы можем затем настроить структуру пространства состояний, описывающую динамику аппарата. Вводите состояния:
x1(t) = v_x(t) Longitudinal velocity [m/s]. x2(t) = v_y(t) Lateral velocity [m/s]. x3(t) = r(t) Yaw rate [rad/s].
пять измеренных или выведенных входных сигналов
u1(t) = s_FL(t) Slip of Front Left tire [ratio]. u2(t) = s_FR(t) Slip of Front Right tire [ratio]. u3(t) = s_RL(t) Slip of Rear Left tire [ratio]. u4(t) = s_RR(t) Slip of Rear Right tire [ratio]. u5(t) = delta(t) Steering angle [rad].
и параметры модели:
m Mass of the vehicle [kg]. a Distance from front axle to COG [m]. b Distance from rear axle to COG [m]. Cx Longitudinal tire stiffness [N]. Cy Lateral tire stiffness [N/rad]. CA Air resistance coefficient [1/m].
Выходами системы являются продольная скорость транспортного средства y1 (t) = x1 (t), боковое ускорение транспортного средства (измеренное акселерометром):
y2(t) = a_y(t) = 1/m*( (F_x,FL(t) + F_x,FR(t))*sin(delta(t))
+ (F_y,FL(t) + F_y,FR(t))*cos(delta(t))
+ F_y,RL(t) + F_y,RR(t))
и скорость рыскания y3 (t) = r (t) (измеренная гироскопом).
В совокупности мы приходим к следующей структуре модели пространства состояний:
d
-- x1(t) = x2(t)*x3(t) + 1/m*( Cx*(u1(t)+u2(t))*cos(u5(t))
dt - 2*Cy*(u5(t)-(x2(t)+a*x3(t))/x1(t))*sin(u5(t))
+ Cx*(u3(t)+u4(t))
- CA*x1(t)^2) d
-- x2(t) = -x1(t)*x3(t) + 1/m*( Cx*(u1(t)+u2(t))*sin(u5(t))
dt + 2*Cy*(u5(t)-(x2(t)+a*x3(t))/x1(t))*cos(u5(t))
+ 2*Cy*(b*x3(t)-x2(t))/x1(t)) d
-- x3(t) = 1/((0.5*(a+b))^2)*m)*( a*( Cx*(u1(t)+u2(t)*sin(u5(t))
dt + 2*Cy*(u5(t) - (x2(t)+a*x3(t))/x1(t))*cos(u5(t)))
- 2*b*Cy*(b*x3(t)-x2(t))/x1(t)) y1(t) = x1(t)
y2(t) = 1/m*( Cx*(u1(t)+u2(t))*sin(u5(t))
+ 2*Cy*(u5(t)-(x2(t)+a*x3(t))/x1(t))*cos(u5(t))
+ 2*Cy*(b*x3(t)-x2(t))/x1(t))
y3(t) = x3(t)Транспортное средство IDNLGREY
В качестве базиса для наших экспериментов по идентификации транспортных средств нам сначала нужно создать файл модели IDNLGREY, описывающий эти уравнения транспортных средств. Здесь мы полагаемся на моделирование C-MEX и создаем файл vehicle_c.c модели, в котором NY установлен на 3. Функции обновления состояния и выхода vehicle_c.c, compute_dx и compute_y, в некоторой степени вовлечены и включают несколько стандартных C-заданных математических функций, таких как cos (.) и sin (.), а также pow (.) для вычисления степени его аргумента.
Функция обновления состояния compute_dx возвращает dx (аргумент 1) и использует 3 входные параметры: вектор состояния x, входной вектор u и шесть скалярных параметров, закодированных в p (t и auxvar файла модели шаблона C-MEX были удалены здесь):
/* State equations. */
void compute_dx(double *dx, double *x, double *u, double **p)
{
/* Retrieve model parameters. */
double *m, *a, *b, *Cx, *Cy, *CA;
m = p[0]; /* Vehicle mass. */
a = p[1]; /* Distance from front axle to COG. */
b = p[2]; /* Distance from rear axle to COG. */
Cx = p[3]; /* Longitudinal tire stiffness. */
Cy = p[4]; /* Lateral tire stiffness. */
CA = p[5]; /* Air resistance coefficient. */ /* x[0]: Longitudinal vehicle velocity. */
/* x[1]: Lateral vehicle velocity. */
/* x[2]: Yaw rate. */
dx[0] = x[1]*x[2]+1/m[0]*(Cx[0]*(u[0]+u[1])*cos(u[4])
-2*Cy[0]*(u[4]-(x[1]+a[0]*x[2])/x[0])*sin(u[4])
+Cx[0]*(u[2]+u[3])-CA[0]*pow(x[0],2));
dx[1] = -x[0]*x[2]+1/m[0]*(Cx[0]*(u[0]+u[1])*sin(u[4])
+2*Cy[0]*(u[4]-(x[1]+a[0]*x[2])/x[0])*cos(u[4])
+2*Cy[0]*(b[0]*x[2]-x[1])/x[0]);
dx[2] = 1/(pow(((a[0]+b[0])/2),2)*m[0])
*(a[0]*(Cx[0]*(u[0]+u[1])*sin(u[4])
+2*Cy[0]*(u[4]-(x[1]+a[0]*x[2])/x[0])*cos(u[4]))
-2*b[0]*Cy[0]*(b[0]*x[2]-x[1])/x[0]);
}Функция compute_y обновления выхода возвращает y (аргумент 1) и использует 3 входные параметры: вектор x состояния, входной вектор u и пять из шести параметров (сопротивление воздуха CA не требуется), закодированных в p:
/* Output equations. */
void compute_y(double *y, double *x, double *u, double **p)
{
/* Retrieve model parameters. */
double *m = p[0]; /* Vehicle mass. */
double *a = p[1]; /* Distance from front axle to COG. */
double *b = p[2]; /* Distance from rear axle to COG. */
double *Cx = p[3]; /* Longitudinal tire stiffness. */
double *Cy = p[4]; /* Lateral tire stiffness. */ /* y[0]: Longitudinal vehicle velocity. */
/* y[1]: Lateral vehicle acceleration. */
/* y[2]: Yaw rate. */
y[0] = x[0];
y[1] = 1/m[0]*(Cx[0]*(u[0]+u[1])*sin(u[4])
+2*Cy[0]*(u[4]-(x[1]+a[0]*x[2])/x[0])*cos(u[4])
+2*Cy[0]*(b[0]*x[2]-x[1])/x[0]);
y[2] = x[2];
}Имея соответствующий файл структуры модели, следующим шагом является создание объекта IDNLGREY, отражающего ситуацию моделирования. Для удобства ведения бухгалтерского учета мы также указываем имена и модули входных и выходных параметров.
FileName = 'vehicle_c'; % File describing the model structure. Order = [3 5 3]; % Model orders [ny nx nu]. Parameters = [1700; 1.5; 1.5; 1.5e5; 4e4; 0.5]; % Initial parameters. InitialStates = [1; 0; 0]; % Initial value of initial states. Ts = 0; % Time-continuous system. nlgr = idnlgrey(FileName, Order, Parameters, InitialStates, Ts, ... 'Name', 'Bicycle vehicle model', 'TimeUnit', 's'); nlgr.InputName = {'Slip on front left tire'; ... % u(1). 'Slip on front right tire'; ... % u(2). 'Slip on rear left tire'; ... % u(3). 'Slip on rear right tire'; ... % u(4). 'Steering angle'}; ... % u(5). nlgr.InputUnit = {'ratio'; 'ratio'; 'ratio'; 'ratio'; 'rad'}; nlgr.OutputName = {'Long. velocity'; ... % y(1); Longitudinal vehicle velocity 'Lat. accel.'; ... % y(2); Lateral vehicle acceleration 'Yaw rate'}; ... % y(3). nlgr.OutputUnit = {'m/s'; 'm/s^2'; 'rad/s'};
Имена и модули (начальных) состояний и параметров модели задаются посредством SETINIT. Мы также используем эту команду, чтобы указать, что первое начальное состояние (продольная скорость) должно быть строго положительным для модели, чтобы она была действительной, и указать, что все параметры модели должны быть строго положительными. Эти ограничения впоследствии будут учитываться при выполнении оценки начального состояния и/или параметра модели.
nlgr = setinit(nlgr, 'Name', {'Longitudinal vehicle velocity' ... % x(1). 'Lateral vehicle velocity' ... % x(2). 'Yaw rate'}); ... % x(3). nlgr = setinit(nlgr, 'Unit', {'m/s'; 'm/s'; 'rad/s'}); nlgr.InitialStates(1).Minimum = eps(0); % Longitudinal velocity > 0 for the model to be valid. nlgr = setpar(nlgr, 'Name', {'Vehicle mass'; ... % m. 'Distance from front axle to COG'; ... % a 'Distance from rear axle to COG'; ... % b. 'Longitudinal tire stiffness'; ... % Cx. 'Lateral tire stiffness'; ... % Cy. 'Air resistance coefficient'}); ... % CA. nlgr = setpar(nlgr, 'Unit', {'kg'; 'm'; 'm'; 'N'; 'N/rad'; '1/m'}); nlgr = setpar(nlgr, 'Minimum', num2cell(eps(0)*ones(6, 1))); % All parameters > 0!
Четыре из шести параметров этой структуры модели могут быть легко получены через лист данных рассматриваемого транспортного средства:
m = 1700 kg a = 1.5 m b = 1.5 m CA = 0.5 or 0.7 1/m (see below)
Поэтому мы не будем оценивать эти параметры:
nlgr.Parameters(1).Fixed = true; nlgr.Parameters(2).Fixed = true; nlgr.Parameters(3).Fixed = true; nlgr.Parameters(6).Fixed = true;
При этом текстовые сводные данные введенной структуры модели IDNLGREY получаются посредством PRESENT следующим образом.
present(nlgr);
nlgr =
Continuous-time nonlinear grey-box model defined by 'vehicle_c' (MEX-file):
dx/dt = F(t, u(t), x(t), p1, ..., p6)
y(t) = H(t, u(t), x(t), p1, ..., p6) + e(t)
with 5 input(s), 3 state(s), 3 output(s), and 2 free parameter(s) (out of 6).
Inputs:
u(1) Slip on front left tire(t) [ratio]
u(2) Slip on front right tire(t) [ratio]
u(3) Slip on rear left tire(t) [ratio]
u(4) Slip on rear right tire(t) [ratio]
u(5) Steering angle(t) [rad]
States: Initial value
x(1) Longitudinal vehicle velocity(t) [m/s] xinit@exp1 1 (fixed) in ]0, Inf]
x(2) Lateral vehicle velocity(t) [m/s] xinit@exp1 0 (fixed) in [-Inf, Inf]
x(3) Yaw rate(t) [rad/s] xinit@exp1 0 (fixed) in [-Inf, Inf]
Outputs:
y(1) Long. velocity(t) [m/s]
y(2) Lat. accel.(t) [m/s^2]
y(3) Yaw rate(t) [rad/s]
Parameters: Value
p1 Vehicle mass [kg] 1700 (fixed) in ]0, Inf]
p2 Distance from front axle to COG [m] 1.5 (fixed) in ]0, Inf]
p3 Distance from rear axle to COG [m] 1.5 (fixed) in ]0, Inf]
p4 Longitudinal tire stiffness [N] 150000 (estimated) in ]0, Inf]
p5 Lateral tire stiffness [N/rad] 40000 (estimated) in ]0, Inf]
p6 Air resistance coefficient [1/m] 0.5 (fixed) in ]0, Inf]
Name: Bicycle vehicle model
Status:
Created by direct construction or transformation. Not estimated.
More information in model's "Report" property.
Входные-выходные данные
На данной точке мы загружаем доступные данные ввода-вывода. Этот файл содержит данные трех различных экспериментов:
A. Simulated data with high stiffness tires [y1 u1]. B. Simulated data with low stiffness tires [y2 u2]. C. Measured data from a Volvo V70 [y3 u3].
Во всех случаях шаг расчета Ts = 0,1 секунды.
load(fullfile(matlabroot, 'toolbox', 'ident', 'iddemos', 'data', 'vehicledata'));
А. Система идентификации с использованием моделируемых данных о высокой жесткости шины
В нашем первом эксперименте по идентификации транспортного средства мы учитываем моделируемые данные о высокой жесткости шин. Сначала создается копия структуры модели nlgr и объекта IDDATA z1, отражающих эту конкретную ситуацию моделирования. 5 входных сигналов сохраняются в u1, а 3 выходных сигнала в y1. Входы скольжения (сгенерированные из сигналов скорости колеса) для передних колес были выбраны синусоидальными с постоянным смещением; скорость рыскания также была синусоидальной, но с другой амплитудой и частотой. В действительности это несколько искусственная ситуация, потому что редко так будоражит транспортное средство в боковом направлении.
nlgr1 = nlgr; nlgr1.Name = 'Bicycle vehicle model with high tire stiffness'; z1 = iddata(y1, u1, 0.1, 'Name', 'Simulated high tire stiffness vehicle data'); z1.InputName = nlgr1.InputName; z1.InputUnit = nlgr1.InputUnit; z1.OutputName = nlgr1.OutputName; z1.OutputUnit = nlgr1.OutputUnit; z1.Tstart = 0; z1.TimeUnit = 's';
Входы и выходы показаны на двух рисунках.
h_gcf = gcf; set(h_gcf,'DefaultLegendLocation','southeast'); h_gcf.Position = [100 100 795 634]; for i = 1:z1.Nu subplot(z1.Nu, 1, i); plot(z1.SamplingInstants, z1.InputData(:,i)); title(['Input #' num2str(i) ': ' z1.InputName{i}]); xlabel(''); axis tight; end xlabel([z1.Domain ' (' z1.TimeUnit ')']);
Фигура 2: Входы в систему транспортного средства с высокой жесткостью шины.
for i = 1:z1.Ny subplot(z1.Ny, 1, i); plot(z1.SamplingInstants, z1.OutputData(:,i)); title(['Output #' num2str(i) ': ' z1.OutputName{i}]); xlabel(''); axis tight; end xlabel([z1.Domain ' (' z1.TimeUnit ')']);
Фигура 3: Выходы системы транспортного средства с высокой жесткостью шины.
Следующим шагом является исследование эффективности начальной модели и для этого мы выполняем симуляцию. Заметьте, что начальное состояние было зафиксировано на ненулевом значении, поскольку первое состояние (продольная скорость транспортного средства) используется в качестве знаменателя в структуре модели. Сравнение между истинным и моделируемым выходами (с начальной моделью) показано в окне графика.
clf compare(z1, nlgr1, [], compareOptions('InitialCondition', 'model'));
Фигура 4: Сравнение истинных выходов и моделируемых выходов исходной модели транспортного средства с высокой жесткостью шины.
В порядок улучшить подгонку модели, далее оцениваются два параметра жесткости шины Cx и Cy, и проводится новая симуляция с предполагаемой моделью.
nlgr1 = nlgreyest(z1, nlgr1);
Сравнение между истинным и моделируемым выходами (с оценочной моделью) показано в окне графика.
compare(z1, nlgr1, [], compareOptions('InitialCondition', 'model'));
Фигура 5: Сравнение истинных выходов и моделируемых выходов предполагаемой модели транспортного средства с высокой жесткостью шины.
Эффективность симуляции предполагаемой модели довольно хорошая. Предполагаемые параметры жесткости также близки к параметрам, используемым в Simulink ®, чтобы сгенерировать истинные выходные данные:
disp(' True Estimated'); fprintf('Longitudinal stiffness: %6.0f %6.0f\n', 2e5, nlgr1.Parameters(4).Value); fprintf('Lateral stiffness : %6.0f %6.0f\n', 5e4, nlgr1.Parameters(5).Value);
True Estimated Longitudinal stiffness: 200000 198517 Lateral stiffness : 50000 53752
B. Система идентификации с использованием моделируемых данных о низкой жесткости шины
Во втором эксперименте мы повторяем моделирование из первого эксперимента, но теперь с имитированными данными о низкой жесткости шин.
nlgr2 = nlgr; nlgr2.Name = 'Bicycle vehicle model with low tire stiffness'; z2 = iddata(y2, u2, 0.1, 'Name', 'Simulated low tire stiffness vehicle data'); z2.InputName = nlgr2.InputName; z2.InputUnit = nlgr2.InputUnit; z2.OutputName = nlgr2.OutputName; z2.OutputUnit = nlgr2.OutputUnit; z2.Tstart = 0; z2.TimeUnit = 's';
Входы и выходы показаны на двух рисунках.
clf for i = 1:z2.Nu subplot(z2.Nu, 1, i); plot(z2.SamplingInstants, z2.InputData(:,i)); title(['Input #' num2str(i) ': ' z2.InputName{i}]); xlabel(''); axis tight; end xlabel([z2.Domain ' (' z2.TimeUnit ')']);
Фигура 6: Входы в систему транспортного средства с низкой жесткостью шины.
clf for i = 1:z2.Ny subplot(z2.Ny, 1, i); plot(z2.SamplingInstants, z2.OutputData(:,i)); title(['Output #' num2str(i) ': ' z2.OutputName{i}]); xlabel(''); axis tight; end xlabel([z2.Domain ' (' z2.TimeUnit ')']);
Фигура 7: Выходы системы транспортного средства с низкой жесткостью шины.
Далее мы исследуем эффективность начальной модели (которая имеет те же параметры, что и исходная модель высокой жесткости шины). Сравнение между истинным и моделируемым выходами (с начальной моделью) показано в окне графика.
clf compare(z2, nlgr2, [], compareOptions('InitialCondition', 'model'));
Фигура 8: Сравнение истинных выходов и моделируемых выходов исходной модели транспортного средства с низкой жесткостью шины.
Далее оцениваются два параметра жесткости.
nlgr2 = nlgreyest(z2, nlgr2);
Сравнение между истинным и моделируемым выходами (с оценочной моделью) показано в окне графика.
compare(z2, nlgr2, [], compareOptions('InitialCondition', 'model'));
Фигура 9: Сравнение истинных выходов и моделируемых выходов предполагаемой модели транспортного средства с низкой жесткостью шины.
Симуляция эффективности предполагаемой модели снова действительно хороша. Даже с той же начальной точкой параметра, что и в случае высокой жесткости шины, предполагаемые параметры жесткости также близки к таковым, используемой в Simulink, чтобы сгенерировать истинные данные о выходе:
disp(' True Estimated'); fprintf('Longitudinal stiffness: %6.0f %6.0f\n', 1e5, nlgr2.Parameters(4).Value); fprintf('Lateral stiffness : %6.0f %6.0f\n', 2.5e4, nlgr2.Parameters(5).Value);
True Estimated Longitudinal stiffness: 100000 99573 Lateral stiffness : 25000 26117
С. Система идентификации с использованием измеренных данных Volvo V70
В заключительном эксперименте мы рассматриваем данные, собранные в V70 Volvo. Как и выше, мы делаем копию типового объекта модели транспортного средства nlgr и создаем новый объект IDDATA, содержащий измеренные данные. Здесь мы также увеличили коэффициент сопротивления воздуха с 0,50 до 0,70, чтобы лучше отразить ситуацию с V70 Volvo.
nlgr3 = nlgr; nlgr3.Name = 'Volvo V70 vehicle model'; nlgr3.Parameters(6).Value = 0.70; % Use another initial CA for the Volvo data. z3 = iddata(y3, u3, 0.1, 'Name', 'Volvo V70 data'); z3.InputName = nlgr3.InputName; z3.InputUnit = nlgr3.InputUnit; z3.OutputName = nlgr3.OutputName; z3.OutputUnit = nlgr3.OutputUnit; z3.Tstart = 0; z3.TimeUnit = 's';
Входы и выходы показаны на двух рисунках. Как видно, измеренные данные довольно зашумлены.
clf for i = 1:z3.Nu subplot(z3.Nu, 1, i); plot(z3.SamplingInstants, z3.InputData(:,i)); title(['Input #' num2str(i) ': ' z3.InputName{i}]); xlabel(''); axis tight; end xlabel([z3.Domain ' (' z3.TimeUnit ')']);
Фигура 10: Измеренные входы от транспортного средства Volvo V70.
clf for i = 1:z3.Ny subplot(z3.Ny, 1, i); plot(z3.SamplingInstants, z3.OutputData(:,i)); title(['Output #' num2str(i) ': ' z3.OutputName{i}]); xlabel(''); axis tight; end xlabel([z3.Domain ' (' z3.TimeUnit ')']);
Фигура 11: Измеренные выходы от транспортного средства Volvo V70.
Затем мы исследуем эффективность начальной модели с оценкой начальных состояний. Сравнение между истинным и моделируемым выходами (с начальной моделью) показано в окне графика.
nlgr3 = setinit(nlgr3, 'Value', {18.7; 0; 0}); % Initial value of initial states. clf compare(z3, nlgr3);
Фигура 12: Сравнение измеренных выходов и моделируемых выходов исходной модели V70 транспортного средства Volvo.
Параметры жесткости шины Cx и Cy оцениваются далее, в этом случае с помощью метода поиска Левенберга-Марквардта, после чего выполняется новая симуляция с оцененной моделью. В сложение мы здесь оцениваем начальное значение продольной скорости, в то время как начальные значения поперечной скорости и скорости рыскания остаются фиксированными.
nlgr3 = setinit(nlgr3, 'Fixed', {false; true; true}); nlgr3 = nlgreyest(z3, nlgr3, nlgreyestOptions('SearchMethod', 'lm'));
Сравнение между истинным и моделируемым выходами (с оценочной моделью) показано в окне графика.
compare(z3, nlgr3);
Фигура 13: Сравнение измеренных выходов и моделируемых выходов первой предполагаемой модели V70 транспортного средства Volvo.
Предполагаемые параметры жесткости конечной модели Volvo V70 разумны, но здесь неизвестно, каковы их реальные значения.
disp(' Estimated'); fprintf('Longitudinal stiffness: %6.0f\n', nlgr3.Parameters(4).Value); fprintf('Lateral stiffness : %6.0f\n', nlgr3.Parameters(5).Value);
Estimated Longitudinal stiffness: 108873 Lateral stiffness : 29964
Дополнительная информация о предполагаемой модели транспортного средства Volvo V70 получена посредством PRESENT. Здесь интересно отметить, что неопределенность, связанная с расчетной жесткостью боковой шины, довольно высока (и значительно выше, чем для продольной жесткости шины). Эта неопределенность частично возникает из-за того, что боковое ускорение изменяется так мало во время тестового привода.
present(nlgr3);
nlgr3 =
Continuous-time nonlinear grey-box model defined by 'vehicle_c' (MEX-file):
dx/dt = F(t, u(t), x(t), p1, ..., p6)
y(t) = H(t, u(t), x(t), p1, ..., p6) + e(t)
with 5 input(s), 3 state(s), 3 output(s), and 2 free parameter(s) (out of 6).
Inputs:
u(1) Slip on front left tire(t) [ratio]
u(2) Slip on front right tire(t) [ratio]
u(3) Slip on rear left tire(t) [ratio]
u(4) Slip on rear right tire(t) [ratio]
u(5) Steering angle(t) [rad]
States: Initial value
x(1) Longitudinal vehicle velocity(t) [m/s] xinit@exp1 17.6049 (estimated) in ]0, Inf]
x(2) Lateral vehicle velocity(t) [m/s] xinit@exp1 0 (fixed) in [-Inf, Inf]
x(3) Yaw rate(t) [rad/s] xinit@exp1 0 (fixed) in [-Inf, Inf]
Outputs:
y(1) Long. velocity(t) [m/s]
y(2) Lat. accel.(t) [m/s^2]
y(3) Yaw rate(t) [rad/s]
Parameters: ValueStandard Deviation
p1 Vehicle mass [kg] 1700 0 (fixed) in ]0, Inf]
p2 Distance from front axle to COG [m] 1.5 0 (fixed) in ]0, Inf]
p3 Distance from rear axle to COG [m] 1.5 0 (fixed) in ]0, Inf]
p4 Longitudinal tire stiffness [N] 108873 26.8501 (estimated) in ]0, Inf]
p5 Lateral tire stiffness [N/rad] 29963.5 217.877 (estimated) in ]0, Inf]
p6 Air resistance coefficient [1/m] 0.7 0 (fixed) in ]0, Inf]
Name: Volvo V70 vehicle model
Status:
Termination condition: Maximum number of iterations reached..
Number of iterations: 20, Number of function evaluations: 41
Estimated using Solver: ode45; Search: lm on time domain data "Volvo V70 data".
Fit to estimation data: [-374.2;29.74;34.46]%
FPE: 2.362e-07, MSE: 0.3106
More information in model's "Report" property.
Заключительные замечания
Оценка параметров жесткости шины на практике является довольно сложной проблемой. Во-первых, приближения, введенные в структуру модели выше, действительны для довольно узкой рабочей области, и данные во время высоких ускорений, торможения и т.д. не могут использоваться. Жесткость также изменяется в зависимости от состояния окружающей среды, например, окружающей температуры, температуры в шинах и условий поверхности дороги, которые не учитываются в используемой структуре модели. Во-вторых, оценка параметров жесткости в большой степени зависит от стиля вождения. Когда в основном идет прямо вперед, как в третьем эксперименте по идентификации, становится трудно оценить параметры жесткости (особенно боковой), или поставить другой путь, неопределенности параметра становятся довольно высокими.