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 $M$ contains complex numbers. Then the conjugate transpose of $M$ is a matrix $M^{H}$ containing the complex conjugates of the matrix transposition of $M$.

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

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

In the next list, $M$ must also be a square matrix.

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

Task

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

Also determine if the matrix is a:

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

See also

• MathWorld entry: conjugate transpose
• MathWorld entry: Hermitian matrix
• MathWorld entry: normal matrix
• MathWorld entry: unitary matrix

11l

Translation of: Nim

-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)

Output:

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

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

ALGOL 68

Uses the same test cases as the Ada sample.

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

Output:

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

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;

matrix
transpose (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;
}

int
isHermitian (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;
}

matrix
multiply (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;
}

int
isNormal (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;
}

int
isUnitary (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;
}


int
main ()
{
  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

C++

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

Output:

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

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));
    }
}

F#

// 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")

Output:

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

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

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)

Output:

true          true          true
true          true          false
false         false         false

Go

package main

import (
    "fmt"
    "math"
    "math/cmplx"
)

// a type to represent matrices
type 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 matricies
func 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 constructors
func 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 stdout
func (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 task
func (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

Haskell

Slow implementation using lists.

import Data.Complex (Complex(..), conjugate)
import Data.List (transpose)

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) -> Bool
isHermitianMatrix = mTest id conjTranspose

isNormalMatrix = mTest mmct (mmul =<< conjTranspose)

isUnitaryMatrix = mTest mmct (ident . length)

mTest
  :: RealFloat a
  => (a2 -> Matrix (Complex a)) -> (a2 -> Matrix (Complex a)) -> a2 -> Bool
mTest f g = (approxEqualMatrix . f) <*> g

mmct
  :: RealFloat a
  => Matrix (Complex a) -> Matrix (Complex a)
mmct = mmul <*> conjTranspose

approxEqualMatrix
  :: (Fractional a, Ord a)
  => 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
  :: Num a
  => 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

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

Conjugate transpose:

   ct =: +@|:                      NB.  Conjugate transpose (ct A is A_ct)

+ when used without a left argument is conjugate, |: is transpose, and @ composes functions.

Examples:

   X         =: +/ . *             NB. Matrix Multiply (x)

   HERMITIAN =:  3 2j1 ,: 2j_1 1  
   (-: ct) HERMITIAN               NB.  A_ct = A
1

   NORMAL    =:  1 1 0 , 0 1 1 ,: 1 0 1
   ((X~ -: X) ct) NORMAL           NB. A_ct x A = A x A_ct
1

   UNITARY   =:  (-:%:2) * 1 1 0 , 0j_1 0j1 0 ,: 0 0 0j1 * %:2
   (ct -: %.)  UNITARY             NB.  A_ct = A^-1
1

Reference (with example matrices):

   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|
+-----+-----+-----+

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

Output:

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

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 complex
def 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 complex
def 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 number
def conjugate:
  if type == "number" then [.,0]
  else [.[0], -(.[1]) ]
  end;

(3) Array operations

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

(4) Matrix operations

# convert a matrix of mixed real/complex entries to all complex entries
def to_complex:
  def toc: if type == "number" then [.,0] else . end;
  map( map(toc) );

# simple matrix pretty-printer
def 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 n
def 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;

Conjugate transposition

# (entries may be real and/or complex)
def conjugate_transpose:
  map( map(conjugate) ) | transpose;

# A Hermitian matrix equals its own conjugate transpose
def is_hermitian:
  to_complex == conjugate_transpose;

# A matrix is normal if it commutes multiplicatively
# with its conjugate transpose
def is_normal:
  . as $M
  | conjugate_transpose as $H
  | matrix_multiply($H; $M) == matrix_multiply($H; $M);

# A unitary matrix (U) has its inverse equal to its conjugate transpose (T)
# i.e. U^-1 == T; NASC is I == UT == TU
def is_unitary:
  . as $M
  | conjugate_transpose as $H
  | complex_identity(length) as $I
  | approximately_equal( $I; matrix_multiply($H;$M); 1e-10)
    and approximately_equal( $I ; matrix_multiply($M;$H); 1e-10)  ;

Examples

def hermitian_example:
  [ [ 3,    [2,1]],
    [[2,-1], 1   ] ];

def normal_example:
  [ [1, 1, 0],
    [0, 1, 1],
    [1, 0, 1] ];

def unitary_example:
  0.707107
  |  [ [ [., 0], [.,  0],   0 ],
       [ [0, -.], [0, .],   0 ],
       [ 0,      0,      [0,1] ] ];

def demo:
  hermitian_example
  | ("Hermitian example:", pp(8)),
    "",
    ("Its conjugate transpose is:",  (to_complex | conjugate_transpose | pp(8))),
    "",
    "Hermitian example: \(hermitian_example | is_hermitian )",
    "",
    "Normal example:    \(normal_example    | is_normal )",
    "",
    "Unitary example:   \(unitary_example   | is_unitary)"
;

demo

Output:

$ jq -r -c -n -f Conjugate_transpose.jq
Hermitian example:

       3   [2,1]
  [2,-1]       1

Conjugate transpose:

  [3,-0]   [2,1]
  [2,-1]  [1,-0]

Hermitian example: true

Normal example:    true

Unitary example:   true

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

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)

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 = Complex
typealias 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 = Conjugate@Transpose@a},a.b === b.a]
UnitaryQ[m_List?MatrixQ] := (Conjugate@Transpose@m.m == IdentityMatrix@Length@m)

m = {{1, 2I, 3}, {3+4I, 5, I}};
m //MatrixForm
->
(1	2I	3
3+4I	5	I)

ConjugateTranspose[m] //MatrixForm
->
(1	3-4I
-2I	5
3	-I)

{HermitianMatrixQ@#, NormalMatrixQ@#, UnitaryQ@#}&@m
-> {False, False, False}

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

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)

Output:

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

