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’а 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.

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