# Conjugate transpose

Conjugate transpose
You are encouraged to solve this task according to the task description, using any language you may know.

Suppose that a matrix ${\displaystyle M}$ contains complex numbers. Then the conjugate transpose of ${\displaystyle M}$ is a matrix ${\displaystyle M^{H}}$ containing the complex conjugates of the matrix transposition of ${\displaystyle M}$.

${\displaystyle (M^{H})_{ji}={\overline {M_{ij}}}}$

This means that row ${\displaystyle j}$, column ${\displaystyle i}$ of the conjugate transpose equals the
complex conjugate of row ${\displaystyle i}$, column ${\displaystyle j}$ of the original matrix.

In the next list, ${\displaystyle M}$ must also be a square matrix.

• A Hermitian matrix equals its own conjugate transpose: ${\displaystyle M^{H}=M}$.
• A normal matrix is commutative in multiplication with its conjugate transpose: ${\displaystyle M^{H}M=MM^{H}}$.
• A unitary matrix has its inverse equal to its conjugate transpose: ${\displaystyle M^{H}=M^{-1}}$.
This is true iff ${\displaystyle M^{H}M=I_{n}}$ and iff ${\displaystyle MM^{H}=I_{n}}$, where ${\displaystyle I_{n}}$ is the identity matrix.

Given some matrix of complex numbers, find its conjugate transpose.

Also determine if the matrix is a:

