Граничные условия путем записи функций

Примечание

ЭТОТ ПЭЙДЖ ОПИСЫВАЕТ УСТАРЕВШИЙ РАБОЧИЙ ПРОЦЕСС. Новые возможности не могут быть совместимы с устаревшим рабочим процессом. Для соответствующего шага в рекомендуемом рабочем процессе смотрите, Задают Граничные условия.

О граничных условиях путем записи функций

Этот раздел показывает, как выразить граничные условия для 2D геометрии с помощью устаревшего синтаксиса функций. Однако рекомендуемый способ выразить граничные условия состоит в том, чтобы использовать, Задают Граничные условия.

Чтобы использовать этот традиционный синтаксис, запишите функции с помощью шаблонов в Граничных условиях для Скалярного УЧП или Граничных условиях для Систем УЧП.

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

Для скалярного УЧП некоторые граничные сегменты могут иметь условия Дирихле, и некоторые граничные сегменты могли обобщить Неймановы условия.

Граничные условия Дирихле

hu = r,

где h и r могут быть функциями x, y, решения u, индекс сегмента ребра, и, для параболических и гиперболических уравнений, время.

Обобщенные Неймановы граничные условия $\vec{n} \cdot (c \nabla u) + q u = g$ на ∂ Ω.

$\vec{n}$ исходящий нормальный модуль. g и q являются функциями, определяемыми на ∂ Ω и могут быть функциями x, y, решения u, индекс сегмента ребра, и, для параболических и гиперболических уравнений, время.

Чтобы записать файл функции, скажем pdebound.m, используют следующий синтаксис:

[qmatrix,gmatrix,hmatrix,rmatrix] = pdebound(p,e,u,time)

Ваша функция возвращает матрицы qmatrix, gmatrix, hmatrix и rmatrix, на основе этих входных параметров:

p Точки в mesh (Данные о Mesh)
e Ребра конечного элемента в mesh, подмножестве всех ребер (Данные о Mesh)
u Решение УЧП
время Время, для параболического или гиперболического УЧП только

Если ваши граничные условия не зависят от u или time, теми входными параметрами является []. Если ваши граничные условия действительно зависят от u или time, то, когда u или time является NaN, гарантируйте, что выходные параметры, такие как qmatrix состоят из матриц NaN правильного размера. Это сигнализирует к решателям, таким как parabolic, использовать зависящий от времени или зависимый решением алгоритм.

Прежде, чем задать граничные условия, необходимо знать граничные метки. Смотрите Идентифицируют Граничные Метки.

Решатель УЧП, такой как assempde или adaptmesh, передает матричный p точек и e ребер. e имеет семь строк и столбцы ne, где вы не обязательно знаете заранее размер ne.

p является 2 _Np матрицей, где p(1,k) является x - координата точки k, и p(2,k) является y - координата точки k.
e является 7 ne матрицей, где
- e(1,k) является индексом первой точки ребра k.
- e(2,k) является индексом второй точки ребра k.
- e(5,k) является меткой ребра геометрии ребра k (см., Идентифицируют Граничные Метки).
e содержит запись для каждого ребра конечного элемента, которое находится на внешнем контуре.

Используйте следующий шаблон для своего массива данных граничных условий.

function [qmatrix,gmatrix,hmatrix,rmatrix] = pdebound(p,e,u,time)

ne = size(e,2); % number of edges
qmatrix = zeros(1,ne);
gmatrix = qmatrix;
hmatrix = zeros(1,2*ne);
rmatrix = hmatrix;

for k = 1:ne
    x1 = p(1,e(1,k)); % x at first point in segment
    x2 = p(1,e(2,k)); % x at second point in segment
    xm = (x1 + x2)/2; % x at segment midpoint
    y1 = p(2,e(1,k)); % y at first point in segment
    y2 = p(2,e(2,k)); % y at second point in segment
    ym = (y1 + y2)/2; % y at segment midpoint
    switch e(5,k)
        case {some_edge_labels}
            % Fill in hmatrix,rmatrix or qmatrix,gmatrix
        case {another_list_of_edge_labels}
            % Fill in hmatrix,rmatrix or qmatrix,gmatrix
        otherwise
            % Fill in hmatrix,rmatrix or qmatrix,gmatrix
    end
