ilu

불완전 LU 분해

구문

[L,U] = ilu(A)

[L,U,P] = ilu(A)

W = ilu(A)

[___] = ilu(A,options)

설명

[L,U] = ilu(A)는 0 채우기로 희소 행렬 A의 불완전 LU 분해를 수행하고, 하부 삼각 행렬 L과 상부 삼각 행렬 U를 반환합니다.

[L,U,P] = ilu(A)는 치환 행렬 P도 반환합니다. 여기서 L과 U는 P*A 또는 A*P의 불완전 인수입니다. 기본적으로 P는 피벗 연산을 수행하지 않은 불완전 LU 분해에 대한 단위 행렬입니다.

W = ilu(A)는 LU 인수의 0이 아닌 요소를 반환합니다. 출력값 W는 L + U - speye(size(A))와 같습니다.

[___] = ilu(A,options)는 구조체 options로 지정된 옵션으로 A에 대한 불완전 LU 분해를 수행합니다.

예를 들어, options의 type 필드를 "ilutp"로 설정하여 피벗 연산으로 불완전 LU 분해를 수행할 수 있습니다. 그런 다음 milu 필드를 "row" 또는 "col"로 설정하여 행의 합계 또는 열의 합계를 유지하는 수정된 불완전 LU 분해를 지정할 수 있습니다. 이 옵션 조합은 치환 행렬 P를 반환합니다. 여기서 "row" 옵션의 경우 L과 U는 A*P의 불완전 인수이거나(이때 U는 열 치환됨), "col" 옵션의 경우 L과 U는 P*A의 불완전 인수입니다(이때 L은 행 치환됨).

예제

다양한 유형의 불완전 LU 분해

ilu 함수는 3가지 유형의 불완전 LU 분해, 즉 0 채우기 분해(ILU(0)), 크라우트 버전(ILUC) 및 임계값 기각과 피벗 연산을 사용한 분해(ILUTP)를 제공합니다.

기본적으로 ilu는 희소 행렬 입력값에 대해 0 채우기 불완전 LU 분해를 수행합니다. 예를 들어, 7840개의 0이 아닌 요소를 갖는 희소 행렬에 대한 완전 및 불완전 행렬 분해를 구합니다. 완전 LU 인수는 126,478개의 0이 아닌 요소를 가지며, 0 채우기가 있는 불완전 LU 인수는 A와 같은 개수인 7840개의 0이 아닌 요소를 가집니다.

A = gallery("neumann",1600) + speye(1600);
n = nnz(A)

n = 7840

n = nnz(lu(A))

n = 126478

n = nnz(ilu(A))

n = 7840

0 채우기 분해는 LU 인수에 입력 행렬의 희소성 패턴을 유지하므로, 분해의 상대 오차는 A의 0이 아닌 요소의 패턴에 대해 본질적으로 0입니다.

[L,U] = ilu(A);
e = norm(A-(L*U).*spones(A),"fro")/norm(A,"fro")

e = 4.8874e-17

그러나 이러한 0 채우기 인수의 곱은 원래 행렬의 적절한 근사가 아닙니다.

e = norm(A-L*U,"fro")/norm(A,"fro")

e = 0.0601

정확도를 높이기 위해 필인(fill-in)이 있는 다른 유형의 불완전 LU 분해를 사용할 수 있습니다. 예를 들어, 1e-6 기각 허용오차와 함께 크라우트 버전을 사용합니다.

options.droptol = 1e-6;
options.type = "crout";
[L,U] = ilu(A,options);

크라우트 버전의 불완전 행렬 분해는 LU 인수에 (완전 행렬 분해보다 적은) 51,482개의 0이 아닌 요소를 갖습니다. 필인을 사용하면 불완전 LU 인수의 곱이 원래 행렬의 더 적절한 근사가 됩니다.

n = nnz(ilu(A,options))

n = 51482

e = norm(A-L*U,"fro")./norm(A,"fro")

e = 9.3040e-07

비교를 위해, 동일한 입력 행렬 A에 대해 임계값 기각과 피벗 연산을 사용한 불완전 행렬 분해를 수행할 경우 크라우트 버전과 유사한 결과가 산출됩니다.

options.type = "ilutp";
[L,U,P] = ilu(A,options);
n = nnz(ilu(A,options))

n = 51541

norm(P*A-L*U,"fro")./norm(A,"fro")

ans = 9.4960e-07

불완전 LU 분해의 기각 허용오차

