24.2 Дифференциально-алгебраические уравнения

Функция daspk может быть использована для решения ДАУ вида

0 = f (x-dot, x, t),    x(t=0) = x_0, x-dot(t=0) = x-dot_0

где x-dot — производная от x. Уравнение решается с помощью решателя ДАУ Petzold’s DASPK.

: [x, xdot, istate, msg] = daspk (fcn, x_0, xdot_0, t, t_crit) ¶

Решить систему дифференциально-алгебраических уравнений.

daspk решает систему уравнений

0 = f (x, xdot, t)

с

x(t_0) = x_0, xdot(t_0) = xdot_0

Решение возвращается в матрицах x и xdot, где каждая строка в результирующих матрицах соответствует одному из элементов в векторе t. Первый элемент t должен быть t_0 и соответствовать начальному состоянию системы x_0 и её производной xdot_0, так что первая строка вывода x равна x_0, а первая строка вывода xdot равна xdot_0.

Первый аргумент, fcn, — это строка, встроенная функция или дескриптор функции, который называет функцию f для вычисления вектора остатков для набора уравнений. Она должна иметь вид

res = f (x, xdot, t)

в которой x, xdot и res — векторы, а t — скаляр.

Если fcn — массив строк из двух элементов или двухэлементный массив ячеек строк, встроенных функций или дескрипторов функций, то первый элемент называет функцию f, описанную выше, а второй элемент называет функцию для вычисления модифицированной якобианной матрицы

df       df
jac = -- + c ------
      dx     d xdot

Функция модифицированной якобианной матрицы должна иметь вид

jac = j (x, xdot, t, c)

Второй и третий аргументы функции daspk задают начальное состояние и его производные, а четвёртый аргумент задаёт вектор выходных времен, в которых требуется решение, включая время, соответствующее начальному состоянию.

Набор начальных состояний и их производных не обязательно должен быть согласован. Если они не согласованы, необходимо использовать функцию daspk_options для предоставления дополнительной информации, чтобы daspk могла вычислить согласованную начальную точку.

Пятый аргумент является необязательным и может быть использован для указания набора времен, до которых решатель ДАУ не должен интегрировать. Это полезно для предотвращения проблем с особенностями и точками, в которых производная имеет разрыв.

После успешного вычисления значение istate будет больше нуля (согласуется с Fortran-версией DASPK).

Если вычисление не успешно, значение istate будет меньше нуля, а msg будет содержать дополнительную информацию.

Вы можете использовать функцию daspk_options для установки необязательных параметров для daspk.

См. также: dassl.

: daspk_options () ¶

: val = daspk_options (opt) ¶

: daspk_options (opt, val) ¶

Запрос или установка параметров для функции daspk.

При вызове без аргументов выводятся имена всех доступных параметров и их текущие значения.

При одном аргументе возвращается значение параметра opt.

При вызове с двумя аргументами daspk_options устанавливает параметр opt в значение val.

Параметры включают

"absolute tolerance"

Абсолюльная погрешность. Может быть вектором или скаляром. Если вектором, он должен соответствовать размеру вектора состояния, и относительная погрешность также должна быть вектором той же длины.

"relative tolerance"

Относительная погрешность. Может быть вектором или скаляром. Если вектором, он должен соответствовать размеру вектора состояния, и абсолютная погрешность также должна быть вектором той же длины.

Тест локальной погрешности, применяемый на каждом шаге интегрирования,

abs (local error in x(i))
       <= rtol(i) * abs (Y(i)) + atol(i)

"compute consistent initial condition"

Обозначая дифференциальные переменные в векторе состояния как ‘Y_d’, а алгебраические переменные как ‘Y_a’, ddaspk может решить одну из двух задач инициализации:

1. При заданных Y_d рассчитать Y_a и Y’_d
2. При заданных Y’ рассчитать Y.

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

Если этот параметр установлен в ненулевое значение, вы также должны установить параметр "algebraic variables" для объявления алгебраических переменных в задаче.

"use initial condition heuristics"

