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

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

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

где x-dot — производная от x. Уравнение решается с помощью решателя ДГУ Петцольда 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"

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

"maximum step size"

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

"step limit"

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

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