희소 행렬을 인수 분해하도록 불완전 LU 분해의 기각 허용오차를 변경합니다.

실수 값의 479×479 비대칭 희소 행렬인 west0479 행렬을 불러옵니다. condest를 사용하여 행렬의 조건수를 추정합니다.

load west0479
A = west0479;
c1 = condest(A)

c1 = 1.4244e+12

equilibrate를 사용하여 행렬의 조건수를 개선합니다. 원래 행렬 A를 치환하고 다시 스케일링하여 새 행렬 B = R*P*A*C를 만듭니다. 이 행렬은 더 나은 조건수를 가지며 대각선 요소가 1과 –1로만 구성되어 있습니다.

[P,R,C] = equilibrate(A);
B = R*P*A*C;
c2 = condest(B)

c2 = 5.1036e+04

B에 대해 행의 합계를 유지하는 임계값 기각과 피벗 연산을 사용한 불완전 LU 분해를 수행하는 옵션을 지정합니다. 비교를 위해, 먼저 필인의 기각 허용오차를 0으로 지정합니다. 이 경우 완전 LU 분해가 생성됩니다.

options.type = "ilutp";
options.milu = "row";
options.droptol = 0;
[L,U,P] = ilu(B,options);

이 분해는 입력 행렬 B를 매우 정확하게 근사하지만, 인수는 B보다 훨씬 더 조밀합니다.

e = norm(B*P-L*U,"fro")/norm(B,"fro")

e = 1.0345e-16

nLU = nnz(L)+nnz(U)-size(B,1)

nLU = 19118

nB = nnz(B)

nB = 1887

기각 허용오차를 변경하여 더 희소하지만 B를 근사하는 데 덜 정확한 불완전 LU 인수를 얻을 수 있습니다. 예를 들어, 다음 플롯은 기각 허용오차에 대해 플로팅된 입력 행렬의 밀도에 대한 불완전 인수의 밀도의 비율과 불완전 행렬 분해의 상대 오차를 보여줍니다.

ntols = 20;
tau = logspace(-6,-2,ntols);
e = zeros(1,ntols);
nLU = zeros(1,ntols);
for k = 1:ntols
    options.droptol = tau(k);
    [L,U,P] = ilu(B,options);
    nLU(k) = nnz(L)+nnz(U)-size(B,1);
    e(k) = norm(B*P-L*U,"fro")/norm(B,"fro");
end
figure
semilogx(tau,nLU./nB,LineWidth=2)
title("Ratio of nonzeros of LU factors with respect to B")
xlabel("drop tolerance")
ylabel("nnz(L)+nnz(U)-size(B,1)/nnz(B)",FontName="FixedWidth")

figure
loglog(tau,e,LineWidth=2)
title("Relative error of the incomplete LU factorization")
xlabel("drop tolerance")
ylabel("$||BP-LU||_F\,/\,||B||_F$",Interpreter="latex")

이 예제에서 임계값 기각을 사용한 불완전 LU 분해의 상대 오차는 기각 허용오차와 동일한 차수(항상 동일하지 않음)입니다.

ilu를 선조건자로 사용하여 선형 시스템 풀기

불완전 LU 분해를 bicgstab에 대한 선조건자로 사용하여 선형 시스템을 풉니다.

479×479 실수 비대칭 희소 행렬인 west0479를 불러옵니다.

load west0479
A = west0479;

$\mathit{Ax}=\mathit{b}$의 참(True) 해가 1로만 구성된 벡터가 되도록 b를 정의합니다.

b = full(sum(A,2));

허용오차와 최대 반복 횟수를 설정합니다.

tol = 1e-12;
maxit = 20;

bicgstab를 사용하여 요청된 허용오차와 반복 횟수로 해를 구합니다. 풀이 과정에 대한 정보를 반환하는 5개의 출력값을 지정합니다.

• x는 A*x = b의 계산된 해입니다.

• fl0은 알고리즘의 수렴 여부를 나타내는 플래그입니다.

• rr0은 계산된 답 x의 상대 잔차입니다.

• it0은 x가 계산된 반복 횟수입니다.

• rv0은 $\Vert \mathit{b}-\mathit{Ax}\Vert$에 대한 각 절반의 반복에서 잔차 내역으로 구성된 벡터입니다.

[x,fl0,rr0,it0,rv0] = bicgstab(A,b,tol,maxit); 
fl0

fl0 = 1

rr0

rr0 = 1

it0

it0 = 0

