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

Настройте и решите УЧП

Задайте функцию для нелинейного f коэффициент в синтаксисе, как дали в f Коэффициенте для specifyCoefficients.

function bcMatrix = myfffun(region,state)

bcMatrix = 1+10*region.z.^2+region.y;

Импортируйте геометрию и исследуйте метки поверхности.

model = createpde;
g = importGeometry(model,'Plate10x10x1.stl');
pdegplot(g,'FaceLabels','on','FaceAlpha',0.5)

Поверхности пронумерованы 1 - 6.

Создайте коэффициенты и граничные условия.

c = 1;
a = 0;
d = 1;
f = @myfffun;
specifyCoefficients(model,'m',0,'d',d,'c',c,'a',a,'f',f);

applyBoundaryCondition(model,'dirichlet','face',1:6,'u',0);

Установите нулевое начальное условие.

setInitialConditions(model,0);

Не создайте mesh со сторонами больше, чем 0,3.

generateMesh(model,'Hmax',0.3);

Установите времена от 0 до 0.2 и решите УЧП.

tlist = 0:0.02:0.2;
results = solvepde(model,tlist);

Постройте срезы через решение

Создайте сетку (x,y,z) точки, где x = 5Y диапазоны от 0 до 10, и z диапазоны от 0 до 1. Интерполируйте решение этих узлов решетки и всех случаев.

yy = 0:0.5:10;
zz = 0:0.25:1;
[YY,ZZ] = meshgrid(yy,zz);
XX = 5*ones(size(YY));
uintrp = interpolateSolution(results,XX,YY,ZZ,1:length(tlist));

Матрица решения uintrp имеет 11 столбцов, один в течение каждого раза в tlist. Возьмите интерполированное решение для второго столбца, который соответствует времени 0.02.

usol = uintrp(:,2);

Элементы usol произойдите из интерполяции решения XX, YY, и ZZ матрицы, которые являются каждым 5 21, соответствуя z-by-y переменные. Измените usol к тому же самому 5 21 размер, и делают объемную поверхностную диаграмму решения. Также сделайте объемные поверхностные диаграммы, соответствующие временам 0.06, 0.10, и 0.20.

figure
usol = reshape(usol,size(XX));
subplot(2,2,1)
surf(usol)
title('t = 0.02')
zlim([0,1.5])
xlim([1,21])
ylim([1,5])

usol = uintrp(:,4);
usol = reshape(usol,size(XX));
subplot(2,2,2)
surf(usol)
title('t = 0.06')
zlim([0,1.5])
xlim([1,21])
ylim([1,5])

usol = uintrp(:,6);
usol = reshape(usol,size(XX));
subplot(2,2,3)
surf(usol)
title('t = 0.10')
zlim([0,1.5])
xlim([1,21])
ylim([1,5])

usol = uintrp(:,11);
usol = reshape(usol,size(XX));
subplot(2,2,4)
surf(usol)
title('t = 0.20')
zlim([0,1.5])
xlim([1,21])
ylim([1,5])

Настройте и решите УЧП

Проблема состоит в том, чтобы решить уравнение Пуассона с нулем граничные условия Дирихле для сложной геометрии. Уравнение Пуассона

Partial Differential Equation Toolbox™ решает уравнения в форме

Таким образом, можно представлять проблему путем установки $$c=1$$ и $$a=0$$. Произвольно набор $$f=10$$.

c = 1;
a = 0;
f = 10;

Первый шаг в решении любой 3-D задачи УЧП должен создать Модель УЧП. Это - контейнер, который содержит количество уравнений, геометрии, mesh и граничных условий для вашего УЧП. Создайте модель, затем импортируйте 'ForearmLink.stl' файл и представление геометрия.

N = 1;
model = createpde(N);
importGeometry(model,'ForearmLink.stl');
pdegplot(model,'FaceAlpha',0.5)
view(-42,24)

Задайте коэффициенты УЧП

Включайте коэффициенты УЧП в model.

specifyCoefficients(model,'m',0,'d',0,'c',c,'a',a,'f',f);

Создайте нуль граничные условия Дирихле на всех поверхностях.

applyBoundaryCondition(model,'dirichlet', ...
                             'Face',1:model.Geometry.NumFaces, ...
                             'u',0);

Создайте mesh и решите УЧП.

generateMesh(model);
result = solvepde(model);

Постройте решение как срезы контура

Поскольку граничные условия $$u=0$$ на всех поверхностях, решении $$u$$ является ненулевым только во внутренней части. Чтобы исследовать внутреннюю часть, возьмите прямоугольную сетку, которая покрывает геометрию интервалом одного модуля в каждом координатном направлении.

[X,Y,Z] = meshgrid(0:135,0:35,0:61);

Для графического вывода и анализа, создайте PDEResults объект из решения. Интерполируйте результат в каждом узле решетки.

