Antworten (2)

Sam Chak
Sam Chak am 1 Sep. 2026 um 18:44
Bearbeitet: Sam Chak am 1 Sep. 2026 um 20:58

3 Stimmen

I know this question is more than a decade old. For anyone trying to implement the top solution in MATLAB, note that there is a slight syntax typo there. The proposed solution I\sqrtm(A) evaluates to just sqrtm(A), but I believe it was intended to be sqrtm(A)\I to achieve the inverse . Else, the literal proposed solution only works in the trivial case where , as shown in Case 1.
With that syntax correction in Case 2, the math holds true under the conjugate transpose . However, it is worth noting that this is not a unique solution. Because matrix multiplication is non-commutative, there are infinitely many valid matrices that satisfy this equation, as demonstrated in Case 3.
While Cholesky is the most common method taught in undergrad linear algebra class (for algorithmic efficiency), I believe there are several unexplored alternative methods to parameterize and generate the infinite solutions, especially if the matrix is close to being singular (See Case 4).
Case 1: when is the identity matrix
% Hermitian positive-definite matrix
A = eye(3)
A = 3×3
1 0 0 0 1 0 0 0 1
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
tof = ishermitian(A) % returns logical 1 if matrix A is legitimately Hermitian
tof = logical
1
% Identity matrix (a.k.a unit matrix)
I = eye(size(A))
I = 3×3
1 0 0 0 1 0 0 0 1
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
% Solution
X = I\sqrtm(A)
X = 3×3
1 0 0 0 1 0 0 0 1
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
% Test: returns 0 if the result is correct
(X')*A*X - I
ans = 3×3
0 0 0 0 0 0 0 0 0
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
Case 2: is any arbitrary Hermitian positive-definite matrix (non-Identity) with
% Hermitian positive-definite matrix
A = [2.5, 2.5, 0.5
2.5, 5.0, 1.5
0.5, 1.5, 2.0];
eig(A)
ans = 3×1
0.7897 1.6913 7.0190
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
tof = ishermitian(A) % returns logical 1 if matrix A is legitimately Hermitian
tof = logical
1
% Identity matrix (a.k.a unit matrix)
I = eye(size(A))
I = 3×3
1 0 0 0 1 0 0 0 1
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
% Solution
X = sqrtm(A)\I % or inv(sqrtm(A))
X = 3×3
0.8565 -0.3013 0.0482 -0.3013 0.6198 -0.1782 0.0482 -0.1782 0.7954
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
% Test: returns 0 (or near double-precision machine epsilon) if the result is numerically accurate
(X')*A*X - I
ans = 3×3
1.0e-14 * -0.0444 0.0134 -0.0291 0.0302 -0.1998 0.0777 -0.0305 0.0749 -0.1332
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
Case 3: using the inverse Cholesky factorization trick
% Cholesky decomposition (A = L*L')
L = chol(A, 'lower'); % same matrix A from Case 2
X_baseline = inv(L');
% generate a random unitary orthogonal matrix Q via QR decomposition
[Q, ~] = qr(rand(3)); % Q'*Q = I
% general solution
X = X_baseline*Q
X = 3×3
-0.4461 0.5473 0.5727 -0.0389 -0.1285 -0.6990 -0.3323 -0.5520 0.5015
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
% Test: returns 0 (or near double-precision machine epsilon) if the result is numerically accurate
(X')*A*X - I
ans = 3×3
1.0e-15 * 0 0 -0.1665 0.1110 0.4441 0 -0.0278 0.1110 0.4441
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
Case 4: when is close to being singular
format longg
lambda_max = 5 - 1e-10; % Largest eigenvalue
lambda_mid = 2 - 1e-11; % Middle eigenvalue
lambda_min = 1e-9; % Smallest eigenvalue
% diagonal eigenvalue matrix (Sigma)
Sigma = diag([lambda_max, lambda_mid, lambda_min]);
% construct a real Hermitian positive-definite matrix using QR decomposition trick
[Q, ~] = qr(rand(3));
A = Q*Sigma*Q'
A = 3×3
3.21244753463304 0.624076812533866 1.91866199897675 0.624076812533866 2.30461855798916 1.23051812127085 1.91866199897675 1.23051812127085 1.4829339082678
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
cond(A)
ans =
4999998354.346
eig(A)
ans = 3×1
9.99999581276109e-10 1.99999999999 4.9999999999
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
tof = ishermitian(A) % returns logical 1 if matrix A is legitimately Hermitian
tof = logical
1
% Cholesky decomposition (A = L*L')
L = chol(A, 'lower'); % same matrix A from Case 2
X_baseline = inv(L');
% generate a random unitary orthogonal matrix Q via QR decomposition
[Q, ~] = qr(rand(3)); % Q'*Q = I
% general solution
X = X_baseline*Q
X = 3×3
6814.00029819071 6905.04323521035 -9809.79099892373 5138.72000419203 5206.83061141719 -7398.70435264685 -13080.89102809 -13254.1598291375 18831.2769607147
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
% Test: returns 0 (or near double-precision machine epsilon) if the result is numerically accurate
(X')*A*X - I
ans = 3×3
3.10737959807739e-08 3.14885255647823e-08 -4.47325874119997e-08 1.79334165295586e-08 1.81726136361249e-08 -2.58169166045263e-08 -8.56152837513946e-08 -8.67557901074179e-08 1.23254721984267e-07
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
Matt J
Matt J am 27 Mai 2014
Bearbeitet: Matt J am 1 Sep. 2026 um 18:49

0 Stimmen

X=sqrtm(A)\I;

Kategorien

Mehr zu Linear Algebra finden Sie in Hilfe-Center und File Exchange

Gefragt:

am 27 Mai 2014

Bearbeitet:

am 1 Sep. 2026 um 20:58

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by