Установите это значение в ненулевое, чтобы использовать начальные условия эвристики, описанные ниже.

"initial condition heuristics"

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

MXNIT

Максимальное количество итераций Ньютона (по умолчанию 5).

MXNJ

Максимальное количество вычислений якобианской матрицы (по умолчанию 6).

MXNH

Максимальное количество значений искусственного параметра шага, которые будут проверены, если параметр "compute consistent initial condition" установлен в 1 (по умолчанию 5).

Обратите внимание, что максимальное общее количество итераций Ньютона, разрешённых, равно MXNIT*MXNJ*MXNH если параметр "compute consistent initial condition" установлен в 1 и MXNIT*MXNJ если он установлен в 2.

LSOFF

Установите это значение в ненулевое, чтобы отключить алгоритм поиска линии (по умолчанию 0).

STPTOL

Минимальный масштабированный шаг в алгоритме поиска линии (по умолчанию eps^(2/3)).

EPINIT

Коэффициент колебания в тесте сходимости итерации Ньютона. Тест применяется к вектору остатков, умноженному слева на приближённую якобианскую матрицу. Для сходимости взвешенная среднеквадратичная норма этого вектора (масштабированная весовыми коэффициентами) должна быть меньше EPINIT*EPCON, где EPCON = 0.33 — аналогичная тестовая константа, используемая в шагах времени. По умолчанию EPINIT = 0.01.

"print initial condition info"

Установите этот параметр в ненулевое значение, чтобы отобразить подробную информацию о расчёте начального состояния (по умолчанию 0).

"exclude algebraic variables from error test"

Установите в ненулевое значение, чтобы исключить алгебраические переменные из теста погрешности. Необходимо также установить параметр "algebraic variables" для объявления алгебраических переменных в задаче (по умолчанию 0).

"algebraic variables"

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

Этот параметр необходим для параметров "compute consistent initial condition" и "exclude algebraic variables from error test".

"enforce inequality constraints"

Установите в одно из следующих значений, чтобы применить ограничения неравенств, заданных параметром "inequality constraint types" (по умолчанию 0).

1. Для проверки ограничений только в расчете начального состояния.
2. Для применения проверки ограничений во время интегрирования.
3. Для применения обоих вариантов 1 и 2.

"inequality constraint types"

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

-2

Меньше нуля.

-1

Меньше или равно нулю.

0

Без ограничений.

1

Больше или равно нулю.

2

Больше нуля.

Этот параметр действует только если параметр "enforce inequality constraints" имеет ненулевое значение.

"initial step size"

Дифференциально-алгебраические задачи иногда страдают от серьёзных проблем со масштабированием на первом шаге. Если вы знаете много о масштабировании вашей задачи, вы можете помочь устранить эту проблему, задав начальный шаг (по умолчанию рассчитывается автоматически).

"maximum order"

Ограничить максимальный порядок метода решения. Этот параметр должен быть от 1 до 5 включительно (по умолчанию 5).

"maximum step size"

Установка максимального шага интегрирования предотвратит пропуск очень больших областей (по умолчанию не задано).

Octave также включает DASSL, более раннюю версию DASPK, и DASRT, которые могут быть использованы для решения ДАУ с ограничениями (условиями остановки).

: [x, xdot, istate, msg] = dassl (fcn, x_0, xdot_0, t, t_crit) ¶

Решить систему дифференциально-алгебраических уравнений.

dassl решает систему уравнений

0 = f (x, xdot, t)

с

x(t_0) = x_0, xdot(t_0) = xdot_0

Решение возвращается в матрицах x и xdot, где каждая строка в результирующих матрицах соответствует одному из элементов в векторе t. Первый элемент t должен быть t_0 и соответствовать начальному состоянию системы x_0 и её производной xdot_0, так что первая строка выходных данных x — x_0, а первая строка выходных данных xdot — xdot_0.

Первый аргумент, fcn, — это строка, встроенная функция или дескриптор функции, который называет функцию f для вычисления вектора остатков для системы уравнений. Она должна иметь вид

