From 2d7f4775d826fd08ebd4d8e99888ec745cc25aad Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Wed, 9 Jun 2021 16:44:39 -0500
Subject: [PATCH 01/15] Created function submodule and interface.

---
 src/stdlib_linalg.fypp               | 18 ++++++++++++++++++
 src/stdlib_linalg_outer_product.fypp | 20 ++++++++++++++++++++
 2 files changed, 38 insertions(+)
 create mode 100644 src/stdlib_linalg_outer_product.fypp

diff --git a/src/stdlib_linalg.fypp b/src/stdlib_linalg.fypp
index 51c2cdd54..fe65df716 100644
--- a/src/stdlib_linalg.fypp
+++ b/src/stdlib_linalg.fypp
@@ -11,6 +11,7 @@ module stdlib_linalg
   public :: diag
   public :: eye
   public :: trace
+  public :: outer_product
 
   interface diag
     !! version: experimental
@@ -52,6 +53,7 @@ module stdlib_linalg
     #:endfor
   end interface
 
+
   ! Matrix trace
   interface trace
     !! version: experimental
@@ -63,6 +65,22 @@ module stdlib_linalg
     #:endfor
   end interface
 
+
+  ! Outer product (of two vectors)
+  interface outer_product
+    !! version: experimental
+    !!
+    !! Computes outer product of two vectors, returning a matrix
+    !! ([Specification](../page/specs/stdlib_linalg.html#description_3))
+    #:for k1, t1 in RCI_KINDS_TYPES
+     module function outer_product_${t1[0]}$${k1}$(u,v) result(res)
+       ${t1}$, intent(in) :: u(:), v(:)
+       ${t1}$ :: res(size(u),size(v))
+       integer :: col
+     end function outer_product_${t1[0]}$${k1}$
+    #:endfor
+  end interface outer_product
+
 contains
 
     function eye(n) result(res)
diff --git a/src/stdlib_linalg_outer_product.fypp b/src/stdlib_linalg_outer_product.fypp
new file mode 100644
index 000000000..6d0803d27
--- /dev/null
+++ b/src/stdlib_linalg_outer_product.fypp
@@ -0,0 +1,20 @@
+#:include "common.fypp"
+#:set RCI_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES + INT_KINDS_TYPES
+submodule (stdlib_linalg) stdlib_linalg_outer_product
+
+  implicit none
+
+contains
+
+    #:for k1, t1 in RCI_KINDS_TYPES
+     module function outer_product_${t1[0]}$${k1}$(u,v) result(res)
+       ${t1}$, intent(in) :: u(:), v(:)
+       ${t1}$ :: res(size(u),size(v))
+       integer :: col
+       do col=1,size(v)
+           res(:,col) = v(col) * u
+       end do
+      end function outer_product_${t1[0]}$${k1}$
+    #:endfor
+
+end submodule

From a85310b30b6aa1d9859f8b63955cec09f2f84b56 Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Wed, 9 Jun 2021 17:11:01 -0500
Subject: [PATCH 02/15] One line addition of submodule file to make files

---
 src/CMakeLists.txt  | 1 +
 src/Makefile.manual | 1 +
 2 files changed, 2 insertions(+)

diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index fa7f6e639..ef1d7963c 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -9,6 +9,7 @@ set(fppFiles
     stdlib_io.fypp
     stdlib_linalg.fypp
     stdlib_linalg_diag.fypp
+    stdlib_linalg_outer_product.fypp
     stdlib_optval.fypp
     stdlib_sorting.fypp
     stdlib_sorting_ord_sort.fypp
diff --git a/src/Makefile.manual b/src/Makefile.manual
index e55bb80fd..f05b25822 100644
--- a/src/Makefile.manual
+++ b/src/Makefile.manual
@@ -6,6 +6,7 @@ SRCFYPP =\
         stdlib_io.fypp \
         stdlib_linalg.fypp \
         stdlib_linalg_diag.fypp \
+	stdlib_linalg_outer_product.fypp \
         stdlib_optval.fypp \
         stdlib_quadrature.fypp \
         stdlib_quadrature_trapz.fypp \

From 2a821864bd943c2ce9376e0a7db959ccd147ed99 Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Wed, 9 Jun 2021 21:34:29 -0500
Subject: [PATCH 03/15] Commit after successful build.

---
 src/tests/linalg/test_linalg.f90 | 20 +++++++++++++++++++-
 1 file changed, 19 insertions(+), 1 deletion(-)

diff --git a/src/tests/linalg/test_linalg.f90 b/src/tests/linalg/test_linalg.f90
index 4ad178d5b..467323099 100644
--- a/src/tests/linalg/test_linalg.f90
+++ b/src/tests/linalg/test_linalg.f90
@@ -2,7 +2,7 @@ program test_linalg
   
   use stdlib_error, only: check
   use stdlib_kinds, only: sp, dp, qp, int8, int16, int32, int64
-  use stdlib_linalg, only: diag, eye, trace
+  use stdlib_linalg, only: diag, eye, trace, outer_product
   
   implicit none
   
@@ -56,6 +56,24 @@ program test_linalg
   call test_trace_int32
   call test_trace_int64
 
+  !
+  ! outer product
+  !
+  !commented until outer_product compiles
+  !call test_outer_product_rsp
+  !call test_outer_product_rsp_k
+  !call test_outer_product_rdp
+  !call test_outer_product_rqp
+
+  !call test_outer_product_csp
+  !call test_outer_product_cdp
+  !call test_outer_product_cqp
+
+  !call test_outer_product_int8
+  !call test_outer_product_int16
+  !call test_outer_product_int32
+  !call test_outer_product_int64
+
 
 contains
 

From 64016f43c1a78035b202e8153fa21e121cf292cb Mon Sep 17 00:00:00 2001
From: Sebastian Ehlert <28669218+awvwgk@users.noreply.github.com>
Date: Sat, 12 Jun 2021 18:30:24 +0200
Subject: [PATCH 04/15] Fix manual makefile build

---
 src/Makefile.manual | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/src/Makefile.manual b/src/Makefile.manual
index f05b25822..f869cbd7f 100644
--- a/src/Makefile.manual
+++ b/src/Makefile.manual
@@ -6,7 +6,7 @@ SRCFYPP =\
         stdlib_io.fypp \
         stdlib_linalg.fypp \
         stdlib_linalg_diag.fypp \
-	stdlib_linalg_outer_product.fypp \
+        stdlib_linalg_outer_product.fypp \
         stdlib_optval.fypp \
         stdlib_quadrature.fypp \
         stdlib_quadrature_trapz.fypp \
@@ -129,3 +129,4 @@ stdlib_stats_distribution_PRNG.o: \
 stdlib_string_type.o: stdlib_ascii.o stdlib_kinds.o
 stdlib_strings.o: stdlib_ascii.o stdlib_string_type.o
 stdlib_math.o: stdlib_kinds.o
+stdlib_linalg_outer_product.o: stdlib_linalg.o

From 89f98efdd1ea2a03e1fb295662759229188969f2 Mon Sep 17 00:00:00 2001
From: Milan Curcic <caomaco@gmail.com>
Date: Sat, 12 Jun 2021 19:21:48 -0700
Subject: [PATCH 05/15] Update src/stdlib_linalg_outer_product.fypp

Co-authored-by: Ivan Pribec <ivan.pribec@gmail.com>
---
 src/stdlib_linalg_outer_product.fypp | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/stdlib_linalg_outer_product.fypp b/src/stdlib_linalg_outer_product.fypp
index 6d0803d27..01f9a5352 100644
--- a/src/stdlib_linalg_outer_product.fypp
+++ b/src/stdlib_linalg_outer_product.fypp
@@ -11,7 +11,7 @@ contains
        ${t1}$, intent(in) :: u(:), v(:)
        ${t1}$ :: res(size(u),size(v))
        integer :: col
-       do col=1,size(v)
+       do col = 1, size(v)
            res(:,col) = v(col) * u
        end do
       end function outer_product_${t1[0]}$${k1}$

From cc336705384e270edecff6a1f32f22e512c37633 Mon Sep 17 00:00:00 2001
From: milancurcic <caomaco@gmail.com>
Date: Sat, 12 Jun 2021 22:37:41 -0400
Subject: [PATCH 06/15] minor style fixes

---
 src/stdlib_linalg.fypp               | 10 +++++-----
 src/stdlib_linalg_outer_product.fypp | 20 ++++++++++----------
 2 files changed, 15 insertions(+), 15 deletions(-)

diff --git a/src/stdlib_linalg.fypp b/src/stdlib_linalg.fypp
index fe65df716..192d68908 100644
--- a/src/stdlib_linalg.fypp
+++ b/src/stdlib_linalg.fypp
@@ -73,11 +73,11 @@ module stdlib_linalg
     !! Computes outer product of two vectors, returning a matrix
     !! ([Specification](../page/specs/stdlib_linalg.html#description_3))
     #:for k1, t1 in RCI_KINDS_TYPES
-     module function outer_product_${t1[0]}$${k1}$(u,v) result(res)
-       ${t1}$, intent(in) :: u(:), v(:)
-       ${t1}$ :: res(size(u),size(v))
-       integer :: col
-     end function outer_product_${t1[0]}$${k1}$
+      module function outer_product_${t1[0]}$${k1}$(u, v) result(res)
+        ${t1}$, intent(in) :: u(:), v(:)
+        ${t1}$ :: res(size(u),size(v))
+        integer :: col
+      end function outer_product_${t1[0]}$${k1}$
     #:endfor
   end interface outer_product
 
diff --git a/src/stdlib_linalg_outer_product.fypp b/src/stdlib_linalg_outer_product.fypp
index 01f9a5352..34c3f32bf 100644
--- a/src/stdlib_linalg_outer_product.fypp
+++ b/src/stdlib_linalg_outer_product.fypp
@@ -6,15 +6,15 @@ submodule (stdlib_linalg) stdlib_linalg_outer_product
 
 contains
 
-    #:for k1, t1 in RCI_KINDS_TYPES
-     module function outer_product_${t1[0]}$${k1}$(u,v) result(res)
-       ${t1}$, intent(in) :: u(:), v(:)
-       ${t1}$ :: res(size(u),size(v))
-       integer :: col
-       do col = 1, size(v)
-           res(:,col) = v(col) * u
-       end do
-      end function outer_product_${t1[0]}$${k1}$
-    #:endfor
+  #:for k1, t1 in RCI_KINDS_TYPES
+    module function outer_product_${t1[0]}$${k1}$(u, v) result(res)
+      ${t1}$, intent(in) :: u(:), v(:)
+      ${t1}$ :: res(size(u),size(v))
+      integer :: col
+      do col = 1, size(v)
+        res(:,col) = v(col) * u
+      end do
+    end function outer_product_${t1[0]}$${k1}$
+  #:endfor
 
 end submodule

From 31c98ae46ace42263e73c367bfddc6d2d14c218a Mon Sep 17 00:00:00 2001
From: milancurcic <caomaco@gmail.com>
Date: Sat, 12 Jun 2021 22:43:04 -0400
Subject: [PATCH 07/15] remove col from the outer_product interface (thanks
 @ivan-pi)

---
 src/stdlib_linalg.fypp | 1 -
 1 file changed, 1 deletion(-)

diff --git a/src/stdlib_linalg.fypp b/src/stdlib_linalg.fypp
index 192d68908..0c1ec48d7 100644
--- a/src/stdlib_linalg.fypp
+++ b/src/stdlib_linalg.fypp
@@ -76,7 +76,6 @@ module stdlib_linalg
       module function outer_product_${t1[0]}$${k1}$(u, v) result(res)
         ${t1}$, intent(in) :: u(:), v(:)
         ${t1}$ :: res(size(u),size(v))
-        integer :: col
       end function outer_product_${t1[0]}$${k1}$
     #:endfor
   end interface outer_product

From 6b0b46e3dcb64314307d69378c8682f5f1ef0e17 Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Sun, 13 Jun 2021 17:51:07 -0500
Subject: [PATCH 08/15] Added tests and confirmed build

---
 src/tests/linalg/test_linalg.f90 | 192 +++++++++++++++++++++++++------
 1 file changed, 156 insertions(+), 36 deletions(-)

diff --git a/src/tests/linalg/test_linalg.f90 b/src/tests/linalg/test_linalg.f90
index 467323099..507ad3b78 100644
--- a/src/tests/linalg/test_linalg.f90
+++ b/src/tests/linalg/test_linalg.f90
@@ -59,18 +59,16 @@ program test_linalg
   !
   ! outer product
   !
-  !commented until outer_product compiles
-  !call test_outer_product_rsp
-  !call test_outer_product_rsp_k
-  !call test_outer_product_rdp
-  !call test_outer_product_rqp
-
-  !call test_outer_product_csp
-  !call test_outer_product_cdp
-  !call test_outer_product_cqp
-
-  !call test_outer_product_int8
-  !call test_outer_product_int16
+  call test_outer_product_rsp
+  call test_outer_product_rdp
+  call test_outer_product_rqp
+
+  call test_outer_product_csp
+  call test_outer_product_cdp
+  call test_outer_product_cqp
+
+  call test_outer_product_int8
+  call test_outer_product_int16
   !call test_outer_product_int32
   !call test_outer_product_int64
 
@@ -93,7 +91,7 @@ subroutine test_eye
     cye = eye(7)
     call check(abs(trace(cye) - cmplx(7.0_sp,0.0_sp,kind=sp)) < sptol, &
           msg="abs(trace(cye) - cmplx(7.0_sp,0.0_sp,kind=sp)) < sptol failed.",warn=warn)
-  end subroutine
+  end subroutine test_eye
 
   subroutine test_diag_rsp
     integer, parameter :: n = 3
@@ -108,7 +106,7 @@ subroutine test_diag_rsp
 
     call check(all(diag(3*a) == 3*v), &
             msg="all(diag(3*a) == 3*v) failed.",warn=warn)
-  end subroutine
+  end subroutine test_diag_rsp
 
   subroutine test_diag_rsp_k
     integer, parameter :: n = 4
@@ -136,7 +134,7 @@ subroutine test_diag_rsp_k
     end do
     call check(size(diag(a,n+1)) == 0, &
         msg="size(diag(a,n+1)) == 0 failed.",warn=warn)
-  end subroutine
+  end subroutine test_diag_rsp_k
 
   subroutine test_diag_rdp
     integer, parameter :: n = 3
@@ -151,7 +149,7 @@ subroutine test_diag_rdp
 
     call check(all(diag(3*a) == 3*v), &
             msg="all(diag(3*a) == 3*v) failed.",warn=warn)
-  end subroutine
+  end subroutine test_diag_rdp
 
   subroutine test_diag_rqp
     integer, parameter :: n = 3
@@ -166,7 +164,7 @@ subroutine test_diag_rqp
 
     call check(all(diag(3*a) == 3*v), &
             msg="all(diag(3*a) == 3*v) failed.", warn=warn)
-  end subroutine
+  end subroutine test_diag_rqp
 
   subroutine test_diag_csp
     integer, parameter :: n = 3
@@ -183,7 +181,7 @@ subroutine test_diag_csp
           msg="all(abs(real(diag(a)) - [(i,i=1,n)]) < sptol)", warn=warn)
     call check(all(abs(aimag(diag(a)) - [(1,i=1,n)]) < sptol), &
           msg="all(abs(aimag(diag(a)) - [(1,i=1,n)]) < sptol)", warn=warn)
-  end subroutine
+  end subroutine test_diag_csp
 
   subroutine test_diag_cdp
     integer, parameter :: n = 3
@@ -193,7 +191,7 @@ subroutine test_diag_cdp
     a = diag([i_],-2) + diag([i_],2)
     call check(a(3,1) == i_ .and. a(1,3) == i_, &
           msg="a(3,1) == i_ .and. a(1,3) == i_ failed.",warn=warn)
-  end subroutine
+  end subroutine test_diag_cdp
 
   subroutine test_diag_cqp
     integer, parameter :: n = 3
@@ -203,7 +201,7 @@ subroutine test_diag_cqp
     a = diag([i_,i_],-1) + diag([i_,i_],1)
     call check(all(diag(a,-1) == i_) .and. all(diag(a,1) == i_), &
           msg="all(diag(a,-1) == i_) .and. all(diag(a,1) == i_) failed.",warn=warn)
-  end subroutine
+  end subroutine test_diag_cqp
 
   subroutine test_diag_int8
     integer, parameter :: n = 3
@@ -217,7 +215,7 @@ subroutine test_diag_int8
           msg="all(diag(a) == pack(a,mask)) failed.", warn=warn)
     call check(all(diag(diag(a)) == merge(a,0_int8,mask)), &
           msg="all(diag(diag(a)) == merge(a,0_int8,mask)) failed.", warn=warn)
-  end subroutine
+  end subroutine test_diag_int8
   subroutine test_diag_int16
     integer, parameter :: n = 4
     integer(int16), allocatable :: a(:,:)
@@ -230,7 +228,7 @@ subroutine test_diag_int16
           msg="all(diag(a) == pack(a,mask))", warn=warn)
     call check(all(diag(diag(a)) == merge(a,0_int16,mask)), &
           msg="all(diag(diag(a)) == merge(a,0_int16,mask)) failed.", warn=warn)
-  end subroutine
+  end subroutine test_diag_int16
   subroutine test_diag_int32
     integer, parameter :: n = 3
     integer(int32) :: a(n,n)
@@ -244,7 +242,7 @@ subroutine test_diag_int32
           msg="all(diag([1,1],-1) == a) failed.", warn=warn)
     call check(all(diag([1,1],1) == transpose(a)), &
           msg="all(diag([1,1],1) == transpose(a)) failed.", warn=warn)
-  end subroutine
+  end subroutine test_diag_int32
   subroutine test_diag_int64
     integer, parameter :: n = 4
     integer(int64) :: a(n,n), c(0:2*n-1)
@@ -275,7 +273,7 @@ subroutine test_diag_int64
     end do
     call check(all(diag(a,-2) == diag(a,2)), &
           msg="all(diag(a,-2) == diag(a,2))", warn=warn)
-  end subroutine
+  end subroutine test_diag_int64
 
 
 
@@ -288,7 +286,7 @@ subroutine test_trace_rsp
     a = reshape([(i,i=1,n**2)],[n,n])
     call check(abs(trace(a) - sum(diag(a))) < sptol, &
           msg="abs(trace(a) - sum(diag(a))) < sptol failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_rsp
 
   subroutine test_trace_rsp_nonsquare
     integer, parameter :: n = 4
@@ -305,7 +303,7 @@ subroutine test_trace_rsp_nonsquare
 
     call check(abs(trace(a) - ans) < sptol, &
           msg="abs(trace(a) - ans) < sptol failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_rsp_nonsquare
 
   subroutine test_trace_rdp
     integer, parameter :: n = 4
@@ -315,7 +313,7 @@ subroutine test_trace_rdp
     a = reshape([(i,i=1,n**2)],[n,n])
     call check(abs(trace(a) - sum(diag(a))) < dptol, &
           msg="abs(trace(a) - sum(diag(a))) < dptol failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_rdp
 
   subroutine test_trace_rdp_nonsquare
     integer, parameter :: n = 4
@@ -332,7 +330,7 @@ subroutine test_trace_rdp_nonsquare
 
     call check(abs(trace(a) - ans) < dptol, &
           msg="abs(trace(a) - ans) < dptol failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_rdp_nonsquare
 
   subroutine test_trace_rqp
     integer, parameter :: n = 3
@@ -342,7 +340,7 @@ subroutine test_trace_rqp
     a = reshape([(i,i=1,n**2)],[n,n])
     call check(abs(trace(a) - sum(diag(a))) < qptol, &
           msg="abs(trace(a) - sum(diag(a))) < qptol failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_rqp
 
 
   subroutine test_trace_csp
@@ -363,7 +361,7 @@ subroutine test_trace_csp
     ! tr(A + B) = tr(A) + tr(B)
     call check(abs(trace(a+b) - (trace(a) + trace(b))) < sptol, &
           msg="abs(trace(a+b) - (trace(a) + trace(b))) < sptol failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_csp
 
   subroutine test_trace_cdp
     integer, parameter :: n = 3
@@ -377,7 +375,7 @@ subroutine test_trace_cdp
 
     call check(abs(trace(a) - ans) < dptol, &
           msg="abs(trace(a) - ans) < dptol failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_cdp
 
   subroutine test_trace_cqp
     integer, parameter :: n = 3
@@ -387,7 +385,7 @@ subroutine test_trace_cqp
     a = 3*eye(n) + 4*eye(n)*i_ ! pythagorean triple
     call check(abs(trace(a)) - 3*5.0_qp < qptol, &
           msg="abs(trace(a)) - 3*5.0_qp < qptol failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_cqp
 
 
   subroutine test_trace_int8
@@ -398,7 +396,7 @@ subroutine test_trace_int8
     a = reshape([(i**2,i=1,n**2)],[n,n])
     call check(trace(a) == (1 + 25 + 81), &
           msg="trace(a) == (1 + 25 + 81) failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_int8
 
   subroutine test_trace_int16
     integer, parameter :: n = 3
@@ -408,7 +406,7 @@ subroutine test_trace_int16
     a = reshape([(i**3,i=1,n**2)],[n,n])
     call check(trace(a) == (1 + 125 + 729), &
           msg="trace(a) == (1 + 125 + 729) failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_int16
 
   subroutine test_trace_int32
     integer, parameter :: n = 3
@@ -418,7 +416,7 @@ subroutine test_trace_int32
     a = reshape([(i**4,i=1,n**2)],[n,n])
     call check(trace(a) == (1 + 625 + 6561), &
           msg="trace(a) == (1 + 625 + 6561) failed.",warn=warn)
-  end subroutine
+  end subroutine test_trace_int32
 
   subroutine test_trace_int64
     integer, parameter :: n = 5
@@ -442,7 +440,129 @@ subroutine test_trace_int64
     call check(trace(h) == sum(c(0:nd:2)), &
           msg="trace(h) == sum(c(0:nd:2)) failed.",warn=warn)
 
-  end subroutine
+  end subroutine test_trace_int64
+
+
+  subroutine test_outer_product_rsp
+    integer, parameter :: n = 2
+    real(sp) :: u(n), v(n), expected(n,n), diff(n,n)
+    write(*,*) "test_outer_product_rsp" 
+    u = [1.,2.]
+    v = [1.,3.]
+    expected = reshape([1.,2.,3.,6.],[n,n])
+    diff = expected - outer_product(u,v)
+    call check(all(abs(diff) < sptol), &
+         msg="all(abs(diff) < sptol) failed.",warn=warn)
+  end subroutine test_outer_product_rsp
+
+  subroutine test_outer_product_rdp
+    integer, parameter :: n = 2
+    real(dp) :: u(n), v(n), expected(n,n), diff(n,n)
+    write(*,*) "test_outer_product_rdp" 
+    u = [1.,2.]
+    v = [1.,3.]
+    expected = reshape([1.,2.,3.,6.],[n,n])
+    diff = expected - outer_product(u,v)
+    call check(all(abs(diff) < dptol), &
+         msg="all(abs(diff) < dptol) failed.",warn=warn)
+  end subroutine test_outer_product_rdp
+
+  subroutine test_outer_product_rqp
+    integer, parameter :: n = 2
+    real(qp) :: u(n), v(n), expected(n,n), diff(n,n)
+    write(*,*) "test_outer_product_rqp" 
+    u = [1.,2.]
+    v = [1.,3.]
+    expected = reshape([1.,2.,3.,6.],[n,n])
+    diff = expected - outer_product(u,v)
+    call check(all(abs(diff) < qptol), &
+         msg="all(abs(diff) < qptol) failed.",warn=warn)
+  end subroutine test_outer_product_rqp
+
+  subroutine test_outer_product_csp
+    integer, parameter :: n = 2
+    complex(sp) :: u(n), v(n), expected(n,n), diff(n,n)
+    write(*,*) "test_outer_product_csp" 
+    u = [cmplx(1.,1.),cmplx(2.,0.)]
+    v = [cmplx(1.,0.),cmplx(3.,1.)]
+    expected = reshape([cmplx(1.,1.),cmplx(2.,0.),cmplx(2.,4.),cmplx(6.,2.)],[n,n])
+    diff = expected - outer_product(u,v)
+    call check(all(abs(diff) < sptol), &
+         msg="all(abs(diff) < sptol) failed.",warn=warn)
+  end subroutine test_outer_product_csp
+
+  subroutine test_outer_product_cdp
+    integer, parameter :: n = 2
+    complex(dp) :: u(n), v(n), expected(n,n), diff(n,n)
+    write(*,*) "test_outer_product_cdp" 
+    u = [cmplx(1.,1.),cmplx(2.,0.)]
+    v = [cmplx(1.,0.),cmplx(3.,1.)]
+    expected = reshape([cmplx(1.,1.),cmplx(2.,0.),cmplx(2.,4.),cmplx(6.,2.)],[n,n])
+    diff = expected - outer_product(u,v)
+    call check(all(abs(diff) < dptol), &
+         msg="all(abs(diff) < dptol) failed.",warn=warn)
+  end subroutine test_outer_product_cdp
+
+  subroutine test_outer_product_cqp
+    integer, parameter :: n = 2
+    complex(qp) :: u(n), v(n), expected(n,n), diff(n,n)
+    write(*,*) "test_outer_product_cqp" 
+    u = [cmplx(1.,1.),cmplx(2.,0.)]
+    v = [cmplx(1.,0.),cmplx(3.,1.)]
+    expected = reshape([cmplx(1.,1.),cmplx(2.,0.),cmplx(2.,4.),cmplx(6.,2.)],[n,n])
+    diff = expected - outer_product(u,v)
+    call check(all(abs(diff) < qptol), &
+         msg="all(abs(diff) < qptol) failed.",warn=warn)
+  end subroutine test_outer_product_cqp
+    
+  subroutine test_outer_product_int8
+    integer, parameter :: n = 2
+    integer(int8) :: u(n), v(n), expected(n,n), diff(n,n)
+    write(*,*) "test_outer_product_int8" 
+    u = [1,2]
+    v = [1,3]
+    expected = reshape([1,2,3,6],[n,n])
+    diff = expected - outer_product(u,v)
+    call check(all(abs(diff) == 0), &
+         msg="all(abs(diff) == 0) failed.",warn=warn)
+  end subroutine test_outer_product_int8
+
+  subroutine test_outer_product_int16
+    integer, parameter :: n = 2
+    integer(int16) :: u(n), v(n), expected(n,n), diff(n,n)
+    write(*,*) "test_outer_product_int16" 
+    u = [1,2]
+    v = [1,3]
+    expected = reshape([1,2,3,6],[n,n])
+    diff = expected - outer_product(u,v)
+    call check(all(abs(diff) == 0), &
+         msg="all(abs(diff) == 0) failed.",warn=warn)
+  end subroutine test_outer_product_int16
+
+  subroutine test_outer_product_int32
+    integer, parameter :: n = 2
+    integer(int32) :: u(n), v(n), expected(n,n), diff(n,n)
+    write(*,*) "test_outer_product_int32" 
+    u = [1,2]
+    v = [1,3]
+    expected = reshape([1,2,3,6],[n,n])
+    diff = expected - outer_product(u,v)
+    call check(all(abs(diff) == 0), &
+         msg="all(abs(diff) == 0) failed.",warn=warn)
+  end subroutine test_outer_product_int32
+
+  subroutine test_outer_product_int64
+    integer, parameter :: n = 2
+    integer(int64) :: u(n), v(n), expected(n,n), diff(n,n)
+    write(*,*) "test_outer_product_int64" 
+    u = [1,2]
+    v = [1,3]
+    expected = reshape([1,2,3,6],[n,n])
+    diff = expected - outer_product(u,v)
+    call check(all(abs(diff) == 0), &
+         msg="all(abs(diff) == 0) failed.",warn=warn)
+  end subroutine test_outer_product_int64
+
 
   pure recursive function catalan_number(n) result(value)
     integer, intent(in) :: n

From dcb5f01781fe19bb5bad2a53760d70c287f8846e Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Sun, 13 Jun 2021 17:57:06 -0500
Subject: [PATCH 09/15] Uncommented two tests and rebuilt

---
 src/tests/linalg/test_linalg.f90 | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/src/tests/linalg/test_linalg.f90 b/src/tests/linalg/test_linalg.f90
index 507ad3b78..cc8d0db68 100644
--- a/src/tests/linalg/test_linalg.f90
+++ b/src/tests/linalg/test_linalg.f90
@@ -69,8 +69,8 @@ program test_linalg
 
   call test_outer_product_int8
   call test_outer_product_int16
-  !call test_outer_product_int32
-  !call test_outer_product_int64
+  call test_outer_product_int32
+  call test_outer_product_int64
 
 
 contains

From 58abf5ff34b963b76c82ebe0671a3b1a4f7d1228 Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Sun, 13 Jun 2021 20:59:05 -0500
Subject: [PATCH 10/15] Fixed grammar in documentation line

---
 src/stdlib_linalg.fypp | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/stdlib_linalg.fypp b/src/stdlib_linalg.fypp
index 0c1ec48d7..53b2ea4bb 100644
--- a/src/stdlib_linalg.fypp
+++ b/src/stdlib_linalg.fypp
@@ -70,7 +70,7 @@ module stdlib_linalg
   interface outer_product
     !! version: experimental
     !!
-    !! Computes outer product of two vectors, returning a matrix
+    !! Computes the outer product of two vectors, returning a matrix
     !! ([Specification](../page/specs/stdlib_linalg.html#description_3))
     #:for k1, t1 in RCI_KINDS_TYPES
       module function outer_product_${t1[0]}$${k1}$(u, v) result(res)

From 03b2dc3ff1078d9e3183f777e4eae8dfb25ae1d7 Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Mon, 14 Jun 2021 00:03:09 -0500
Subject: [PATCH 11/15] Added specs documentation.

---
 doc/specs/stdlib_linalg.md | 38 ++++++++++++++++++++++++++++++++++++++
 1 file changed, 38 insertions(+)

diff --git a/doc/specs/stdlib_linalg.md b/doc/specs/stdlib_linalg.md
index 6aa76fa88..902d54cd4 100644
--- a/doc/specs/stdlib_linalg.md
+++ b/doc/specs/stdlib_linalg.md
@@ -168,3 +168,41 @@ program demo_trace
     print *, trace(A) ! 1 + 5 + 9
 end program demo_trace
 ```
+
+## `outer_product` - Computes the outer product of two vectors
+
+### Status
+
+Experimental
+
+### Description
+
+Create a diagonal array or extract the diagonal elements of an array
+
+### Syntax
+
+`d = [[stdlib_linalg(module):outer_product(interface)]](u, v)`
+
+### Arguments
+
+`u`: Shall be a rank-1 array
+
+`v`: Shall be a rank-1 array
+
+### Return value
+
+Returns a rank-2 array equal to `u v^T` (assuming `u, v` are column vectors). The dimensions of the returned array are `[size(u), size(v)]`
+
+### Example
+
+```fortran
+program demo_outer_product
+    use stdlib_linalg, only: outer_product
+    implicit none
+    real, allocatable :: A(:,:), u(:), v(:)
+    u = [1.,2.,3.]
+    v = [3., 4.]
+    A=outer_product(u,v)
+    !A = transpose(reshape([3., 4., 6., 8., 9., 12.], [3,2]))
+end program demo_outer_product
+```
\ No newline at end of file

From 45e9458a6512354fe589cd9344d786f5b0826060 Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Mon, 14 Jun 2021 09:05:40 -0500
Subject: [PATCH 12/15] Fixed errors and corrected style in specs doc

---
 doc/specs/stdlib_linalg.md | 10 +++++-----
 1 file changed, 5 insertions(+), 5 deletions(-)

diff --git a/doc/specs/stdlib_linalg.md b/doc/specs/stdlib_linalg.md
index 902d54cd4..c755dfb38 100644
--- a/doc/specs/stdlib_linalg.md
+++ b/doc/specs/stdlib_linalg.md
@@ -177,7 +177,7 @@ Experimental
 
 ### Description
 
-Create a diagonal array or extract the diagonal elements of an array
+Computes the outer product of two vectors
 
 ### Syntax
 
@@ -191,7 +191,7 @@ Create a diagonal array or extract the diagonal elements of an array
 
 ### Return value
 
-Returns a rank-2 array equal to `u v^T` (assuming `u, v` are column vectors). The dimensions of the returned array are `[size(u), size(v)]`
+Returns a rank-2 array equal to `u v^T` (assuming `u, v` are column vectors). The shape of the returned array is `[size(u), size(v)]`.
 
 ### Example
 
@@ -200,9 +200,9 @@ program demo_outer_product
     use stdlib_linalg, only: outer_product
     implicit none
     real, allocatable :: A(:,:), u(:), v(:)
-    u = [1.,2.,3.]
+    u = [1. ,2. ,3. ]
     v = [3., 4.]
-    A=outer_product(u,v)
-    !A = transpose(reshape([3., 4., 6., 8., 9., 12.], [3,2]))
+    A = outer_product(u,v)
+    !A = reshape([3., 6., 9., 4., 8., 12.], [3,2])
 end program demo_outer_product
 ```
\ No newline at end of file

From bdc0ed4dc52ee09ae12430076664168bbbadce3e Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Mon, 21 Jun 2021 09:04:11 -0500
Subject: [PATCH 13/15] Fixed spacing in specs docs

Co-authored-by: Jeremie Vandenplas <jeremie.vandenplas@gmail.com>
---
 doc/specs/stdlib_linalg.md | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/doc/specs/stdlib_linalg.md b/doc/specs/stdlib_linalg.md
index c755dfb38..0430ab1e3 100644
--- a/doc/specs/stdlib_linalg.md
+++ b/doc/specs/stdlib_linalg.md
@@ -200,9 +200,9 @@ program demo_outer_product
     use stdlib_linalg, only: outer_product
     implicit none
     real, allocatable :: A(:,:), u(:), v(:)
-    u = [1. ,2. ,3. ]
+    u = [1., 2., 3. ]
     v = [3., 4.]
     A = outer_product(u,v)
     !A = reshape([3., 6., 9., 4., 8., 12.], [3,2])
 end program demo_outer_product
-```
\ No newline at end of file
+```

From 6aac9cade787e825a15030108c91819ef5cf384a Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Mon, 21 Jun 2021 09:04:56 -0500
Subject: [PATCH 14/15] matix -> rank-2 array in specs docs

Co-authored-by: Jeremie Vandenplas <jeremie.vandenplas@gmail.com>
---
 src/stdlib_linalg.fypp | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/stdlib_linalg.fypp b/src/stdlib_linalg.fypp
index 53b2ea4bb..b258dc084 100644
--- a/src/stdlib_linalg.fypp
+++ b/src/stdlib_linalg.fypp
@@ -70,7 +70,7 @@ module stdlib_linalg
   interface outer_product
     !! version: experimental
     !!
-    !! Computes the outer product of two vectors, returning a matrix
+    !! Computes the outer product of two vectors, returning a rank-2 array
     !! ([Specification](../page/specs/stdlib_linalg.html#description_3))
     #:for k1, t1 in RCI_KINDS_TYPES
       module function outer_product_${t1[0]}$${k1}$(u, v) result(res)

From e852bb55ebf803e35fbe68bbafa36f494d2c64f6 Mon Sep 17 00:00:00 2001
From: GHBrown <gabriel.h.brown@gmail.com>
Date: Wed, 23 Jun 2021 20:40:47 -0500
Subject: [PATCH 15/15] Fixed docs wording and made pure function

---
 doc/specs/stdlib_linalg.md           | 2 +-
 src/stdlib_linalg.fypp               | 2 +-
 src/stdlib_linalg_outer_product.fypp | 2 +-
 3 files changed, 3 insertions(+), 3 deletions(-)

diff --git a/doc/specs/stdlib_linalg.md b/doc/specs/stdlib_linalg.md
index 0430ab1e3..cab16279c 100644
--- a/doc/specs/stdlib_linalg.md
+++ b/doc/specs/stdlib_linalg.md
@@ -191,7 +191,7 @@ Computes the outer product of two vectors
 
 ### Return value
 
-Returns a rank-2 array equal to `u v^T` (assuming `u, v` are column vectors). The shape of the returned array is `[size(u), size(v)]`.
+Returns a rank-2 array equal to `u v^T` (where `u, v` are considered column vectors). The shape of the returned array is `[size(u), size(v)]`.
 
 ### Example
 
diff --git a/src/stdlib_linalg.fypp b/src/stdlib_linalg.fypp
index b258dc084..5e0388c0b 100644
--- a/src/stdlib_linalg.fypp
+++ b/src/stdlib_linalg.fypp
@@ -73,7 +73,7 @@ module stdlib_linalg
     !! Computes the outer product of two vectors, returning a rank-2 array
     !! ([Specification](../page/specs/stdlib_linalg.html#description_3))
     #:for k1, t1 in RCI_KINDS_TYPES
-      module function outer_product_${t1[0]}$${k1}$(u, v) result(res)
+      pure module function outer_product_${t1[0]}$${k1}$(u, v) result(res)
         ${t1}$, intent(in) :: u(:), v(:)
         ${t1}$ :: res(size(u),size(v))
       end function outer_product_${t1[0]}$${k1}$
diff --git a/src/stdlib_linalg_outer_product.fypp b/src/stdlib_linalg_outer_product.fypp
index 34c3f32bf..26c726435 100644
--- a/src/stdlib_linalg_outer_product.fypp
+++ b/src/stdlib_linalg_outer_product.fypp
@@ -7,7 +7,7 @@ submodule (stdlib_linalg) stdlib_linalg_outer_product
 contains
 
   #:for k1, t1 in RCI_KINDS_TYPES
-    module function outer_product_${t1[0]}$${k1}$(u, v) result(res)
+    pure module function outer_product_${t1[0]}$${k1}$(u, v) result(res)
       ${t1}$, intent(in) :: u(:), v(:)
       ${t1}$ :: res(size(u),size(v))
       integer :: col