exponenta event banner

Решение BVP с несколькими граничными условиями

В этом примере показано, как решить задачу многоточечного граничного значения, где интересующее решение удовлетворяет условиям внутри интервала интегрирования.

Для x в [0, λ] рассмотрим уравнения

v′=C-1n,

C′=vC-min (x, 1) start.

Известные параметры проблемы - n, κ, λ> 1, и η =λ2n ⋅κ2.

Член min (x, 1) в уравнении для C (x) не является гладким при x = 1, поэтому задача не может быть решена напрямую. Вместо этого можно разбить задачу на два: один набор в интервале [0,1], а другой набор в интервале [1, λ]. Связь между двумя областями заключается в том, что решения должны быть непрерывными при x = 1. Решение должно также удовлетворять граничным условиям

v (0) = 0,

C (λ) = 1.

Уравнения для каждой области:

Регион 1: 0≤x≤1

v′=C-1n,

C′=vC-xη.

Регион 2: 1≤x≤λ

v′=C-1n,

C′=vC-1η.

Точка интерфейса x = 1 включена в обе области. В этот момент решатель создает как левое, так и правое решения, которые должны быть равными для обеспечения непрерывности решения.

Чтобы решить эту систему уравнений в MATLAB, необходимо кодировать уравнения, граничные условия и начальное предположение, прежде чем вызывать решатель задач граничных значений bvp5c. Требуемые функции можно либо включить в качестве локальных функций в конце файла (как здесь сделано), либо сохранить их как отдельные именованные файлы в каталоге по пути MATLAB.

Уравнения кода

Уравнения для v (x) и C (x) зависят от решаемой области. Для многоточечных граничных задач производная функция должна принимать третий входной аргументregion, который используется для идентификации области, где оценивается производная. Решатель нумерует области слева направо, начиная с 1.

Создание функции для кодирования уравнений с сигнатурой dydx = f(x,y,region,p), где:

  • x является независимой переменной.

  • y является зависимой переменной.

  • dydx(1) дает уравнение для v (x), иdydx(2) уравнение для C (x).

  • region - номер области, в которой вычисляется производная (в этой двухобластной задаче, region является 1 или 2).

  • p вектор, содержащий ценности постоянных параметров [n, κ, λ, η].

Оператор switch используется для возврата различных уравнений в зависимости от области решения. Функция:

function dydx = f(x,y,region,p) % equations being solved
n = p(1);
eta = p(4);

dydx = zeros(2,1);
dydx(1) = (y(2) - 1)/n;

switch region
    case 1    % x in [0 1]
        dydx(2) = (y(1)*y(2) - x)/eta;
    case 2    % x in [1 lambda]
        dydx(2) = (y(1)*y(2) - 1)/eta;
end
end

Примечание.В конце примера все функции включаются как локальные.

Граничные условия кода

Решение двух дифференциальных уравнений первого порядка в двух областях требует четырех граничных условий. Два из этих условий вытекают из первоначальной проблемы:

v (0) = 0,

C (λ) -1 = 0.

Два других условия обеспечивают непрерывность левого и правого решений в точке интерфейса x = 1:

vL (1) -vR (1) = 0,

CL (1) -CR (1) = 0.

Для многоточечных BVP, аргументы граничных условий функции YL и YR стать матрицами. В частности, kВ-й колонне YL(:,k) - решение на левой границе kВ-й регион. Аналогично, YR(:,k) - решение на правой границе kВ-й регион.

В этой задаче y (0) аппроксимируется YL(:,1), в то время как y (λ) аппроксимируется YR(:,end). Непрерывность решения при x = 1 требует, чтобыYR(:,1) = YL(:,2).

Функция, которая кодирует остаточное значение четырех граничных условий,

function res = bc(YL,YR)
res = [YL(1,1)               % v(0) = 0
       YR(1,1) - YL(1,2)     % Continuity of v(x) at x=1
       YR(2,1) - YL(2,2)     % Continuity of C(x) at x=1
       YR(2,end) - 1];       % C(lambda) = 1
end

Начальное предположение формы

Для многоточечных BVP граничные условия автоматически применяются в начале и конце интервала интегрирования. Однако необходимо указать двойные записи в xmesh для других точек интерфейса. Простое предположение, удовлетворяющее граничным условиям, является постоянным предположением y = [1; 1].

xc = 1;
xmesh = [0 0.25 0.5 0.75 xc xc 1.25 1.5 1.75 2];
yinit = [1; 1];
sol = bvpinit(xmesh,yinit);

Уравнение решения

Определите значения постоянных параметров и поместите их в вектор p. Предоставление функции bvp5c с синтаксисом @(x,y,r) f(x,y,r,p) для обеспечения вектора параметров.

Вычислите решение для нескольких значений, используя каждое решение в качестве начального предположения для следующего. Вычислите для каждого значения, равного, значение осмолярности Os = 1v (λ). Для каждой итерации цикла сравните вычисленное значение с приблизительным аналитическим решением.

lambda = 2;
n = 5e-2;
for kappa = 2:5
   eta = lambda^2/(n*kappa^2);
   p = [n kappa lambda eta];
   sol = bvp5c(@(x,y,r) f(x,y,r,p), @bc, sol);
   K2 = lambda*sinh(kappa/lambda)/(kappa*cosh(kappa));
   approx = 1/(1 - K2);
   computed = 1/sol.y(1,end);
   fprintf('  %2i    %10.3f     %10.3f \n',kappa,computed,approx);
end
   2         1.462          1.454 
   3         1.172          1.164 
   4         1.078          1.071 
   5         1.039          1.034 

Решение для построения графика

Постройте график компонентов решения для v (x) и C (x) и вертикальной линии в точке сопряжения x = 1. Отображаемое решение для δ = 5 является результатом конечной итерации цикла.

plot(sol.x,sol.y(1,:),'--o',sol.x,sol.y(2,:),'--o')
line([1 1], [0 2], 'Color', 'k')
legend('v(x)','C(x)')
title('A Three-Point BVP Solved with bvp5c')
xlabel({'x', '\lambda = 2, \kappa = 5'})
ylabel('v(x) and C(x)')

Figure contains an axes. The axes with title A Three-Point BVP Solved with bvp5c contains 3 objects of type line. These objects represent v(x), C(x).

Локальные функции

Здесь перечислены локальные вспомогательные функции, которые решатель BVP bvp5c вызывает для вычисления решения. Кроме того, эти функции можно сохранить в виде собственных файлов в каталоге по пути MATLAB.

function dydx = f(x,y,region,p) % equations being solved
n = p(1);
eta = p(4);

dydx = zeros(2,1);
dydx(1) = (y(2) - 1)/n;

switch region
    case 1    % x in [0 1]
        dydx(2) = (y(1)*y(2) - x)/eta;
    case 2    % x in [1 lambda]
        dydx(2) = (y(1)*y(2) - 1)/eta;
end
end
%-------------------------------------------
function res = bc(YL,YR) % boundary conditions
res = [YL(1,1)               % v(0) = 0
       YR(1,1) - YL(1,2)     % Continuity of v(x) at x=1
       YR(2,1) - YL(2,2)     % Continuity of C(x) at x=1
       YR(2,end) - 1];       % C(lambda) = 1
end
%-------------------------------------------

См. также

| |

Связанные темы