Оценка градиентов решений PDE в произвольных точках
[___] = evaluateGradient( возвращает интерполированные значения градиентов в точках, указанных в results,querypoints)querypoints.
[___] = evaluateGradient(___, возвращает интерполированные значения градиентов для системы уравнений для индексов (компонентов) уравнений iU)iU. При решении системы эллиптических ПДЭ указать iU после входных аргументов в любом из предыдущих синтаксисов.
Первое измерение gradx, grady, и, в 3-D случае, gradz соответствует точкам запроса. Вторая размерность соответствует индексам уравнений iU.
[___] = evaluateGradient(___, возвращает интерполированные значения градиентов для зависящего от времени уравнения или системы зависящих от времени уравнений в моменты времени iT)iT. При вычислении градиента для PDE, зависящего от времени, укажите iT после входных аргументов в любом из предыдущих синтаксисов. Для системы зависимых от времени уравнений укажите оба временных индекса iT и индексы уравнений (компоненты) iU.
Первое измерение gradx, grady, и, в 3-D случае, gradz соответствует точкам запроса. Для одного PDE, зависящего от времени, второй размер соответствует шагам времени. iT. Для системы зависящих от времени PDE вторая размерность соответствует индексам уравнений iU, и третье измерение соответствует шагам времени iT.
Оцените градиенты решения скалярной эллиптической задачи вдоль прямой. Постройте график результатов.
Создайте решение задачи 1 на L-образной мембране с нулевыми граничными условиями Дирихле.
model = createpde; geometryFromEdges(model,@lshapeg); applyBoundaryCondition(model,'dirichlet','Edge',1:model.Geometry.NumEdges,'u',0); specifyCoefficients(model,'m',0,... 'd',0,... 'c',1,... 'a',0,... 'f',1); generateMesh(model,'Hmax',0.05); results = solvepde(model);
Оцените градиенты решения вдоль прямой линии из (x,y)=(-1,-1) кому (1,1). Постройте график результатов в виде графика quiver с помощью quiver.
xq = linspace(-1,1,101); yq = xq; [gradx,grady] = evaluateGradient(results,xq,yq); gradx = reshape(gradx,size(xq)); grady = reshape(grady,size(yq)); quiver(xq,yq,gradx,grady) xlabel('x') ylabel('y')
Рассчитайте градиенты для среднего времени выхода броуновской частицы из области, которая содержит поглощающие (бегущие) границы и отражающие границы. Используйте уравнение Пуассона с постоянными коэффициентами и 3-D геометрию прямоугольного блока для моделирования этой задачи.
Создайте решение для этой проблемы.
model = createpde; importGeometry(model,'Block.stl'); applyBoundaryCondition(model,'dirichlet','Face',[1,2,5],'u',0); specifyCoefficients(model,'m',0,... 'd',0,... 'c',1,... 'a',0,... 'f',2); generateMesh(model); results = solvepde(model);
Создайте сетку и интерполируйте градиенты решения в сетку.
[X,Y,Z] = meshgrid(1:16:100,1:6:20,1:7:50); [gradx,grady,gradz] = evaluateGradient(results,X,Y,Z);
Измените форму градиентов на форму сетки и постройте график градиентов.
gradx = reshape(gradx,size(X)); grady = reshape(grady,size(Y)); gradz = reshape(gradz,size(Z)); quiver3(X,Y,Z,gradx,grady,gradz) axis equal xlabel('x') ylabel('y') zlabel('z')
Решение скалярной эллиптической задачи и интерполяция градиентов решения в плотную сетку. Используйте матрицу запроса для задания сетки.
Создайте решение задачи -1 на L-образной мембране с нулевыми граничными условиями Дирихле.
model = createpde; geometryFromEdges(model,@lshapeg); applyBoundaryCondition(model,'dirichlet','Edge',1:model.Geometry.NumEdges,'u',0); specifyCoefficients(model,'m',0,... 'd',0,... 'c',1,... 'a',0,... 'f',1); generateMesh(model,'Hmax',0.05); results = solvepde(model);
Интерполяция градиентов решения в сетку от -1 до 1 в каждом направлении. Постройте график результата с помощью quiver функция печати.
v = linspace(-1,1,101); [X,Y] = meshgrid(v); querypoints = [X(:),Y(:)]'; [gradx,grady] = evaluateGradient(results,querypoints); quiver(X(:),Y(:),gradx,grady) xlabel('x') ylabel('y')
Увеличьте изображение определенной части графика, чтобы просмотреть дополнительные сведения. Например, ограничьте диапазон печати 0,2 в каждом направлении.
axis([-0.2 0.2 -0.2 0.2])
Оцените градиенты решения для двухкомпонентной эллиптической системы и постройте график результатов.
Создайте модель PDE для двух компонентов.
model = createpde(2);
Создайте геометрию 2-D в виде прямоугольника с круглым отверстием в центре. Дополнительные сведения о создании геометрии см. в примере в разделе «Решение PDE с постоянными граничными условиями».
R1 = [3,4,-0.3,0.3,0.3,-0.3,-0.3,-0.3,0.3,0.3]'; C1 = [1,0,0,0.1]'; C1 = [C1;zeros(length(R1)-length(C1),1)]; geom = [R1,C1]; ns = (char('R1','C1'))'; sf = 'R1 - C1'; g = decsg(geom,sf,ns);
Включите геометрию в модель и просмотрите геометрию.
geometryFromEdges(model,g); pdegplot(model,'EdgeLabels','on') axis equal axis([-0.4,0.4,-0.4,0.4])
Задайте граничные условия и коэффициенты.
specifyCoefficients(model,'m',0,... 'd',0,... 'c',1,... 'a',0,... 'f',[2; -2]); applyBoundaryCondition(model,'dirichlet','Edge',3,'u',[-1,1]); applyBoundaryCondition(model,'dirichlet','Edge',1,'u',[1,-1]); applyBoundaryCondition(model,'neumann','Edge',[2,4:8],'g',[0,0]);
Создайте сетку и решите проблему.
generateMesh(model,'Hmax',0.1);
results = solvepde(model);Интерполяция градиентов решения в сетку от -0,3 до 0,3 в каждом направлении для каждого из двух компонентов.
v = linspace(-0.3,0.3,15); [X,Y] = meshgrid(v); [gradx,grady] = evaluateGradient(results,X,Y,[1,2]);
Постройте график градиентов для первого компонента.
figure gradx1 = gradx(:,1); grady1 = grady(:,1); quiver(X(:),Y(:),gradx1,grady1) title('Component 1') axis equal xlim([-0.3,0.3])
Постройте график градиентов для второго компонента.
figure gradx2 = gradx(:,2); grady2 = grady(:,2); quiver(X(:),Y(:),gradx2,grady2) title('Component 2') axis equal xlim([-0.3,0.3])
Решите систему гиперболических PDE и оцените градиенты.
Импорт геометрии перекрытия для 3-D проблемы с тремя компонентами решения. Постройте график геометрии.
model = createpde(3); importGeometry(model,'Plate10x10x1.stl'); pdegplot(model,'FaceLabels','on','FaceAlpha',0.5)
Задайте граничные условия так, чтобы грань 2 была фиксированной (нулевое отклонение в любом направлении) и грань 5 имела нагрузку 1e3 в положительном z-направление. Эта нагрузка приводит к изгибу перекрытия вверх. Установите начальное условие, что решение равно нулю, и его производная по времени также равна нулю.
applyBoundaryCondition(model,'dirichlet','Face',2,'u',[0,0,0]); applyBoundaryCondition(model,'neumann','Face',5,'g',[0,0,1e3]); setInitialConditions(model,0,0);
Создайте коэффициенты PDE для уравнений линейной упругости. Задайте свойства материала, аналогичные свойствам стали. См. раздел Уравнения линейной упругости.
E = 200e9; nu = 0.3; specifyCoefficients(model,'m',1,... 'd',0,... 'c',elasticityC3D(E,nu),... 'a',0,... 'f',[0;0;0]);
Создание сетки, настройка Hmax кому 1.
generateMesh(model,'Hmax',1);Решение проблемы на время 0 через 5e-3 в шагах 1e-4. Возможно, вам придется подождать решения несколько минут.
tlist = 0:5e-4:5e-3; results = solvepde(model,tlist);
Оценка градиентов решения при фиксированном x- и z- координаты в центрах их диапазонов, 5 и 0,5 соответственно. Оценка для y от 0 до 10 с шагом 0,2. Получить только компонент 3, z-компонент.
yy = 0:0.2:10; zz = 0.5*ones(size(yy)); xx = 10*zz; component = 3; [gradx,grady,gradz] = evaluateGradient(results,xx,yy,zz,component,1:length(tlist));
Три проекции градиентов решения представляют собой массивы 51 на 1 на 51. Использовать squeeze для удаления одиночного размера. Удаление одиночного размера преобразует эти массивы в матрицы 51 на 51, что упрощает индексирование в них.
gradx = squeeze(gradx); grady = squeeze(grady); gradz = squeeze(gradz);
Печать интерполированного градиентного компонента grady вдоль y ось для следующих шести значений из интервала времени tlist.
figure t = [1:2:11]; for i = t p(i) = plot(yy,grady(:,i),'DisplayName', strcat('t=', num2str(tlist(i)))); hold on end legend(p(t)) xlabel('y') ylabel('grady') ylim([-5e-6, 20e-6])
results объект содержит решение и его градиент, рассчитанный в узловых точках треугольной или тетраэдрической сетки. Можно получить доступ к решению и трем компонентам градиента в узловых точках с помощью точечной нотации.
interpolateSolution и evaluateGradient разрешить интерполяцию решения и его градиента в пользовательскую сетку, например, заданную meshgrid.
contour | evaluateCGradient | interpolateSolution | PDEModel | quiver | quiver3 | StationaryResults | TimeDependentResults