shuffle complex vectors

shuffle.f90

 
 MODULE shuffle
   use globals
+  use utils
   implicit none
 
 CONTAINS
   
   FUNCTION append_to_each_row(a, m) result(res)
     ! appends element a to each row of m
-    integer :: a, m(:,:), i
-    integer :: res(size(m,1),size(m,2)+1)
+    complex(kind=prec) :: a, m(:,:)
+    integer :: i
+    complex(kind=prec) :: res(size(m,1),size(m,2)+1)
     do i=1,size(m,1)
       res(i,:) = [a,m(i,:)]
     end do
@@ -16,31 +18,25 @@ CONTAINS
 
   FUNCTION stack_matrices_vertically(m1, m2) result(res)
     ! appends to matrix m1 the rows of matrix m2
-    integer :: m1(:,:), m2(:,:)
-    integer :: res(size(m1,1)+size(m2,1), size(m1,2))
+    complex(kind=prec) :: m1(:,:), m2(:,:)
+    complex(kind=prec) :: res(size(m1,1)+size(m2,1), size(m1,2))
     res(1:size(m1,1), :) = m1
     res(size(m1,1)+1:size(res,1),:) = m2 
   END FUNCTION stack_matrices_vertically
-  
-  RECURSIVE FUNCTION factorial(n) result(res)
-    integer, intent(in) :: n
-    integer :: res
-    res = merge(1,n*factorial(n-1),n==0)
-  END FUNCTION factorial
 
   RECURSIVE FUNCTION shuffle_product(v1, v2) result(res)
-    integer :: v1(:), v2(:)
+    complex(kind=prec) :: v1(:), v2(:)
     integer :: i
-    integer :: res(product((/(i,i=1,size(v1)+size(v2))/))/  & 
+    complex(kind=prec) :: res(product((/(i,i=1,size(v1)+size(v2))/))/  & 
       (product((/(i,i=1,size(v1))/))*product((/(i,i=1,size(v2))/))), & 
       size(v1) + size(v2))
-    integer :: p1(product((/(i,i=1,size(v1)+size(v2)-1)/))/  & 
+    complex(kind=prec) :: p1(product((/(i,i=1,size(v1)+size(v2)-1)/))/  & 
       (product((/(i,i=1,size(v1)-1)/))*product((/(i,i=1,size(v2))/))), & 
       size(v1) + size(v2) - 1)
-    integer :: p2(product((/(i,i=1,size(v1)+size(v2)-1)/))/  & 
+    complex(kind=prec) :: p2(product((/(i,i=1,size(v1)+size(v2)-1)/))/  & 
       (product((/(i,i=1,size(v1))/))*product((/(i,i=1,size(v2)-1)/))), & 
       size(v1) + size(v2) - 1)
-    integer :: alpha, beta, w1(size(v1)-1), w2(size(v2)-1)
+    complex(kind=prec) :: alpha, beta, w1(size(v1)-1), w2(size(v2)-1)
 
     res = 0
     if(size(v1) == 0) then 
@@ -59,22 +55,23 @@ CONTAINS
     res = stack_matrices_vertically( &
       append_to_each_row(alpha, shuffle_product(w1, v2)), & 
       append_to_each_row(beta, shuffle_product(v1, w2)) )
+
   END FUNCTION shuffle_product
 
 END MODULE shuffle
 
-! PROGRAM test
-!   use shuffle_algebra
-!   implicit none
+PROGRAM test
+  use utils
+  use shuffle
+  implicit none
 
-!   integer :: v1(3), v2(3)
-!   integer :: amount_shuffles
-!   integer :: res(3,3)
+  complex(kind=prec) :: v1(3), v2(2)
+  integer :: amount_shuffles
 
-!   v1 = (/1,2,3/)
-!   v2 = (/-1,-2,-3/)
+  v1 = cmplx((/1,2,3/))
+  v2 = cmplx((/4,5/))
 
-!   call print_as_matrix(shuffle_product(v1, v2))
+  call print_matrix(shuffle_product(v1,v2))
 
-! END PROGRAM test
+END PROGRAM test
 

utils.f90

@@ -8,26 +8,12 @@
 MODULE utils
   use globals
   implicit none
+
   ! logical :: print_enabled = .true.
   ! logical :: warnings_enabled = .true.
 
-  INTERFACE print_matrix
-    ! prints 2d array as matrix. For complex it takes absolutes
-    
-    SUBROUTINE print_integer_matrix(m) 
-      integer :: m(:,:)
-      integer :: s(2)
-    END SUBROUTINE print_integer_matrix
-
-  SUBROUTINE print_complex_matrix(m) 
-    complex :: m(:,:)
-    integer :: s(2)
-  END SUBROUTINE print_complex_matrix
-
-  END INTERFACE print_matrix
-
 CONTAINS
-
+  
   FUNCTION  get_condensed_m(z) result(m)
     ! returns condensed m where the ones not needed are filled with 0
     complex(kind=prec), intent(in) :: z(:)
@@ -101,6 +87,12 @@ CONTAINS
     complex(kind=prec) :: res(n)
     res = 0
   END FUNCTION zero_array
+  
+  RECURSIVE FUNCTION factorial(n) result(res)
+    integer, intent(in) :: n
+    integer :: res
+    res = merge(1,n*factorial(n-1),n==0)
+  END FUNCTION factorial
 
   FUNCTION shuffle_with_zero(a) result(res)
     ! rows of result are shuffles of a with 0
@@ -118,6 +110,15 @@ CONTAINS
     end do
   END FUNCTION shuffle_with_zero
 
+  SUBROUTINE print_matrix(m) 
+    complex(kind=prec) :: m(:,:)
+    integer :: s(2), i
+    s = shape(m)
+    do i = 1,s(1)
+      print*, abs(m(i,:))
+    end do
+  END SUBROUTINE print_matrix
+
   ! subroutine print(s1,s2,s3,s4,s5)
   !   character(len = *), intent(in), optional :: s1, s2, s3, s4, s5
   !   if(print_enabled) then
@@ -135,24 +136,6 @@ CONTAINS
 END MODULE utils
 
 
-SUBROUTINE print_integer_matrix(m) 
-  integer :: m(:,:)
-  integer :: s(2), i
-  s = shape(m)
-  do i = 1,s(1)
-    print*, m(i,:)
-  end do
-END SUBROUTINE print_integer_matrix
-
-SUBROUTINE print_complex_matrix(m) 
-  complex :: m(:,:)
-  integer :: s(2), i
-  s = shape(m)
-  do i = 1,s(1)
-    print*, abs(m(i,:))
-  end do
-END SUBROUTINE print_complex_matrix
-
 ! PROGRAM test
 !   use  utils
 !   implicit none
@@ -172,5 +155,4 @@ END SUBROUTINE print_complex_matrix
 !     condensed_size = find_first_zero(m_prime)-1 
 !   end if
 !   print*, condensed_size
-
 ! END  PROGRAM test