res = f (x, xdot, t)

в котором x, xdot и res — векторы, а t — скаляр.

Если fcn — массив строк из двух элементов или двумерный массив ячеек из строк, встроенных функций или дескрипторов функций, то первый элемент называет функцию f, описанную выше, а второй элемент называет функцию для вычисления модифицированной якобиана

df       df
jac = -- + c ------
      dx     d xdot

Функция модифицированного якобиана должна иметь вид

jac = j (x, xdot, t, c)

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

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

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

После успешного вычисления значение istate будет больше нуля (в соответствии с Fortran-версией DASSL).

Если вычисление не удалось, значение istate будет меньше нуля, и msg будет содержать дополнительную информацию.

Вы можете использовать функцию dassl_options для задания необязательных параметров для dassl.

См. также: daspk, dasrt, lsode.

: dassl_options () ¶

: val = dassl_options (opt) ¶

: dassl_options (opt, val) ¶

Запрос или установка параметров для функции dassl.

При вызове без аргументов отображаются имена всех доступных параметров и их текущие значения.

При задании одного аргумента возвращается значение параметра opt.

При вызове с двумя аргументами dassl_options устанавливает параметр opt в значение val.

Параметры включают

"absolute tolerance"

Абсолютная погрешность. Может быть вектором или скаляром. Если это вектор, он должен соответствовать размерности вектора состояния, а относительная погрешность также должна быть вектором той же длины.

"relative tolerance"

Относительная погрешность. Может быть вектором или скаляром. Если это вектор, он должен соответствовать размерности вектора состояния, а абсолютная погрешность также должна быть вектором той же длины.

Тест локальной погрешности, применяемый на каждом шаге интегрирования,

abs (local error in x(i))
       <= rtol(i) * abs (Y(i)) + atol(i)

"compute consistent initial condition"

Если не равно нулю, dassl попытается вычислить согласованный набор начальных условий. Это, как правило, не надёжно, поэтому лучше предоставить согласованный набор и оставить этот параметр равным нулю.

"enforce nonnegativity constraints"

Если известно, что решения ваших уравнений всегда неотрицательны, это может помочь установить этот параметр в ненулевое значение. Однако, вероятно, лучше сначала попробовать оставить этот параметр равным нулю и установить его в ненулевое значение только в том случае, если это не работает достаточно хорошо.

"initial step size"

Дифференциально-алгебраические задачи иногда могут страдать от серьёзных проблем масштабирования на первом шаге. Если вы знаете много о масштабировании вашей задачи, вы можете помочь решить эту проблему, указав начальный шаг.

"maximum order"

Ограничить максимальный порядок метода решения. Этот параметр должен быть в диапазоне от 1 до 5 включительно.

"maximum step size"

Установление максимального шага позволит избежать прохождения очень больших областей (значение по умолчанию не указано).

"step limit"

Максимальное количество шагов интегрирования, которые необходимо предпринять при единственном вызове подлежащего Fortran-кода.

: [x, xdot, t_out, istat, msg] = dasrt (fcn, g, x_0, xdot_0, t) ¶

: … = dasrt (fcn, g, x_0, xdot_0, t, t_crit) ¶

: … = dasrt (fcn, x_0, xdot_0, t) ¶

: … = dasrt (fcn, x_0, xdot_0, t, t_crit) ¶

Решает систему дифференциально-алгебраических уравнений.

dasrt решает систему уравнений

0 = f (x, xdot, t)

с

x(t_0) = x_0, xdot(t_0) = xdot_0

с функциональными критериями остановки (решение корней).

Решение возвращается в матрицах x и xdot, где каждая строка в результирующих матрицах соответствует одному из элементов в векторе t_out. Первый элемент t должен быть t_0 и соответствовать начальному состоянию системы x_0 и её производной xdot_0, так что первая строка выходных данных x равна x_0, а первая строка выходных данных xdot равна xdot_0.

