Conjugate transpose: Difference between revisions

m
m (→‎{{header|Sidef}}: Fix link: Perl 6 --> Raku)
m (→‎{{header|Wren}}: Minor tidy)
 
(26 intermediate revisions by 15 users not shown)
Line 34:
* MathWorld entry: [http://mathworld.wolfram.com/UnitaryMatrix.html unitary matrix]
<br><br>
 
=={{header|11l}}==
{{trans|Nim}}
 
<syntaxhighlight lang="11l">-V eps = 1e-10
 
F to_str(m)
V r = ‘’
L(row) m
V i = L.index
r ‘’= I i == 0 {‘[’} E ‘ ’
L(val) row
V j = L.index
I j != 0
r ‘’= ‘ ’
r ‘’= ‘(#2.4, #2.4)’.format(val.real, val.imag)
r ‘’= I i == m.len - 1 {‘]’} E "\n"
R r
 
F conjugateTransposed(m)
V r = [[0i] * m.len] * m.len
L(i) 0 .< m.len
L(j) 0 .< m.len
r[j][i] = conjugate(m[i][j])
R r
 
F mmul(m1, m2)
V r = [[0i] * m1.len] * m1.len
L(i) 0 .< m1.len
L(j) 0 .< m1.len
L(k) 0 .< m1.len
r[i][j] += m1[i][k] * m2[k][j]
R r
 
F isHermitian(m)
L(i) 0 .< m.len
L(j) 0 .< m.len
I m[i][j] != conjugate(m[j][i])
R 0B
R 1B
 
F isEqual(m1, m2)
L(i) 0 .< m1.len
L(j) 0 .< m1.len
I m1[i][j] != m2[i][j]
R 0B
R 1B
 
F isNormal(m)
V h = conjugateTransposed(m)
R isEqual(mmul(m, h), mmul(h, m))
 
F isIdentity(m)
L(i) 0 .< m.len
L(j) 0 .< m.len
I i == j
I abs(m[i][j] - 1.0) > :eps
R 0B
E
I abs(m[i][j]) > :eps
R 0B
R 1B
 
F isUnitary(m)
V h = conjugateTransposed(m)
R isIdentity(mmul(m, h)) & isIdentity(mmul(h, m))
 
F test(m)
print(‘Matrix’)
print(‘------’)
print(to_str(m))
print(‘’)
print(‘Conjugate transposed’)
print(‘--------------------’)
print(to_str(conjugateTransposed(m)))
print(‘’)
print(‘Hermitian: ’(I isHermitian(m) {‘true’} E ‘false’))
print(‘Normal: ’(I isNormal(m) {‘true’} E ‘false’))
print(‘Unitary: ’(I isUnitary(m) {‘true’} E ‘false’))
 
V M2 = [[3.0 + 0.0i, 2.0 + 1.0i],
[2.0 - 1.0i, 1.0 + 0.0i]]
 
V M3 = [[1.0 + 0.0i, 1.0 + 0.0i, 0.0 + 0.0i],
[0.0 + 0.0i, 1.0 + 0.0i, 1.0 + 0.0i],
[1.0 + 0.0i, 0.0 + 0.0i, 1.0 + 0.0i]]
 
V SR2 = 1 / sqrt(2.0)
V SR2i = SR2 * 1i
V M4 = [[SR2 + 0.0i, SR2 + 0.0i, 0.0 + 0.0i],
[0.0 + SR2i, 0.0 - SR2i, 0.0 + 0.0i],
[0.0 + 0.0i, 0.0 + 0.0i, 0.0 + 1.0i]]
 
test(M2)
print("\n")
test(M3)
print("\n")
test(M4)</syntaxhighlight>
 
{{out}}
<pre>
Matrix
------
[( 3.0000, 0.0000) ( 2.0000, 1.0000)
( 2.0000, -1.0000) ( 1.0000, 0.0000)]
 
Conjugate transposed
--------------------
[( 3.0000, 0.0000) ( 2.0000, 1.0000)
( 2.0000, -1.0000) ( 1.0000, 0.0000)]
 
Hermitian: true
Normal: true
Unitary: false
 
 
Matrix
------
[( 1.0000, 0.0000) ( 1.0000, 0.0000) ( 0.0000, 0.0000)
( 0.0000, 0.0000) ( 1.0000, 0.0000) ( 1.0000, 0.0000)
( 1.0000, 0.0000) ( 0.0000, 0.0000) ( 1.0000, 0.0000)]
 
Conjugate transposed
--------------------
[( 1.0000, 0.0000) ( 0.0000, 0.0000) ( 1.0000, 0.0000)
( 1.0000, 0.0000) ( 1.0000, 0.0000) ( 0.0000, 0.0000)
( 0.0000, 0.0000) ( 1.0000, 0.0000) ( 1.0000, 0.0000)]
 
Hermitian: false
Normal: true
Unitary: false
 
 
Matrix
------
[( 0.7071, 0.0000) ( 0.7071, 0.0000) ( 0.0000, 0.0000)
( 0.0000, 0.7071) ( 0.0000, -0.7071) ( 0.0000, 0.0000)
( 0.0000, 0.0000) ( 0.0000, 0.0000) ( 0.0000, 1.0000)]
 
Conjugate transposed
--------------------
[( 0.7071, 0.0000) ( 0.0000, -0.7071) ( 0.0000, 0.0000)
( 0.7071, 0.0000) ( 0.0000, 0.7071) ( 0.0000, 0.0000)
( 0.0000, 0.0000) ( 0.0000, 0.0000) ( 0.0000, -1.0000)]
 
Hermitian: false
Normal: true
Unitary: true
</pre>
 
=={{header|Ada}}==
<langsyntaxhighlight Adalang="ada">with Ada.Text_IO; use Ada.Text_IO;
with Ada.Complex_Text_IO; use Ada.Complex_Text_IO;
with Ada.Numerics.Complex_Types; use Ada.Numerics.Complex_Types;
Line 82 ⟶ 231:
Put_Line("nmat:"); Examine(nmat); New_Line;
Put_Line("umat:"); Examine(umat);
end ConTrans;</langsyntaxhighlight>
{{out}}
<pre>hmat:
Line 117 ⟶ 266:
Normal?: TRUE
Unitary?: TRUE</pre>
 
=={{header|ALGOL 68}}==
Uses the same test cases as the Ada sample.
<syntaxhighlight lang="algol68">BEGIN # find and classify the complex conjugate transpose of a complex matrix #
# returns the conjugate transpose of m #
OP CONJUGATETRANSPOSE = ( [,]COMPL m )[,]COMPL:
BEGIN
[ 2 LWB m : 2 UPB m, 1 LWB m : 1 UPB m ]COMPL result;
FOR i FROM 1 LWB m TO 1 UPB m DO
FOR j FROM 2 LWB m TO 2 UPB m DO
result[ j, i ] := CONJ m[ i, j ]
OD
OD;
result
END # CONJUGATETRANSPOSE # ;
# returns TRUE if m is an identity matrix, FALSE otherwise #
OP ISIDENTITY = ( [,]COMPL m )BOOL:
IF 1 LWB m /= 2 LWB m OR 1 UPB m /= 2 UPB m THEN
# non-square matrix #
FALSE
ELSE
# the matrix is square #
# returns TRUE IF v - e is nearly 0, FALSE Otherwise #
PROC nearly equal = ( COMPL v, REAL e )BOOL: ABS re OF v - e < 1e-14 AND ABS im OF v < 1e-14;
BOOL result := TRUE;
FOR i FROM 1 LWB m TO 1 UPB m WHILE result DO
IF result := nearly equal( m[ i, i ], 1 ) THEN
# the diagonal element is 1 - test the non-diagonals #
FOR j FROM 1 LWB m TO 1 UPB m WHILE result DO
IF i /= j THEN result := nearly equal( m[ i, j ], 0 ) FI
OD
FI
OD;
result
FI # ISIDENTITY # ;
# returns m multiplied by n #
PRIO X = 7;
OP X = ( [,]COMPL m, n )[,]COMPL:
BEGIN
[ 1 : 1 UPB m, 1 : 2 UPB n ]COMPL r;
FOR i FROM 1 LWB m TO 1 UPB m DO
FOR j FROM 2 LWB n TO 2 UPB n DO
r[ i, j ] := 0;
FOR k TO 2 UPB n DO
r[ i, j ] +:= m[ i, k ] * n[ k, j ]
OD
OD
OD;
r
END # X # ;
# prints the complex matris m #
PROC show matrix = ( [,]COMPL m )VOID:
FOR i FROM 1 LWB m TO 1 UPB m DO
print( ( " " ) );
FOR j FROM 2 LWB m TO 2 UPB m DO
print( ( "( ", fixed( re OF m[ i, j ], -8, 4 )
, ", ", fixed( im OF m[ i, j ], -8, 4 )
, "i )"
)
)
OD;
print( ( newline ) )
OD # show matrix # ;
# display the matrix m, its conjugate transpose and whether it is Hermitian, Normal and Unitary #
PROC show = ( [,]COMPL m )VOID:
BEGIN
[,]COMPL c = CONJUGATETRANSPOSE m;
[,]COMPL cm = c X m;
[,]COMPL mc = m X c;
print( ( "Matrix:", newline ) );
show matrix( m );
print( ( "Conjugate Transpose:", newline ) );
show matrix( c );
BOOL is normal = cm = mc;
BOOL is unitary = IF NOT is normal THEN FALSE
ELSE ISIDENTITY mc
FI;
print( ( IF c = m THEN "" ELSE "not " FI, "Hermitian; "
, IF is normal THEN "" ELSE "not " FI, "Normal; "
, IF is unitary THEN "" ELSE "not " FI, "Unitary"
, newline
)
);
print( ( newline ) )
END # show # ;
# test some matrices for Hermitian, Normal and Unitary #
show( ( ( ( 3.0000 I 0.0000 ), ( 2.0000 I 1.0000 ) )
, ( ( 2.0000 I -1.0000 ), ( 1.0000 I 0.0000 ) )
)
);
show( ( ( ( 1.0000 I 0.0000 ), ( 1.0000 I 0.0000 ), ( 0.0000 I 0.0000 ) )
, ( ( 0.0000 I 0.0000 ), ( 1.0000 I 0.0000 ), ( 1.0000 I 0.0000 ) )
, ( ( 1.0000 I 0.0000 ), ( 0.0000 I 0.0000 ), ( 1.0000 I 0.0000 ) )
)
);
REAL rh = sqrt( 0.5 );
show( ( ( ( rh I 0.0000 ), ( rh I 0.0000 ), ( 0.0000 I 0.0000 ) )
, ( ( 0.0000 I rh ), ( 0.0000 I - rh ), ( 0.0000 I 0.0000 ) )
, ( ( 0.0000 I 0.0000 ), ( 0.0000 I 0.0000 ), ( 0.0000 I 1.0000 ) )
)
)
END</syntaxhighlight>
{{out}}
<pre>
Matrix:
( 3.0000, 0.0000i )( 2.0000, 1.0000i )
( 2.0000, -1.0000i )( 1.0000, 0.0000i )
Conjugate Transpose:
( 3.0000, 0.0000i )( 2.0000, 1.0000i )
( 2.0000, -1.0000i )( 1.0000, 0.0000i )
Hermitian; Normal; not Unitary
 
Matrix:
( 1.0000, 0.0000i )( 1.0000, 0.0000i )( 0.0000, 0.0000i )
( 0.0000, 0.0000i )( 1.0000, 0.0000i )( 1.0000, 0.0000i )
( 1.0000, 0.0000i )( 0.0000, 0.0000i )( 1.0000, 0.0000i )
Conjugate Transpose:
( 1.0000, 0.0000i )( 0.0000, 0.0000i )( 1.0000, 0.0000i )
( 1.0000, 0.0000i )( 1.0000, 0.0000i )( 0.0000, 0.0000i )
( 0.0000, 0.0000i )( 1.0000, 0.0000i )( 1.0000, 0.0000i )
not Hermitian; Normal; not Unitary
 
Matrix:
( 0.7071, 0.0000i )( 0.7071, 0.0000i )( 0.0000, 0.0000i )
( 0.0000, 0.7071i )( 0.0000, -0.7071i )( 0.0000, 0.0000i )
( 0.0000, 0.0000i )( 0.0000, 0.0000i )( 0.0000, 1.0000i )
Conjugate Transpose:
( 0.7071, 0.0000i )( 0.0000, -0.7071i )( 0.0000, 0.0000i )
( 0.7071, 0.0000i )( 0.0000, 0.7071i )( 0.0000, 0.0000i )
( 0.0000, 0.0000i )( 0.0000, 0.0000i )( 0.0000, -1.0000i )
not Hermitian; Normal; Unitary</pre>
 
=={{header|C}}==
<langsyntaxhighlight lang="c">/* Uses C99 specified complex.h, complex datatype has to be defined and operation provided if used on non-C99 compilers */
 
#include<stdlib.h>
Line 299 ⟶ 579:
 
return 0;
}</langsyntaxhighlight>
{{out}}
<pre>
Line 318 ⟶ 598:
 
Complex Matrix A is not normal</pre>
 
=={{header|C++}}==
<syntaxhighlight lang="cpp">#include <cassert>
#include <cmath>
#include <complex>
#include <iomanip>
#include <iostream>
#include <sstream>
#include <vector>
 
template <typename scalar_type> class complex_matrix {
public:
using element_type = std::complex<scalar_type>;
 
complex_matrix(size_t rows, size_t columns)
: rows_(rows), columns_(columns), elements_(rows * columns) {}
 
complex_matrix(size_t rows, size_t columns, element_type value)
: rows_(rows), columns_(columns), elements_(rows * columns, value) {}
 
complex_matrix(size_t rows, size_t columns,
const std::initializer_list<std::initializer_list<element_type>>& values)
: rows_(rows), columns_(columns), elements_(rows * columns) {
assert(values.size() <= rows_);
size_t i = 0;
for (const auto& row : values) {
assert(row.size() <= columns_);
std::copy(begin(row), end(row), &elements_[i]);
i += columns_;
}
}
 
size_t rows() const { return rows_; }
size_t columns() const { return columns_; }
 
const element_type& operator()(size_t row, size_t column) const {
assert(row < rows_);
assert(column < columns_);
return elements_[row * columns_ + column];
}
element_type& operator()(size_t row, size_t column) {
assert(row < rows_);
assert(column < columns_);
return elements_[row * columns_ + column];
}
 
friend bool operator==(const complex_matrix& a, const complex_matrix& b) {
return a.rows_ == b.rows_ && a.columns_ == b.columns_ &&
a.elements_ == b.elements_;
}
 
private:
size_t rows_;
size_t columns_;
std::vector<element_type> elements_;
};
 
template <typename scalar_type>
complex_matrix<scalar_type> product(const complex_matrix<scalar_type>& a,
const complex_matrix<scalar_type>& b) {
assert(a.columns() == b.rows());
size_t arows = a.rows();
size_t bcolumns = b.columns();
size_t n = a.columns();
complex_matrix<scalar_type> c(arows, bcolumns);
for (size_t i = 0; i < arows; ++i) {
for (size_t j = 0; j < n; ++j) {
for (size_t k = 0; k < bcolumns; ++k)
c(i, k) += a(i, j) * b(j, k);
}
}
return c;
}
 
template <typename scalar_type>
complex_matrix<scalar_type>
conjugate_transpose(const complex_matrix<scalar_type>& a) {
size_t rows = a.rows(), columns = a.columns();
complex_matrix<scalar_type> b(columns, rows);
for (size_t i = 0; i < columns; i++) {
for (size_t j = 0; j < rows; j++) {
b(i, j) = std::conj(a(j, i));
}
}
return b;
}
 
template <typename scalar_type>
std::string to_string(const std::complex<scalar_type>& c) {
std::ostringstream out;
const int precision = 6;
out << std::fixed << std::setprecision(precision);
out << std::setw(precision + 3) << c.real();
if (c.imag() > 0)
out << " + " << std::setw(precision + 2) << c.imag() << 'i';
else if (c.imag() == 0)
out << " + " << std::setw(precision + 2) << 0.0 << 'i';
else
out << " - " << std::setw(precision + 2) << -c.imag() << 'i';
return out.str();
}
 
template <typename scalar_type>
void print(std::ostream& out, const complex_matrix<scalar_type>& a) {
size_t rows = a.rows(), columns = a.columns();
for (size_t row = 0; row < rows; ++row) {
for (size_t column = 0; column < columns; ++column) {
if (column > 0)
out << ' ';
out << to_string(a(row, column));
}
out << '\n';
}
}
 
template <typename scalar_type>
bool is_hermitian_matrix(const complex_matrix<scalar_type>& matrix) {
if (matrix.rows() != matrix.columns())
return false;
return matrix == conjugate_transpose(matrix);
}
 
template <typename scalar_type>
bool is_normal_matrix(const complex_matrix<scalar_type>& matrix) {
if (matrix.rows() != matrix.columns())
return false;
auto c = conjugate_transpose(matrix);
return product(c, matrix) == product(matrix, c);
}
 
bool is_equal(const std::complex<double>& a, double b) {
constexpr double e = 1e-15;
return std::abs(a.imag()) < e && std::abs(a.real() - b) < e;
}
 
template <typename scalar_type>
bool is_identity_matrix(const complex_matrix<scalar_type>& matrix) {
if (matrix.rows() != matrix.columns())
return false;
size_t rows = matrix.rows();
for (size_t i = 0; i < rows; ++i) {
for (size_t j = 0; j < rows; ++j) {
if (!is_equal(matrix(i, j), scalar_type(i == j ? 1 : 0)))
return false;
}
}
return true;
}
 
template <typename scalar_type>
bool is_unitary_matrix(const complex_matrix<scalar_type>& matrix) {
if (matrix.rows() != matrix.columns())
return false;
auto c = conjugate_transpose(matrix);
auto p = product(c, matrix);
return is_identity_matrix(p) && p == product(matrix, c);
}
 
template <typename scalar_type>
void test(const complex_matrix<scalar_type>& matrix) {
std::cout << "Matrix:\n";
print(std::cout, matrix);
std::cout << "Conjugate transpose:\n";
print(std::cout, conjugate_transpose(matrix));
std::cout << std::boolalpha;
std::cout << "Hermitian: " << is_hermitian_matrix(matrix) << '\n';
std::cout << "Normal: " << is_normal_matrix(matrix) << '\n';
std::cout << "Unitary: " << is_unitary_matrix(matrix) << '\n';
}
 
int main() {
using matrix = complex_matrix<double>;
 
matrix matrix1(3, 3, {{{2, 0}, {2, 1}, {4, 0}},
{{2, -1}, {3, 0}, {0, 1}},
{{4, 0}, {0, -1}, {1, 0}}});
 
double n = std::sqrt(0.5);
matrix matrix2(3, 3, {{{n, 0}, {n, 0}, {0, 0}},
{{0, -n}, {0, n}, {0, 0}},
{{0, 0}, {0, 0}, {0, 1}}});
 
matrix matrix3(3, 3, {{{2, 2}, {3, 1}, {-3, 5}},
{{2, -1}, {4, 1}, {0, 0}},
{{7, -5}, {1, -4}, {1, 0}}});
 
test(matrix1);
std::cout << '\n';
test(matrix2);
std::cout << '\n';
test(matrix3);
return 0;
}</syntaxhighlight>
 
{{out}}
<pre>
Matrix:
2.000000 + 0.000000i 2.000000 + 1.000000i 4.000000 + 0.000000i
2.000000 - 1.000000i 3.000000 + 0.000000i 0.000000 + 1.000000i
4.000000 + 0.000000i 0.000000 - 1.000000i 1.000000 + 0.000000i
Conjugate transpose:
2.000000 + 0.000000i 2.000000 + 1.000000i 4.000000 + 0.000000i
2.000000 - 1.000000i 3.000000 + 0.000000i 0.000000 + 1.000000i
4.000000 + 0.000000i 0.000000 - 1.000000i 1.000000 + 0.000000i
Hermitian: true
Normal: true
Unitary: false
 
Matrix:
0.707107 + 0.000000i 0.707107 + 0.000000i 0.000000 + 0.000000i
0.000000 - 0.707107i 0.000000 + 0.707107i 0.000000 + 0.000000i
0.000000 + 0.000000i 0.000000 + 0.000000i 0.000000 + 1.000000i
Conjugate transpose:
0.707107 + 0.000000i 0.000000 + 0.707107i 0.000000 + 0.000000i
0.707107 + 0.000000i 0.000000 - 0.707107i 0.000000 + 0.000000i
0.000000 + 0.000000i 0.000000 + 0.000000i 0.000000 - 1.000000i
Hermitian: false
Normal: true
Unitary: true
 
Matrix:
2.000000 + 2.000000i 3.000000 + 1.000000i -3.000000 + 5.000000i
2.000000 - 1.000000i 4.000000 + 1.000000i 0.000000 + 0.000000i
7.000000 - 5.000000i 1.000000 - 4.000000i 1.000000 + 0.000000i
Conjugate transpose:
2.000000 - 2.000000i 2.000000 + 1.000000i 7.000000 + 5.000000i
3.000000 - 1.000000i 4.000000 - 1.000000i 1.000000 + 4.000000i
-3.000000 - 5.000000i 0.000000 + 0.000000i 1.000000 + 0.000000i
Hermitian: false
Normal: false
Unitary: false
</pre>
 
=={{header|Common Lisp}}==
<syntaxhighlight lang="lisp">
<lang Lisp>
(defun matrix-multiply (m1 m2)
(mapcar
Line 351 ⟶ 863:
(defun unitary-p (m)
(identity-p (matrix-multiply m (conjugate-transpose m))) )
</syntaxhighlight>
</lang>
 
{{out}}
Line 374 ⟶ 886:
=={{header|D}}==
{{trans|Python}} A well typed and mostly imperative version:
<langsyntaxhighlight lang="d">import std.stdio, std.complex, std.math, std.range, std.algorithm,
std.numeric;
 
Line 480 ⟶ 992:
writefln("Unitary? %s.\n", isUnitary(mat, ct));
}
}</langsyntaxhighlight>
{{out}}
<pre>Matrix:
Line 519 ⟶ 1,031:
===Alternative Version===
A more functional version that contains some typing problems (same output).
<langsyntaxhighlight lang="d">import std.stdio, std.complex, std.math, std.range, std.algorithm,
std.numeric, std.exception, std.traits;
 
Line 608 ⟶ 1,120:
writefln("Unitary? %s.\n", isUnitary(mat, ct));
}
}</langsyntaxhighlight>
 
