Этот пример показывает, как решить уравнение Пуассона с помощью предварительно обусловленного итерационного решателя и распределенных массивов. При помощи распределенных массивов можно масштабировать вычисление с помощью памяти кластера машин, а не только памяти одной машины. Вы можете масштабировать, не меняя свой код.
Этот пример продолжает темы, рассмотренные в «Использование распределенных массивов для решения систем линейных уравнений с итерационными методами». Основываясь на [1], пример моделирует распределение тепла в комнате с помощью уравнения Пуассона в форме, известной как однородное установившееся тепловое уравнение. Устойчивое состояние означает, что тепло не изменяется со временем и однородное означает, что нет внешнего источника тепла.
В уравнении, представляет температуру в каждой точке в комнате. Чтобы решить уравнение, вы сначала аппроксимируете его системой линейных уравнений с помощью метода дискретизации конечного различия. Затем вы используете предварительно обусловленные сопряженные градиенты (pcg) метод для решения системы. Предварительное кондиционирование преобразует задачу, чтобы улучшить эффективность числового решателя. При помощи распределенных массивов можно использовать объединенную память кластера машин и допускать более мелкие дискретизации.
Узнать, как:
Настройка дискретных 3-D сетки и граничных условий.
Задайте дискретизацию и мультигрид-прекондиционер.
Примените предварительно обусловленный численный решатель, чтобы решить уравнение тепла через объем 3-D.
В этом примере кубик боковых 1 моделирует комнату. Первым шагом является дискретизация его с помощью 3-D сетки.
Метод предварительного обусловления в этом примере использует несколько сеток с различными уровнями зернистости. Каждый уровень сокращает сетку в множитель 2 в каждой размерности. Определите количество уровней multigrid.
multigridLevels = 2;
Определите число точек в каждой размерности, X, Y, и Z, из лучших сетки. Метод предварительного обусловления требует, чтобы число точек в этой сетке было делено на 2^multigridLevels. В этом случае число точек должно быть делено на 4, потому что количество уровней multigrid 2.
numPoints.X = 32; numPoints.Y = 32; numPoints.Z = 32;
Дискретизируйте пространственные размерности с помощью 3-D сетки при помощи meshgrid функция. Разделите каждую размерность равномерно в соответствии с числом точек при помощи linspace. Обратите внимание, что для включения контуров куба необходимо добавить две дополнительные точки.
[X,Y,Z] = meshgrid(linspace(0,1,numPoints.X+2), ... linspace(0,1,numPoints.Y+2), ... linspace(0,1,numPoints.Z+2));
Предположим, в комнате есть окно и дверь. Стенки и потолок имеют постоянную температуру 0 степеней, окно имеет постоянную температуру 16 степеней, а дверь имеет постоянную температуру 15 степеней. Пол на 0,5 степени. Цель состоит в том, чтобы определить распределение температуры по внутренней части помещения.
Определите координаты пола, окна и двери с помощью реляционных операторов и определите температуру на этих краевых элементах. Контуры являются гранями куба и, следовательно, одной из X, Y, или Z должен быть 0 или 1. Установите остальную часть контура и внутреннюю часть куба равной 0.
floor = (0.0 <= X & X <= 1.0) & (0.0 <= Y & Y <= 1) & (Z == 0.0); window = (X == 1) & (0.2 <= Y & Y <= 0.8) & (0.4 <= Z & Z <= 0.6); door = (0.4 <= X & X <= 0.6) & (Y == 1.0) & (0.0 <= Z & Z <= 0.6); u = zeros(size(X)); u(floor) = 0.5; u(window) = 16; u(door) = 15;
Эти граничные условия задают постоянные значения, которые решение должно принимать вдоль контура области. Этот тип граничного условия известен как граничное условие Дирихле.
Визуализируйте граничные условия с помощью slice функция. Используйте срезы, расположенные на контуре куба, которые показывают ненулевые граничные условия.
xSlices = 1; ySlices = 1; zSlices = 0; f = slice(X,Y,Z,u,xSlices,ySlices,zSlices,'nearest'); title('Constant nonzero boundary conditions'), xlabel('x'), ylabel('y'), zlabel('z'); colorbar, colormap cool; shading interp; set(f,'EdgeColor',[0 0 0]);
Этот пример дискретизирует дифференциальное уравнение в линейную систему с помощью метода приближения конечных различий и использует мультигрид-предкондиционер, чтобы улучшить эффективность итерационного решателя. В данном примере используйте дискретизацию и предварительный кондиционер в вспомогательных функциях discretizePoissonEquation и multigridPreconditioner. Для других задач выберите дискретизацию и предварительное условие, которые подходят для вашего приложения.
В этом примере discretizePoissonEquation дискретизирует уравнение Пуассона методом семиточечных конечных различий в несколько сеток с различными уровнями гранулярности. Функция создает мультигридную структуру дискретизаций с предварительно вычисленными треугольными факторизациями и операторами, которые отображают между грубым и мелким уровнями. Предварительный кондиционер использует эту мультигридную информацию, чтобы эффективно аппроксимировать решение.
Среди других методов, этот предварительный кондиционер применяет сглаживание, чтобы минимизировать ошибки с серией приближений. Определите количество шагов сглаживания. Использование большего количества шагов делает приближения более точными, но и более интенсивными в вычислительном отношении. Затем дискретизируйте дифференциальное уравнение и настройте предварительный кондиционер.
numberOfSmootherSteps = 1; [A,b,multigridData] = discretizePoissonEquation(numPoints,multigridLevels,numberOfSmootherSteps,u);
Level 0: The problem is of dimension 32768 with 223232 nonzeros. Level 1: The problem is of dimension 4096 with 27136 nonzeros. Level 2: The problem is of dimension 512 with 3200 nonzeros.
preconditioner = setupPreconditioner(multigridData);
Решить линейную систему можно используя предварительно обусловленные сопряженные градиенты.
tol = 1e-12; maxit = numel(A); pcg(A,b,tol,maxit,preconditioner);
pcg converged at iteration 45 to a solution with relative residual 5.4e-13.
Если вам нужно больше вычислительных ресурсов, таких как память, можно масштабировать с помощью распределенных массивов, не меняя код. Распределенные массивы распределяют ваши данные между несколькими работниками, и они могут использовать вычислительную эффективность и память кластера машин.
Создайте пул параллельных рабочих. По умолчанию parpool использует ваш кластер по умолчанию. Проверьте профиль кластера по умолчанию на вкладке MATLAB Home, в области окружение, в Parallel > Select a Default Cluster.
parpool;
Starting parallel pool (parpool) using the 'MyCluster' profile ... Connected to the parallel pool (number of workers: 12).
Распределите температурную переменную u по памяти работников вашего кластера при помощи distributed функция.
distU = distributed(u);
Можно использовать тот же код, что и раньше; изменения не требуются, потому что функции дискретизации и предварительной подготовки создают распределенные массивы, если вход является распределенным массивом. Многие функции MATLAB усилены для распределенных массивов, поэтому с ними можно работать так же, как с массивами в памяти.
Обратите внимание, что discretizePoissonEquation возвращает структуру, содержащую распределенные данные. Чтобы использовать распределенные данные внутри структуры распределенным способом, необходимо создать структуру внутри spmd блок. Вы также должны вызвать любую функцию, которая использует ее внутри spmd блок.
spmd [A,b,multigridData] = discretizePoissonEquation(numPoints,multigridLevels,numberOfSmootherSteps,distU); preconditioner = setupPreconditioner(multigridData); end
Analyzing and transferring files to the workers ...done. Lab 1: Level 0: The problem is of dimension 32768 with 223232 nonzeros. Level 1: The problem is of dimension 4096 with 27136 nonzeros. Level 2: The problem is of dimension 512 with 3200 nonzeros.
Использование pcg внутри spmd блок для решения линейной системы распределенным способом.
spmd x = pcg(A,b,tol,maxit,preconditioner); end
Lab 1: pcg converged at iteration 45 to a solution with relative residual 5.4e-13.
Решением от решателя является вектор, который помещается в памяти. Отправка данных от работников клиенту при помощи gather. Измените форму данных обратно в трехмерный массив и переупорядочите размерности, чтобы получить окончательное решение. Установите внутреннюю часть u к этому решению. Внешняя часть, контур, уже содержит значение граничных условий.
x3D = reshape(gather(x),numPoints.X,numPoints.Y,numPoints.Z); u(2:end-1,2:end-1,2:end-1) = permute(x3D, [2, 1, 3]);
Визуализируйте решение с помощью slice функция. Добавьте дополнительные срезы для построения графика температуры внутри куба. Можно использовать Rotate инструмент, или изменить положение срезов, чтобы исследовать решение.
xSlices = [.5,1]; ySlices = [.5,1]; zSlices = [0,.5]; f = slice(X,Y,Z,u,xSlices,ySlices,zSlices,'nearest'); title('Heat distribution'), xlabel('x'), ylabel('y'), zlabel('z'); colorbar, colormap cool; shading interp;
Можно попробовать различные значения numPoints в этом примере для проверки различных уровней дискретизации. Использование большего значения увеличивает разрешение, но требует большей памяти. В сложение, большая multigridLevels Таким образом, чем эффективнее память, тем выше предварительное условие. Однако большая multigridLevels подразумевает менее точный предварительный кондиционер, так как грубость снижает точность на каждом уровне. В результате решателю может потребоваться больше итераций, чтобы достичь того же уровня точности.
Задайте мультигридный предварительный кондиционер для использования с предварительно обусловленным методом сопряженных градиентов. Этот тип предкондиционера использует несколько дискретизационных сеток с различными уровнями гранулярности, чтобы более эффективно аппроксимировать решение системы линейных уравнений. Метод предварительного кондиционирования в этом примере основан на [2] и следует этим основным этапам:
Пресмут с использованием метода приближения Гаусса-Зайделя.
Вычислите остаточное решение на более грубом уровне.
Рекурсивно предусловие на более грубом уровне или решите непосредственно, если на самом грубом уровне.
Обновите решение с помощью более грубой сетки.
Постсмут с использованием метода приближения Гаусса-Зайделя.
function x = multigridPreconditioner(mgData,r,level) if(level < mgData(level).MaxLevel) x = zeros(size(r),'like',r); % Presmooth using Gauss-Seidel for i=1:mgData(level).NumberOfSmootherSteps x = mgData(level).Matrices.LD \ (-mgData(level).Matrices.U*x + r); x = mgData(level).Matrices.DU \ (-mgData(level).Matrices.L*x + r); end % Compute residual on a coarser level Axf = mgData(level).Matrices.A*x; rc = r(mgData(level).Fine2CoarseOperator)- Axf(mgData(level).Fine2CoarseOperator); % Recursive call until coarsest level is reached xc = multigridPreconditioner(mgData, rc, level+1); % Update solution with prolonged coarse grid solution x(mgData(level).Fine2CoarseOperator) = x(mgData(level).Fine2CoarseOperator)+xc; % Postsmooth using Gauss-Seidel for i = 1:mgData(level).NumberOfSmootherSteps x = mgData(level).Matrices.LD \ (-mgData(level).Matrices.U*x + r); x = mgData(level).Matrices.DU \ (-mgData(level).Matrices.L*x + r); end else % Obtain exact solution on the coarsest level x = mgData(level).Matrices.A \ r; end end
Создайте функцию, которая принимает мультигридные данные и возвращает указатель на функцию, который применяет предварительный кондиционер к входным данным. В этом примере этот указатель на функцию является входом предкондиционера в pcg. Вы должны создать эту функцию, потому что невозможно задать анонимные функции внутри spmd блоки.
function preconditioner = setupPreconditioner(multigridData) if ~isempty(multigridData) preconditioner = @(x,varargin) multigridPreconditioner(multigridData,x,1); else preconditioner = []; end end
[1] Донгарра, Дж., М. А. Херу и П. Лющек. HPCG Benchmark: A New Metric for Ranking High Performance Computing Systems (неопр.) (недоступная ссылка). Ноксвилл, TN: Университет Теннесси, 2015.
[2] Elman, H. C., D. J. Silvester, and A. J. Wathen. Конечные элементы и быстрые итеративные решатели: с приложениями в несжимаемой динамике жидкости. Оксфорд, Великобритания: Oxford University Press, 2005, Section 2.5.
distributed | gather | pcg | spmd