PARI/GP

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

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 warnings;
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->row($ARG)->as_list) . q{]} }
                          (1 .. $n_rows);
    return join("\n", @row_strings);
}

sub identity {
    my $N = shift;
    my $m = Math::MatrixReal->new($N, $N);
    $m->one();
    return $m;
}

sub example1 {
    my $m = Math::MatrixReal->new(2, 2);
    $m->assign(1, 1, cplx(3, 0));
    $m->assign(1, 2, cplx(2, 1));
    $m->assign(2, 1, cplx(2, -1));
    $m->assign(2, 2, cplx(1, 0));
    return $m;
}

sub example2 {
    my $m = Math::MatrixReal->new(3, 3);
    $m->assign(1, 1, cplx(1, 0));
    $m->assign(1, 2, cplx(1, 0));
    $m->assign(1, 3, cplx(0, 0));
    $m->assign(2, 1, cplx(0, 0));
    $m->assign(2, 2, cplx(1, 0));
    $m->assign(2, 3, cplx(1, 0));
    $m->assign(3, 1, cplx(1, 0));
    $m->assign(3, 2, cplx(0, 0));
    $m->assign(3, 3, cplx(1, 0));
    return $m;
}

sub example3 {
    my $m = Math::MatrixReal->new(3, 3);
    $m->assign(1, 1, cplx(0.70710677, 0));
    $m->assign(1, 2, cplx(0.70710677, 0));
    $m->assign(1, 3, cplx(0, 0));
    $m->assign(2, 1, cplx(0, -0.70710677));
    $m->assign(2, 2, cplx(0, 0.70710677));
    $m->assign(2, 3, cplx(0, 0));
    $m->assign(3, 1, cplx(0, 0));
    $m->assign(3, 2, cplx(0, 0));
    $m->assign(3, 3, cplx(0, 1));
    return $m;
}

Output:

Starting matrix:
[3, 2+i]
[2-i, 1]
Its conjugate transpose:
[3, 2+i]
[2-i, 1]
Is Hermitian? TRUE
Is normal? TRUE
Is unitary? FALSE

Starting matrix:
[1, 1, 0]
[0, 1, 1]
[1, 0, 1]
Its conjugate transpose:
[1, 0, 1]
[1, 1, 0]
[0, 1, 1]
Is Hermitian? FALSE
Is normal? TRUE
Is unitary? FALSE

Starting matrix:
[0.70710677, 0.70710677, 0]
[-0.70710677i, 0.70710677i, 0]
[0, 0, i]
Its conjugate transpose:
[0.70710677, 0.70710677i, 0]
[0.70710677, -0.70710677i, 0]
[0, 0, -i]
Is Hermitian? FALSE
Is normal? TRUE
Is unitary? TRUE

Phix

Note this code has no testing for non-square matrices.

with javascript_semantics
include complex.e
 
procedure m_print(sequence a)
    a = deep_copy(a)
    integer l = length(a)
    for i=1 to l do
        for j=1 to l do
            a[i][j] = complex_sprint(a[i][j])
        end for
        a[i] = "["&join(a[i],",")&"]"
    end for
    puts(1,join(a,"\n")&"\n")
end procedure
 
function conjugate_transpose(sequence a)
    sequence res = deep_copy(a)
    integer l = length(a)
    for i=1 to l do
        for j=1 to l do
          res[i][j] = complex_conjugate(a[j][i])
        end for
    end for
    return res
end function
 
function m_unitary(sequence act)
-- note: a was normal and act = a*ct already
    integer l = length(act)
    for i=1 to l do
        for j=1 to l do
            atom {re,im} = act[i,j]
            -- round to nearest billionth
            -- (powers of 2 help the FPU out)
            re = round(re,1024*1024*1024)
            im = round(im,1024*1024*1024)
            if im!=0
            or (i=j and re!=1)
            or (i!=j and re!=0) then
                return 0
            end if
        end for
    end for
    return 1
end function
 
function m_mul(sequence a, sequence b)
    sequence res = sq_mul(a,0)
    integer l = length(a)
    for i=1 to l do
        for j=1 to l do
            for k=1 to l do
                res[i][j] = complex_add(res[i][j],complex_mul(a[i][k],b[k][j]))
            end for
        end for
    end for
    return res
end function
 
procedure test(sequence a)
    sequence ct = conjugate_transpose(a)
    printf(1,"Original matrix:\n")
    m_print(a)
    printf(1,"Conjugate transpose:\n")
    m_print(ct)
    -- 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, 
    --       but as it is, all tests pass with the builtin = operator.
    printf(1,"Hermitian?: %t\n",a=ct) -- (this one)
    sequence act = m_mul(a,ct), cta = m_mul(ct,a)
    bool normal = act=cta                                   -- (&this one)
    printf(1,"Normal?: %t\n",normal)
    printf(1,"Unitary?: %t\n\n",normal and m_unitary(act))
end procedure
 
constant x = sqrt(2)/2,
         tests = {{{{3, 0},{2,1}},
                   {{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}}}}
 
papply(tests,test)

Output:

Original matrix:
[3,2+i]
[2-i,1]
Conjugate transpose:
[3,2+i]
[2-i,1]
Hermitian?: true
Normal?: true
Unitary?: false

Original matrix:
[1,1+i,2i]
[1-i,5,-3]
[-2i,-3,0]
Conjugate transpose:
[1,1+i,2i]
[1-i,5,-3]
[-2i,-3,0]
Hermitian?: true
Normal?: true
Unitary?: false

Original matrix:
[0.5+0.5i,0.5-0.5i]
[0.5-0.5i,0.5+0.5i]
Conjugate transpose:
[0.5-0.5i,0.5+0.5i]
[0.5+0.5i,0.5-0.5i]
Hermitian?: false
Normal?: true
Unitary?: true

Original matrix:
[1,1,0]
[0,1,1]
[1,0,1]
Conjugate transpose:
[1,0,1]
[1,1,0]
[0,1,1]
Hermitian?: false
Normal?: true
Unitary?: false

