Изучение якобианской структуры решателя модели

Открыть модель

В этом примере показано, как использовать Simulink ® для изучения шаблона разреженности якобиана решателя и связи между шаблоном разреженности якобиана решателя и зависимостью между компонентами физической системы. Используется модель Simulink, моделирующая синхронизацию трех метроном, размещенных на свободном движущемся основании.

Якобианский шаблон решателя

Как правило, непрерывная часть модели Simulink может быть записана следующим образом:

$
\left\{ {\begin{array}{*{20}c}
 {\dot x = f(x,u)} \\
 {y = g(x,u)} \\
\end{array}} \right.
$

где $$ x $$ - непрерывные состояния и $$ u $$ - входы. $$ y $$ являются выходами.

Матрица:

$J_x = \frac{{\partial f}}{{\partial x}}$

называется матрицей якобиана решателя системы. Если для решения системных уравнений используется неявный решатель ОДУ $$ J_x $$ , его необходимо вычислить при необходимости. Например, матрица Якобиана набора ОДУ первого порядка

$
\left\{ {\begin{array}{*{20}c}
 {\dot x_1 = x_1 x_2 } \\
 {\dot x_2 = x_2 } \\
\end{array}} \right.
$

является

$
J_x = \left[ {\begin{array}{*{20}c}
 {\frac{{\partial f_1 }}{{\partial x_1 }}} & {\frac{{\partial f_1 }}{{\partial x_2 }}} \\
 {\frac{{\partial f_2 }}{{\partial x_1 }}} & {\frac{{\partial f_2 }}{{\partial x_2 }}} \\
\end{array}} \right] = \left[ {\begin{array}{*{20}c}
 {x_2 } & {x_1 } \\
 0 & 1 \\
\end{array}} \right]
$

Можно преобразовать матрицу якобиана решателя в логическую разреженную матрицу, представив каждый ненулевой элемент как 1, а каждый элемент, который всегда равен нулю (жесткий ноль), как 0. Например, булевой матрицей, соответствующей вышеуказанному якобиану $$ J_x $$ , является:

$
J_{xp} = \left[ {\begin{array}{*{20}c}
 1 & 1 \\
 0 & 1 \\
\end{array}} \right],
$

где $J_{xp}$ называется решателем Jacobian матрицы массива.

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

$
J_{xp} (i,j) = \left\{ {\begin{array}{*{20}c}
 {1~{\rm }if~{\rm }x_j~{\rm }is~{\rm }explicitly~{\rm }used~{\rm }to~{\rm }compute~{\rm }\dot x_i } \\
 {0~{\rm }otherwise} \\
\end{array}} \right.
$

Как видно, решающая матрица якобиана фактически представляет зависимость между переменными состояния и их производными. Если для вычисления $\dot x_i$ требуется значение, $$ x_j $$ то существует зависимость и $\dot x_i \leftarrow x_j$ . $J_{xp} (i,j) = 1$ Эти зависимости определяются физической природой системы, и, таким образом, изучая решатель Jacobian матрицы, можно исследовать физическую структуру физической системы, представленной моделью. Simulink предоставляет API для пользователя, чтобы получить решатель Jacobian матрицы шаблона. Ниже показано, как получить доступ к шаблону якобиана решателя и использовать его для исследования модели.

Шаблон и зависимость: синхронизация метрономов

Синхронизация определяется как корректировка ритмов колеблющихся объектов из-за их слабого взаимодействия [1]. Об одном из первых научно задокументированных наблюдений синхронизации сообщил голландский ученый Кристиан Гюйгенс, изобретатель маятниковых часов [2]. Гюйгенс обнаружил, что два маятника, прикрепленные к одной балке, поддерживаемой стульями, будут качать в прямо противоположном направлении через некоторое время. Аналогичная конфигурация, используемая в этом примере, показана на фиг.1.

Рис. 1: Настройка, используемая в этом примере: три метронома на движущейся базе

Моделирование системы

Модель этой физической системы можно разделить на две части:

Механизм метрономов
Движущееся основание

Механизм метрономов

Как показано на фиг.1, динамические уравнения одного метронома на движущейся основе могут быть получены как [3]:

$
\underbrace {\ddot \theta }_{Angular~acceleration} +
\underbrace {\frac{{mrg}}{J}\sin \theta }_{Gravitational~term} +
\underbrace {\frac{b}{J}\left[ {\left( {\frac{\theta }{{\theta _r }}} \right)^2 - 1}
\right]\dot \theta }_{Escapement~and~damping~term} + \underbrace
{\left( {\frac{{rm\cos \theta }}{J}} \right)\ddot
x}_{Coupled~inertial~force~term} = 0~~~~~~~~~~(eq.1)
$

Первые два термина описывают простой маятник без трения. Третий термин используется для аппроксимации выхода * и любого демпфирования метронома. Этот термин увеличивает угловую скорость для $~\theta < \theta_{0}$ и уменьшает ее для $\theta > \theta _{0}$ [3]. Последний член - это связанный эффект от движущегося основания, в терминах инерционной силы.

Движущееся основание

Предполагая, что движение основания безфрикционно, центр масс системы не изменится и можно вывести следующее:

$
x = - \frac{m}{{M + 3m}}r(\sin \theta _1 + \sin \theta _2 + \sin
\theta _3 ),
$

где $$ M $$ - масса основания и $$ m $$ - масса маятника.

Затем eq.1 можно изменить на:

$
\ddot x = - \frac{m}{{M + 3m}}r\left[ {\ddot \theta _1 \cos \theta _1 -
\sin \theta _1 (\dot \theta _1 )^2 + \ddot \theta _2 \cos \theta _2 -
\sin \theta _2 (\dot \theta _2 )^2 + \ddot \theta _3 \cos \theta _3 -
\sin \theta _3 (\dot \theta _3 )^2 } \right]~~~~~~~~~~(eq.2)
$

Ниже показана вышеописанная система, реализованная с использованием Simulink. Модель содержит три метрономные подсистемы и движущуюся базу.

Рис. 2: Модель системы метрономов Simulink

Моделирование этой системы показывает интересный феномен: Синхронизация. Это показывает, что все три метронома с разным начальным фазовым углом в конечном итоге синхронизируются друг с другом. На фиг.3 показаны результаты моделирования. Основной причиной синхронизации является движущаяся база, связывающая все эти метрономы вместе. Это физическое соединение можно наблюдать из динамического уравнения каждого метронома.

Кроме того, это физическое соединение также может наблюдаться из якобианского шаблона решателя этой модели. Следующий код MATLAB ® показывает, как получить шаблон якобиана решателя.

Рис. 3: Синхронизированные метрономы

Действия по получению шаблона якобиана решателя

% 1. Set the solver to be any implicit solver

   set_param(mdl, 'Solver', 'ode15s');

% 2. Set the solver Jacobian method to be Sparse perturbation *

   set_param(mdl, 'SolverJacobianMethodControl', 'SparsePerturbation');

% 3. Get the solver Jacobian object.

   J = Simulink.Solver.getSlvrJacobianPattern(mdl);
   disp('J = ');
   disp(J);

% 4. Show the pattern in a figure. use the method J.show

   J.show;

% 5. Explore the pattern with the given state name and other information

   stateNames = J.stateNames;
   disp('stateNames = ');
   disp(stateNames);

Вы увидите следующие результаты:

J = 
  SlvrJacobianPattern with properties:

        Jpattern: [8x8 double]
     numColGroup: 6
        colGroup: [8x1 double]
      stateNames: {8x1 cell}
    blockHandles: [8x1 double]

stateNames = 
    {'sldemo_metro/Moving base/Integrator1'              }
    {'(sldemo_metro/Metronome1/Integrator2).(Theta1)'    }
    {'(sldemo_metro/Metronome2/Integrator2).(Theta2)'    }
    {'(sldemo_metro/Metronome3/Integrator2).(Theta3)'    }
    {'(sldemo_metro/Metronome3/Integrator1).(Thetadot_3)'}
    {'(sldemo_metro/Metronome2/Integrator1).(Thetadot_2)'}
    {'(sldemo_metro/Metronome1/Integrator1).(Thetadot_1)'}
    {'sldemo_metro/Moving base/Integrator'               }

Рис. 4. Якобианский шаблон решателя

Свойства объекта Jacobian Pattern решателя

Решатель якобианский массив J является объектом. Он содержит следующие свойства:

     Jpattern: A sparse mxArray which is the Jacobian pattern.

  numColGroup: Number of column groups.

     colgroup: A column partition of the sparse pattern matrix.

   stateNames: A cell array containing the state name of each state.

blockHandles: Block handles of the owner of each state.

Изучение шаблона якобиана решателя

Как показано на рис. 4. Во-первых, решатель Якобиана этой модели разрежен, а число ненулевых элементов равно 28. Во-вторых, число групп столбцов равно 6; меньше числа состояний 8.

Строка 1-4 указывает следующие зависимости:

$
 \begin{array}{l}
 \dot x_1 \leftarrow x_8 (\dot {baseposition} \leftarrow {basepseed}) \\
 \dot x_2 \leftarrow x_7 (\dot {Theta1} \leftarrow {Thetadot1} ) \\
 \dot x_3 \leftarrow x_6 (\dot {Theta2} \leftarrow {Thetadot2} ) \\
 \dot x_4 \leftarrow x_5 (\dot {Theta3} \leftarrow {Thetadot3}) \\
 \end{array}
$

Эти отношения ясны, поскольку скорость является производной позиции. В строке 5-8 показаны следующие уравнения:

$
 \begin{array}{l}
 \dot x_5 \leftarrow x_{2 - 7} \\
 \dot x_6 \leftarrow x_{2 - 7} \\
 \dot x_7 \leftarrow x_{2 - 7} \\
 \dot x_8 \leftarrow x_{2 - 7} \\
 \end{array}
$

Эти соотношения показывают, что для вычисления углового ускорения метрономов или ускорения движущегося основания необходимы угловое положение и угловая скорость метрономов, но не положение и скорость основания. Эти отношения можно найти, изучая экв. (1) и экв. (2) непосредственно.

Заключение

Якобианский шаблон решателя - это инструмент для изучения зависимостей данных между производными переменных состояния и переменными состояния. Эти соотношения обычно отражают определенные физические связи в физической системе. Используя предоставленные инструменты, можно обнаружить эти уравнения, связанные с моделью Simulink, даже без исходных динамических уравнений физической модели.

Ссылки

[1] Аркадий Пикосвкий, Майкл Розенблюм и Юрген Куртс. Синхронизация. Издательство Кембриджского университета, 2001 год.

[2] Уорд Т. Уд, Дизайн и экспериментальные результаты синхронизации метрономов, вдохновленные Кристианом Гюйгенсом, магистерская диссертация, Эйндховенский технологический университет, 2006.

[3] Панталеоне, Джеймс, Американский журнал физики, том 70, выпуск 10, стр. 992-1000, 2002.

Escapement - это набор механизмов, который управляет метрономом. Для получения дополнительной информации см. [2].

Документация