bicgstab가 요청된 반복 횟수인 20회 이내에 요청된 허용오차 1e-12로 수렴하지 않으므로 fl0은 1이 됩니다. 사실상 bicgstab의 동작이 좋지 않기 때문에 초기 추측값 x0 = zeros(size(A,2),1)이 최적의 해가 되고 it0 = 0에서 볼 수 있듯이 이 해가 반환됩니다.

이러한 느린 수렴을 개선할 수 있도록 선조건자 행렬을 지정할 수 있습니다. A는 비대칭 행렬이므로 ilu를 사용하여 선조건자 $\mathit{M}=\mathit{L}\mathit{U}$를 생성합니다. 대각선상에 있지 않으면서 값이 1e-6보다 작은 요소는 무시하도록 기각 허용오차를 지정합니다. L과 U를 bicgstab에 대한 입력값으로 지정하여 선조건이 적용된 시스템 ${\mathit{A}\mathit{M}}^{-1}\mathit{M}\mathit{x}=\mathit{b}$를 풉니다.

options = struct("type","ilutp","droptol",1e-6);
[L,U] = ilu(A,options);
[x_precond,fl1,rr1,it1,rv1] = bicgstab(A,b,tol,maxit,L,U);
fl1

fl1 = 0

rr1

rr1 = 5.9892e-14

it1

it1 = 3

ilu 선조건자를 사용한 결과 3번째 반복에서 요청된 허용오차 1e-12보다 작은 상대 잔차 rr1이 생성됩니다. 출력값 rv1(1)은 norm(b)가 되고, 출력값 rv1(end)는 norm(b-A*x1)이 됩니다.

각 절반의 반복에서 상대 잔차를 플로팅하면 bicgstab의 진행률을 추적할 수 있습니다. 각 해의 잔차 내역을 플로팅하고, 지정된 허용오차에 대한 선을 추가합니다.

semilogy((0:numel(rv0)-1)/2,rv0/norm(b),"-o")
hold on
semilogy((0:numel(rv1)-1)/2,rv1/norm(b),"-o")
yline(tol,"r--");
legend("No preconditioner","ILU preconditioner","Tolerance",Location="East")
xlabel("Iteration number")
ylabel("Relative residual")

입력 인수

A — 입력 행렬

희소 정사각 행렬

입력 행렬로, 희소 정사각 행렬로 지정됩니다.

options — LU 분해 옵션

구조체

LU 분해 옵션으로, 다음과 같은 최대 5개의 필드를 가진 구조체로 지정됩니다.

필드 이름	요약	설명
type	불완전 LU 분해의 유형	type에 유효한 값은 다음과 같습니다. "nofill"(디폴트 값) — ILU(0)으로 알려진 0 레벨 필인의 ILU 분해를 수행합니다. type이 "nofill"로 설정된 상태에서 milu 옵션만 사용됩니다. 다른 모든 필드는 무시됩니다. "crout" — ILUC로 알려진 ILU 분해의 크라우트 버전을 수행합니다. type이 "crout"로 설정되면 droptol 옵션과 milu 옵션만 사용됩니다. 다른 모든 필드는 무시됩니다. "ilutp" — 임계값 및 피벗 연산을 사용한 ILU 분해(ILUTP라고 함)를 수행합니다. 피벗 연산은 이 유형에서만 수행됩니다. 디폴트 값은 "nofill"입니다.
droptol	불완전 LU 분해의 기각 허용오차	droptol에 유효한 값은 음이 아닌 스칼라입니다. U의 0이 아닌 요소는 abs(U(i,j)) >= droptol*norm(A(:,j))를 충족합니다. 이 기준의 충족 여부와 관계없이 유지되는 대각선 요소는 예외입니다. L의 요소는 피벗에 의해 스케일링되기 전에 국소 기각 허용오차에 대해 테스트됩니다. 따라서 L의 0이 아닌 요소는 abs(L(i,j)) >= droptol*norm(A(:,j))/U(j,j)를 충족합니다. 디폴트 값은 완전 LU 분해를 생성하는 0입니다.
milu	수정된 불완전 LU 분해의 유형	milu에 유효한 값은 다음과 같습니다. "off"(디폴트 값) — 수정된 불완전 LU 분해가 생성되지 않습니다. "row" — 행의 합계에 대해 수정된 불완전 LU 분해를 생성합니다. 원래 행렬의 행 합계가 유지되도록 상부 삼각 인수 U의 대각선 요소가 보정됩니다. 즉, A*e = L*U*e입니다. 여기서 e는 1로 구성된 열 벡터입니다. "col" — 열의 합계에 대해 수정된 불완전 LU 분해를 생성합니다. 원래 행렬의 열 합계가 유지되도록 상부 삼각 인수 U의 대각선 요소가 보상됩니다. 즉, e'*A = e'*L*U입니다. 디폴트 값은 "off"입니다.
udiag	U의 0 대각선 요소를 바꿀지 여부	udiag에 유효한 값은 0과 1입니다. udiag가 1인 경우, 상부 삼각 인수 U의 모든 0 대각선 요소는 특이값 인수를 피하기 위해 국소 기각 허용오차에 의해 대체됩니다. 디폴트 값은 0입니다.
thresh	피벗 임계값	유효한 값은 0(이때 알고리즘은 대각선 피벗 선택)과 1(이때 알고리즘은 피벗 연산이 수행될 열의 최대 크기 항목 선택) 사이입니다. 피벗 연산은 열에 있는 대각선 요소의 크기가 그 열에 있는 임의의 하부대각선 요소의 크기에 thresh를 곱한 값보다 작은 경우 발생합니다. 디폴트 값은 1입니다.