Original matrix:
[0.707107,0.707107,0]
[-0.707107i,0.707107i,0]
[0,0,i]
Conjugate transpose:
[0.707107,0.707107i,0]
[0.707107,-0.707107i,0]
[0,0,-i]
Hermitian?: false
Normal?: true
Unitary?: true

Original matrix:
[2+7i,9-5i]
[3+4i,8-6i]
Conjugate transpose:
[2-7i,3-4i]
[9+5i,8+6i]
Hermitian?: false
Normal?: false
Unitary?: false

PL/I

test: procedure options (main); /* 1 October 2012 */
   declare n fixed binary;

   put ('Conjugate a complex square matrix.');
   put skip list ('What is the order of the matrix?:');
   get (n);
   begin;
      declare (M, MH, MM, MM_MMH, MM_MHM, IDENTITY)(n,n) fixed complex;
      declare i fixed binary;

      IDENTITY = 0; do i = 1 to n; IDENTITY(I,I) = 1; end;
      put skip list ('Please type the matrix:');
      get list (M);
      do i = 1 to n;
         put skip list (M(i,*));
      end;
      do i = 1 to n;
         MH(i,*) = conjg(M(*,i));
      end;
      put skip list ('The conjugate transpose is:');
      do i = 1 to n;
         put skip list (MH(i,*));
      end;
      if all(M=MH) then
         put skip list ('Matrix is Hermitian');
      call MMULT(M, MH, MM_MMH);
      call MMULT(MH, M, MM_MHM);

      if all(MM_MMH = MM_MHM) then
         put skip list ('Matrix is Normal');

      if all(ABS(MM_MMH - IDENTITY) < 0.0001) then
         put skip list ('Matrix is unitary');
      if all(ABS(MM_MHM - IDENTITY) < 0.0001) then
         put skip list ('Matrix is unitary');
   end;

MMULT: procedure (M, MH, MM);
   declare (M, MH, MM)(*,*) fixed complex;
   declare (i, j, n) fixed binary;

   n = hbound(M,1);
   do i = 1 to n;
      do j = 1 to n;
         MM(i,j) = sum(M(i,*) * MH(*,j) );
      end;
   end;
end MMULT;
end test;

Outputs from separate runs:

Please type the matrix: 

       1+0I                    1+0I                    1+0I       
       1+0I                    1+0I                    1+0I       
       1+0I                    1+0I                    1+0I       
The conjugate transpose is: 
       1-0I                    1-0I                    1-0I       
       1-0I                    1-0I                    1-0I       
       1-0I                    1-0I                    1-0I       
Matrix is Hermitian 
Matrix is Normal 

       1+0I                    1+0I                    0+0I
       0+0I                    1+0I                    1+0I       
       1+0I                    0+0I                    1+0I       
The conjugate transpose is: 
       1-0I                    0-0I                    1-0I       
       1-0I                    1-0I                    0-0I       
       0-0I                    1-0I                    1-0I       
Matrix is Normal 

Next test performed with declaration of matrixes changed to decimal precision (10,5).

Please type the matrix:

      0.70710+0.00000I        0.70710+0.00000I        0.00000+0.00000I
      0.00000+0.70710I        0.00000-0.70710I        0.00000+0.00000I
      0.00000+0.00000I        0.00000+0.00000I        0.00000+1.00000I
    
The conjugate transpose is: 
      0.70710-0.00000I        0.00000-0.70710I        0.00000-0.00000I
      0.70710-0.00000I        0.00000+0.70710I        0.00000-0.00000I
      0.00000-0.00000I        0.00000-0.00000I        0.00000-1.00000I

Matrix is Normal 
Matrix is unitary 
Matrix is unitary

PowerShell

function conjugate-transpose($a) {
    $arr = @()
    if($a) { 
        $n = $a.count - 1 
        if(0 -lt $n) { 
            $m = ($a | foreach {$_.count} | measure-object -Minimum).Minimum - 1
            if( 0 -le $m) {
                if (0 -lt $m) {
                    $arr =@(0)*($m+1)
                    foreach($i in 0..$m) {
                        $arr[$i] = foreach($j in 0..$n) {@([System.Numerics.complex]::Conjugate($a[$j][$i]))}    
                    }
                } else {$arr = foreach($row in $a) {[System.Numerics.complex]::Conjugate($row[0])}}
            }
        } else {$arr = foreach($row in $a) {[System.Numerics.complex]::Conjugate($row[0])}}
    }
    $arr
}

function multarrays-complex($a, $b) {
    $c = @()
    if($a -and $b) {
        $n = $a.count - 1
        $m = $b[0].count - 1
        $c = @([System.Numerics.complex]::new(0,0))*($n+1)
        foreach ($i in 0..$n) {    
            $c[$i] = foreach ($j in 0..$m) { 
                [System.Numerics.complex]$sum = [System.Numerics.complex]::new(0,0)
                foreach ($k in 0..$n){$sum = [System.Numerics.complex]::Add($sum, ([System.Numerics.complex]::Multiply($a[$i][$k],$b[$k][$j])))}
                $sum
            }
        }
    }
    $c
}

function identity-complex($n) {
    if(0 -lt $n) { 
        $array = @(0) * $n
        foreach ($i in 0..($n-1)) {
            $array[$i] = @([System.Numerics.complex]::new(0,0)) * $n
            $array[$i][$i] = [System.Numerics.complex]::new(1,0)
        }  
        $array 
    } else { @() }
}

function are-eq ($a,$b) { -not (Compare-Object $a $b -SyncWindow 0)}

function show($a) {
    if($a) { 
        0..($a.Count - 1) | foreach{ if($a[$_]){"$($a[$_])"}else{""} }
    }
}
function complex($a,$b) {[System.Numerics.complex]::new($a,$b)}
 
