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'a 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) для получения дополнительной информации об реализации DASSL.