end

Для каждого столбца k в e запись k rmatrix является значением rmatrix в первой точке в ребре и записью ne + , k является значением во второй точке в ребре. Например, если r = x ² + y ⁴, то запишите эти строки:

rmatrix(k) = x1^2 + y1^4;
rmatrix(k+ne) = x2^2 + y2^4;

Синтаксис для hmatrix идентичен: запись k hmatrix является значением r в первой точке в ребре и записью k + ne, является значением во второй точке в ребре.

Для каждого столбца k в e запись k qmatrix является значением qmatrix в средней точке в ребре. Например, если q = x ² + y ⁴, то запишите эти строки:

qmatrix(k) = xm^2 + ym^4;

Синтаксис для gmatrix идентичен: запись k gmatrix является значением gmatrix в средней точке в ребре.

Если коэффициенты зависят от решения u, используйте элемент u(e(1,k)) в качестве значения решения в первой точке ребра k и u(e(2,k)) как значение решения во второй точке ребра k.

Например, рассмотрите следующую геометрию, прямоугольник с круговой дырой.

% Rectangle is code 3, 4 sides,
% followed by x-coordinates and then y-coordinates
R1 = [3,4,-1,1,1,-1,-.4,-.4,.4,.4]';
% Circle is code 1, center (.5,0), radius .2
C1 = [1,.5,0,.2]';
% Pad C1 with zeros to enable concatenation with R1
C1 = [C1;zeros(length(R1)-length(C1),1)];
geom = [R1,C1];

% Names for the two geometric objects
ns = (char('R1','C1'))';

% Set formula
sf = 'R1-C1';

% Create geometry
gd = decsg(geom,sf,ns);

% View geometry
pdegplot(gd,'EdgeLabels','on')
xlim([-1.1 1.1])
axis equal

Предположим, что граничными условиями на внешней границе (сегменты 1 - 4) является Дирихле, со значением u (x, y) = t (x – y), где t время. Предположим, что круговой контур (сегменты 5 - 8) имеет обобщенное Нейманово условие с q = 1 и g = x ² + y ².

Запишите следующий файл данных граничных условий, чтобы представлять граничные условия:

function [qmatrix,gmatrix,hmatrix,rmatrix] = pdebound(p,e,u,time)

ne = size(e,2); % number of edges
qmatrix = zeros(1,ne);
gmatrix = qmatrix;
hmatrix = zeros(1,2*ne);
rmatrix = hmatrix;

for k = 1:ne
    x1 = p(1,e(1,k)); % x at first point in segment
    x2 = p(1,e(2,k)); % x at second point in segment
    xm = (x1 + x2)/2; % x at segment midpoint
    y1 = p(2,e(1,k)); % y at first point in segment
    y2 = p(2,e(2,k)); % y at second point in segment
    ym = (y1 + y2)/2; % y at segment midpoint
    switch e(5,k)
        case {1,2,3,4} % rectangle boundaries
            hmatrix(k) = 1;
            hmatrix(k+ne) = 1;
            rmatrix(k) = time*(x1 - y1);
            rmatrix(k+ne) = time*(x2 - y2);
        otherwise % same as case {5,6,7,8}, circle boundaries
            qmatrix(k) = 1;
            gmatrix(k) = xm^2 + ym^2;
    end
end

Граничные условия для систем УЧП

Общие смешанные граничные условия для систем УЧП уравнений N (см. уравнения, которые Можно Решить Используя Устаревшие Функции)