$id2 = identity-complex 2
$m = @(@((complex 2 7), (complex 9 -5)),@((complex 3 4), (complex 8 -6)))
$hm = conjugate-transpose $m
$mhm = multarrays-complex $m $hm
$hmm = multarrays-complex $hm $m
"`$m ="
show $m
""
"`$hm = conjugate-transpose `$m ="
show $hm
""
"`$m * `$hm ="
show $mhm
""
"`$hm * `$m ="
show $hmm
""
"Hermitian? `$m = $(are-eq $m $hm)"
"Normal? `$m = $(are-eq $mhm $hmm)"
"Unitary? `$m = $((are-eq $id2 $hmm) -and (are-eq $id2 $mhm))"

Output:

$m =
(2, 7) (9, -5)
(3, 4) (8, -6)

$hm = conjugate-transpose $m =
(2, -7) (3, -4)
(9, 5) (8, 6)

$m * $hm =
(159, 0) (136, 27)
(136, -27) (125, 0)

$hm * $m =
(78, 0) (-17, -123)
(-17, 123) (206, 0)

Hermitian? $m = False
Normal? $m = False
Unitary? $m = False

Python

Internally, matrices must be represented as rectangular tuples of tuples of complex numbers.

def conjugate_transpose(m):
    return tuple(tuple(n.conjugate() for n in row) for row in zip(*m))

def mmul( ma, mb):
    return tuple(tuple(sum( ea*eb for ea,eb in zip(a,b)) for b in zip(*mb)) for a in ma)

def mi(size):
    'Complex Identity matrix'
    sz = range(size)
    m = [[0 + 0j for i in sz] for j in sz]
    for i in range(size):
        m[i][i] = 1 + 0j
    return tuple(tuple(row) for row in m)

def __allsame(vector):
    first, rest = vector[0], vector[1:]
    return all(i == first for i in rest)

def __allnearsame(vector, eps=1e-14):
    first, rest = vector[0], vector[1:]
    return all(abs(first.real - i.real) < eps and abs(first.imag - i.imag) < eps
               for i in rest)

def isequal(matrices, eps=1e-14):
    'Check any number of matrices for equality within eps'
    x = [len(m) for m in matrices]
    if not __allsame(x): return False
    y = [len(m[0]) for m in matrices]
    if not __allsame(y): return False
    for s in range(x[0]):
        for t in range(y[0]):
            if not __allnearsame([m[s][t] for m in matrices], eps): return False
    return True
    

def ishermitian(m, ct):
    return isequal([m, ct])

def isnormal(m, ct):
    return isequal([mmul(m, ct), mmul(ct, m)])

def isunitary(m, ct):
    mct, ctm = mmul(m, ct), mmul(ct, m)
    mctx, mcty, cmx, ctmy = len(mct), len(mct[0]), len(ctm), len(ctm[0])
    ident = mi(mctx)
    return isequal([mct, ctm, ident])

def printm(comment, m):
    print(comment)
    fields = [['%g%+gj' % (f.real, f.imag) for f in row] for row in m]
    width = max(max(len(f) for f in row) for row in fields)
    lines = (', '.join('%*s' % (width, f) for f in row) for row in fields)
    print('\n'.join(lines))

if __name__ == '__main__':
    for matrix in [
            ((( 3.000+0.000j), (+2.000+1.000j)), 
            (( 2.000-1.000j), (+1.000+0.000j))),

            ((( 1.000+0.000j), (+1.000+0.000j), (+0.000+0.000j)), 
            (( 0.000+0.000j), (+1.000+0.000j), (+1.000+0.000j)), 
            (( 1.000+0.000j), (+0.000+0.000j), (+1.000+0.000j))),

            ((( 2**0.5/2+0.000j), (+2**0.5/2+0.000j), (+0.000+0.000j)), 
            (( 0.000+2**0.5/2j), (+0.000-2**0.5/2j), (+0.000+0.000j)), 
            (( 0.000+0.000j), (+0.000+0.000j), (+0.000+1.000j)))]:
        printm('\nMatrix:', matrix)
        ct = conjugate_transpose(matrix)
        printm('Its conjugate transpose:', ct)
        print('Hermitian? %s.' % ishermitian(matrix, ct))
        print('Normal?    %s.' % isnormal(matrix, ct))
        print('Unitary?   %s.' % isunitary(matrix, ct))

Output:

Matrix:
3+0j, 2+1j
2-1j, 1+0j
Its conjugate transpose:
3-0j, 2+1j
2-1j, 1-0j
Hermitian? True.
Normal?    True.
Unitary?   False.

Matrix:
1+0j, 1+0j, 0+0j
0+0j, 1+0j, 1+0j
1+0j, 0+0j, 1+0j
Its conjugate transpose:
1-0j, 0-0j, 1-0j
1-0j, 1-0j, 0-0j
0-0j, 1-0j, 1-0j
Hermitian? False.
Normal?    True.
Unitary?   False.

Matrix:
0.707107+0j, 0.707107+0j,        0+0j
0-0.707107j, 0+0.707107j,        0+0j
       0+0j,        0+0j,        0+1j
Its conjugate transpose:
0.707107-0j, 0+0.707107j,        0-0j
0.707107-0j, 0-0.707107j,        0-0j
       0-0j,        0-0j,        0-1j
Hermitian? False.
Normal?    True.
Unitary?   True.

Racket

#lang racket
(require math)
(define H matrix-hermitian)

(define (normal? M)
  (define MH (H M))
  (equal? (matrix* MH M) 
          (matrix* M MH)))

(define (unitary? M)
  (define MH (H M))
  (and (matrix-identity? (matrix* MH M))
       (matrix-identity? (matrix* M MH))))

(define (hermitian? M)
  (equal? (H M) M))

Test:

(define M (matrix [[3.000+0.000i +2.000+1.000i]
                   [2.000-1.000i +1.000+0.000i]]))
(H M)
(normal? M)
(unitary? M)
(hermitian? M)

Output:

