Этот пример показывает, как использовать продолжение, чтобы постепенно расширить решение BVP больших интервалов.
Краевые задачи Falkner-Skan [1] являются результатом решений для подобия вязкого, несжимаемого, ламинарного течения по плоской пластине. Уравнение в качестве примера
.
Проблема создана на бесконечном интервале с Согласно граничным условиям
,
,
.
BVP не может быть решен на бесконечном интервале, и это невозможно решить BVP в очень большом конечном интервале. Вместо этого этот пример решает последовательность проблем, созданных на меньшем интервале проверять, что решение имеет сопоставимое поведение как . Эта практика разбивания проблемы в более простые проблемы, с решением каждой проблемы, возвращающейся в как исходное предположение для следующей проблемы, называется продолжением.
Чтобы решить эту систему уравнений в MATLAB, необходимо закодировать уравнения, граничные условия и опции прежде, чем вызвать решатель для краевой задачи bvp4c
. Вы любой может включать необходимые функции как локальные функции в конце файла (как сделано здесь) или сохранить их как отдельные, именованные файлы в директории на пути MATLAB.
Создайте функцию, чтобы закодировать уравнения. Эта функция должна иметь подпись dfdx = fsode(x,f)
, где:
x является независимой переменной.
f является зависимой переменной.
Можно переписать уравнение третьего порядка как систему уравнений первого порядка с помощью замен , , и . Уравнения становятся
,
,
.
Соответствующая функция
function dfdeta = fsode(x,f) b = 0.5; dfdeta = [ f(2) f(3) -f(1)*f(3) - b*(1 - f(2)^2) ]; end
Примечание: Все функции включены как локальные функции в конце примера.
Теперь, запишите функцию, которая возвращает остаточное значение граничных условий в граничных точках. Эта функция должна иметь подпись res = fsbc(f0,finf)
, где:
f0
является значением граничного условия в начале интервала.
finf
является значением граничного условия в конце интервала.
Чтобы вычислить остаточные значения, необходимо поместить граничные условия в форму . В этой форме граничные условия
,
,
.
Соответствующая функция
function res = fsbc(f0,finf) res = [f0(1) f0(2) finf(2) - 1]; end
Наконец, необходимо обеспечить исходное предположение для решения. Грубая сетка пяти точек и постоянного предположения, которое удовлетворяет граничные условия, достаточно хороша, чтобы получить сходимость на интервале . Переменная infinity
обозначает правый предел интервала интегрирования. Когда значение infinity
увеличивается на последующих итерациях от 3 до его максимального значения 6, решение от каждой предыдущей итерации действия как исходное предположение для следующей итерации.
infinity = 3; maxinfinity = 6; solinit = bvpinit(linspace(0,infinity,5),[0 0 1]);
Решите проблему в начальном интервале . Постройте решение и сравните значение к аналитическому значению [1].
sol = bvp4c(@fsode,@fsbc,solinit); x = sol.x; f = sol.y; plot(x,f(2,:),x(end),f(2,end),'o'); axis([0 maxinfinity 0 1.4]); title('Falkner-Skan Equation, Positive Wall Shear, \beta = 0.5.') xlabel('x') ylabel('df/dx') hold on
fprintf('Cebeci & Keller report that f''''(0) = 0.92768.\n')
Cebeci & Keller report that f''(0) = 0.92768.
fprintf('Value computed using infinity = %g is %7.5f.\n', ... infinity,f(3,1))
Value computed using infinity = 3 is 0.92915.
Теперь, решите проблему на прогрессивно больших интервалах путем увеличения значения infinity
для каждой итерации. Функция bvpinit
экстраполирует каждое решение нового интервала, чтобы действовать как исходное предположение для следующего значения infinity
. Каждая итерация распечатывает расчетное значение и накладывает график решения по предыдущим решениям. Когда infinity = 6
, сопоставимое поведение решения становится очевидным и значение очень близко к ожидаемому значению.
for Bnew = infinity+1:maxinfinity solinit = bvpinit(sol,[0 Bnew]); % Extend solution to new interval sol = bvp4c(@fsode,@fsbc,solinit); x = sol.x; f = sol.y; fprintf('Value computed using infinity = %g is %7.5f.\n', ... Bnew,f(3,1)) plot(x,f(2,:),x(end),f(2,end),'o'); drawnow end
Value computed using infinity = 4 is 0.92774.
Value computed using infinity = 5 is 0.92770.
Value computed using infinity = 6 is 0.92770.
hold off
Перечисленный здесь локальные функции помощника, которые решатель BVP bvp4c
вызывает, чтобы вычислить решение. Также можно сохранить эти функции как их собственные файлы в директории на пути MATLAB.
function dfdeta = fsode(x,f) % equation being solved dfdeta = [ f(2) f(3) -f(1)*f(3) - 0.5*(1 - f(2)^2) ]; end %------------------------------------------- function res = fsbc(f0,finf) % boundary conditions res = [f0(1) f0(2) finf(2) - 1]; end %-------------------------------------------
[1] Cebeci, T. и Х. Б. Келлер. "Стрельба и Параллельные Методы Стрельбы для Решения Уравнения пограничного слоя Falkner-Skan". J. Аккомпанируйте Физике, Изданию 7, 1971, стр 289-300.