물리정보 신경망을 사용하여 단위 원판에서 푸아송 방정식 풀기

다음 제품이 필요합니다.

이 예제에서는 물리정보 신경망(PINN)을 사용하여 디리클레 경계 조건을 갖는 푸아송 방정식을 푸는 방법을 보여줍니다. PDE 모델 설정을 사용하여 PINN 훈련에 필요한 데이터를 생성합니다.

단위 원판에서 디리클레 경계 조건이 0인 푸아송 방정식은 $Ω$ 에서 $- Δ u = 1$ ( $δ Ω$ 에서는 $u = 0$ )로 작성할 수 있습니다. 여기서 $Ω$ 는 단위 원판입니다. 엄밀해는 $u (x, y) = (1 - x^{2} - y^{2}) / 4$ 입니다.

PDE 모델 설정

PDE 모델을 만들고 지오메트리를 포함시킵니다.

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

지오메트리를 플로팅하고 경계 조건 정의에 사용할 모서리 레이블을 표시합니다.

figure 
pdegplot(model,EdgeLabels="on"); 
axis equal

Figure contains an axes object. The axes object contains 5 objects of type line, text.

모든 모서리에 0 디리클레 경계 조건을 지정합니다.

applyBoundaryCondition(model,"dirichlet", ...
    Edge=1:model.Geometry.NumEdges,u=0);

계수로 구성된 구조체 배열을 만듭니다. PDE 모델에 대한 계수를 지정합니다.

pdeCoeffs.m = 0;
pdeCoeffs.d = 0;
pdeCoeffs.c = 1;
pdeCoeffs.a = 0;
pdeCoeffs.f = 1;
specifyCoefficients(model,m=pdeCoeffs.m,d=pdeCoeffs.d, ...
    c=pdeCoeffs.c,a=pdeCoeffs.a,f=pdeCoeffs.f);

경계에 많은 수의 노드가 있는 메시를 생성하고 플로팅합니다.

msh = generateMesh(model,Hmax=0.05,Hgrad=2, ...
    Hedge={1:model.Geometry.NumEdges,0.005});
figure 
pdemesh(model); 
axis equal

Figure contains an axes object. The axes object contains 2 objects of type line.

PINN 훈련을 위한 공간 데이터 생성하기

PINN을 훈련시키려면 도메인과 경계의 연결 점(collocation point)에서 손실을 모델링합니다. 이 예제에서 연결 점은 메시 노드입니다.

boundaryNodes = findNodes(msh,"region", ...
                          Edge=1:model.Geometry.NumEdges);
domainNodes = setdiff(1:size(msh.Nodes,2),boundaryNodes);
domainCollocationPoints = msh.Nodes(:,domainNodes)';

신경망 아키텍처 정의하기

각각 50개의 은닉 뉴런을 갖는 4개의 완전 연결 연산을 포함하는 다층 퍼셉트론 아키텍처를 정의합니다. 첫 번째 완전 연결 연산에는 입력값 x 및 y에 대응하는 두 개의 입력 채널이 있습니다. 마지막 완전 연결 연산에는 u(x,y)에 대응하는 출력이 하나 있습니다.

numNeurons = 50;
layers = [
    featureInputLayer(2,Name="featureinput")
    fullyConnectedLayer(numNeurons,Name="fc1")
    tanhLayer(Name="tanh_1")
    fullyConnectedLayer(numNeurons,Name="fc2")
    tanhLayer(Name="tanh_2")
    fullyConnectedLayer(numNeurons,Name="fc3")
    tanhLayer(Name="tanh_3")
    fullyConnectedLayer(1,Name="fc4")
    ];

pinn = dlnetwork(layers);

훈련 옵션 지정하기

Epoch 횟수, 미니 배치 크기, 초기 학습률, 학습률 감쇠를 지정합니다.

numEpochs = 50;
miniBatchSize = 500;
initialLearnRate = 0.01;
learnRateDecay = 0.005;

훈련 데이터를 dlarray (Deep Learning Toolbox) 객체로 변환합니다.

ds = arrayDatastore(domainCollocationPoints);
mbq = minibatchqueue(ds,MiniBatchSize=miniBatchSize, ...
                        MiniBatchFormat="BC", ...
                        OutputEnvironment="cpu");

평균 기울기와 제곱 평균 기울기를 초기화합니다.

averageGrad = [];
averageSqGrad = [];

훈련 진행 상황 모니터의 총 반복 횟수를 계산합니다.

numIterations = numEpochs* ...
    ceil(size(domainCollocationPoints,1)/miniBatchSize);

trainingProgressMonitor (Deep Learning Toolbox) 객체를 초기화합니다.

monitor = trainingProgressMonitor(Metrics="Loss", ...
                                  Info="Epoch", ...
                                  XLabel="Iteration");

Training Progess Window. The figure is set to track the loss with respect to the learnable parameters in pinn.

PINN 훈련시키기

사용자 지정 훈련 루프를 사용하여 모델을 훈련시킵니다. adamupdate (Deep Learning Toolbox) 함수를 사용하여 신경망 파라미터를 업데이트합니다. 각 반복의 끝에서 훈련 진행 상황을 표시합니다.

이 훈련 코드는 modelLoss 헬퍼 함수를 사용합니다. 자세한 내용은 모델 손실 함수 항목을 참조하십시오.