(array #[#[3.0-0.0i 2.0+1.0i] #[2.0-1.0i 1.0-0.0i]])
#t
#f
#f

Raku

(formerly Perl 6)

Works with: Rakudo version 2015-12-13

for [ # Test Matrices
       [   1, 1+i, 2i],
       [ 1-i,   5, -3],
       [0-2i,  -3,  0]
    ],
    [
       [1, 1, 0],
       [0, 1, 1],
       [1, 0, 1]
    ],
    [
       [0.707 ,    0.707,  0],
       [0.707i, 0-0.707i,  0],
       [0     ,        0,  i]
    ]
    -> @m {
        say "\nMatrix:";
        @m.&say-it;
        my @t = @m».conj.&mat-trans;
        say "\nTranspose:";
        @t.&say-it;
        say "Is Hermitian?\t{is-Hermitian(@m, @t)}";
        say "Is Normal?\t{is-Normal(@m, @t)}";
        say "Is Unitary?\t{is-Unitary(@m, @t)}";
    }

sub is-Hermitian (@m, @t, --> Bool) {
    so @m».Complex eqv @t».Complex
 }

sub is-Normal (@m, @t, --> Bool) {
    so mat-mult(@m, @t)».Complex eqv mat-mult(@t, @m)».Complex
}

sub is-Unitary (@m, @t, --> Bool) {
    so mat-mult(@m, @t, 1e-3)».Complex eqv mat-ident(+@m)».Complex;
}

sub mat-trans (@m) { map { [ @m[*;$_] ] }, ^@m[0] }

sub mat-ident ($n) { [ map { [ flat 0 xx $_, 1, 0 xx $n - 1 - $_ ] }, ^$n ] }

sub mat-mult (@a, @b, \ε = 1e-15) {
    my @p;
    for ^@a X ^@b[0] -> ($r, $c) {
        @p[$r][$c] += @a[$r][$_] * @b[$_][$c] for ^@b;
        @p[$r][$c].=round(ε); # avoid floating point math errors
    }
    @p
}

sub say-it (@array) { $_».fmt("%9s").say for @array }

Output:

Matrix:
[        1      1+1i      0+2i]
[     1-1i         5        -3]
[     0-2i        -3         0]

Transpose:
[        1      1+1i      0+2i]
[     1-1i         5        -3]
[     0-2i        -3         0]
Is Hermitian?	True
Is Normal?	True
Is Unitary?	False

Matrix:
[        1         1         0]
[        0         1         1]
[        1         0         1]

Transpose:
[        1         0         1]
[        1         1         0]
[        0         1         1]
Is Hermitian?	False
Is Normal?	True
Is Unitary?	False

Matrix:
[    0.707     0.707         0]
[ 0+0.707i  0-0.707i         0]
[        0         0      0+1i]

Transpose:
[    0.707  0-0.707i         0]
[    0.707  0+0.707i         0]
[        0         0      0-1i]
Is Hermitian?	False
Is Normal?	True
Is Unitary?	True

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*/
                  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              /*    "     "   "      "        "      */
say '  M is Normal:   '           word('no yes', isHermitian('MMH', "MHM", N) + 1)
say '  M is Unary:    '           word('no yes',     isUnary('M',   N) + 1)
say 'MMH is Unary:    '           word('no yes',     isUnary('MMH', N) + 1)
say 'MHM is Unary:    '           word('no yes',     isUnary('MHM', N) + 1)
exit                                             /*stick a fork in it,  we're all done. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
cP:       procedure; arg ',' c;        return word( strip( translate(c, , 'IJ') )  0, 1)
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*/;                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*/;                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*/;                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*/;              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*/;                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*/
             say _
             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

output   when using the default input:

───────────────────────────────────matrix M────────────────────────────────────
        1                1                1
        1                1                1
        1                1                1


───────────────────────────────────matrix MH───────────────────────────────────
        1                1                1
        1                1                1
        1                1                1

M is Hermitian:   yes
  M is Normal:    yes
  M is Unary:     no
MMH is Unary:     no
MHM is Unary:     no

output   when using the input of:     3   .7071   .7071   0   0,.7071   0,-.7071   0   0   0   0,1

───────────────────────────────────matrix M────────────────────────────────────
        0.7071           0.7071           0
        0,0.7071j        0,-0.7071        0
        0                0                0,1j


───────────────────────────────────matrix MH───────────────────────────────────
        0.7071           0,-0.7071        0
        0.7071           0,0.7071j        0
        0                0                0,-1j

M is Hermitian:   no
  M is Normal:    yes
  M is Unary:     no
MMH is Unary:     yes
MHM is Unary:     yes

RPL

Although basic, RPL's matrix handling capabilities help to keep the code compact but still readable.

Works with: Halcyon Calc version 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

Output:

3: "Hermitian. Normal. "
2: "Normal. "
1: "Normal. Unitary. "

Ruby

Works with: Ruby version 2.0

require 'matrix'

# Start with some matrix.
i = Complex::I
matrix = Matrix[[i, 0, 0],
                [0, i, 0],
                [0, 0, i]]

# Find the conjugate transpose.
#   Matrix#conjugate appeared in Ruby 1.9.2.
conjt = matrix.conj.t           # aliases for matrix.conjugate.tranpose
print 'conjugate tranpose: '; puts conjt

if matrix.square?
  # These predicates appeared in Ruby 1.9.3.
  print 'Hermitian? '; puts matrix.hermitian?
  print '   normal? '; puts matrix.normal?
  print '  unitary? '; puts matrix.unitary?
else
  # Matrix is not square. These predicates would
  # raise ExceptionForMatrix::ErrDimensionMismatch.
  print 'Hermitian? false'
  print '   normal? false'
  print '  unitary? false'
end

Note: Ruby 1.9 had a bug in the Matrix#hermitian? method. It's fixed in 2.0.

Rust

Uses external crate 'num', version 0.1.34

extern crate num; // crate for complex numbers

use num::complex::Complex;
use std::ops::Mul;
use std::fmt;


#[derive(Debug, PartialEq)]
struct Matrix<f32> {
    grid: [[Complex<f32>; 2]; 2], // used to represent matrix
}


impl Matrix<f32> { // implements a method call for calculating the conjugate transpose
    fn conjugate_transpose(&self) -> Matrix<f32> {
        Matrix {grid: [[self.grid[0][0].conj(), self.grid[1][0].conj()],
        [self.grid[0][1].conj(), self.grid[1][1].conj()]]}
    }
}

impl Mul for Matrix<f32> { // implements '*' (multiplication) for the matrix
    type Output = Matrix<f32>;

    fn mul(self, other: Matrix<f32>) -> Matrix<f32> {
        Matrix {grid: [[self.grid[0][0]*other.grid[0][0] + self.grid[0][1]*other.grid[1][0],
            self.grid[0][0]*other.grid[0][1] + self.grid[0][1]*other.grid[1][1]],
            [self.grid[1][0]*other.grid[0][0] + self.grid[1][1]*other.grid[1][0],
            self.grid[1][0]*other.grid[1][0] + self.grid[1][1]*other.grid[1][1]]]}
    }
}

impl Copy for Matrix<f32> {} // implemented to prevent 'moved value' errors in if statements below
impl Clone for Matrix<f32> {
    fn clone(&self) -> Matrix<f32> {
        *self
    }
}

impl fmt::Display for Matrix<f32> { // implemented to make output nicer
    fn fmt(&self, f: &mut fmt::Formatter) -> fmt::Result {
        write!(f, "({}, {})\n({}, {})", self.grid[0][0], self.grid[0][1], self.grid[1][0], self.grid[1][1])
    }
}

fn main() {
    let a = Matrix {grid: [[Complex::new(3.0, 0.0), Complex::new(2.0, 1.0)],
        [Complex::new(2.0, -1.0), Complex::new(1.0, 0.0)]]};

    let b = Matrix {grid: [[Complex::new(0.5, 0.5), Complex::new(0.5, -0.5)],
        [Complex::new(0.5, -0.5), Complex::new(0.5, 0.5)]]};

    test_type(a);
    test_type(b);
}

fn test_type(mat: Matrix<f32>) {
    let identity = Matrix {grid: [[Complex::new(1.0, 0.0), Complex::new(0.0, 0.0)],
        [Complex::new(0.0, 0.0), Complex::new(1.0, 0.0)]]};
    let mat_conj = mat.conjugate_transpose();

    println!("Matrix: \n{}\nConjugate transpose: \n{}", mat, mat_conj);

    if mat == mat_conj {
        println!("Hermitian?: TRUE");
    } else {
        println!("Hermitian?: FALSE");
    }

    if mat*mat_conj == mat_conj*mat {
        println!("Normal?: TRUE");
    } else {
        println!("Normal?: FALSE");
    }

    if mat*mat_conj == identity {
        println!("Unitary?: TRUE");
    } else {
        println!("Unitary?: FALSE");
    }
}

Output:

Matrix:
(3+0i, 2+1i)
(2-1i, 1+0i)
Conjugate transpose:
(3+0i, 2+1i)
(2-1i, 1+0i)
Hermitian?: TRUE
Normal?: TRUE
Unitary?: FALSE
Matrix:
(0.5+0.5i, 0.5-0.5i)
(0.5-0.5i, 0.5+0.5i)
Conjugate transpose:
(0.5-0.5i, 0.5+0.5i)
(0.5+0.5i, 0.5-0.5i)
Hermitian?: FALSE
Normal?: TRUE
Unitary?: TRUE

Scala

object ConjugateTranspose {
  
  case class Complex(re: Double, im: Double) {
    def conjugate(): Complex = Complex(re, -im)
    def +(other: Complex) = Complex(re + other.re, im + other.im)
    def *(other: Complex) = Complex(re * other.re - im * other.im, re * other.im + im * other.re)
    override def toString(): String = {
      if (im < 0) {
        s"${re}${im}i"
      } else {
        s"${re}+${im}i"
      }
    }
  }
  
  case class Matrix(val entries: Vector[Vector[Complex]]) {
    
    def *(other: Matrix): Matrix = {
      new Matrix(
        Vector.tabulate(entries.size, other.entries(0).size)((r, c) => {
          val rightRow = entries(r)
          val leftCol = other.entries.map(_(c))
          rightRow.zip(leftCol)
            .map{ case (x, y) => x * y } // multiply pair-wise
            .foldLeft(new Complex(0,0)){ case (x, y) => x + y } // sum over all
        })
      )
    }
    
    def conjugateTranspose(): Matrix = {
      new Matrix(
        Vector.tabulate(entries(0).size, entries.size)((r, c) => entries(c)(r).conjugate)
      )
    }
    
    def isHermitian(): Boolean = {
      this == conjugateTranspose()
    }
    
    def isNormal(): Boolean = {
      val ct = conjugateTranspose()
      this * ct == ct * this
    }
    
    def isIdentity(): Boolean = {
      val entriesWithIndexes = for (r <- 0 until entries.size; c <- 0 until entries(r).size) yield (r, c, entries(r)(c))
      entriesWithIndexes.forall { case (r, c, x) =>
        if (r == c) {
          x == Complex(1.0, 0.0)
        } else {
          x == Complex(0.0, 0.0)
        }
      }
    }
    
    def isUnitary(): Boolean = {
      (this * conjugateTranspose()).isIdentity()
    }
    
    override def toString(): String = {
      entries.map("  " + _.mkString("[", ",", "]")).mkString("[\n", "\n", "\n]")
    }
    
  }
  
  def main(args: Array[String]): Unit = {
    val m = new Matrix(
      Vector.fill(3, 3)(new Complex(Math.random() * 2 - 1.0, Math.random() * 2 - 1.0))
    )
    println("Matrix: " + m)
    println("Conjugate Transpose: " + m.conjugateTranspose())
    println("Hermitian: " + m.isHermitian())
    println("Normal: " + m.isNormal())
    println("Unitary: " + m.isUnitary())
  }
  
}

Output:

Matrix: [
  [-0.7679977131543951-0.439979346567841i,-0.6011221529373452+0.510336881376179i,-0.22458301626795674-0.2036390034398219i]
  [-0.29309032295973036+0.3034337168992096i,-0.06392399629070344-0.8178102917845342i,0.06006452944412022-0.6141208421036348i]
  [0.34841978725201117+0.3778314407778909i,0.6768867572228499+0.7323625144544055i,-0.8246879334889017-0.009443253424316733i]
]
Conjugate Transpose: [
  [-0.7679977131543951+0.439979346567841i,-0.29309032295973036-0.3034337168992096i,0.34841978725201117-0.3778314407778909i]
  [-0.6011221529373452-0.510336881376179i,-0.06392399629070344+0.8178102917845342i,0.6768867572228499-0.7323625144544055i]
  [-0.22458301626795674+0.2036390034398219i,0.06006452944412022+0.6141208421036348i,-0.8246879334889017+0.009443253424316733i]
]
Hermitian: false
Normal: false
Unitary: false

Sidef

Translation of: Raku

func is_Hermitian (Array m, Array t) -> Bool { m == t }

func mat_mult (Array a, Array b, Number ε = -3) {
    var p = []
    for r, c in (^a ~X ^b[0]) {
        for k in (^b) {
            p[r][c] := 0 += (a[r][k] * b[k][c]) -> round!(ε)
        }
    }
    return p
}

func mat_trans (Array m) {
    var r = []
    for i,j in (^m ~X ^m[0]) {
        r[j][i] = m[i][j]
    }
    return r
}

func mat_ident (Number n) {
    ^n -> map {|i|
        [i.of(0)..., 1, (n - i - 1).of(0)...]
    }
}

func is_Normal (Array m, Array t) -> Bool {
    mat_mult(m, t) == mat_mult(t, m)
}

func is_Unitary (Array m, Array t) -> Bool {
    mat_mult(m, t) == mat_ident(m.len)
}

func say_it (Array a) {
    a.each {|b|
        b.map { "%9s" % _ }.join(' ').say
    }
}

[
    [
       [   1, 1+1i, 2i],
       [1-1i,    5, -3],
       [0-2i,   -3,  0]
    ],
    [
       [1, 1, 0],
       [0, 1, 1],
       [1, 0, 1]
    ],
    [
       [0.707 ,   0.707,  0],
       [0.707i, -0.707i,  0],
       [0     ,       0,  1i]
    ]
].each { |m|
    say "\nMatrix:"
    say_it(m)
    var t = mat_trans(m.map{.map{.conj}})
    say "\nTranspose:"
    say_it(t)
    say "Is Hermitian?\t#{is_Hermitian(m, t)}"
    say "Is Normal?\t#{is_Normal(m, t)}"
    say "Is Unitary?\t#{is_Unitary(m, t)}"
}

Output:

Matrix:
        1       1+i        2i
      1-i         5        -3
      -2i        -3         0

Transpose:
        1       1+i        2i
      1-i         5        -3
      -2i        -3         0
Is Hermitian?	true
Is Normal?	true
Is Unitary?	false

Matrix:
        1         1         0
        0         1         1
        1         0         1

Transpose:
        1         0         1
        1         1         0
        0         1         1
Is Hermitian?	false
Is Normal?	true
Is Unitary?	false

Matrix:
    0.707     0.707         0
   0.707i   -0.707i         0
        0         0         i

Transpose:
    0.707   -0.707i         0
    0.707    0.707i         0
        0         0        -i
Is Hermitian?	false
Is Normal?	true
Is Unitary?	true

Sparkling

Sparkling has support for basic complex algebraic operations, but complex matrix operations are not in the standard library.

# Computes conjugate transpose of M
let conjTransp = function conjTransp(M) {
	return map(range(sizeof M[0]), function(row) {
		return map(range(sizeof M), function(col) {
			return cplx_conj(M[col][row]);
		});
	});
};

# Helper for cplxMatMul
let cplxVecScalarMul = function cplxVecScalarMul(A, B, row, col) {
	var M = { "re": 0.0, "im": 0.0 };
	let N = sizeof A;
	for (var i = 0; i < N; i++) {
		let P = cplx_mul(A[row][i], B[i][col]);
		M = cplx_add(M, P);
	}
	return M;
};

# Multiplies matrices A and B
# A and B are assumed to be square and of the same size,
# this condition is not checked.
let cplxMatMul = function cplxMatMul(A, B) {
	var R = {};
	let N = sizeof A;
	for (var row = 0; row < N; row++) {
		R[row] = {};
		for (var col = 0; col < N; col++) {
			R[row][col] = cplxVecScalarMul(A, B, row, col);
		}
	}
	return R;
};

# Helper for creating an array representing a complex number
# given its textual representation
let _ = function makeComplex(str) {
	let sep = indexof(str, "+", 1);
	if sep < 0 {
		sep = indexof(str, "-", 1);
	}
	let reStr = substrto(str, sep);
	let imStr = substrfrom(str, sep);
	return { "re": tofloat(reStr), "im": tofloat(imStr) };
};

# Formats a complex matrix
let printCplxMat = function printCplxMat(M) {
	foreach(M, function(i, row) {
		foreach(row, function(j, elem) {
			printf("    %.2f%+.2fi", elem.re, elem.im);
		});
		print();
	});
};

# A Hermitian matrix
let H = {
	{ _("3+0i"), _("2+1i") },
	{ _("2-1i"), _("0+0i") }
};

# A normal matrix
let N = {
	{ _("1+0i"), _("1+0i"), _("0+0i") },
	{ _("0+0i"), _("1+0i"), _("1+0i") },
	{ _("1+0i"), _("0+0i"), _("1+0i") }
};

# A unitary matrix
let U = {
	{ _("0.70710678118+0i"), _("0.70710678118+0i"), _("0+0i") },
	{ _("0-0.70710678118i"), _("0+0.70710678118i"), _("0+0i") },
	{ _("0+0i"),             _("0+0i"),             _("0+1i") }
};


print("Hermitian matrix:\nH = ");
printCplxMat(H);
print("H* = ");
printCplxMat(conjTransp(H));
print();

print("Normal matrix:\nN = ");
printCplxMat(N);
print("N* = ");
printCplxMat(conjTransp(N));
print("N* x N = ");
printCplxMat(cplxMatMul(conjTransp(N), N));
print("N x N* = ");
printCplxMat(cplxMatMul(N, conjTransp(N)));
print();

print("Unitary matrix:\nU = ");
printCplxMat(U);
print("U* = ");
printCplxMat(conjTransp(U));
print("U x U* = ");
printCplxMat(cplxMatMul(U, conjTransp(U)));
print();

Stata

In Mata, the ' operator is always the conjugate transpose. To get only the matrix transpose without complex conjugate, use the [ transposeonly] function.

: a=1,2i\3i,4

: a
        1    2
    +-----------+
  1 |   1   2i  |
  2 |  3i    4  |
    +-----------+

: a'
         1     2
    +-------------+
  1 |    1   -3i  |
  2 |  -2i     4  |
    +-------------+

: transposeonly(a)
        1    2
    +-----------+
  1 |   1   3i  |
  2 |  2i    4  |
    +-----------+

: a*a'==a'*a
  0

: a'==a
  0

: a'*a==I(rows(a))
  0

Tcl

Tcl's matrixes (in Tcllib) do not assume that the contents are numeric at all. As such, they do not provide mathematical operations over them and this considerably increases the complexity of the code below. Note the use of lambda terms to simplify access to the complex number package.

Library: Tcllib (Package: math::complexnumbers)

Library: Tcllib (Package: struct::matrix)

package require struct::matrix
package require math::complexnumbers

proc complexMatrix.equal {m1 m2 {epsilon 1e-14}} {
    if {[$m1 rows] != [$m2 rows] || [$m1 columns] != [$m2 columns]} {
	return 0
    }
    # Compute the magnitude of the difference between two complex numbers
    set ceq [list apply {{epsilon a b} {
	expr {[mod [- $a $b]] < $epsilon}
    } ::math::complexnumbers} $epsilon]
    for {set i 0} {$i<[$m1 columns]} {incr i} {
	for {set j 0} {$j<[$m1 rows]} {incr j} {
	    if {![{*}$ceq [$m1 get cell $i $j] [$m2 get cell $i $j]]} {
		return 0
	    }
	}
    }
    return 1
}

proc complexMatrix.multiply {a b} {
    if {[$a columns] != [$b rows]} {
        error "incompatible sizes"
    }
    # Simplest to use a lambda in the complex NS
    set cpm {{sum a b} {
	+ $sum [* $a $b]
    } ::math::complexnumbers}
    set c0 [math::complexnumbers::complex 0.0 0.0];   # Complex zero
    set c [struct::matrix]
    $c add columns [$b columns]
    $c add rows [$a rows]
    for {set i 0} {$i < [$a rows]} {incr i} {
        for {set j 0} {$j < [$b columns]} {incr j} {
            set sum $c0
	    foreach rv [$a get row $i] cv [$b get column $j] {
		set sum [apply $cpm $sum $rv $cv]
            }
	    $c set cell $j $i $sum
        }
    }
    return $c
}

proc complexMatrix.conjugateTranspose {matrix} {
    set mat [struct::matrix]
    $mat = $matrix
    $mat transpose
    for {set c 0} {$c < [$mat columns]} {incr c} {
	for {set r 0} {$r < [$mat rows]} {incr r} {
	    set val [$mat get cell $c $r]
	    $mat set cell $c $r [math::complexnumbers::conj $val]
	}
    }
    return $mat
}

Using these tools to test for the properties described in the task:

proc isHermitian {matrix {epsilon 1e-14}} {
    if {[$matrix rows] != [$matrix columns]} {
	# Must be square!
	return 0
    }
    set cc [complexMatrix.conjugateTranspose $matrix]
    set result [complexMatrix.equal $matrix $cc $epsilon]
    $cc destroy
    return $result
}

proc isNormal {matrix {epsilon 1e-14}} {
    if {[$matrix rows] != [$matrix columns]} {
	# Must be square!
	return 0
    }
    set mh [complexMatrix.conjugateTranspose $matrix]
    set mhm [complexMatrix.multiply $mh $matrix]
    set mmh [complexMatrix.multiply $matrix $mh]
    $mh destroy
    set result [complexMatrix.equal $mhm $mmh $epsilon]
    $mhm destroy
    $mmh destroy
    return $result
}

proc isUnitary {matrix {epsilon 1e-14}} {
    if {[$matrix rows] != [$matrix columns]} {
	# Must be square!
	return 0
    }
    set mh [complexMatrix.conjugateTranspose $matrix]
    set mhm [complexMatrix.multiply $mh $matrix]
    set mmh [complexMatrix.multiply $matrix $mh]
    $mh destroy
    set result [complexMatrix.equal $mhm $mmh $epsilon]
    $mhm destroy
    if {$result} {
	set id [struct::matrix]
	$id = $matrix;   # Just for its dimensions
	for {set c 0} {$c < [$id columns]} {incr c} {
	    for {set r 0} {$r < [$id rows]} {incr r} {
		$id set cell $c $r \
		    [math::complexnumbers::complex [expr {$c==$r}] 0]
	    }
	}
	set result [complexMatrix.equal $mmh $id $epsilon]
	$id destroy
    }
    $mmh destroy
    return $result
}

Wren

Library: Wren-fmt

Library: 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.

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))")

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

For the final example if we use a tolerance of 1e-14:
Unitary   : true