Уравнение Пуассона с источником точек и адаптивным уточнением сетки

Этот пример показывает, как решить уравнение Пуассона с источником точек дельта-функции на единичном диске, используя adaptmesh функция.

В частности, решите уравнение Пуассона

на единичном диске с нулевыми граничными условиями Дирихле. Точное решение, выраженное в полярных координатах,

который сингулярен в источник.

Используя адаптивное уточнение сетки, Partial Equation Toolbox™ может точно найти решение везде, вдали от источника.

Следующие переменные определяют задачу:

c = 1;
a = 0;
f = @circlef;

Создайте модель УЧП с одной зависимой переменной.

numberOfPDE = 1;
model = createpde(numberOfPDE);

Создайте геометрию и включите ее в модель.

g = @circleg;
geometryFromEdges(model,g);

Постройте график геометрии и отобразите метки ребра.

figure; 
pdegplot(model,'EdgeLabels','on'); 
axis equal
title 'Geometry With Edge Labels Displayed';

Задайте нулевое решение на всех четырех внешних краях окружности.

applyBoundaryCondition(model,'dirichlet','Edge',(1:4),'u',0);

adaptmesh решает эллиптические PDE с помощью генерации адаптивной сетки. The tripick позволяет вам задать функцию, которая возвращает треугольники, которые будут уточнены в следующей итерации. circlepick возвращает треугольники с вычисленными оценками ошибок, больше заданный допуск. Допуск предоставляется circlepick использование параметра 'par'.

[u,p,e,t] = adaptmesh(g,model,c,a,f,'tripick','circlepick','maxt',2000,'par',1e-3);

Number of triangles: 258
Number of triangles: 515
Number of triangles: 747
Number of triangles: 1003
Number of triangles: 1243
Number of triangles: 1481
Number of triangles: 1705
Number of triangles: 1943
Number of triangles: 2155

Maximum number of triangles obtained.

Постройте график самого мелкого mesh.

figure; 
pdemesh(p,e,t); 
axis equal

Постройте график значений ошибок.

x = p(1,:)';
y = p(2,:)';
r = sqrt(x.^2+y.^2);
uu = -log(r)/2/pi;
figure;
pdeplot(p,e,t,'XYData',u-uu,'ZData',u-uu,'Mesh','off');

Постройте график решения КЭМ на самом мелком mesh.

figure;
pdeplot(p,e,t,'XYData',u,'ZData',u,'Mesh','off');