출력 인수

L — 하부 삼각 인수

희소 정사각 행렬

하부 삼각 인수로, 희소 정사각 행렬로 반환됩니다.

• options.type가 "nofill"(디폴트 값) 또는 "crout"로 지정된 경우 L은 단위 하부 삼각 행렬(즉, 주대각선 요소가 1인 하부 삼각 행렬)로 반환됩니다.

• options.type이 "ilutp"로 지정되고 options.milu가 "col"로 지정된 경우 L은 행 치환된 단위 하부 삼각 행렬로 반환됩니다.

U — 상부 삼각 인수

희소 정사각 행렬

상부 삼각 인수로, 희소 정사각 행렬로 반환됩니다.

• options.type이 "nofill"(디폴트 값) 또는 "crout"로 지정된 경우 U는 상부 삼각 행렬로 반환됩니다.

• options.type이 "ilutp"로 지정되고 options.milu가 "row"로 지정된 경우 U는 열 치환된 상부 삼각 행렬로 반환됩니다.

P — 치환 행렬

희소 정사각 행렬

치환 행렬로, 희소 정사각 행렬로 반환됩니다.

options.type이 "nofill"(디폴트 값) 또는 "crout"로 지정된 경우 P는 A와 동일한 크기의 단위 행렬로 반환됩니다(이들 모두 피벗 연산을 수행하지 않기 때문).

W — LU 인수의 0이 아닌 요소

희소 정사각 행렬

LU 인수의 0이 아닌 요소로, 희소 정사각 행렬로 반환됩니다(여기서 W = L + U - speye(size(A))).

팁

• 이 함수에 의해 반환되는 불완전 행렬 분해는 bicg, bicgstab, gmres와 같은 반복법을 사용하여 풀고 있는 선형 연립방정식을 위한 선조건자로 유용할 수 있습니다.

참고 문헌

[1] Saad, Y. "Preconditioning Techniques", chap. 10 in Iterative Methods for Sparse Linear Systems. The PWS Series in Computer Science. Boston: PWS Pub. Co, 1996.

확장 기능

스레드 기반 환경

MATLAB®의 backgroundPool을 사용해 백그라운드에서 코드를 실행하거나 Parallel Computing Toolbox™의 ThreadPool을 사용해 코드 실행 속도를 높일 수 있습니다.

이 함수는 스레드 기반 환경을 완전히 지원합니다. 자세한 내용은 스레드 기반 환경에서 MATLAB 함수 실행하기 항목을 참조하십시오.

분산 배열

Parallel Computing Toolbox™를 사용하여 대규모 배열을 클러스터의 결합된 메모리에 걸쳐 분할할 수 있습니다.

사용법 관련 참고 및 제한 사항:

• 구조체형 배열 options에 필드 type을 포함한 경우 이를 'nofill'로 설정해야 합니다.

자세한 내용은 분산 배열을 사용하여 MATLAB 함수 실행 (Parallel Computing Toolbox) 항목을 참조하십시오.

버전 내역

R2007a에 개발됨

참고 항목

lu | ichol | bicg | bicgstab | gmres

도움말 항목

• 희소 행렬 연산
• 선형 시스템을 위한 반복법

×

예제 열기

이 예제의 수정된 버전이 있습니다. 사용자가 편집한 내용을 반영하여 이 예제를 여시겠습니까?

아니요, 수정된 버전을 덮어씁니다. 예