Иногда дифференциальные уравнения определяются с граничными условиями не только на концах интервала, но и с дополнительным условием в некоторой промежуточной точке расчетного интервала. Чаще всего такие задачи содержат данные о негладких в некоторой внутренней точке интервала решениях. Для них имеется встроенная функция
bvaifit, также реализующая алгоритм стрельбы.
- z1 — вектор, присваивающий недостающим начальным условиям на левой границе интервала начальные значения.
- z2 — вектор того же размера, присваивающий недостающим начальным условиям на правой границе интервала начальные значения.
- х0 — левая граница расчетного интервала.
- x1 — правая граница расчетного интервала.
- xf — точка внутри интервала.
- Dо(х,у) — векторная функция, описывающая систему N ОДУ размера Nx1 и двух аргументов — скалярного х и векторного у. При этом у — это неизвестная векторная функция аргумента х того же размера Nx1.
- load1(x0,z) — векторная функция размера Nx1 левых граничных условий, причем недостающие начальные условия поименовываются соответствующими компонентами векторного аргумента z.
- load2(x1,z) — векторная функция размера Nx1 левых граничных условий, причем недостающие начальные условия поименовываются соответствующими компонентами векторного аргумента z.
- score (xf, у) — векторная функция размера Nx1, выражающая внутреннее условие для векторной функции у в точке xf.
Система уравнений и левое краевое условие вводится так же, как и
в предыдущем листинге для функции sbval. Обратите внимание, что таким же образом записано и правое краевое условие. Для того чтобы ввести условие отражения на правой границе, пришлось определить еще один неизвестный пристрелочный параметр z20. Строка листинга, в которой определена функция
score, задает условие стрельбы — сшивку двух решений в точке xf. В самой последней строке листинга выдан ответ — определенные численным методом значения обоих пристрелочных параметров, которые объединены в вектор а (мы применили в предпоследней строке операцию транспонирования, чтобы результат получился в форме вектора, а не матрицы-строки). Для корректного построения графика решения лучше составить его из двух частей — решения задачи Коши на интервале (x0,xf) и другой задачи Коши для интервала (xf,x1). Реализация этого способа приведена в листинге 10.4, который является продолжением листинга 10.3. В последней строке листинга 10.4 выведено значение второй искомой функции на правой границе интервала. Всегда полезно проконтролировать, что оно совпадает с соответствующим пристрелочным параметром (выведенным в последней строке листинга 10.3).
Листинг 10.4. Решение краевой задачи (продолжение
листинга 10.3)
Решение краевой задачи приведено на рис. 10.5. С физической точки зрения естественно, что интенсивность света уменьшается быстрее по мере распространения в более плотной среде в правой половине расчетного интервала. В средней точке
xf=0.5, как и ожидалось, производные обоих решений имеют разрыв.
ПРИМЕЧАНИЕ
Еще один пример решения краевых задач с разрывными коэффициентами ОДУ приведен в справочной системе Mathcad.
Рис. 10.5. Решение краевой задачи с разрывом в xf=0.5 (продолжение листингов 10.3—10.4)
Ради справедливости необходимо заметить, что разобранную краевую задачу легко решить и с помощью функции
sbval, заменив в листинге 10.2 зависимость а (х) на третью строку листинга 10.3. В этом случае листинг 10.2 даст в
точности тот же ответ, что показан на рис. 10.5. Однако в определенных случаях (в том числе из соображений быстроты расчетов) удобнее использовать функцию
bvaifit, т. е. вести пристрелку с обеих границ интервала.
СОВЕТ
Если вы имеете дело с подобными уравнениями, попробуйте сначала решить их как обычную краевую задачу с помощью более надежной и легкой в применении функции
sbval.