V = interpolateSolution(result,X,Y,Z);
V = reshape(V,size(X));

Постройте срезы контура для различных значений $$z$$.

figure
colormap jet
contourslice(X,Y,Z,V,[],[],0:5:60)
xlabel('x')
ylabel('y')
zlabel('z')
colorbar
view(-11,14)
axis equal

Постройте срезы контура для различных значений $$y$$.

figure
colormap jet
contourslice(X,Y,Z,V,[],1:6:31,[])
xlabel('x')
ylabel('y')
zlabel('z')
colorbar
view(-62,34)
axis equal

Сохраните память путем оценки по мере необходимости

Для больших проблем у можно закончиться память при создании прекрасной 3-D сетки. Кроме того, это может быть длительно, чтобы оценить решение на полной сетке. Чтобы сохранить память и время, оцените только в точках, которые вы строите. Можно также использовать этот метод, чтобы интерполировать к наклоненным сеткам, или на другие поверхности.

Например, интерполируйте решение сетки на наклоненной плоскости $$0\le x\le 135$$, $$0\le y\le 35$$, и $$z=x/10+y/2$$. Постройте оба контура, и окрасил поверхностные данные. Используйте прекрасную сетку с разрядкой 0.2.

[X,Y] = meshgrid(0:0.2:135,0:0.2:35);
Z = X/10 + Y/2;
V = interpolateSolution(result,X,Y,Z);
V = reshape(V,size(X));
figure
subplot(2,1,1)
contour(X,Y,V);
axis equal
title('Contour Plot on Tilted Plane')
xlabel('x')
ylabel('y')
colorbar
subplot(2,1,2)
surf(X,Y,V,'LineStyle','none');
axis equal
view(0,90)
title('Colored Plot on Tilted Plane')
xlabel('x')
ylabel('y')
colorbar

Импортируйте и просмотрите геометрию

model = createpde;
importGeometry(model,'Block.stl');
pdegplot(model,'FaceLabels','on','FaceAlpha',0.5)
view(-42,24)

Условия границы множества

Установите поверхности 1, 2, и 5 быть местами, где частица может выйти. На этих поверхностях, решении $$u=0$$. Сохраните значение по умолчанию, отражающее граничные условия на поверхностях 3, 4, и 6.

applyBoundaryCondition(model,'dirichlet','Face',[1,2,5],'u',0);

Создайте коэффициенты УЧП

УЧП

В синтаксисе Partial Differential Equation Toolbox™,

Это уравнение переводит в коэффициенты c = 1, a = 0, и f = 2. Введите коэффициенты.

c = 1;
a = 0;
f = 2;
specifyCoefficients(model,'m',0,'d',0,'c',c','a',a,'f',f);

Создайте Mesh и решите УЧП

Инициализируйте mesh.

generateMesh(model);

Решите УЧП.

results = solvepde(model);

Исследуйте решение в графике среза контура

Создайте сетку и интерполируйте решение сетки.

[X,Y,Z] = meshgrid(0:135,0:35,0:61);
V = interpolateSolution(results,X,Y,Z);
V = reshape(V,size(X));

Создайте график среза контура для пяти фиксированных значений y- координата.

figure
colormap jet
contourslice(X,Y,Z,V,[],0:4:16,[])
xlabel('x')
ylabel('y')
zlabel('z')
xlim([0,100])
ylim([0,20])
zlim([0,50])
axis equal
view(-50,22)
colorbar

Частица имеет самое большое среднее выходное время около точки $$(x,y,z)=(100,0,0)$$.

Используйте градиенты для дрожи и оптимизируйте графики

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

[X,Y,Z] = meshgrid(1:9:99,1:3:20,1:6:50);
[gradx,grady,gradz] = evaluateGradient(results,X,Y,Z);

Постройте векторы градиента. Сначала измените аппроксимированные градиенты к форме mesh.

gradx = reshape(gradx,size(X));
grady = reshape(grady,size(Y));
gradz = reshape(gradz,size(Z));

figure
quiver3(X,Y,Z,gradx,grady,gradz)
axis equal
xlabel 'x'
ylabel 'y'
zlabel 'z'
title('Quiver Plot of Estimated Gradient of Solution')

Постройте потоки аппроксимированного градиента. Запустите потоки с более разреженного набора начальных точек.

hold on
[sx,sy,sz] = meshgrid([1,46],1:6:20,1:12:50);
streamline(X,Y,Z,gradx,grady,gradz,sx,sy,sz)
title('Quiver Plot with Streamlines')
hold off

Потоки показывают что маленькие значения y и z дайте большие средние выходные времена. Они также показывают что x- координата оказывает значительное влияние на u когда x мал, но когда x больше 40, большие значения оказывают мало влияния на u. Точно так же, когда z меньше 20, его значения оказывают мало влияния на u.