iteration = 0;
epoch = 0;
learningRate = initialLearnRate;
while epoch < numEpochs && ~monitor.Stop
    epoch = epoch + 1;
    reset(mbq);
    while hasdata(mbq) && ~monitor.Stop
        iteration = iteration + 1;
        XY = next(mbq);
        % Evaluate the model loss and gradients using dlfeval.
        [loss,gradients] = dlfeval(@modelLoss,model,pinn,XY,pdeCoeffs);
        % Update the network parameters using the adamupdate function.
        [pinn,averageGrad,averageSqGrad] = ...
            adamupdate(pinn,gradients,averageGrad, ...
                       averageSqGrad,iteration,learningRate);
    end
    % Update learning rate.
    learningRate = initialLearnRate / (1+learnRateDecay*iteration);
    % Update the training progress monitor.
    recordMetrics(monitor,iteration,Loss=loss);
    updateInfo(monitor,Epoch=epoch + " of " + numEpochs);
    monitor.Progress = 100 * iteration/numIterations;
end

PINN 테스트하기

메시 노드에서 참(True) 해를 구해 플로팅합니다.

trueSolution = @(msh) (1-msh.Nodes(1,:).^2-msh.Nodes(2,:).^2)/4;
Utrue = trueSolution(msh);

figure;
pdeplot(model,XYData=Utrue);
xlabel('$x$',interpreter='latex')
ylabel('$y$',interpreter='latex')
zlabel('$u(x,y)$',interpreter='latex')
title('True Solution: $u(x,y) = (1-x^2-y^2)/4$', ...
                      interpreter='latex')
axis equal

Figure contains an axes object. The axes object with title True Solution: u leftParenthesis x , y rightParenthesis equals leftParenthesis 1 minus x Squared baseline minus y Squared baseline rightParenthesis slash 4, xlabel $x$, ylabel $y$ contains an object of type patch.

이제 PINN이 예측한 해를 구해 플로팅합니다.

nodesDLarry = dlarray(msh.Nodes,"CB");
Upinn = gather(extractdata(predict(pinn,nodesDLarry)));

figure;
pdeplot(model,XYData=Upinn);
xlabel('$x$',interpreter='latex')
ylabel('$y$',interpreter='latex')
zlabel('$u(x,y)$',interpreter='latex')
title(sprintf(['PINN Predicted Solution: ' ...
    'L2 Error = %0.1e'],norm(Upinn-Utrue)))
axis equal

Figure contains an axes object. The axes object with title PINN Predicted Solution: L2 Error = 1.7e-01, xlabel $x$, ylabel $y$ contains an object of type patch.

모델 손실 함수

modelLoss 헬퍼 함수는 dlnetwork 객체 pinn과 입력 데이터의 미니 배치 XY를 입력값으로 받습니다. 이 함수는 손실과 pinn의 학습 가능한 파라미터에 대한 손실 기울기를 반환합니다. 기울기를 자동으로 계산하기 위해 dlgradient (Deep Learning Toolbox) 함수를 사용합니다. 입력값(x,y)이 주어졌을 때, 신경망의 출력값 u(x,y)가 푸아송 방정식과 경계 조건을 충족하도록 모델을 훈련시킵니다.

function [loss,gradients] = modelLoss(model,net,XY,pdeCoeffs)
dim = 2;
[U,~] = forward(net,XY);

% Compute gradients of U and Laplacian of U.
gradU = dlgradient(sum(U,"all"),XY,EnableHigherDerivatives=true);
Laplacian = 0;
for i=1:dim
    gradU2 = dlgradient(sum(pdeCoeffs.c.*gradU(i,:),"all"), ...
                          XY,EnableHigherDerivatives=true);
    Laplacian = gradU2(i,:)+Laplacian;
end

% Enforce PDE. Calculate lossF.
res = -pdeCoeffs.f - Laplacian + pdeCoeffs.a.*U;
zeroTarget = zeros(size(res),"like",res);
lossF = mse(res,zeroTarget);

% Enforce boundary conditions. Calculate lossU.
actualBC = []; % Boundary information
BC_XY = []; % Boundary coordinates
% Loop over the boundary edges and find boundary
% coordinates and actual BC assigned to PDE model.
numBoundaries = model.Geometry.NumEdges;
for i=1:numBoundaries
    BCi = findBoundaryConditions(model.BoundaryConditions,Edge=i);
    BCiNodes = findNodes(model.Mesh,"region",Edge=i);
    BC_XY = [BC_XY,model.Mesh.Nodes(:,BCiNodes)]; %#ok<AGROW> 
    actualBCi = ones(1,numel(BCiNodes))*BCi.u;
    actualBC = [actualBC actualBCi]; %#ok<AGROW>
end
BC_XY = dlarray(BC_XY,"CB"); % Format the coordinates
[predictedBC,~] = forward(net,BC_XY);
lossU = mse(predictedBC,actualBC);

% Combine the losses. Use a weighted linear combination 
% of loss (lossF and lossU) terms.
lambda = 0.4; % Weighting factor
loss = lambda*lossF + (1-lambda)*lossU;

% Calculate gradients with respect to the learnable parameters.
gradients = dlgradient(loss,net.Learnables); 
end

참고 항목

dlarray (Deep Learning Toolbox) | dlfeval (Deep Learning Toolbox) | dlgradient (Deep Learning Toolbox) | adamupdate (Deep Learning Toolbox)

도움말 항목

Solve PDE Using Fourier Neural Operator (Deep Learning Toolbox)
물리정보 신경망을 사용하여 PDE 해 구하기 (Deep Learning Toolbox)
Solve Heat Equation Using Graph Neural Network