$\begin{array}{l} h u = r \\ n \cdot (c \otimes \nabla u) + q u = g + h^{'} μ \end{array}$

Обозначение $n \cdot (c \otimes \nabla u)$ означает N-by-1 матрица с (i, 1) - компонент

$\sum_{j = 1}^{N} (потому что (α) c_{i, j, 1, 1} \frac{\partial}{\partial x} + потому что (α) c_{i, j, 1, 2} \frac{\partial}{\partial y} + \sin (α) c_{i, j, 2, 1} \frac{\partial}{\partial x} + \sin (α) c_{i, j, 2, 2} \frac{\partial}{\partial y}) u_{j}$

где исходящий вектор нормали контура $n = (потому что (α), \sin (α))$ . Для каждого сегмента ребра существует M, условиями Дирихле и h-матрицей является M-by-N, M ≥ 0. Обобщенное Нейманово условие содержит источник $h^{'} μ$ где решатель вычисляет множители Лагранжа µ, таким образом, что условия Дирихле удовлетворены.

Чтобы записать файл функции, скажем pdebound.m, используют следующий синтаксис:

[qmatrix,gmatrix,hmatrix,rmatrix] = pdebound(p,e,u,time)

Ваша функция возвращает матрицы qmatrix, gmatrix, hmatrix и rmatrix, на основе этих входных параметров:

p Точки в mesh (Данные о Mesh)
e Ребра конечного элемента в mesh, подмножестве всех ребер (Данные о Mesh)
u Решение УЧП
время Время, для параболического или гиперболического УЧП только

p является 2 _Np матрицей, где p(1,k) является x - координата точки k, и p(2,k) является y - координата точки k.
e является 7 ne матрицей, где
- e(1,k) является индексом первой точки ребра k.
- e(2,k) является индексом второй точки ребра k.
- e(5,k) является меткой ребра геометрии ребра k (см., Идентифицируют Граничные Метки).
e содержит запись для каждого ребра конечного элемента, которое находится на внешнем контуре.

Позвольте N быть размерностью системы УЧП; смотрите уравнения, которые Можно Решить Используя Устаревшие Функции. Используйте следующий шаблон для своего массива данных граничных условий.

function [qmatrix,gmatrix,hmatrix,rmatrix] = pdebound(p,e,u,time)

N = 3; % Set N = the number of equations
ne = size(e,2); % number of edges
qmatrix = zeros(N^2,ne);
gmatrix = zeros(N,ne);
hmatrix = zeros(N^2,2*ne);
rmatrix = zeros(N,2*ne);

for k = 1:ne
    x1 = p(1,e(1,k)); % x at first point in segment
    x2 = p(1,e(2,k)); % x at second point in segment
    xm = (x1 + x2)/2; % x at segment midpoint
    y1 = p(2,e(1,k)); % y at first point in segment
    y2 = p(2,e(2,k)); % y at second point in segment
    ym = (y1 + y2)/2; % y at segment midpoint
    switch e(5,k)
        case {some_edge_labels}
            % Fill in hmatrix,rmatrix or qmatrix,gmatrix
        case {another_list_of_edge_labels}
            % Fill in hmatrix,rmatrix or qmatrix,gmatrix
        otherwise
            % Fill in hmatrix,rmatrix or qmatrix,gmatrix
        
    end
end

Для массива данных граничных условий вы представляете матрицу h для каждого сегмента ребра как вектор, беря матрицу по столбцам, как hmatrix(:). Столбец k hmatrix соответствует матрице в первой точке ребра e(1,k) и столбец k + ne, соответствует матрице во второй точке ребра e(2,k).

Точно так же вы представляете каждый вектор r для ребра как столбец в матричном rmatrix. Столбец k соответствует вектору в первой точке ребра e(1,k) и столбец k + ne, соответствует вектору во второй точке ребра e(2,k).

Представляйте записи для матрицы q для каждого сегмента ребра как вектор, qmatrix(:), подобный матричному hmatrix(:). Точно так же представляйте g для каждого сегмента ребра, вектор-столбец в матричном gmatrix. В отличие от h и r, которые имеют два столбца для каждого сегмента, q и g имеют всего один столбец для каждого сегмента, который является значением функции в средней точке сегмента ребра.

Например, рассмотрите следующую геометрию, прямоугольник с круговой дырой.

% Rectangle is code 3, 4 sides,
% followed by x-coordinates and then y-coordinates
R1 = [3,4,-1,1,1,-1,-.4,-.4,.4,.4]';
% Circle is code 1, center (.5,0), radius .2
C1 = [1,.5,0,.2]';
% Pad C1 with zeros to enable concatenation with R1
C1 = [C1;zeros(length(R1)-length(C1),1)];
geom = [R1,C1];

% Names for the two geometric objects
ns = (char('R1','C1'))';

% Set formula
sf = 'R1-C1';

% Create geometry
gd = decsg(geom,sf,ns);

% View geometry
pdegplot(gd,'EdgeLabels','on')
xlim([-1.1 1.1])
axis equal

Предположим N = 3. Предположим, что граничные условия смешаны. Существует M = 1 условие Дирихле:

Первый компонент u = 0 на прямоугольных сегментах (числа 1-4). Так h (1,1) = 1 и r (1) = 0 для тех сегментов.
Вторые компоненты u = 0 на круговых сегментах (числа 5-8). Так h (2,2) = 1 и r (2) = 0 для тех сегментов.
На прямоугольных сегментах (числа 1-4),
$q = (\begin{matrix} 011000110 \end{matrix})$
и
$g = (\begin{matrix} 1 + x^{2} \\ 0 \\ 1 + y^{2} \end{matrix})$
На круговых сегментах (числа 5-8),
$q = (\begin{matrix} 0 & 1 + x^{2} & 2 + y^{2} \\ 000 \\ 1 + x^{4} & 1 + y^{4} & 0 \end{matrix})$
и
$g = (\begin{matrix} потому что (π x) \\ 0 \\ \tanh (x + y) \end{matrix})$

Запишите следующий файл данных граничных условий, чтобы представлять граничные условия:

function [qmatrix,gmatrix,hmatrix,rmatrix] = pdebound(p,e,u,time)

N = 3;
ne = size(e,2); % number of edges
qmatrix = zeros(N^2,ne);
gmatrix = zeros(N,ne);
hmatrix = zeros(N^2,2*ne);
rmatrix = zeros(N,2*ne);

for k = 1:ne
    x1 = p(1,e(1,k)); % x at first point in segment
    x2 = p(1,e(2,k)); % x at second point in segment
    xm = (x1 + x2)/2; % x at segment midpoint
    y1 = p(2,e(1,k)); % y at first point in segment
    y2 = p(2,e(2,k)); % y at second point in segment
    ym = (y1 + y2)/2; % y at segment midpoint
    switch e(5,k)
        case {1,2,3,4}
            hk = zeros(N);
            hk(1,1) = 1;
            hk = hk(:);
            hmatrix(:,k) = hk;
            hmatrix(:,k+ne) = hk;
            
            rk = zeros(N,1); % Not strictly necessary
            rmatrix(:,k) = rk; % These are already 0
            rmatrix(:,k+ne) = rk;
            
            qk = zeros(N);
            qk(1,2) = 1;
            qk(1,3) = 1;
            qk(3,1) = 1;
            qk(3,2) = 1;
            qk = qk(:);
            qmatrix(:,k) = qk;
            
            gk = zeros(N,1);
            gk(1) = 1+xm^2;
            gk(3) = 1+ym^2;
            gmatrix(:,k) = gk;
            
        case {5,6,7,8}
            hk = zeros(N);
            hk(2,2) = 1;
            hk = hk(:);
            hmatrix(:,k) = hk;
            hmatrix(:,k+ne) = hk;
            
            rk = zeros(N,1); % Not strictly necessary
            rmatrix(:,k) = rk; % These are already 0
            rmatrix(:,k+ne) = rk;

            qk = zeros(N);
            qk(1,2) = 1+xm^2;
            qk(1,3) = 2+ym^2;
            qk(3,1) = 1+xm^4;
            qk(3,2) = 1+ym^4;
            qk = qk(:);
            qmatrix(:,k) = qk;
            
            gk = zeros(N,1);
            gk(1) = cos(pi*xm);
            gk(3) = tanh(xm*ym);
            gmatrix(:,k) = gk;
            
        end
end

Документация