=={{header|F_Sharp|F#}}==
<syntaxhighlight lang="fsharp">
// Conjugate transpose. Nigel Galloway: January 10th., 2022
let fN g=let g=g|>List.map(List.map(fun(n,g)->System.Numerics.Complex(n,g)))|>MathNet.Numerics.LinearAlgebra.MatrixExtensions.matrix in (g,g.ConjugateTranspose())
let fG n g=(MathNet.Numerics.LinearAlgebra.Matrix.inverse n-g)|>MathNet.Numerics.LinearAlgebra.Matrix.forall(fun(n:System.Numerics.Complex)->abs n.Real<1e-14&&abs n.Imaginary<1e-14)
let test=[fN [[(3.0,0.0);(2.0,1.0)];[(2.0,-1.0);(1.0,0.0)]];fN [[(1.0,0.0);(1.0,0.0);(0.0,0.0)];[(0.0,0.0);(1.0,0.0);(1.0,0.0)];[(1.0,0.0);(0.0,0.0);(1.0,0.0)]];fN [[(1.0/sqrt 2.0,0.0);(1.0/sqrt 2.0,0.0);(0.0,0.0)];[(0.0,1.0/sqrt 2.0);(0.0,-1.0/sqrt 2.0);(0.0,0.0)];[(0.0,0.0);(0.0,0.0);(0.0,1.0)]]]
test|>List.iter(fun(n,g)->printfn $"Matrix\n------\n%A{n}\nConjugate transposed\n--------------------\n%A{g}\nIs hermitian: %A{n.IsHermitian()}\nIs normal: %A{n*g=g*n}\nIs unitary: %A{fG n g}\n")
</syntaxhighlight>
{{out}}
<pre>
Matrix
------
DenseMatrix 2x2-Complex
(3, 0) (2, 1)
(2, -1) (1, 0)
 
Conjugate transposed
--------------------
DenseMatrix 2x2-Complex
(3, -0) (2, 1)
(2, -1) (1, -0)
 
Is hermitian: true
Is normal: true
Is unitary: false
 
Matrix
------
DenseMatrix 3x3-Complex
(1, 0) (1, 0) (0, 0)
(0, 0) (1, 0) (1, 0)
(1, 0) (0, 0) (1, 0)
 
Conjugate transposed
--------------------
DenseMatrix 3x3-Complex
(1, -0) (0, -0) (1, -0)
(1, -0) (1, -0) (0, -0)
(0, -0) (1, -0) (1, -0)
 
Is hermitian: false
Is normal: true
Is unitary: false
 
Matrix
------
DenseMatrix 3x3-Complex
(0.707107, 0) (0.707107, 0) (0, 0)
(0, 0.707107) (0, -0.707107) (0, 0)
(0, 0) (0, 0) (0, 1)
 
Conjugate transposed
--------------------
DenseMatrix 3x3-Complex
(0.707107, -0) (0, -0.707107) (0, -0)
(0.707107, -0) (0, 0.707107) (0, -0)
(0, -0) (0, -0) (0, -1)
 
Is hermitian: false
Is normal: true
Is unitary: true
</pre>
=={{header|Factor}}==
Before the fix to [https://github.com/slavapestov/factor/issues/484 Factor bug #484], <code>m.</code> gave the wrong answer and this code failed. Factor 0.94 is too old to work.
 
{{works with|Factor|development (future 0.95)}}
<langsyntaxhighlight lang="factor">USING: kernel math.functions math.matrices sequences ;
IN: rosetta.hermitian
 
Line 627 ⟶ 1,201:
 
: unitary-matrix? ( matrix -- ? )
[ dup conj-t m. ] [ length identity-matrix ] bi = ;</langsyntaxhighlight>
 
Usage:
Line 645 ⟶ 1,219:
The examples and algorithms are taken from the j solution, except for UnitaryQ. The j solution uses the matrix inverse verb. Compilation on linux, assuming the program is file f.f08 :<pre>
gfortran -std=f2008 -Wall -fopenmp -ffree-form -fall-intrinsics -fimplicit-none f.f08 -o f</pre>
<syntaxhighlight lang="fortran">
<lang FORTRAN>
program conjugate_transpose
 
Line 710 ⟶ 1,284:
 
end program conjugate_transpose
</syntaxhighlight>
</lang>
<pre>
-*- mode: compilation; default-directory: "/tmp/" -*-
Line 755 ⟶ 1,329:
Compilation finished at Fri Jun 7 16:31:38
</pre>
 
=={{header|FreeBASIC}}==
<syntaxhighlight lang="freebasic">'complex type and operators for it
type complex
real as double
imag as double
end type
 
operator + ( a as complex, b as complex ) as complex
dim as complex r
r.real = a.real + b.real
r.imag = a.imag + b.imag
return r
end operator
 
operator * ( a as complex, b as complex ) as complex
dim as complex r
r.real = a.real*b.real - a.imag*b.imag
r.imag = a.real*b.imag + b.real*a.imag
return r
end operator
 
operator = ( a as complex, b as complex ) as boolean
if not a.real = b.real then return false
if not a.imag = b.imag then return false
return true
end operator
 
function complex_conjugate( a as complex ) as complex
dim as complex r
r.real = a.real
r.imag = -a.imag
return r
end function
 
'matrix type and operations for it
'reuses code from the matrix multiplication task
type Matrix
dim as complex m( any , any )
declare constructor ( )
declare constructor ( byval x as uinteger )
end type
 
constructor Matrix ( )
end constructor
 
constructor Matrix ( byval x as uinteger )
redim this.m( x - 1 , x - 1 )
end constructor
 
operator * ( byref a as Matrix , byref b as Matrix ) as Matrix
dim as Matrix ret
dim as uinteger i, j, k
redim ret.m( ubound( a.m , 1 ) , ubound( a.m , 1 ) )
for i = 0 to ubound( a.m , 1 )
for j = 0 to ubound( b.m , 2 )
for k = 0 to ubound( b.m , 1 )
ret.m( i , j ) += a.m( i , k ) * b.m( k , j )
next k
next j
next i
return ret
end operator
 
function conjugate_transpose( byref a as Matrix ) as Matrix
dim as Matrix ret
dim as uinteger i, j
redim ret.m( ubound( a.m , 1 ) , ubound( a.m , 1 ) )
for i = 0 to ubound( a.m , 1 )
for j = 0 to ubound( a.m , 2 )
ret.m( i, j ) = complex_conjugate(a.m( j, i ))
next j
next i
return ret
end function
 
'tests if matrices are unitary, hermitian, or normal
 
operator = (byref a as Matrix, byref b as matrix) as boolean
dim as integer i, j
if ubound(a.m, 1) <> ubound(b.m, 1) then return false
for i = 0 to ubound( a.m , 1 )
for j = 0 to ubound( a.m , 2 )
if not a.m(i,j)=b.m(i,j) then return false
next j
next i
return true
end operator
 
function is_identity( byref a as Matrix ) as boolean
dim as integer i, j
for i = 0 to ubound( a.m , 1 )
for j = 0 to ubound( a.m , 2 )
if i = j and ( not a.m(i,j).real = 1.0 or not a.m(i,j).imag = 0.0 ) then return false
if i <> j and ( not a.m(i,j).real = 0.0 or not a.m(i,j).imag = 0.0 ) then return false
next j
next i
return true
end function
 
function is_hermitian( byref a as Matrix ) as boolean
if a = conjugate_transpose(a) then return true
return false
end function
 
function is_normal( byref a as Matrix ) as boolean
dim as Matrix aa = conjugate_transpose(a)
if a*aa = aa*a then return true else return false
end function
 
function is_unitary( byref a as Matrix ) as boolean
dim as Matrix aa = conjugate_transpose(a)
if not is_identity( a*aa ) or not is_identity( aa*a ) then return false
return true
end function
 
'''now some example matrices
dim as Matrix A = Matrix(2) 'an identity matrix
A.m(0,0).real = 1.0 : A.m(0,0).imag = 0.0 : A.m(0,1).real = 0.0 : A.m(0,1).imag = 0.0
A.m(1,0).real = 0.0 : A.m(1,0).imag = 0.0 : A.m(1,1).real = 1.0 : A.m(1,1).imag = 0.0
 
dim as Matrix B = Matrix(2) 'a hermitian matrix
B.m(0,0).real = 1.0 : B.m(0,0).imag = 0.0 : B.m(0,1).real = 1.0 : B.m(0,1).imag = -1.0
B.m(1,0).real = 1.0 : B.m(1,0).imag = 1.0 : B.m(1,1).real = 1.0 : B.m(1,1).imag = 0.0
 
dim as Matrix C = Matrix(2) 'a random matrix
C.m(0,0).real = rnd : C.m(0,0).imag = rnd : C.m(0,1).real = rnd : C.m(0,1).imag = rnd
C.m(1,0).real = rnd : C.m(1,0).imag = rnd : C.m(1,1).real = rnd : C.m(1,1).imag = rnd
 
print is_hermitian(A), is_normal(A), is_unitary(A)
print is_hermitian(B), is_normal(B), is_unitary(B)
print is_hermitian(C), is_normal(C), is_unitary(C)</syntaxhighlight>
{{out}}
<pre>true true true
true true false
false false false</pre>
 
=={{header|Go}}==
<langsyntaxhighlight lang="go">package main
 
import (
Line 868 ⟶ 1,578:
}
return m3
}</langsyntaxhighlight>
Output:
<pre>
Line 911 ⟶ 1,621:
=={{header|Haskell}}==
Slow implementation using lists.
<langsyntaxhighlight lang="haskell">import Data.ListComplex (transposeComplex(..), conjugate)
import Data.ComplexList (transpose)
 
type Matrix a = [[a]]
Line 918 ⟶ 1,628:
main :: IO ()
main =
mapM_ (\a -> do
(\a -> do
putStrLn "\nMatrix:"
putStrLn mapM_ print a"\nMatrix:"
mapM_ putStrLnprint "Conjugate Transpose:"a
putStrLn "Conjugate Transpose:"
mapM_ print (conjTranspose a)
mapM_ print (conjTranspose a)
putStrLn $ "Hermitian? " ++ show (isHermitianMatrix a)
putStrLn $ "NormalHermitian? " ++ show (isNormalMatrixisHermitianMatrix a)
putStrLn $ "UnitaryNormal? " ++ show (isUnitaryMatrixisNormalMatrix a))
putStrLn ([[[3,$ "Unitary? " ++ show (isUnitaryMatrix 2:+1],a))
([ [[3, 2 :+ 1], [2 :+ (-1), 1 ]],
, [[1, 1, 0], [0, 1, 1], [1, 0, 1]]
, [ [sqrt 2 / 2 :+ 0, sqrt 2 / 2 :+ 0, 0]
, [0 :+ sqrt 2 / 2, 0 :+ (-sqrt 2 / 2), 0]
, [0, 0, 0 :+ 1]
]
] :: [Matrix (Complex Double)])
 
isHermitianMatrix, isNormalMatrix, isUnitaryMatrix
[[1, 1, 0],
:: RealFloat a
[0, 1, 1],
=> Matrix (Complex a) -> Bool
[1, 0, 1]],
isHermitianMatrix = mTest id conjTranspose
 
isNormalMatrix = mTest mmct (mmul =<< conjTranspose)
[[sqrt 2/2:+0, sqrt 2/2:+0, 0 ],
[0:+sqrt 2/2, 0:+ (-sqrt 2/2), 0 ],
[0, 0, 0:+1]]] :: [Matrix (Complex Double)])
 
isHermitianMatrix, isNormalMatrix, isUnitaryMatrix ::= RealFloatmTest ammct =>(ident Matrix. (Complex alength) -> Bool
isHermitianMatrix a = a `approxEqualMatrix` conjTranspose a
isNormalMatrix a = (a `mmul` conjTranspose a) `approxEqualMatrix` (conjTranspose a `mmul` a)
isUnitaryMatrix a = (a `mmul` conjTranspose a) `approxEqualMatrix` ident (length a)
 
mTest
approxEqualMatrix :: (Fractional a, Ord a) => Matrix (Complex a) -> Matrix (Complex a) -> Bool
:: RealFloat a
approxEqualMatrix a b = length a == length b && length (head a) == length (head b) &&
=> (a2 -> Matrix (Complex a)) -> (a2 -> Matrix (Complex a)) -> a2 -> Bool
and (zipWith approxEqualComplex (concat a) (concat b))
mTest f g = (approxEqualMatrix . f) <*> g
where approxEqualComplex (rx :+ ix) (ry :+ iy) = abs (rx - ry) < eps && abs (ix - iy) < eps
eps = 1e-14
 
mmct
mmul :: Num a => Matrix a -> Matrix a -> Matrix a
:: RealFloat a
mmul a b = [[sum (zipWith (*) row column) | column <- transpose b] | row <- a]
=> Matrix (Complex a) -> Matrix (Complex a)
mmct = mmul <*> conjTranspose
 
approxEqualMatrix
ident :: Num a => Int -> Matrix a
:: (Fractional a, Ord a)
ident size = [[fromIntegral $ div a b * div b a | a <- [1..size]] | b <- [1..size]]
=> Matrix (Complex a) -> Matrix (Complex a) -> Bool
approxEqualMatrix a b =
length a == length b &&
length (head a) == length (head b) &&
and (zipWith approxEqualComplex (concat a) (concat b))
where
approxEqualComplex (rx :+ ix) (ry :+ iy) =
abs (rx - ry) < eps && abs (ix - iy) < eps
eps = 1e-14
 
mmul
conjTranspose :: Num a => Matrix (Complex a) -> Matrix (Complex a)
:: Num a
conjTranspose = map (map conjugate) . transpose</lang>
=> Matrix a -> Matrix a -> Matrix a
mmul a b =
[ [ sum (zipWith (*) row column)
| column <- transpose b ]
| row <- a ]
 
ident
:: Num a
=> Int -> Matrix a
ident size =
[ [ fromIntegral $ div a b * div b a
| a <- [1 .. size] ]
| b <- [1 .. size] ]
 
conjTranspose
:: Num a
=> Matrix (Complex a) -> Matrix (Complex a)
conjTranspose = map (map conjugate) . transpose</syntaxhighlight>
Output:
<pre>
Line 995 ⟶ 1,735:
=={{header|J}}==
 
'''SolutionConjugate transpose:''': <langsyntaxhighlight lang="j"> ct =: +@|: NB. Conjugate transpose (ct A is A_ct)</langsyntaxhighlight>
 
'''Examples''': <lang j> X =: +/ . * NB. Matrix Multiply (x)
<code>+</code> when used without a left argument is conjugate, <code>|:</code> is transpose, and <code>@</code> composes functions.
 
'''Examples''': <syntaxhighlight lang="j"> X =: +/ . * NB. Matrix Multiply (x)
 
HERMITIAN =: 3 2j1 ,: 2j_1 1
Line 1,008 ⟶ 1,751:
UNITARY =: (-:%:2) * 1 1 0 , 0j_1 0j1 0 ,: 0 0 0j1 * %:2
(ct -: %.) UNITARY NB. A_ct = A^-1
1</langsyntaxhighlight>
 
'''Reference''' (with example matrices for other langs to use):<syntaxhighlight lang ="j"> HERMITIAN;NORMAL;UNITARY
+--------+-----+--------------------------+
| 3 2j1|1 1 0| 0.707107 0.707107 0|
Line 1,025 ⟶ 1,768:
+-----+-----+-----+
|1 1 0|0 1 0|0 1 1|
+-----+-----+-----+</langsyntaxhighlight>
 
=={{header|Java}}==
<syntaxhighlight lang ="java">
 
import java.util.Arrays;
import java.util.List;
 
public final class ConjugateTranspose {
 
public static void main(String[] aArgs) {
ComplexMatrix one = new ComplexMatrix( new Complex[][] { { new Complex(0, 4), new Complex(-1, 1) },
{ new Complex(1, -1), new Complex(0, 4) } } );
 
ComplexMatrix two = new ComplexMatrix(
new Complex[][] { { new Complex(1, 0), new Complex(1, 1), new Complex(0, 2) },
{ new Complex(1, -1), new Complex(5, 0), new Complex(-3, 0) },
{ new Complex(0, -2), new Complex(-3, 0), new Complex(0, 0) } } );
 
final double term = 1.0 / Math.sqrt(2.0);
ComplexMatrix three = new ComplexMatrix( new Complex[][] { { new Complex(term, 0), new Complex(term, 0) },
{ new Complex(0, term), new Complex(0, -term) } } );
List<ComplexMatrix> matricies = List.of( one, two, three );
for ( ComplexMatrix matrix : matricies ) {
System.out.println("Matrix:");
matrix.display();
System.out.println("Conjugate transpose:");
matrix.conjugateTranspose().display();
System.out.println("Hermitian: " + matrix.isHermitian());
System.out.println("Normal: " + matrix.isNormal());
System.out.println("Unitary: " + matrix.isUnitary() + System.lineSeparator());
}
}
 
}
 
final class ComplexMatrix {
public ComplexMatrix(Complex[][] aData) {
rowCount = aData.length;
colCount = aData[0].length;
data = Arrays.stream(aData).map( row -> Arrays.copyOf(row, row.length) ).toArray(Complex[][]::new);
}
public ComplexMatrix multiply(ComplexMatrix aOther) {
if ( colCount != aOther.rowCount ) {
throw new RuntimeException("Incompatible matrix dimensions.");
}
Complex[][] newData = new Complex[rowCount][aOther.colCount];
Arrays.stream(newData).forEach( row -> Arrays.fill(row, new Complex(0, 0)) );
for ( int row = 0; row < rowCount; row++ ) {
for ( int col = 0; col < aOther.colCount; col++ ) {
for ( int k = 0; k < colCount; k++ ) {
newData[row][col] = newData[row][col].add(data[row][k].multiply(aOther.data[k][col]));
}
}
}
return new ComplexMatrix(newData);
}
 
public ComplexMatrix conjugateTranspose() {
if ( rowCount != colCount ) {
throw new IllegalArgumentException("Only applicable to a square matrix");
}
Complex[][] newData = new Complex[colCount][rowCount];
for ( int row = 0; row < rowCount; row++ ) {
for ( int col = 0; col < colCount; col++ ) {
newData[col][row] = data[row][col].conjugate();
}
}
return new ComplexMatrix(newData);
}
public static ComplexMatrix identity(int aSize) {
Complex[][] data = new Complex[aSize][aSize];
for ( int row = 0; row < aSize; row++ ) {
for ( int col = 0; col < aSize; col++ ) {
data[row][col] = ( row == col ) ? new Complex(1, 0) : new Complex(0, 0);
}
}
return new ComplexMatrix(data);
}
public boolean equals(ComplexMatrix aOther) {
if ( aOther.rowCount != rowCount || aOther.colCount != colCount ) {
return false;
}
for ( int row = 0; row < rowCount; row++ ) {
for ( int col = 0; col < colCount; col++ ) {
if ( data[row][col].subtract(aOther.data[row][col]).modulus() > EPSILON ) {
return false;
}
}
}
return true;
}
public void display() {
for ( int row = 0; row < rowCount; row++ ) {
System.out.print("[");
for ( int col = 0; col < colCount - 1; col++ ) {
System.out.print(data[row][col] + ", ");
}
System.out.println(data[row][colCount - 1] + " ]");
}
}
public boolean isHermitian() {
return equals(conjugateTranspose());
}
public boolean isNormal() {
ComplexMatrix conjugateTranspose = conjugateTranspose();
return multiply(conjugateTranspose).equals(conjugateTranspose.multiply(this));
}
public boolean isUnitary() {
ComplexMatrix conjugateTranspose = conjugateTranspose();
return multiply(conjugateTranspose).equals(identity(rowCount)) &&
conjugateTranspose.multiply(this).equals(identity(rowCount));
}
private final int rowCount;
private final int colCount;
private final Complex[][] data;
private static final double EPSILON = 0.000_000_000_001;
}
 
final class Complex {
public Complex(double aReal, double aImag) {
real = aReal;
imag = aImag;
}
public Complex add(Complex aOther) {
return new Complex(real + aOther.real, imag + aOther.imag);
}
public Complex multiply(Complex aOther) {
return new Complex(real * aOther.real - imag * aOther.imag, real * aOther.imag + imag * aOther.real);
}
public Complex negate() {
return new Complex(-real, -imag);
}
public Complex subtract(Complex aOther) {
return this.add(aOther.negate());
}
public Complex conjugate() {
return new Complex(real, -imag);
}
public double modulus() {
return Math.hypot(real, imag);
}
public boolean equals(Complex aOther) {
return real == aOther.real && imag == aOther.imag;
}
@Override
public String toString() {
String prefix = ( real < 0.0 ) ? "" : " ";
String realPart = prefix + String.format("%.3f", real);
String sign = ( imag < 0.0 ) ? " - " : " + ";
return realPart + sign + String.format("%.3f", Math.abs(imag)) + "i";
}
private final double real;
private final double imag;
}
</syntaxhighlight>
{{ out }}
<pre>
Matrix:
[ 0.000 + 4.000i, -1.000 + 1.000i ]
[ 1.000 - 1.000i, 0.000 + 4.000i ]
Conjugate transpose:
[ 0.000 - 4.000i, 1.000 + 1.000i ]
[-1.000 - 1.000i, 0.000 - 4.000i ]
Hermitian: false
Normal: true
Unitary: false
 
Matrix:
[ 1.000 + 0.000i, 1.000 + 1.000i, 0.000 + 2.000i ]
[ 1.000 - 1.000i, 5.000 + 0.000i, -3.000 + 0.000i ]
[ 0.000 - 2.000i, -3.000 + 0.000i, 0.000 + 0.000i ]
Conjugate transpose:
[ 1.000 + 0.000i, 1.000 + 1.000i, 0.000 + 2.000i ]
[ 1.000 - 1.000i, 5.000 + 0.000i, -3.000 + 0.000i ]
[ 0.000 - 2.000i, -3.000 + 0.000i, 0.000 + 0.000i ]
Hermitian: true
Normal: true
Unitary: false
 
Matrix:
[ 0.707 + 0.000i, 0.707 + 0.000i ]
[ 0.000 + 0.707i, 0.000 - 0.707i ]
Conjugate transpose:
[ 0.707 + 0.000i, 0.000 - 0.707i ]
[ 0.707 + 0.000i, 0.000 + 0.707i ]
Hermitian: false
Normal: true
Unitary: true
</pre>
 
=={{header|jq}}==
Line 1,037 ⟶ 1,992:
 
If your jq does not have "transpose" then the following may be used:
<langsyntaxhighlight lang="jq"># transpose/0 expects its input to be a rectangular matrix
# (an array of equal-length arrays):
def transpose:
if (.[0] | length) == 0 then []
else [map(.[0])] + (map(.[1:]) | transpose)
end ;</langsyntaxhighlight>
'''(2) Operations on real/complex numbers'''
<langsyntaxhighlight lang="jq"># x must be real or complex, and ditto for y;
# always return complex
def plus(x; y):
Line 1,070 ⟶ 2,025:
if type == "number" then [.,0]
else [.[0], -(.[1]) ]
end;</langsyntaxhighlight>
'''(3) Array operations'''
<langsyntaxhighlight lang="jq"># a and b are arrays of real/complex numbers
def dot_product(a; b):
a as $a | b as $b
| reduce range(0;$a|length) as $i
(0; . as $s | plus($s; multiply($a[$i]; $b[$i]) ));</langsyntaxhighlight>
'''(4) Matrix operations'''
<langsyntaxhighlight lang="jq"># convert a matrix of mixed real/complex entries to all complex entries
def to_complex:
def toc: if type == "number" then [.,0] else . end;
Line 1,115 ⟶ 2,070:
reduce range(0;M|length) as $i
(0; reduce range(0; M[0]|length) as $j
(.; 0 + sqdiff( M[$i][$j]; N[$i][$j] ) ) ) <= epsilon;</langsyntaxhighlight>
====Conjugate transposition====
<langsyntaxhighlight lang="jq"># (entries may be real and/or complex)
def conjugate_transpose:
map( map(conjugate) ) | transpose;
Line 1,139 ⟶ 2,094:
| complex_identity(length) as $I
| approximately_equal( $I; matrix_multiply($H;$M); 1e-10)
and approximately_equal( $I ; matrix_multiply($M;$H); 1e-10) ; </langsyntaxhighlight>
 
====Examples====
<langsyntaxhighlight lang="jq">def hermitian_example:
[ [ 3, [2,1]],
[[2,-1], 1 ] ];
Line 1,170 ⟶ 2,125:
;
 
demo</langsyntaxhighlight>
{{out}}
<langsyntaxhighlight lang="sh">$ jq -r -c -n -f Conjugate_transpose.jq
Hermitian example:
 
Line 1,187 ⟶ 2,142:
Normal example: true
 
Unitary example: true</langsyntaxhighlight>
 
=={{header|Julia}}==
Julia has a built-in matrix type, and the conjugate-transpose of a complex matrix <code>A</code> is simply:
<syntaxhighlight lang ="julia">A'</langsyntaxhighlight>
(similar to Matlab). You can check whether <code>A</code> is Hermitian via the built-in function
<syntaxhighlight lang ="julia">ishermitian(A)</langsyntaxhighlight>
Ignoring the possibility of roundoff errors for floating-point matrices (like most of the examples in the other languages), you can check whether a matrix is normal or unitary by the following functions
<langsyntaxhighlight lang="julia">eye(A) = A^0
isnormal(A) = size(A,1) == size(A,2) && A'*A == A*A'
isunitary(A) = size(A,1) == size(A,2) && A'*A == eye(A)</langsyntaxhighlight>
 
=={{header|Kotlin}}==
As Kotlin doesn't have built in classes for complex numbers or matrices, some basic functionality needs to be coded in order to tackle this task:
<langsyntaxhighlight lang="scala">// version 1.1.3
 
typealias C = Complex
Line 1,325 ⟶ 2,280:
println("Unitary? ${mct.isUnitary()}\n")
}
}</langsyntaxhighlight>
 
{{out}}
Line 1,372 ⟶ 2,327:
=={{header|Maple}}==
The commands <code>HermitianTranspose</code> and <code>IsUnitary</code> are provided by the <code>LinearAlgebra</code> package.
<langsyntaxhighlight Maplelang="maple">M:=<<3|2+I>,<2-I|1>>:
 
with(LinearAlgebra):
Line 1,381 ⟶ 2,336:
type(M,'Matrix'(hermitian));
IsNormal(M);
IsUnitary(M);</langsyntaxhighlight>
Output:
<pre> [ 3 2 + I]
Line 1,398 ⟶ 2,353:
 
=={{header|Mathematica}} / {{header|Wolfram Language}}==
<langsyntaxhighlight Mathematicalang="mathematica">NormalMatrixQ[a_List?MatrixQ] := Module[{b = Conjugate@Transpose@a},a.b === b.a]
UnitaryQ[m_List?MatrixQ] := (Conjugate@Transpose@m.m == IdentityMatrix@Length@m)
 
Line 1,414 ⟶ 2,369:
 
{HermitianMatrixQ@#, NormalMatrixQ@#, UnitaryQ@#}&@m
-> {False, False, False}</langsyntaxhighlight>
 
=={{header|Nim}}==
 
The complex type is defined as generic regarding the type of real an imaginary part. We have chosen to use Complex[float] and make only our Matrix type generic regarding the dimensions. Thus, a Matrix has a two dimensions M and N which are static, i.e. known at compile time. We have enforced the condition M = N for square matrices (also at compile time).
 
<syntaxhighlight lang="nim">import complex, strformat
 
type Matrix[M, N: static Positive] = array[M, array[N, Complex[float]]]
 
const Eps = 1e-10 # Tolerance used for float comparisons.
 
 
####################################################################################################
# Templates.
 
template `[]`(m: Matrix; i, j: Natural): Complex[float] =
## Allow to get value of an element using m[i, j] syntax.
m[i][j]
 
template `[]=`(m: var Matrix; i, j: Natural; val: Complex[float]) =
## Allow to set value of an element using m[i, j] syntax.
m[i][j] = val
 
 
####################################################################################################
# General operations.
 
func `$`(m: Matrix): string =
## Return the string representation of a matrix using one line per row.
 
for i, row in m:
result.add(if i == 0: '[' else: ' ')
for j, val in row:
if j != 0: result.add(' ')
result.add(&"({val.re:7.4f}, {val.im:7.4f})")
result.add(if i == m.high: ']' else: '\n')
 
#---------------------------------------------------------------------------------------------------
 
func conjugateTransposed[M, N: static int](m: Matrix[M, N]): Matrix[N, M] =
## Return the conjugate transpose of a matrix.
 
for i in 0..<m.M:
for j in 0..<m.N:
result[j, i] = m[i, j].conjugate()
 
#---------------------------------------------------------------------------------------------------
 
func `*`[M, K, N: static int](m1: Matrix[M, K]; m2: Matrix[K, N]): Matrix[M, N] =
# Compute the product of two matrices.
 
for i in 0..<M:
for j in 0..<N:
for k in 0..<K:
result[i, j] = result[i, j] + m1[i, k] * m2[k, j]
 
 
####################################################################################################
# Properties.
 
func isHermitian(m: Matrix): bool =
## Check if a matrix is hermitian.
 
when m.M != m.N:
{.error: "hermitian test only allowed for square matrices".}
else:
for i in 0..<m.M:
for j in i..<m.N:
if m[i, j] != m[j, i].conjugate:
return false
result = true
 
#---------------------------------------------------------------------------------------------------
 
func isNormal(m: Matrix): bool =
## Check if a matrix is normal.
 
when m.M != m.N:
{.error: "normal test only allowed for square matrices".}
else:
let h = m.conjugateTransposed
result = m * h == h * m
 
#---------------------------------------------------------------------------------------------------
 
func isIdentity(m: Matrix): bool =
## Check if a matrix is the identity matrix.
 
when m.M != m.N:
{.error: "identity test only allowed for square matrices".}
else:
for i in 0..<m.M:
for j in 0..<m.N:
if i == j:
if abs(m[i, j] - 1.0) > Eps:
return false
else:
if abs(m[i, j]) > Eps:
return false
result = true
 
#---------------------------------------------------------------------------------------------------
 
func isUnitary(m: Matrix): bool =
## Check if a matrix is unitary.
 
when m.M != m.N:
{.error: "unitary test only allowed for square matrices".}
else:
let h = m.conjugateTransposed
result = (m * h).isIdentity and (h * m).isIdentity
 
#———————————————————————————————————————————————————————————————————————————————————————————————————
 
when isMainModule:
 
import math
 
proc test(m: Matrix) =
echo "\n"
echo "Matrix"
echo "------"
echo m
echo ""
echo "Conjugate transposed"
echo "--------------------"
echo m.conjugateTransposed
 
when m.M == m.N:
# Only for squares matrices.
echo ""
echo "Hermitian: ", m.isHermitian
echo "Normal: ", m.isNormal
echo "Unitary: ", m.isUnitary
 
#-------------------------------------------------------------------------------------------------
 
# Non square matrix.
const M1: Matrix[2, 3] = [[1.0 + im 2.0, 3.0 + im 0.0, 2.0 + im 5.0],
[3.0 - im 1.0, 2.0 + im 0.0, 0.0 + im 3.0]]
 
# Square matrices.
const M2: Matrix[2, 2] = [[3.0 + im 0.0, 2.0 + im 1.0],
[2.0 - im 1.0, 1.0 + im 0.0]]
 
const M3: Matrix[3, 3] = [[1.0 + im 0.0, 1.0 + im 0.0, 0.0 + im 0.0],
[0.0 + im 0.0, 1.0 + im 0.0, 1.0 + im 0.0],
[1.0 + im 0.0, 0.0 + im 0.0, 1.0 + im 0.0]]
 
const SR2 = 1 / sqrt(2.0)
const M4: Matrix[3, 3] = [[SR2 + im 0.0, SR2 + im 0.0, 0.0 + im 0.0],
[0.0 + im SR2, 0.0 - im SR2, 0.0 + im 0.0],
[0.0 + im 0.0, 0.0 + im 0.0, 0.0 + im 1.0]]
 
test(M1)
test(M2)
test(M3)
test(M4)</syntaxhighlight>
 
{{out}}
<pre>Matrix
------
[( 1.0000, 2.0000) ( 3.0000, 0.0000) ( 2.0000, 5.0000)
( 3.0000, -1.0000) ( 2.0000, 0.0000) ( 0.0000, 3.0000)]
 
Conjugate transposed
--------------------
[( 1.0000, -2.0000) ( 3.0000, 1.0000)
( 3.0000, -0.0000) ( 2.0000, -0.0000)
( 2.0000, -5.0000) ( 0.0000, -3.0000)]
 
 
Matrix
------
[( 3.0000, 0.0000) ( 2.0000, 1.0000)
( 2.0000, -1.0000) ( 1.0000, 0.0000)]
 
Conjugate transposed
--------------------
[( 3.0000, -0.0000) ( 2.0000, 1.0000)
( 2.0000, -1.0000) ( 1.0000, -0.0000)]
 
Hermitian: true
Normal: true
Unitary: false
 
 
Matrix
------
[( 1.0000, 0.0000) ( 1.0000, 0.0000) ( 0.0000, 0.0000)
( 0.0000, 0.0000) ( 1.0000, 0.0000) ( 1.0000, 0.0000)
( 1.0000, 0.0000) ( 0.0000, 0.0000) ( 1.0000, 0.0000)]
 
Conjugate transposed
--------------------
[( 1.0000, -0.0000) ( 0.0000, -0.0000) ( 1.0000, -0.0000)
( 1.0000, -0.0000) ( 1.0000, -0.0000) ( 0.0000, -0.0000)
( 0.0000, -0.0000) ( 1.0000, -0.0000) ( 1.0000, -0.0000)]
 
Hermitian: false
Normal: true
Unitary: false
 
 
Matrix
------
[( 0.7071, 0.0000) ( 0.7071, 0.0000) ( 0.0000, 0.0000)
( 0.0000, 0.7071) ( 0.0000, -0.7071) ( 0.0000, 0.0000)
( 0.0000, 0.0000) ( 0.0000, 0.0000) ( 0.0000, 1.0000)]
 
Conjugate transposed
--------------------
[( 0.7071, -0.0000) ( 0.0000, -0.7071) ( 0.0000, -0.0000)
( 0.7071, -0.0000) ( 0.0000, 0.7071) ( 0.0000, -0.0000)
( 0.0000, -0.0000) ( 0.0000, -0.0000) ( 0.0000, -1.0000)]
 
Hermitian: false
Normal: true
Unitary: true</pre>
 
=={{header|PARI/GP}}==
<syntaxhighlight lang="text">conjtranspose(M)=conj(M~)
isHermitian(M)=M==conj(M~)
isnormal(M)=my(H=conj(M~));H*M==M*H
isunitary(M)=M*conj(M~)==1</langsyntaxhighlight>
 
=={{header|Perl}}==
In general, using two or more modules which overload operators can be problematic. For this task, using both Math::Complex and Math::MatrixReal gives us the behavior we want for everything except matrix I/O, i.e. parsing and stringification.
<langsyntaxhighlight lang="perl">use strict;
use warnings;
use English;
use Math::Complex;
Line 1,467 ⟶ 2,642:
sub identity {
my $N = shift;
my $m = new Math::MatrixReal->new($N, $N);
$m->one();
return $m;
Line 1,473 ⟶ 2,648:
 
sub example1 {
my $m = new Math::MatrixReal->new(2, 2);
$m->assign(1, 1, cplx(3, 0));
$m->assign(1, 2, cplx(2, 1));
Line 1,482 ⟶ 2,657:
 
sub example2 {
my $m = new Math::MatrixReal->new(3, 3);
$m->assign(1, 1, cplx(1, 0));
$m->assign(1, 2, cplx(1, 0));
Line 1,496 ⟶ 2,671:
 
sub example3 {
my $m = new Math::MatrixReal->new(3, 3);
$m->assign(1, 1, cplx(0.70710677, 0));
$m->assign(1, 2, cplx(0.70710677, 0));
Line 1,507 ⟶ 2,682:
$m->assign(3, 3, cplx(0, 1));
return $m;
}</langsyntaxhighlight>
{{out}}
<pre>
Line 1,546 ⟶ 2,721:
 
=={{header|Phix}}==
Phix has no support for complex numbers, so roll our own, ditto matrix maths. Note this code has no testing for non-square matrices.
<!--<syntaxhighlight lang="phix">(phixonline)-->
<lang Phix>enum REAL, IMAG
<span style="color: #008080;">with</span> <span style="color: #008080;">javascript_semantics</span>
<span style="color: #008080;">include</span> <span style="color: #004080;">complex</span><span style="color: #0000FF;">.</span><span style="color: #000000;">e</span>
type complex(sequence s)
return length(s)=2 and atom(s[REAL]) and atom(s[IMAG])
<span style="color: #008080;">procedure</span> <span style="color: #000000;">m_print</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
end type
<span style="color: #000000;">a</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">deep_copy</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
<span style="color: #004080;">integer</span> <span style="color: #000000;">l</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
function c_add(complex a, complex b)
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">l</span> <span style="color: #008080;">do</span>
return sq_add(a,b)
<span style="color: #008080;">for</span> <span style="color: #000000;">j</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">l</span> <span style="color: #008080;">do</span>
end function
<span style="color: #000000;">a</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">complex_sprint</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">])</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
function c_mul(complex a, complex b)
<span style="color: #000000;">a</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #008000;">"["</span><span style="color: #0000FF;">&</span><span style="color: #7060A8;">join</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">],</span><span style="color: #008000;">","</span><span style="color: #0000FF;">)&</span><span style="color: #008000;">"]"</span>
return {a[REAL] * b[REAL] - a[IMAG] * b[IMAG],
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
a[REAL] * b[IMAG] + a[IMAG] * b[REAL]}
<span style="color: #7060A8;">puts</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #7060A8;">join</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"\n"</span><span style="color: #0000FF;">)&</span><span style="color: #008000;">"\n"</span><span style="color: #0000FF;">)</span>
end function
<span style="color: #008080;">end</span> <span style="color: #008080;">procedure</span>
 
function c_conj(complex a)
<span style="color: #008080;">function</span> <span style="color: #000000;">conjugate_transpose</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
return {a[REAL],-a[IMAG]}
<span style="color: #004080;">sequence</span> <span style="color: #000000;">res</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">deep_copy</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
end function
<span style="color: #004080;">integer</span> <span style="color: #000000;">l</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
 
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">l</span> <span style="color: #008080;">do</span>
function c_print(complex a)
<span style="color: #008080;">for</span> <span style="color: #000000;">j</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">l</span> <span style="color: #008080;">do</span>
if a[IMAG]=0 then return sprintf("%g",a[REAL]) end if
<span style="color: #000000;">res</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">complex_conjugate</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">[</span><span style="color: #000000;">j</span><span style="color: #0000FF;">][</span><span style="color: #000000;">i</span><span style="color: #0000FF;">])</span>
return sprintf("%g%+gi",a)
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
end function
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
 
<span style="color: #008080;">return</span> <span style="color: #000000;">res</span>
procedure m_print(sequence a)
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
integer l = length(a)
for i=1 to l do
<span style="color: #008080;">function</span> <span style="color: #000000;">m_unitary</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">act</span><span style="color: #0000FF;">)</span>
for j=1 to l do
<span style="color: #000080;font-style:italic;">-- note: a was normal and act = a*ct already</span>
a[i][j] = c_print(a[i][j])
<span style="color: #004080;">integer</span> <span style="color: #000000;">l</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">act</span><span style="color: #0000FF;">)</span>
end for
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">l</span> <span style="color: #008080;">do</span>
a[i] = "["&join(a[i],",")&"]"
<span style="color: #008080;">for</span> <span style="color: #000000;">j</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">l</span> <span style="color: #008080;">do</span>
end for
<span style="color: #004080;">atom</span> <span style="color: #0000FF;">{</span><span style="color: #000000;">re</span><span style="color: #0000FF;">,</span><span style="color: #000000;">im</span><span style="color: #0000FF;">}</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">act</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">,</span><span style="color: #000000;">j</span><span style="color: #0000FF;">]</span>
puts(1,join(a,"\n")&"\n")
<span style="color: #000080;font-style:italic;">-- round to nearest billionth
end procedure
-- (powers of 2 help the FPU out)</span>
 
<span style="color: #000000;">re</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">round</span><span style="color: #0000FF;">(</span><span style="color: #000000;">re</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1024</span><span style="color: #0000FF;">*</span><span style="color: #000000;">1024</span><span style="color: #0000FF;">*</span><span style="color: #000000;">1024</span><span style="color: #0000FF;">)</span>
 
<span style="color: #000000;">im</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">round</span><span style="color: #0000FF;">(</span><span style="color: #000000;">im</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1024</span><span style="color: #0000FF;">*</span><span style="color: #000000;">1024</span><span style="color: #0000FF;">*</span><span style="color: #000000;">1024</span><span style="color: #0000FF;">)</span>
function conjugate_transpose(sequence a)
<span style="color: #008080;">if</span> <span style="color: #000000;">im</span><span style="color: #0000FF;">!=</span><span style="color: #000000;">0</span>
sequence res = a
<span style="color: #008080;">or</span> <span style="color: #0000FF;">(</span><span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">j</span> <span style="color: #008080;">and</span> <span style="color: #000000;">re</span><span style="color: #0000FF;">!=</span><span style="color: #000000;">1</span><span style="color: #0000FF;">)</span>
integer l = length(a)
<span style="color: #008080;">or</span> <span style="color: #0000FF;">(</span><span style="color: #000000;">i</span><span style="color: #0000FF;">!=</span><span style="color: #000000;">j</span> <span style="color: #008080;">and</span> <span style="color: #000000;">re</span><span style="color: #0000FF;">!=</span><span style="color: #000000;">0</span><span style="color: #0000FF;">)</span> <span style="color: #008080;">then</span>
for i=1 to l do
<span style="color: #008080;">return</span> <span style="color: #000000;">0</span>
for j=1 to l do
<span style="color: #008080;">end</span> <span style="color: #008080;">if</span>
res[i][j] = c_conj(a[j][i])
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
end for
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
end for
<span style="color: #008080;">return</span> <span style="color: #000000;">1</span>
return res
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
end function
 
<span style="color: #008080;">function</span> <span style="color: #000000;">m_mul</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">,</span> <span style="color: #004080;">sequence</span> <span style="color: #000000;">b</span><span style="color: #0000FF;">)</span>
function m_unitary(sequence act)
<span style="color: #004080;">sequence</span> <span style="color: #000000;">res</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">sq_mul</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">,</span><span style="color: #000000;">0</span><span style="color: #0000FF;">)</span>
-- note: a was normal and act = a*ct already
<span style="color: #004080;">integer</span> <span style="color: #000000;">l</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">length</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
integer l = length(act)
<span style="color: #008080;">for</span> <span style="color: #000000;">i</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">l</span> <span style="color: #008080;">do</span>
for i=1 to l do
<span style="color: #008080;">for</span> <span style="color: #000000;">j</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">l</span> <span style="color: #008080;">do</span>
for j=1 to l do
<span style="color: #008080;">for</span> <span style="color: #000000;">k</span><span style="color: #0000FF;">=</span><span style="color: #000000;">1</span> <span style="color: #008080;">to</span> <span style="color: #000000;">l</span> <span style="color: #008080;">do</span>
atom {re,im} = act[i,j]
<span style="color: #000000;">res</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">]</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">complex_add</span><span style="color: #0000FF;">(</span><span style="color: #000000;">res</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">],</span><span style="color: #7060A8;">complex_mul</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">[</span><span style="color: #000000;">i</span><span style="color: #0000FF;">][</span><span style="color: #000000;">k</span><span style="color: #0000FF;">],</span><span style="color: #000000;">b</span><span style="color: #0000FF;">[</span><span style="color: #000000;">k</span><span style="color: #0000FF;">][</span><span style="color: #000000;">j</span><span style="color: #0000FF;">]))</span>
-- round to nearest billionth
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
-- (powers of 2 help the FPU out)
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
re = round(re,1024*1024*1024)
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
im = round(im,1024*1024*1024)
<span style="color: #008080;">return</span> <span style="color: #000000;">res</span>
if im!=0
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
or (i=j and re!=1)
or (i!=j and re!=0) then
<span style="color: #008080;">procedure</span> <span style="color: #000000;">test</span><span style="color: #0000FF;">(</span><span style="color: #004080;">sequence</span> <span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
return 0
<span style="color: #004080;">sequence</span> <span style="color: #000000;">ct</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">conjugate_transpose</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
end if
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Original matrix:\n"</span><span style="color: #0000FF;">)</span>
end for
<span style="color: #000000;">m_print</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
end for
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Conjugate transpose:\n"</span><span style="color: #0000FF;">)</span>
return 1
<span style="color: #000000;">m_print</span><span style="color: #0000FF;">(</span><span style="color: #000000;">ct</span><span style="color: #0000FF;">)</span>
end function
<span style="color: #000080;font-style:italic;">-- note: rounding similar to that in m_unitary may be rqd (in a similar
 
-- loop in a new m_equal function) on these two equality tests,
function m_mul(sequence a, sequence b)
-- but as it is, all tests pass with the builtin = operator.</span>
sequence res = sq_mul(a,0)
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Hermitian?: %t\n"</span><span style="color: #0000FF;">,</span><span style="color: #000000;">a</span><span style="color: #0000FF;">=</span><span style="color: #000000;">ct</span><span style="color: #0000FF;">)</span> <span style="color: #000080;font-style:italic;">-- (this one)</span>
integer l = length(a)
<span style="color: #004080;">sequence</span> <span style="color: #000000;">act</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">m_mul</span><span style="color: #0000FF;">(</span><span style="color: #000000;">a</span><span style="color: #0000FF;">,</span><span style="color: #000000;">ct</span><span style="color: #0000FF;">),</span> <span style="color: #000000;">cta</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">m_mul</span><span style="color: #0000FF;">(</span><span style="color: #000000;">ct</span><span style="color: #0000FF;">,</span><span style="color: #000000;">a</span><span style="color: #0000FF;">)</span>
for i=1 to l do
<span style="color: #004080;">bool</span> <span style="color: #000000;">normal</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">act</span><span style="color: #0000FF;">=</span><span style="color: #000000;">cta</span> <span style="color: #000080;font-style:italic;">-- (&this one)</span>
for j=1 to l do
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Normal?: %t\n"</span><span style="color: #0000FF;">,</span><span style="color: #000000;">normal</span><span style="color: #0000FF;">)</span>
for k=1 to l do
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Unitary?: %t\n\n"</span><span style="color: #0000FF;">,</span><span style="color: #000000;">normal</span> <span style="color: #008080;">and</span> <span style="color: #000000;">m_unitary</span><span style="color: #0000FF;">(</span><span style="color: #000000;">act</span><span style="color: #0000FF;">))</span>
res[i][j] = c_add(res[i][j],c_mul(a[i][k],b[k][j]))
<span style="color: #008080;">end</span> <span style="color: #008080;">procedure</span>
end for
end for
<span style="color: #008080;">constant</span> <span style="color: #000000;">x</span> <span style="color: #0000FF;">=</span> <span style="color: #7060A8;">sqrt</span><span style="color: #0000FF;">(</span><span style="color: #000000;">2</span><span style="color: #0000FF;">)/</span><span style="color: #000000;">2</span><span style="color: #0000FF;">,</span>
end for
<span style="color: #000000;">tests</span> <span style="color: #0000FF;">=</span> <span style="color: #0000FF;">{{{{</span><span style="color: #000000;">3</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">2</span><span style="color: #0000FF;">,</span><span style="color: #000000;">1</span><span style="color: #0000FF;">}},</span>
return res
<span style="color: #0000FF;">{{</span><span style="color: #000000;">2</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #000000;">0</span><span style="color: #0000FF;">}}},</span>
end function
 
<span style="color: #0000FF;">{{{</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">},{</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">2</span><span style="color: #0000FF;">}},</span>
procedure test(sequence a)
<span style="color: #0000FF;">{{</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">1</span><span style="color: #0000FF;">},{</span> <span style="color: #000000;">5</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{-</span><span style="color: #000000;">3</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">}},</span>
sequence ct = conjugate_transpose(a)
<span style="color: #0000FF;">{{</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">2</span><span style="color: #0000FF;">},{-</span><span style="color: #000000;">3</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">}}},</span>
printf(1,"Original matrix:\n")
m_print(a)
<span style="color: #0000FF;">{{{</span><span style="color: #000000;">0.5</span><span style="color: #0000FF;">,+</span><span style="color: #000000;">0.5</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">0.5</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">0.5</span><span style="color: #0000FF;">}},</span>
printf(1,"Conjugate transpose:\n")
<span style="color: #0000FF;">{{</span><span style="color: #000000;">0.5</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">0.5</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">0.5</span><span style="color: #0000FF;">,+</span><span style="color: #000000;">0.5</span><span style="color: #0000FF;">}}},</span>
m_print(ct)
-- note: rounding similar to that in m_unitary may be rqd (in a similar
<span style="color: #0000FF;">{{{</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">}},</span>
-- loop in a new m_equal function) on these two equality tests,
<span style="color: #0000FF;">{{</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">}},</span>
-- but as it is, all tests pass with the builtin = operator.
<span style="color: #0000FF;">{{</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">}}},</span>
printf(1,"Hermitian?: %s\n",{iff(a=ct?"TRUE":"FALSE")}) -- (this one)
sequence act = m_mul(a,ct), cta = m_mul(ct,a)
<span style="color: #0000FF;">{{{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">x</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">}},</span>
bool normal = act=cta -- (&this one)
<span style="color: #0000FF;">{{</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">x</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">x</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">}},</span>
printf(1,"Normal?: %s\n",{iff(normal?"TRUE":"FALSE")})
<span style="color: #0000FF;">{{</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">0</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">0</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">1</span><span style="color: #0000FF;">}}},</span>
printf(1,"Unitary?: %s\n\n",{iff(normal and m_unitary(act)?"TRUE":"FALSE")})
end procedure
<span style="color: #0000FF;">{{{</span><span style="color: #000000;">2</span><span style="color: #0000FF;">,</span><span style="color: #000000;">7</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">9</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">5</span><span style="color: #0000FF;">}},</span>
 
<span style="color: #0000FF;">{{</span><span style="color: #000000;">3</span><span style="color: #0000FF;">,</span><span style="color: #000000;">4</span><span style="color: #0000FF;">},{</span><span style="color: #000000;">8</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">6</span><span style="color: #0000FF;">}}}}</span>
constant x = sqrt(2)/2
 
<span style="color: #7060A8;">papply</span><span style="color: #0000FF;">(</span><span style="color: #000000;">tests</span><span style="color: #0000FF;">,</span><span style="color: #000000;">test</span><span style="color: #0000FF;">)</span>
constant tests = {{{{3, 0},{2,1}},
<!--</syntaxhighlight>-->
{{2,-1},{1,0}}},
 
{{{ 1, 0},{ 1, 1},{ 0, 2}},
{{ 1,-1},{ 5, 0},{-3, 0}},
{{ 0,-2},{-3, 0},{ 0, 0}}},
 
{{{0.5,+0.5},{0.5,-0.5}},
{{0.5,-0.5},{0.5,+0.5}}},
 
{{{ 1, 0},{ 1, 0},{ 0, 0}},
{{ 0, 0},{ 1, 0},{ 1, 0}},
{{ 1, 0},{ 0, 0},{ 1, 0}}},
 
{{{x, 0},{x, 0},{0, 0}},
{{0,-x},{0, x},{0, 0}},
{{0, 0},{0, 0},{0, 1}}},
 
{{{2,7},{9,-5}},
{{3,4},{8,-6}}}}
 
for i=1 to length(tests) do test(tests[i]) end for</lang>
{{out}}
<pre>
Original matrix:
[3,2+1ii]
[2-1ii,1]
Conjugate transpose:
[3,2+1ii]
[2-1ii,1]
Hermitian?: TRUEtrue
Normal?: TRUEtrue
Unitary?: FALSEfalse
 
Original matrix:
[1,1+1ii,0+2i]
[1-1ii,5,-3]
[0-2i,-3,0]
Conjugate transpose:
[1,1+1ii,0+2i]
[1-1ii,5,-3]
[0-2i,-3,0]
Hermitian?: TRUEtrue
Normal?: TRUEtrue
Unitary?: FALSEfalse
 
Original matrix:
Line 1,697 ⟶ 2,852:
[0.5-0.5i,0.5+0.5i]
[0.5+0.5i,0.5-0.5i]
Hermitian?: FALSEfalse
Normal?: TRUEtrue
Unitary?: TRUEtrue
 
Original matrix:
Line 1,709 ⟶ 2,864:
[1,1,0]
[0,1,1]
Hermitian?: FALSEfalse
Normal?: TRUEtrue
Unitary?: FALSEfalse
 
Original matrix:
[0.707107,0.707107,0]
[0-0.707107i,0+0.707107i,0]
[0,0,0+1ii]
Conjugate transpose:
[0.707107,0+0.707107i,0]
[0.707107,0-0.707107i,0]
[0,0,0-1ii]
Hermitian?: FALSEfalse
Normal?: TRUEtrue
Unitary?: TRUEtrue
 
Original matrix:
Line 1,731 ⟶ 2,886:
[2-7i,3-4i]
[9+5i,8+6i]
Hermitian?: FALSEfalse
Normal?: FALSEfalse
Unitary?: FALSEfalse
</pre>
 
=={{header|PL/I}}==
<syntaxhighlight lang="pl/i">
<lang PL/I>
test: procedure options (main); /* 1 October 2012 */
declare n fixed binary;
Line 1,787 ⟶ 2,942:
end MMULT;
end test;
</syntaxhighlight>
</lang>
Outputs from separate runs:
<pre>
Line 1,831 ⟶ 2,986:
 
=={{header|PowerShell}}==
<syntaxhighlight lang="powershell">
<lang PowerShell>
function conjugate-transpose($a) {
$arr = @()
Line 1,908 ⟶ 3,063:
"Normal? `$m = $(are-eq $mhm $hmm)"
"Unitary? `$m = $((are-eq $id2 $hmm) -and (are-eq $id2 $mhm))"
</syntaxhighlight>
</lang>
<b>Output:</b>
<pre>
Line 1,934 ⟶ 3,089:
=={{header|Python}}==
Internally, matrices must be represented as rectangular tuples of tuples of complex numbers.
<langsyntaxhighlight lang="python">def conjugate_transpose(m):
return tuple(tuple(n.conjugate() for n in row) for row in zip(*m))
 
Line 2,005 ⟶ 3,160:
print('Hermitian? %s.' % ishermitian(matrix, ct))
print('Normal? %s.' % isnormal(matrix, ct))
print('Unitary? %s.' % isunitary(matrix, ct))</langsyntaxhighlight>
 
{{out}}
Line 2,043 ⟶ 3,198:
 
=={{header|Racket}}==
<langsyntaxhighlight lang="racket">
#lang racket
(require math)
Line 2,060 ⟶ 3,215:
(define (hermitian? M)
(equal? (H M) M))
</syntaxhighlight>
</lang>
Test:
<langsyntaxhighlight lang="racket">
(define M (matrix [[3.000+0.000i +2.000+1.000i]
[2.000-1.000i +1.000+0.000i]]))
Line 2,069 ⟶ 3,224:
(unitary? M)
(hermitian? M)
</syntaxhighlight>
</lang>
Output:
<langsyntaxhighlight lang="racket">
(array #[#[3.0-0.0i 2.0+1.0i] #[2.0-1.0i 1.0-0.0i]])
#t
#f
#f
</syntaxhighlight>
</lang>
 
=={{header|Raku}}==
(formerly Perl 6)
{{works with|Rakudo|2015-12-13}}
<syntaxhighlight lang="raku" perl6line>for [ # Test Matrices
[ 1, 1+i, 2i],
[ 1-i, 5, -3],
Line 2,132 ⟶ 3,287:
}
 
sub say-it (@array) { $_».fmt("%9s").say for @array }</langsyntaxhighlight>
{{out}}
<pre>Matrix:
Line 2,175 ⟶ 3,330:
 
=={{header|REXX}}==
<langsyntaxhighlight lang="rexx">/*REXX program performs a conjugate transpose on a complex square matrix. */
parse arg N elements; if N==''|N=="," then N=3 /*Not specified? Then use the default.*/
k= 0; do r=1 for N
do c=1 for N; k= k+1; M.r.c= word( word(elements, k) 1,1); end /*c*/1)
end /*rc*/
end /*r*/
call showCmat 'M' ,N /*display a nicely formatted matrix. */
identity.= 0; do d=1 for N; identity.d.d= 1; end /*d*/
call conjCmat 'MH', "M" ,N /*conjugate the M matrix ───► MH */
call showCmat 'MH' ,N /*display a nicely formatted matrix. */
say 'M is Hermitian: ' word('no yes', isHermitian('M', "MH", N) + 1)
call multCmat 'M', 'MH', 'MMH', N /*multiple the two matrices together. */
call multCmat 'MH', 'M', 'MHM', N /* " " " " " */
Line 2,196 ⟶ 3,352:
rP: procedure; parse arg r ','; return word( r 0, 1) /*◄──maybe return a 0 ↑ */
/*──────────────────────────────────────────────────────────────────────────────────────*/
conjCmat: parse arg matX,matY,rows 1 cols; call normCmat matY, rows
do r=1 for rows; _=
do c=1 for cols; v= value(matY'.'r"."c)
rP= rP(v); cP= -cP(v); call value matX'.'c"."r, rP','cP
end /*c*/
end /*r*/; end /*r*/ return
return
/*──────────────────────────────────────────────────────────────────────────────────────*/
isHermitian: parse arg matX,matY,rows 1 cols; call normCmat matX, rows
call normCmat matY, rows
do r=1 for rows; _=
do c=1 for cols
if value(matX'.'r"."c) \= value(matY'.'r"."c) then return 0
end /*c*/
end /*r*/; end /*r*/ return 1
return 1
/*──────────────────────────────────────────────────────────────────────────────────────*/
isUnary: parse arg matX,rows 1 cols
do r=1 for rows; _=
do c=1 for cols; z= value(matX'.'r"."c); rP= rP(z); cP= cP(z)
if abs( sqrt( rP(z) **2 + cP(z)**2) - (r==c)) >= .0001 then return 0
end /*c*/
end /*r*/; end /*r*/ return 1
return 1
/*──────────────────────────────────────────────────────────────────────────────────────*/
multCmat: parse arg matA,matB,matT,rows 1 cols; call value matT'.', 0
do r=1 for rows; _=
do c=1 for cols
do k=1 for cols; T= value(matT'.'r"."c); Tr= rP(T); Tc= cP(T)
A= value(matA'.'r"."k); Ar= rP(A); Ac= cP(A)
B= value(matB'.'k"."c); Br= rP(B); Bc= cP(B)
Pr= Ar*Br - Ac*Bc; Pc= Ac*Br + Ar*Bc; Tr= Tr+Pr; Tc= Tc+Pc
call value matT'.'r"."c,Tr','Tc
end /*k*/
end /*c*/
end /*r*/; end /*r*/ return
return
/*──────────────────────────────────────────────────────────────────────────────────────*/
normCmat: parse arg matN,rows 1 cols
do r=1 to rows; _=
do c=1 to cols; v= translate( value(matN'.'r"."c), , "IiJj")
parse upper var v real ',' cplx
if real\=='' then real= real / 1
if cplx\=='' then cplx= cplx / 1; if cplx=0 then cplx=
if cplx\=='' then cplx= cplx"j"
call value matN'.'r"."c, strip(real','cplx, "T", ',')
end /*c*/
end /*r*/; end /*r*/ return
return
/*──────────────────────────────────────────────────────────────────────────────────────*/
showCmat: parse arg matX,rows,cols; if cols=='' then cols= rows; @@= left('', 6)
say; say center('matrix' matX, 79, '─'); call normCmat matX, rows, cols
do r=1 to rows; _=
do c=1 to cols; _= _ @@ left( value(matX'.'r"."c), 9); end /*c*/
end say _/*c*/
say end /*r*/_
end /*r*/; say; return
/*──────────────────────────────────────────────────────────────────────────────────────*/
sqrt: procedure; parse arg x; if x=0 then return 0; d=digits(); numeric form; h=d+6
numeric digits; parse value format(x,2,1,,0) 'E0' with g 'E' _ .; g=g *.5'e'_ % 2
m.=9; do j=0 while h>9; m.j=h; h=h%2+1; end /*j*/
do k=j+5 to 0 by -1; numeric digits m.k; g=(g+x/g)*.5; end /*k*/; return g</langsyntaxhighlight>
'''{{out|output''' |text=&nbsp; when using the default input:}}
<pre>
───────────────────────────────────matrix M────────────────────────────────────
Line 2,277 ⟶ 3,428:
MHM is Unary: no
</pre>
 
'''output''' &nbsp; when using the input of: &nbsp; <tt> 3 &nbsp; .7071 &nbsp; .7071 &nbsp; 0 &nbsp; 0,.7071 &nbsp; 0,-.7071 &nbsp; 0 &nbsp; 0 &nbsp; 0 &nbsp; 0,1 </tt>
{{out|output|text=&nbsp; when using the input of: &nbsp; &nbsp; <tt> 3 &nbsp; .7071 &nbsp; .7071 &nbsp; 0 &nbsp; 0,.7071 &nbsp; 0,-.7071 &nbsp; 0 &nbsp; 0 &nbsp; 0 &nbsp; 0,1 </tt>}}
<pre>
───────────────────────────────────matrix M────────────────────────────────────
Line 2,295 ⟶ 3,447:
MMH is Unary: yes
MHM is Unary: yes
</pre>
 
=={{header|RPL}}==
Although basic, RPL's matrix handling capabilities help to keep the code compact but still readable.
{{works with|Halcyon Calc|4.2.7}}
≪ - ABS 1E-10 < ≫
´SAME?´ STO
≪ DUP TRN → m mh
≪ m mh SAME? "Hermitian. " "" IFTE
m mh * mh m * SAME? "Normal. " "" IFTE +
m INV mh SAME? "Unitary. " "" IFTE +
´CNJTRN’ STO
[[(3,0) (2,1)][(2,-1) (1,0)]]
'Hm' STO
[[1 1 0][0 1 1][1 0 1]]
'Nm' STO
≪ [[1 1 0][(0,1) (0,-1) 0][0 0 0]] 2 √ 2 / * {3 3} (0,1) PUT ≫
'Um' STO
Hm CNJTRN
Nm CNJTRN
Um CNJTRN
{{out}}
<pre>
3: "Hermitian. Normal. "
2: "Normal. "
1: "Normal. Unitary. "
</pre>
 
=={{header|Ruby}}==
{{works with|Ruby|2.0}}
<langsyntaxhighlight lang="ruby">require 'matrix'
 
# Start with some matrix.
Line 2,323 ⟶ 3,505:
print ' normal? false'
print ' unitary? false'
end</langsyntaxhighlight>
Note: Ruby 1.9 had a bug in the Matrix#hermitian? method. It's fixed in 2.0.
 
=={{header|Rust}}==
Uses external crate 'num', version 0.1.34
<langsyntaxhighlight lang="rust">
extern crate num; // crate for complex numbers
 
Line 2,408 ⟶ 3,590:
println!("Unitary?: FALSE");
}
}</langsyntaxhighlight>
Output:
<pre>
Line 2,432 ⟶ 3,614:
 
=={{header|Scala}}==
<langsyntaxhighlight Scalalang="scala">object ConjugateTranspose {
case class Complex(re: Double, im: Double) {
Line 2,508 ⟶ 3,690:
}
}</langsyntaxhighlight>
{{out}}
<pre>
Line 2,528 ⟶ 3,710:
=={{header|Sidef}}==
{{trans|Raku}}
<langsyntaxhighlight lang="ruby">func is_Hermitian (Array m, Array t) -> Bool { m == t }
 
func mat_mult (Array a, Array b, Number ε = -3) {
Line 2,593 ⟶ 3,775:
say "Is Normal?\t#{is_Normal(m, t)}"
say "Is Unitary?\t#{is_Unitary(m, t)}"
}</langsyntaxhighlight>
{{out}}
<pre>
Line 2,639 ⟶ 3,821:
Sparkling has support for basic complex algebraic operations, but complex matrix operations are not in the standard library.
 
<langsyntaxhighlight lang="sparkling"># Computes conjugate transpose of M
let conjTransp = function conjTransp(M) {
return map(range(sizeof M[0]), function(row) {
Line 2,739 ⟶ 3,921:
print("U x U* = ");
printCplxMat(cplxMatMul(U, conjTransp(U)));
print();</langsyntaxhighlight>
 
=={{header|Stata}}==
In Mata, the ' operator is always the conjugate transpose. To get only the matrix transpose without complex conjugate, use the [ transposeonly] function.
 
<langsyntaxhighlight lang="stata">
: a=1,2i\3i,4
 
Line 2,776 ⟶ 3,958:
: a'*a==I(rows(a))
0
</syntaxhighlight>
</lang>
 
=={{header|Tcl}}==
Line 2,782 ⟶ 3,964:
{{tcllib|math::complexnumbers}}
{{tcllib|struct::matrix}}
<langsyntaxhighlight lang="tcl">package require struct::matrix
package require math::complexnumbers
 
Line 2,838 ⟶ 4,020:
}
return $mat
}</langsyntaxhighlight>
Using these tools to test for the properties described in the task:
<langsyntaxhighlight lang="tcl">proc isHermitian {matrix {epsilon 1e-14}} {
if {[$matrix rows] != [$matrix columns]} {
# Must be square!
Line 2,891 ⟶ 4,073:
$mmh destroy
return $result
}</langsyntaxhighlight>
<!-- Wot, no demonstration? -->
 
=={{header|Wren}}==
{{libheader|Wren-fmt}}
{{libheader|Wren-complex}}
Although the third example is in fact a unitary matrix, the ''isUnitary'' method of the above module returns false.
 
This is because the methods in the module work as accurately as they can within the confines of 64-bit floating point arithmetic and don't therefore allow for the small rounding error that occurs due to the use of the irrational number, sqrt(2).
 
However, if we use the ''almostEquals'' method with the default tolerance of 1.0e-14, then we do get a ''true'' result.
<syntaxhighlight lang="wren">import "./complex" for Complex, CMatrix
import "./fmt" for Fmt
 
var cm1 = CMatrix.new(
[
[Complex.new(3), Complex.new(2, 1)],
[Complex.new(2, -1), Complex.one ]
]
)
var cm2 = CMatrix.fromReals([ [1, 1, 0], [0, 1, 1], [1, 0, 1] ])
var x = 2.sqrt/2
var cm3 = CMatrix.new(
[
[Complex.new(x), Complex.new(x), Complex.zero],
[Complex.new(0, -x), Complex.new(0, x), Complex.zero],
[Complex.zero, Complex.zero, Complex.imagOne]
]
)
 
for (cm in [cm1, cm2, cm3]) {
System.print("Matrix:")
Fmt.mprint(cm, 5, 3)
System.print("\nConjugate transpose:")
Fmt.mprint(cm.conjTranspose, 5, 3)
System.print("\nHermitian : %(cm.isHermitian)")
System.print("Normal : %(cm.isNormal)")
System.print("Unitary : %(cm.isUnitary)")
System.print()
}
 
System.print("For the final example if we use a tolerance of 1e-14:")
var cm4 = cm3 * cm3.conjTranspose
var id = CMatrix.identity(3)
System.print("Unitary : %(cm4.almostEquals(id))")</syntaxhighlight>
 
{{out}}
<pre>
Matrix:
|3.000 + 0.000i 2.000 + 1.000i|
|2.000 - 1.000i 1.000 + 0.000i|
 
Conjugate transpose:
|3.000 + 0.000i 2.000 + 1.000i|
|2.000 - 1.000i 1.000 + 0.000i|
 
Hermitian : true
Normal : true
Unitary : false
 
Matrix:
|1.000 + 0.000i 1.000 + 0.000i 0.000 + 0.000i|
|0.000 + 0.000i 1.000 + 0.000i 1.000 + 0.000i|
|1.000 + 0.000i 0.000 + 0.000i 1.000 + 0.000i|
 
Conjugate transpose:
|1.000 + 0.000i 0.000 + 0.000i 1.000 + 0.000i|
|1.000 + 0.000i 1.000 + 0.000i 0.000 + 0.000i|
|0.000 + 0.000i 1.000 + 0.000i 1.000 + 0.000i|
 
Hermitian : false
Normal : true
Unitary : false
 
Matrix:
|0.707 + 0.000i 0.707 + 0.000i 0.000 + 0.000i|
|0.000 - 0.707i 0.000 + 0.707i 0.000 + 0.000i|
|0.000 + 0.000i 0.000 + 0.000i 0.000 + 1.000i|
 
Conjugate transpose:
|0.707 + 0.000i 0.000 + 0.707i 0.000 + 0.000i|
|0.707 + 0.000i 0.000 - 0.707i 0.000 + 0.000i|
|0.000 + 0.000i 0.000 + 0.000i 0.000 - 1.000i|
 
Hermitian : false
Normal : true
Unitary : false
 
For the final example if we use a tolerance of 1e-14:
Unitary : true
</pre>
9,476

edits