• Hermitian matrix,
• normal matrix, or
• unitary matrix.

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;with Ada.Numerics.Complex_Arrays; use Ada.Numerics.Complex_Arrays;procedure ConTrans is   subtype CM is Complex_Matrix;   S2O2 : constant Float := 0.7071067811865;    procedure Print (mat : CM) is begin      for row in mat'Range(1) loop for col in mat'Range(2) loop         Put(mat(row,col), Exp=>0, Aft=>4);      end loop; New_Line; end loop;   end Print;    function almostzero(mat : CM; tol : Float) return Boolean is begin      for row in mat'Range(1) loop for col in mat'Range(2) loop         if abs(mat(row,col)) > tol then return False; end if;      end loop; end loop;      return True;   end almostzero;    procedure Examine (mat : CM) is      CT : CM := Conjugate (Transpose(mat));      isherm, isnorm, isunit : Boolean;   begin      isherm := almostzero(mat-CT, 1.0e-6);      isnorm := almostzero(mat*CT-CT*mat, 1.0e-6);      isunit := almostzero(CT-Inverse(mat), 1.0e-6);      Print(mat);      Put_Line("Conjugate transpose:"); Print(CT);      Put_Line("Hermitian?: " & isherm'Img);      Put_Line("Normal?: " & isnorm'Img);      Put_Line("Unitary?: " & isunit'Img);   end Examine;    hmat : CM := ((3.0+0.0*i, 2.0+1.0*i), (2.0-1.0*i, 1.0+0.0*i));   nmat : CM := ((1.0+0.0*i, 1.0+0.0*i, 0.0+0.0*i),                 (0.0+0.0*i, 1.0+0.0*i, 1.0+0.0*i),                 (1.0+0.0*i, 0.0+0.0*i, 1.0+0.0*i));   umat : CM := ((S2O2+0.0*i, S2O2+0.0*i, 0.0+0.0*i),                 (0.0+S2O2*i, 0.0-S2O2*i, 0.0+0.0*i),                 (0.0+0.0*i, 0.0+0.0*i, 0.0+1.0*i));begin   Put_Line("hmat:"); Examine(hmat); New_Line;   Put_Line("nmat:"); Examine(nmat); New_Line;   Put_Line("umat:"); Examine(umat);end ConTrans;
Output:
hmat:
( 3.0000, 0.0000)( 2.0000, 1.0000)
( 2.0000,-1.0000)( 1.0000, 0.0000)
Conjugate transpose:
( 3.0000,-0.0000)( 2.0000, 1.0000)
( 2.0000,-1.0000)( 1.0000,-0.0000)
Hermitian?: TRUE
Normal?: TRUE
Unitary?: FALSE

nmat:
( 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 transpose:
( 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

umat:
( 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 transpose:
( 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

## C

/* Uses C99 specified complex.h, complex datatype has to be defined and operation provided if used on non-C99 compilers */ #include<stdlib.h>#include<stdio.h>#include<complex.h> typedef struct{  int rows, cols;  complex **z;} matrix; matrixtranspose (matrix a){  int i, j;  matrix b;   b.rows = a.cols;  b.cols = a.rows;   b.z = malloc (b.rows * sizeof (complex *));   for (i = 0; i < b.rows; i++)    {      b.z[i] = malloc (b.cols * sizeof (complex));      for (j = 0; j < b.cols; j++)        {          b.z[i][j] = conj (a.z[j][i]);        }    }   return b;} intisHermitian (matrix a){  int i, j;  matrix b = transpose (a);   if (b.rows == a.rows && b.cols == a.cols)    {      for (i = 0; i < b.rows; i++)        {          for (j = 0; j < b.cols; j++)            {              if (b.z[i][j] != a.z[i][j])                return 0;            }        }    }   else    return 0;   return 1;} matrixmultiply (matrix a, matrix b){  matrix c;  int i, j;   if (a.cols == b.rows)    {      c.rows = a.rows;      c.cols = b.cols;       c.z = malloc (c.rows * (sizeof (complex *)));       for (i = 0; i < c.rows; i++)        {          c.z[i] = malloc (c.cols * sizeof (complex));          c.z[i][j] = 0 + 0 * I;          for (j = 0; j < b.cols; j++)            {              c.z[i][j] += a.z[i][j] * b.z[j][i];            }        }     }   return c;} intisNormal (matrix a){  int i, j;  matrix a_ah, ah_a;   if (a.rows != a.cols)    return 0;   a_ah = multiply (a, transpose (a));  ah_a = multiply (transpose (a), a);   for (i = 0; i < a.rows; i++)    {      for (j = 0; j < a.cols; j++)        {          if (a_ah.z[i][j] != ah_a.z[i][j])            return 0;        }    }   return 1;} intisUnitary (matrix a){  matrix b;  int i, j;  if (isNormal (a) == 1)    {      b = multiply (a, transpose(a));       for (i = 0; i < b.rows; i++)        {          for (j = 0; j < b.cols; j++)            {              if ((i == j && b.z[i][j] != 1) || (i != j && b.z[i][j] != 0))                return 0;            }        }      return 1;    }  return 0;}  intmain (){  complex z = 3 + 4 * I;  matrix a, aT;  int i, j;  printf ("Enter rows and columns :");  scanf ("%d%d", &a.rows, &a.cols);   a.z = malloc (a.rows * sizeof (complex *));  printf ("Randomly Generated Complex Matrix A is : ");  for (i = 0; i < a.rows; i++)    {      printf ("\n");      a.z[i] = malloc (a.cols * sizeof (complex));      for (j = 0; j < a.cols; j++)        {          a.z[i][j] = rand () % 10 + rand () % 10 * I;          printf ("\t%f + %fi", creal (a.z[i][j]), cimag (a.z[i][j]));        }    }   aT = transpose (a);   printf ("\n\nTranspose of Complex Matrix A is : ");  for (i = 0; i < aT.rows; i++)    {      printf ("\n");      aT.z[i] = malloc (aT.cols * sizeof (complex));      for (j = 0; j < aT.cols; j++)        {          aT.z[i][j] = rand () % 10 + rand () % 10 * I;          printf ("\t%f + %fi", creal (aT.z[i][j]), cimag (aT.z[i][j]));        }    }   printf ("\n\nComplex Matrix A %s hermitian",          isHermitian (a) == 1 ? "is" : "is not");  printf ("\n\nComplex Matrix A %s unitary",          isUnitary (a) == 1 ? "is" : "is not");  printf ("\n\nComplex Matrix A %s normal",          isNormal (a) == 1 ? "is" : "is not");     return 0;}
Output:
Enter rows and columns :3 3
Randomly Generated Complex Matrix A is :
3.000000 + 6.000000i    7.000000 + 5.000000i    3.000000 + 5.000000i
6.000000 + 2.000000i    9.000000 + 1.000000i    2.000000 + 7.000000i
0.000000 + 9.000000i    3.000000 + 6.000000i    0.000000 + 6.000000i

Transpose of Complex Matrix A is :
2.000000 + 6.000000i    1.000000 + 8.000000i    7.000000 + 9.000000i
2.000000 + 0.000000i    2.000000 + 3.000000i    7.000000 + 5.000000i
9.000000 + 2.000000i    2.000000 + 8.000000i    9.000000 + 7.000000i

Complex Matrix A is not hermitian

Complex Matrix A is not unitary

Complex Matrix A is not normal

## Common Lisp

 (defun matrix-multiply (m1 m2) (mapcar  (lambda (row)   (apply #'mapcar    (lambda (&rest column)     (apply #'+ (mapcar #'* row column))) m2)) m1)) (defun identity-p (m &optional (tolerance 1e-6)) "Is m an identity matrix?"  (loop for row in m    for r = 1 then (1+ r) do      (loop for col in row        for c = 1 then (1+ c) do          (if (eql r c)            (unless (< (abs (- col 1)) tolerance) (return-from identity-p nil))            (unless (< (abs col) tolerance) (return-from identity-p nil)) )))  T ) (defun conjugate-transpose (m)  (apply #'mapcar #'list (mapcar #'(lambda (r) (mapcar #'conjugate r)) m)) ) (defun hermitian-p (m)  (equalp m (conjugate-transpose m))) (defun normal-p (m)  (let ((m* (conjugate-transpose m)))    (equalp (matrix-multiply m m*) (matrix-multiply m* m)) )) (defun unitary-p (m)  (identity-p (matrix-multiply m (conjugate-transpose m))) )
Output:
(hermitian-p
'((3        #C(2 1))
(#C(2 -1) 1) ))
=> T

(normal-p
'((#C(0 1) 0)
(0       #C(3 -5)) ))
==> T

(unitary-p
'((0.70710677        0.70710677       0)
(#C(0 -0.70710677) #C(0 0.70710677) 0)
(0                 0                1) ))
==> T


## D

Translation of: Python
A well typed and mostly imperative version:
import std.stdio, std.complex, std.math, std.range, std.algorithm,       std.numeric; T[][] conjugateTranspose(T)(in T[][] m) pure nothrow @safe {    auto r = new typeof(return)(m[0].length, m.length);    foreach (immutable nr, const row; m)        foreach (immutable nc, immutable c; row)            r[nc][nr] = c.conj;    return r;} bool isRectangular(T)(in T[][] M) pure nothrow @safe @nogc {    return M.all!(row => row.length == M[0].length);} T[][] matMul(T)(in T[][] A, in T[][] B) pure nothrow /*@safe*/in {    assert(A.isRectangular && B.isRectangular &&           !A.empty && !B.empty && A[0].length == B.length);} body {    auto result = new T[][](A.length, B[0].length);    auto aux = new T[B.length];     foreach (immutable j; 0 .. B[0].length) {        foreach (immutable k, const row; B)            aux[k] = row[j];        foreach (immutable i, const ai; A)            result[i][j] = dotProduct(ai, aux);    }     return result;} /// Check any number of complex matrices for equality within/// some bits of mantissa.bool areEqual(T)(in Complex!T[][][] matrices, in size_t nBits=20)pure nothrow /*@safe*/ {    static bool allSame(U)(in U[] v) pure nothrow @nogc {        return v[1 .. $].all!(c => c == v[0]); } bool allNearSame(in Complex!T[] v) pure nothrow @nogc { auto v0 = v[0].Complex!T; // To avoid another cast. return v[1 ..$].all!(c => feqrel(v0.re, c.re) >= nBits &&                                   feqrel(v0.im, c.im) >= nBits);    }     immutable x = matrices.map!(m => m.length).array;    if (!allSame(x))        return false;    immutable y = matrices.map!(m => m[0].length).array;    if (!allSame(y))        return false;    foreach (immutable s; 0 .. x[0])        foreach (immutable t; 0 .. y[0])            if (!allNearSame(matrices.map!(m => m[s][t]).array))                return false;    return true;} bool isHermitian(T)(in Complex!T[][] m, in Complex!T[][] ct)pure nothrow /*@safe*/ {    return [m, ct].areEqual;} bool isNormal(T)(in Complex!T[][] m, in Complex!T[][] ct)pure nothrow /*@safe*/ {    return [matMul(m, ct), matMul(ct, m)].areEqual;} auto complexIdentitymatrix(in size_t side) pure nothrow /*@safe*/ {    return side.iota.map!(r => side.iota.map!(c => complex(r == c)).array).array;} bool isUnitary(T)(in Complex!T[][] m, in Complex!T[][] ct)pure nothrow /*@safe*/ {    immutable mct = matMul(m, ct);    immutable ident = mct.length.complexIdentitymatrix;    return [mct, matMul(ct, m), ident].areEqual;} void main() /*@safe*/ {    alias C = complex;    immutable x = 2 ^^ 0.5 / 2;     immutable data = [[[C(3.0,  0.0), C(2.0, 1.0)],                       [C(2.0, -1.0), C(1.0, 0.0)]],                       [[C(1.0, 0.0), C(1.0, 0.0), C(0.0, 0.0)],                       [C(0.0, 0.0), C(1.0, 0.0), C(1.0, 0.0)],                       [C(1.0, 0.0), C(0.0, 0.0), C(1.0, 0.0)]],                       [[C(x,    0.0), C(x,   0.0), C(0.0, 0.0)],                       [C(0.0, -x),   C(0.0, x),   C(0.0, 0.0)],                       [C(0.0,  0.0), C(0.0, 0.0), C(0.0, 1.0)]]];     foreach (immutable mat; data) {        enum mFormat = "[%([%(%1.3f, %)],\n %)]]";        writefln("Matrix:\n" ~ mFormat, mat);        immutable ct = conjugateTranspose(mat);        "Its conjugate transpose:".writeln;        writefln(mFormat, ct);        writefln("Hermitian? %s.", isHermitian(mat, ct));        writefln("Normal?    %s.", isNormal(mat, ct));        writefln("Unitary?   %s.\n", isUnitary(mat, ct));    }}
Output:
Matrix:
[[3.000+0.000i, 2.000+1.000i],
[2.000-1.000i, 1.000+0.000i]]
Its 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]]
Its 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]]
Its 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?   true.


### Alternative Version

A more functional version that contains some typing problems (same output).

import std.stdio, std.complex, std.math, std.range, std.algorithm,       std.numeric, std.exception, std.traits; // alias CM(T) = Complex!T[][]; // Not yet useful. auto conjugateTranspose(T)(in Complex!T[][] m) pure nothrow /*@safe*/if (!hasIndirections!T) {    return iota(m[0].length).map!(i => m.transversal(i).map!conj.array).array;} T[][] matMul(T)(immutable T[][] A, immutable T[][] B) pure nothrow /*@safe*/ {    immutable Bt = B[0].length.iota.map!(i => B.transversal(i).array).array;    return A.map!(a => Bt.map!(b => a.dotProduct(b)).array).array;} /// Check any number of complex matrices for equality within/// some bits of mantissa.bool areEqual(T)(in Complex!T[][][] matrices, in size_t nBits=20)pure nothrow /*@safe*/ {    static bool allSame(U)(in U[] v) pure nothrow @nogc @safe {        return v[1 .. $].all!(c => c == v[0]); } bool allNearSame(in Complex!T[] v) pure nothrow @nogc @safe { auto v0 = v[0].Complex!T; // To avoid another cast. return v[1 ..$].all!(c => feqrel(v0.re, c.re) >= nBits &&                                   feqrel(v0.im, c.im) >= nBits);    }     immutable x = matrices.map!(m => m.length).array;    if (!allSame(x))        return false;    immutable y = matrices.map!(m => m[0].length).array;    if (!allSame(y))        return false;    foreach (immutable s; 0 .. x[0])        foreach (immutable t; 0 .. y[0])            if (!allNearSame(matrices.map!(m => m[s][t]).array))                return false;    return true;} bool isHermitian(T)(in Complex!T[][] m, in Complex!T[][] ct)pure nothrow /*@safe*/ {    return [m, ct].areEqual;} bool isNormal(T)(immutable Complex!T[][] m, immutable Complex!T[][] ct)pure nothrow /*@safe*/ {    return [matMul(m, ct), matMul(ct, m)].areEqual;} auto complexIdentitymatrix(in size_t side) pure nothrow /*@safe*/ {    return side.iota.map!(r => side.iota.map!(c => complex(r == c)).array).array;} bool isUnitary(T)(immutable Complex!T[][] m, immutable Complex!T[][] ct)pure nothrow /*@safe*/ {    immutable mct = matMul(m, ct);    immutable ident = mct.length.complexIdentitymatrix;    return [mct, matMul(ct, m), ident].areEqual;} void main() {    alias C = complex;    immutable x = 2 ^^ 0.5 / 2;     foreach (/*immutable*/ const matrix;        [[[C(3.0,  0.0), C(2.0, 1.0)],          [C(2.0, -1.0), C(1.0, 0.0)]],          [[C(1.0, 0.0), C(1.0, 0.0), C(0.0, 0.0)],          [C(0.0, 0.0), C(1.0, 0.0), C(1.0, 0.0)],          [C(1.0, 0.0), C(0.0, 0.0), C(1.0, 0.0)]],          [[C(x,    0.0), C(x,   0.0), C(0.0, 0.0)],          [C(0.0, -x),   C(0.0, x),   C(0.0, 0.0)],          [C(0.0,  0.0), C(0.0, 0.0), C(0.0, 1.0)]]]) {        immutable mat = matrix.assumeUnique; //*         enum mFormat = "[%([%(%1.3f, %)],\n %)]]";        writefln("Matrix:\n" ~ mFormat, mat);        immutable ct = conjugateTranspose(mat);        "Its conjugate transpose:".writeln;        writefln(mFormat, ct);        writefln("Hermitian? %s.", isHermitian(mat, ct));        writefln("Normal?    %s.", isNormal(mat, ct));        writefln("Unitary?   %s.\n", isUnitary(mat, ct));    }}

## Factor

Before the fix to Factor bug #484, m. gave the wrong answer and this code failed. Factor 0.94 is too old to work.

Works with: Factor version development (future 0.95)
USING: kernel math.functions math.matrices sequences ;IN: rosetta.hermitian : conj-t ( matrix -- conjugate-transpose )    flip [ [ conjugate ] map ] map ; : hermitian-matrix? ( matrix -- ? )    dup conj-t = ; : normal-matrix? ( matrix -- ? )    dup conj-t [ m. ] [ swap m. ] 2bi = ; : unitary-matrix? ( matrix -- ? )    [ dup conj-t m. ] [ length identity-matrix ] bi = ;

Usage:

USE: rosetta.hermitian
IN: scratchpad { { C{ 1 2 } 0 }
{ 0 C{ 3 4 } } }
[ hermitian-matrix? . ]
[ normal-matrix? . ]
[ unitary-matrix? . ] tri
f
t
f


## Fortran

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 :
gfortran -std=f2008 -Wall -fopenmp -ffree-form -fall-intrinsics -fimplicit-none f.f08 -o f
 program conjugate_transpose   complex, dimension(3, 3) :: a  integer :: i  a = reshape((/ (i, i=1,9) /), shape(a))  call characterize(a)  a(:2,:2) = reshape((/cmplx(3,0),cmplx(2,-1),cmplx(2,1),cmplx(1,0)/),(/2,2/))  call characterize(a(:2,:2))  call characterize(cmplx(reshape((/1,0,1,1,1,0,0,1,1/),(/3,3/)),0))  a(3,:) = (/cmplx(0,0), cmplx(0,0), cmplx(0,1)/)*sqrt(2.0)  a(2,:) = (/cmplx(0,-1),cmplx(0,1),cmplx(0,0)/)  a(1,:) = (/1,1,0/)  a = a * sqrt(2.0)/2.0  call characterize(a) contains   subroutine characterize(a)    complex, dimension(:,:), intent(in) :: a    integer :: i, j    do i=1, size(a,1)       print *,(a(i, j), j=1,size(a,1))    end do    print *,'Is Hermitian?  ',HermitianQ(a)    print *,'Is normal?  ',NormalQ(a)    print *,'Unitary?  ',UnitaryQ(a)    print '(/)'  end subroutine characterize   function ct(a) result(b) ! return the conjugate transpose of a matrix    complex, dimension(:,:), intent(in) :: a    complex, dimension(size(a,1),size(a,1)) :: b    b = conjg(transpose(a))  end function ct   function identity(n) result(b) ! return identity matrix    integer, intent(in) :: n    real, dimension(n,n) :: b    integer :: i    b = 0    do i=1, n       b(i,i) = 1    end do  end function identity   logical function HermitianQ(a)    complex, dimension(:,:), intent(in) :: a    HermitianQ = all(a .eq. ct(a))  end function HermitianQ   logical function NormalQ(a)    complex, dimension(:,:), intent(in) :: a    NormalQ = all(matmul(ct(a),a) .eq. matmul(a,ct(a)))  end function NormalQ   logical function UnitaryQ(a)    ! if  A inverse equals A star    ! then multiplying each side by A should result in the identity matrix    ! Thus show that  A times A star  is sufficiently close to  I .    complex, dimension(:,:), intent(in) :: a    UnitaryQ = all(abs(matmul(a,ct(a)) - identity(size(a,1))) .lt. 1e-6)  end function UnitaryQ end program conjugate_transpose
-*- mode: compilation; default-directory: "/tmp/" -*-
Compilation started at Fri Jun  7 16:31:38

a=./f && make $a && time$a
gfortran -std=f2008 -Wall -fopenmp -ffree-form -fall-intrinsics -fimplicit-none f.f08 -o f
(  1.00000000    ,  0.00000000    ) (  4.00000000    ,  0.00000000    ) (  7.00000000    ,  0.00000000    )
(  2.00000000    ,  0.00000000    ) (  5.00000000    ,  0.00000000    ) (  8.00000000    ,  0.00000000    )
(  3.00000000    ,  0.00000000    ) (  6.00000000    ,  0.00000000    ) (  9.00000000    ,  0.00000000    )
Is Hermitian?   F
Is normal?   F
Unitary?   F

(  3.00000000    ,  0.00000000    ) (  2.00000000    ,  1.00000000    )
(  2.00000000    , -1.00000000    ) (  1.00000000    ,  0.00000000    )
Is Hermitian?   T
Is normal?   T
Unitary?   F

(  1.00000000    ,  0.00000000    ) (  1.00000000    ,  0.00000000    ) (  0.00000000    ,  0.00000000    )
(  0.00000000    ,  0.00000000    ) (  1.00000000    ,  0.00000000    ) (  1.00000000    ,  0.00000000    )
(  1.00000000    ,  0.00000000    ) (  0.00000000    ,  0.00000000    ) (  1.00000000    ,  0.00000000    )
Is Hermitian?   F
Is normal?   T
Unitary?   F

( 0.707106769    ,  0.00000000    ) ( 0.707106769    ,  0.00000000    ) (  0.00000000    ,  0.00000000    )
(  0.00000000    ,-0.707106769    ) (  0.00000000    , 0.707106769    ) (  0.00000000    ,  0.00000000    )
(  0.00000000    ,  0.00000000    ) (  0.00000000    ,  0.00000000    ) (  0.00000000    , 0.999999940    )
Is Hermitian?   F
Is normal?   T
Unitary?   T

real	0m0.002s
user	0m0.000s
sys	0m0.000s

Compilation finished at Fri Jun  7 16:31:38


## Go

package main import (    "fmt"    "math"    "math/cmplx") // a type to represent matricestype matrix struct {    ele  []complex128    cols int} // conjugate transpose, implemented here as a method on the matrix type.func (m *matrix) conjTranspose() *matrix {    r := &matrix{make([]complex128, len(m.ele)), len(m.ele) / m.cols}    rx := 0    for _, e := range m.ele {        r.ele[rx] = cmplx.Conj(e)        rx += r.cols        if rx >= len(r.ele) {            rx -= len(r.ele) - 1        }    }    return r} // program to demonstrate capabilites on example matriciesfunc main() {    show("h", matrixFromRows([][]complex128{        {3, 2 + 1i},        {2 - 1i, 1}}))     show("n", matrixFromRows([][]complex128{        {1, 1, 0},        {0, 1, 1},        {1, 0, 1}}))     show("u", matrixFromRows([][]complex128{        {math.Sqrt2 / 2, math.Sqrt2 / 2, 0},        {math.Sqrt2 / -2i, math.Sqrt2 / 2i, 0},        {0, 0, 1i}}))} func show(name string, m *matrix) {    m.print(name)    ct := m.conjTranspose()    ct.print(name + "_ct")     fmt.Println("Hermitian:", m.equal(ct, 1e-14))     mct := m.mult(ct)    ctm := ct.mult(m)    fmt.Println("Normal:", mct.equal(ctm, 1e-14))     i := eye(m.cols)    fmt.Println("Unitary:", mct.equal(i, 1e-14) && ctm.equal(i, 1e-14))} // two constructorsfunc matrixFromRows(rows [][]complex128) *matrix {    m := &matrix{make([]complex128, len(rows)*len(rows[0])), len(rows[0])}    for rx, row := range rows {        copy(m.ele[rx*m.cols:(rx+1)*m.cols], row)    }    return m} func eye(n int) *matrix {    r := &matrix{make([]complex128, n*n), n}    n++    for x := 0; x < len(r.ele); x += n {        r.ele[x] = 1    }    return r} // print method outputs matrix to stdoutfunc (m *matrix) print(heading string) {    fmt.Print("\n", heading, "\n")    for e := 0; e < len(m.ele); e += m.cols {        fmt.Printf("%6.3f ", m.ele[e:e+m.cols])        fmt.Println()    }} // equal method uses Îµ to allow for floating point error.func (a *matrix) equal(b *matrix, Îµ float64) bool {    for x, aEle := range a.ele {        if math.Abs(real(aEle)-real(b.ele[x])) > math.Abs(real(aEle))*Îµ ||            math.Abs(imag(aEle)-imag(b.ele[x])) > math.Abs(imag(aEle))*Îµ {            return false        }    }    return true} // mult method taken from matrix multiply taskfunc (m1 *matrix) mult(m2 *matrix) (m3 *matrix) {    m3 = &matrix{make([]complex128, (len(m1.ele)/m1.cols)*m2.cols), m2.cols}    for m1c0, m3x := 0, 0; m1c0 < len(m1.ele); m1c0 += m1.cols {        for m2r0 := 0; m2r0 < m2.cols; m2r0++ {            for m1x, m2x := m1c0, m2r0; m2x < len(m2.ele); m2x += m2.cols {                m3.ele[m3x] += m1.ele[m1x] * m2.ele[m2x]                m1x++            }            m3x++        }    }    return m3}

Output:

h
[( 3.000+0.000i) (+2.000+1.000i)]
[( 2.000-1.000i) (+1.000+0.000i)]

h_ct
[( 3.000-0.000i) (+2.000+1.000i)]
[( 2.000-1.000i) (+1.000-0.000i)]
Hermitian: true
Normal: true
Unitary: false

n
[( 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)]

n_ct
[( 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

u
[( 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)]

u_ct
[( 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: true


Slow implementation using lists.

import Data.List (transpose)import Data.Complex type Matrix a = [[a]] main :: IO ()main =    mapM_ (\a -> do        putStrLn "\nMatrix:"        mapM_ print a        putStrLn "Conjugate Transpose:"        mapM_ print (conjTranspose a)        putStrLn $"Hermitian? " ++ show (isHermitianMatrix a) putStrLn$ "Normal? " ++ show (isNormalMatrix a)        putStrLn $"Unitary? " ++ show (isUnitaryMatrix 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 :: RealFloat a => Matrix (Complex a) -> BoolisHermitianMatrix a = a approxEqualMatrix conjTranspose aisNormalMatrix a = (a mmul conjTranspose a) approxEqualMatrix (conjTranspose a mmul a)isUnitaryMatrix a = (a mmul conjTranspose a) approxEqualMatrix ident (length a) approxEqualMatrix :: (Fractional a, Ord a) => Matrix (Complex a) -> Matrix (Complex a) -> BoolapproxEqualMatrix 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 :: Num a => Matrix a -> Matrix a -> Matrix ammul a b = [[sum (zipWith (*) row column) | column <- transpose b] | row <- a] ident :: Num a => Int -> Matrix aident 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

Output:

Matrix:
[3.0 :+ 0.0,2.0 :+ 1.0]
[2.0 :+ (-1.0),1.0 :+ 0.0]
Conjugate Transpose:
[3.0 :+ (-0.0),2.0 :+ 1.0]
[2.0 :+ (-1.0),1.0 :+ (-0.0)]
Hermitian? True
Normal? True
Unitary? False

Matrix:
[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]
Conjugate Transpose:
[1.0 :+ (-0.0),0.0 :+ (-0.0),1.0 :+ (-0.0)]
[1.0 :+ (-0.0),1.0 :+ (-0.0),0.0 :+ (-0.0)]
[0.0 :+ (-0.0),1.0 :+ (-0.0),1.0 :+ (-0.0)]
Hermitian? False
Normal? True
Unitary? False

Matrix:
[0.7071067811865476 :+ 0.0,0.7071067811865476 :+ 0.0,0.0 :+ 0.0]
[0.0 :+ 0.7071067811865476,0.0 :+ (-0.7071067811865476),0.0 :+ 0.0]
[0.0 :+ 0.0,0.0 :+ 0.0,0.0 :+ 1.0]
Conjugate Transpose:
[0.7071067811865476 :+ (-0.0),0.0 :+ (-0.7071067811865476),0.0 :+ (-0.0)]
[0.7071067811865476 :+ (-0.0),0.0 :+ 0.7071067811865476,0.0 :+ (-0.0)]
[0.0 :+ (-0.0),0.0 :+ (-0.0),0.0 :+ (-1.0)]
Hermitian? False
Normal? True
Unitary? True


## J

Solution:
   ct =: [email protected]|:                      NB.  Conjugate transpose (ct A is A_ct)
Examples:
   X         =: +/ . *             NB. Matrix Multiply (x)    HERMITIAN =:  3 2j1 ,: 2j_1 1     (-: ct) HERMITIAN               NB.  A_ct = A1    NORMAL    =:  1 1 0 , 0 1 1 ,: 1 0 1   ((X~ -: X) ct) NORMAL           NB. A_ct x A = A x A_ct1    UNITARY   =:  (-:%:2) * 1 1 0 , 0j_1 0j1 0 ,: 0 0 0j1 * %:2   (ct -: %.)  UNITARY             NB.  A_ct = A^-11
Reference (example matrices for other langs to use):
   HERMITIAN;NORMAL;UNITARY+--------+-----+--------------------------+|   3 2j1|1 1 0|   0.707107   0.707107   0||2j_1   1|0 1 1|0j_0.707107 0j0.707107   0||        |1 0 1|          0          0 0j1|+--------+-----+--------------------------+   NB. In J, PjQ is P + Q*i and the 0.7071... is sqrt(2)    hermitian=: -: ct   normal =: (X~ -: X) ct   unitary=: ct -: %.    (hermitian,normal,unitary)&.>HERMITIAN;NORMAL;UNITARY+-----+-----+-----+|1 1 0|0 1 0|0 1 1|+-----+-----+-----+

## jq

Works with: jq version 1.4

In the following, we use the array [x,y] to represent the complex number x + iy, but the following functions also accept a number wherever a complex number is acceptable.

#### Infrastructure

(1) transpose/0:

If your jq does not have "transpose" then the following may be used:

# 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 ;

(2) Operations on real/complex numbers

# x must be real or complex, and ditto for y;# always return complexdef plus(x; y):    if (x|type) == "number" then       if  (y|type) == "number" then [ x+y, 0 ]       else [ x + y[0], y[1]]       end    elif (y|type) == "number" then plus(y;x)    else [ x[0] + y[0], x[1] + y[1] ]    end; # x must be real or complex, and ditto for y;# always return complexdef multiply(x; y):    if (x|type) == "number" then       if  (y|type) == "number" then [ x*y, 0 ]       else [x * y[0], x * y[1]]       end    elif (y|type) == "number" then multiply(y;x)    else [ x[0] * y[0] - x[1] * y[1],  x[0] * y[1] + x[1] * y[0]]    end; # conjugate of a real or complex numberdef conjugate:  if type == "number" then [.,0]  else [.[0], -(.[1]) ]  end;

(3) Array operations

# a and b are arrays of real/complex numbersdef 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]) ));

(4) Matrix operations

# convert a matrix of mixed real/complex entries to all complex entriesdef to_complex:  def toc: if type == "number" then [.,0] else . end;  map( map(toc) ); # simple matrix pretty-printerdef pp(wide):  def pad: tostring | (wide - length) * " " + .;  def row: reduce .[] as $x (""; . + ($x|pad));  reduce .[] as $row (""; . + "\n\($row|row)"); # Matrix multiplication# A and B should both be real/complex matrices,# A being m by n, and B being n by p.def matrix_multiply(A; B):  A as $A | B as$B  | ($B[0]|length) as$p  | ($B|transpose) as$BT  | reduce range(0; $A|length) as$i       ([]; reduce range(0; $p) as$j          (.; .[$i][$j] = dot_product( $A[$i]; $BT[$j] ) )) ; # Complex identity matrix of dimension ndef complex_identity(n):  def indicator(i;n):  [range(0;n)] | map( [0,0]) | .[i] = [1,0];  reduce range(0; n) as $i ([]; . + [indicator($i; n )] ); # Approximate equality of two matrices# Are two real/complex matrices essentially equal# in the sense that the sum of the squared element-wise differences# is less than or equal to epsilon?# The two matrices must be conformal.def approximately_equal(M; N; epsilon):  def norm: multiply(. ; conjugate ) | .[0];  def sqdiff( x; y): plus(x; multiply(y; -1)) | norm;  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;

## Julia

Julia has a built-in matrix type, and the conjugate-transpose of a complex matrix A is simply:

A'

(similar to Matlab). You can check whether A is Hermitian via the built-in function

ishermitian(A)

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

isnormal(A) = size(A,1) == size(A,2) && A'*A == A*A'isunitary(A) = size(A,1) == size(A,2) && A'*A == eye(A)

## 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:

// version 1.1.3 typealias C = Complextypealias Vector = Array<C>typealias Matrix = Array<Vector> class Complex(val real: Double, val imag: Double) {     operator fun plus(other: Complex) =        Complex(this.real + other.real, this.imag + other.imag)     operator fun times(other: Complex) =        Complex(this.real * other.real - this.imag * other.imag,                this.real * other.imag + this.imag * other.real)     fun conj() = Complex(this.real, -this.imag)     /* tolerable equality allowing for rounding of Doubles */    infix fun teq(other: Complex) =        Math.abs(this.real - other.real) <= 1e-14 &&        Math.abs(this.imag - other.imag) <= 1e-14     override fun toString() = "${"%.3f".format(real)} " + when { imag > 0.0 -> "+${"%.3f".format(imag)}i"        imag == 0.0  -> "+ 0.000i"        else         -> "- ${"%.3f".format(-imag)}i" }} fun Matrix.conjTranspose(): Matrix { val rows = this.size val cols = this[0].size return Matrix(cols) { i -> Vector(rows) { j -> this[j][i].conj() } }} operator fun Matrix.times(other: Matrix): Matrix { val rows1 = this.size val cols1 = this[0].size val rows2 = other.size val cols2 = other[0].size require(cols1 == rows2) val result = Matrix(rows1) { Vector(cols2) { C(0.0, 0.0) } } for (i in 0 until rows1) { for (j in 0 until cols2) { for (k in 0 until rows2) { result[i][j] += this[i][k] * other[k][j] } } } return result} /* tolerable matrix equality using the same concept as for complex numbers */infix fun Matrix.teq(other: Matrix): Boolean { if (this.size != other.size || this[0].size != other[0].size) return false for (i in 0 until this.size) { for (j in 0 until this[0].size) if (!(this[i][j] teq other[i][j])) return false } return true} fun Matrix.isHermitian() = this teq this.conjTranspose() fun Matrix.isNormal(): Boolean { val ct = this.conjTranspose() return (this * ct) teq (ct * this)} fun Matrix.isUnitary(): Boolean { val ct = this.conjTranspose() val prod = this * ct val ident = identityMatrix(prod.size) val prod2 = ct * this return (prod teq ident) && (prod2 teq ident)} fun Matrix.print() { val rows = this.size val cols = this[0].size for (i in 0 until rows) { for (j in 0 until cols) { print(this[i][j]) print(if(j < cols - 1) ", " else "\n") } } println()} fun identityMatrix(n: Int): Matrix { require(n >= 1) val ident = Matrix(n) { Vector(n) { C(0.0, 0.0) } } for (i in 0 until n) ident[i][i] = C(1.0, 0.0) return ident} fun main(args: Array<String>) { val x = Math.sqrt(2.0) / 2.0 val matrices = arrayOf( arrayOf( arrayOf(C(3.0, 0.0), C(2.0, 1.0)), arrayOf(C(2.0, -1.0), C(1.0, 0.0)) ), arrayOf( arrayOf(C(1.0, 0.0), C(1.0, 0.0), C(0.0, 0.0)), arrayOf(C(0.0, 0.0), C(1.0, 0.0), C(1.0, 0.0)), arrayOf(C(1.0, 0.0), C(0.0, 0.0), C(1.0, 0.0)) ), arrayOf( arrayOf(C(x, 0.0), C(x, 0.0), C(0.0, 0.0)), arrayOf(C(0.0, -x), C(0.0, x), C(0.0, 0.0)), arrayOf(C(0.0, 0.0), C(0.0, 0.0), C(0.0, 1.0)) ) ) for (m in matrices) { println("Matrix:") m.print() val mct = m.conjTranspose() println("Conjugate transpose:") mct.print() println("Hermitian?${mct.isHermitian()}")        println("Normal?    ${mct.isNormal()}") println("Unitary?${mct.isUnitary()}\n")    }}
Output:
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?   true


## Maple

The commands HermitianTranspose and IsUnitary are provided by the LinearAlgebra package.

M:=<<3|2+I>,<2-I|1>>: with(LinearAlgebra):IsNormal:=A->EqualEntries(A^%H.A,A.A^%H): M^%H;HermitianTranspose(M);type(M,'Matrix'(hermitian));IsNormal(M);IsUnitary(M);

Output:

                               [  3    2 + I]
[            ]
[2 - I    1  ]

[  3    2 + I]
[            ]
[2 - I    1  ]

true

true

false

## Mathematica / Wolfram Language

NormalMatrixQ[a_List?MatrixQ] := Module[{b = [email protected]@a},a.b === b.a]UnitaryQ[m_List?MatrixQ] := ([email protected]@m.m == [email protected]@m) m = {{1, 2I, 3}, {3+4I, 5, I}};m //MatrixForm->(1	2I	33+4I	5	I) ConjugateTranspose[m] //MatrixForm->(1	3-4I-2I	53	-I) {[email protected]#, [email protected]#, [email protected]#}&@m-> {False, False, False}

## PARI/GP

conjtranspose(M)=conj(M~)isHermitian(M)=M==conj(M~)isnormal(M)=my(H=conj(M~));H*M==M*Hisunitary(M)=M*conj(M~)==1

## 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.

use strict;use English;use Math::Complex;use Math::MatrixReal; my @examples = (example1(), example2(), example3());foreach my $m (@examples) { print "Starting matrix:\n", cmat_as_string($m), "\n";    my $m_ct = conjugate_transpose($m);    print "Its conjugate transpose:\n", cmat_as_string($m_ct), "\n"; print "Is Hermitian? ", (cmats_are_equal($m, $m_ct) ? 'TRUE' : 'FALSE'), "\n"; my$product = $m_ct *$m;    print "Is normal? ", (cmats_are_equal($product,$m * $m_ct) ? 'TRUE' : 'FALSE'), "\n"; my$I = identity(($m->dim())[0]); print "Is unitary? ", (cmats_are_equal($product, $I) ? 'TRUE' : 'FALSE'), "\n"; print "\n";}exit 0; sub cmats_are_equal { my ($m1, $m2) = @ARG; my$max_norm = 1.0e-7;    return abs($m1 -$m2) < $max_norm; # Math::MatrixReal overloads abs().} # Note that Math::Complex and Math::MatrixReal both overload '~', for# complex conjugates and matrix transpositions respectively.sub conjugate_transpose { my$m_T = ~ shift;    my $result =$m_T->each(sub {~ $ARG[0]}); return$result;} sub cmat_as_string {    my $m = shift; my$n_rows = ($m->dim())[0]; my @row_strings = map { q{[} . join(q{, },$m-