Вектор t задаёт верхнюю границу длины интегрирования. Если условие остановки выполняется, вектор t_out будет короче t, а последний элемент t_out будет точкой, в которой выполнилось условие остановки, и может не соответствовать ни одному элементу вектора t.

Первый аргумент, fcn, является строкой, встроенной функцией или дескриптором функции, который называет функцию f для вычисления вектора остатков для системы уравнений. Она должна иметь вид

res = f (x, xdot, t)

в котором x, xdot и res являются векторами, а t является скаляром.

Если fcn является двумерным массивом строк или двумерным массивом ячеек строк, встроенных функций или дескрипторов функций, первый элемент называет функцию f, описанную выше, а второй элемент называет функцию для вычисления модифицированной якобианной матрицы

df       df
jac = -- + c ------
      dx     d xdot

Функция модифицированной якобианной матрицы должна иметь вид

jac = j (x, xdot, t, c)

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

g_out = g (x, t)

и возвращать вектор значений функций ограничений. Если значение какой-либо из функций ограничений меняет знак, DASRT попытается остановить интегрирование в точке изменения знака.

Если имя функции ограничения опущено, dasrt решает ту же задачу, что и daspk или dassl.

Обратите внимание, что из-за числовых ошибок в функциях ограничений, обусловленных округлением и ошибками интегрирования, DASRT может возвращать ложные корни или возвращать один и тот же корень в двух или более близких значениях T. Если такие ложные корни вызывают подозрение, пользователь должен рассмотреть меньшие погрешности или более высокую точность в оценке функций ограничений.

Если корень какой-либо функции ограничения определяет конец задачи, вход в DASRT должен тем не менее позволять интегрирование в точку немного за этим корнем, чтобы DASRT мог найти корень с помощью интерполяции.

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

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

Шестой аргумент является необязательным и может использоваться для задания набора времён, до которых решатель DAE не должен интегрировать. Это полезно для избежания трудностей с особенностями и точками, где есть разрыв в производной.

После успешного вычисления значение istate будет больше нуля (согласовано с Fortran-версией DASSL).

Если вычисление не завершилось успешно, значение istate будет меньше нуля, а msg будет содержать дополнительную информацию.

Вы можете использовать функцию dasrt_options для установки необязательных параметров для dasrt.

См. также: dasrt_options, daspk, dasrt, lsode.

: dasrt_options () ¶

: val = dasrt_options (opt) ¶

: dasrt_options (opt, val) ¶

Запрос или установка параметров для функции dasrt.

При вызове без аргументов выводятся имена всех доступных параметров и их текущие значения.

При вызове с одним аргументом возвращается значение параметра opt.

При вызове с двумя аргументами dasrt_options устанавливает параметр opt в значение val.

Параметры включают

"absolute tolerance"

Абсолютная погрешность. Может быть вектором или скаляром. Если вектор, он должен соответствовать размеру вектора состояния, и относительная погрешность также должна быть вектором той же длины.

"relative tolerance"

Относительная погрешность. Может быть вектором или скаляром. Если вектор, он должен соответствовать размеру вектора состояния, и абсолютная погрешность также должна быть вектором той же длины.

Локальный тест на погрешность, применяемый на каждом шаге интегрирования, имеет вид

abs (local error in x(i)) <= ...
      rtol(i) * abs (Y(i)) + atol(i)

"initial step size"

Дифференциально-алгебраические задачи иногда могут испытывать серьёзные трудности со масштабированием на первом шаге. Если вам известно много о масштабировании вашей задачи, вы можете помочь уменьшить эту проблему, указав начальный шаг.

"maximum order"

Ограничение максимального порядка метода решения.

"maximum step size"

Установка максимального шага позволит избежать прохождения очень больших участков.

"step limit"

Максимальное количество шагов интегрирования, которые нужно попробовать при одном вызове базового Fortran-кода.

См. K. E. Brenan и др., Численное решение задач Коши в дифференциально-алгебраических уравнениях, North-Holland (1989), DOI: https://doi.org/10.1137/1.9781611971224, для получения дополнительной информации об реализации DASSL.