diff --git a/example/linalg/example_repmat.f90 b/example/linalg/example_repmat.f90
new file mode 100644
index 000000000..c907ec879
--- /dev/null
+++ b/example/linalg/example_repmat.f90
@@ -0,0 +1,12 @@
+program example_repmat
+ use stdlib_linalg, only: repmat
+ implicit none
+ real, allocatable :: a(:, :),b(:,:)
+ a = reshape([1., 2., 3., 4.],[2, 2],order=[2, 1])
+ b = repmat(a, 2, 2)
+! A = reshape([1., 2., 1., 2.,&
+! 3., 4., 3., 4.,&
+! 1., 2., 1., 2.,&
+! 3., 4., 3., 4.,&
+! ], [4, 4],order=[2, 1])
+end program example_repmat
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index 0fb95a2d3..a056e791b 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -22,6 +22,7 @@ set(fppFiles
stdlib_linalg.fypp
stdlib_linalg_diag.fypp
stdlib_linalg_outer_product.fypp
+ stdlib_linalg_repmat.fypp
stdlib_optval.fypp
stdlib_selection.fypp
stdlib_sorting.fypp
diff --git a/src/stdlib_linalg.fypp b/src/stdlib_linalg.fypp
index bc1017f0a..d064c948e 100644
--- a/src/stdlib_linalg.fypp
+++ b/src/stdlib_linalg.fypp
@@ -14,6 +14,7 @@ module stdlib_linalg
public :: eye
public :: trace
public :: outer_product
+ public :: repmat
public :: is_square
public :: is_diagonal
public :: is_symmetric
@@ -92,6 +93,21 @@ module stdlib_linalg
#:endfor
end interface outer_product
+ ! repeat matrix (of 2d array)
+ interface repmat
+ !! version: experimental
+ !!
+ !! Creates large matrices from a small array, `repmat()` repeats the given values of the array to create the large matrix.
+ !! ([Specification](../page/specs/stdlib_linalg.html#
+ !! repmat-creates-large-matrices-from-a-small-array))
+ #:for k1, t1 in RCI_KINDS_TYPES
+ pure module function repmat_${t1[0]}$${k1}$(a,m,n) result(res)
+ ${t1}$, intent(in) :: a(:,:)
+ integer, intent(in) :: m,n
+ ${t1}$ :: res(m*size(a,1),n*size(a,2))
+ end function repmat_${t1[0]}$${k1}$
+ #:endfor
+ end interface repmat
! Check for squareness
interface is_square
diff --git a/src/stdlib_linalg_repmat.fypp b/src/stdlib_linalg_repmat.fypp
new file mode 100644
index 000000000..3ea47a064
--- /dev/null
+++ b/src/stdlib_linalg_repmat.fypp
@@ -0,0 +1,30 @@
+#:include "common.fypp"
+#:set RCI_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES + INT_KINDS_TYPES
+submodule (stdlib_linalg) stdlib_linalg_repmat
+
+ implicit none
+
+contains
+
+ #:for k1, t1 in RCI_KINDS_TYPES
+ pure module function repmat_${t1[0]}$${k1}$(a, m, n) result(res)
+ ${t1}$, intent(in) :: a(:,:)
+ integer, intent(in) :: m,n
+ ${t1}$ :: res(m*size(a,1),n*size(a,2))
+ integer ::i,j,k,l
+ associate(ma=>size(a,1),na=>size(a,2))
+ do j=1,n
+ do l=1,na
+ do i=1,m
+ do k=1,ma
+ res((i-1)*ma+k,(j-1)*na+l)=a(k,l)
+ end do
+ end do
+ end do
+ end do
+ end associate
+ end function repmat_${t1[0]}$${k1}$
+
+ #:endfor
+
+end submodule
diff --git a/test/linalg/test_linalg.fypp b/test/linalg/test_linalg.fypp
index f74cbff6b..f278345da 100644
--- a/test/linalg/test_linalg.fypp
+++ b/test/linalg/test_linalg.fypp
@@ -3,7 +3,7 @@
module test_linalg
use testdrive, only : new_unittest, unittest_type, error_type, check, skip_test
use stdlib_kinds, only: sp, dp, xdp, qp, int8, int16, int32, int64
- use stdlib_linalg, only: diag, eye, trace, outer_product
+ use stdlib_linalg, only: diag, eye, trace, outer_product, repmat
implicit none
@@ -57,7 +57,17 @@ contains
new_unittest("outer_product_int8", test_outer_product_int8), &
new_unittest("outer_product_int16", test_outer_product_int16), &
new_unittest("outer_product_int32", test_outer_product_int32), &
- new_unittest("outer_product_int64", test_outer_product_int64) &
+ new_unittest("outer_product_int64", test_outer_product_int64), &
+ new_unittest("test_repmat_rsp",test_repmat_rsp), &
+ new_unittest("test_repmat_rdp",test_repmat_rdp), &
+ new_unittest("test_repmat_rqp",test_repmat_rqp), &
+ new_unittest("test_repmat_csp",test_repmat_csp), &
+ new_unittest("test_repmat_cdp",test_repmat_cdp), &
+ new_unittest("test_repmat_cqp",test_repmat_cqp), &
+ new_unittest("test_repmat_int8",test_repmat_int8), &
+ new_unittest("test_repmat_int16",test_repmat_int16), &
+ new_unittest("test_repmat_int32",test_repmat_int32), &
+ new_unittest("test_repmat_int64",test_repmat_int64) &
]
end subroutine collect_linalg
@@ -703,6 +713,197 @@ contains
end subroutine test_outer_product_int64
+ subroutine test_repmat_rqp(error)
+ !> Error handling
+ type(error_type), allocatable, intent(out) :: error
+#:if WITH_QP
+ integer, parameter :: n = 2
+ real(qp) :: a(n,n)
+ real(qp),allocatable:: expected(:,:),diff(:,:)
+
+ a=reshape([1,2,3,4],shape(a),order=[2,1])
+ expected = reshape([&
+ 1,2,1,2,&
+ 3,4,3,4,&
+ 1,2,1,2,&
+ 3,4,3,4 ],[4,4],order=[2,1])
+ diff = expected - repmat(a, 2, 2)
+ call check(error, all(abs(diff) == 0), &
+ "all(abs(diff) == 0) failed.")
+#:else
+ call skip_test(error, "Quadruple precision is not enabled")
+#:endif
+ end subroutine test_repmat_rqp
+ subroutine test_repmat_cqp(error)
+ !> Error handling
+ type(error_type), allocatable, intent(out) :: error
+#:if WITH_QP
+ integer, parameter :: n = 2
+ complex(qp) :: a(n,n)
+ complex(qp),allocatable:: expected(:,:),diff(:,:)
+
+ a=reshape([1,2,3,4],shape(a),order=[2,1])
+ expected = reshape([&
+ 1,2,1,2,&
+ 3,4,3,4,&
+ 1,2,1,2,&
+ 3,4,3,4 ],[4,4],order=[2,1])
+ diff = expected - repmat(a, 2, 2)
+ call check(error, all(abs(diff) == 0), &
+ "all(abs(diff) == 0) failed.")
+#:else
+ call skip_test(error, "Quadruple precision is not enabled")
+#:endif
+ end subroutine test_repmat_cqp
+
+
+ subroutine test_repmat_rsp(error)
+ !> Error handling
+ type(error_type), allocatable, intent(out) :: error
+ integer, parameter :: n = 2
+ real(sp) :: a(n,n)
+ real(sp),allocatable:: expected(:,:),diff(:,:)
+
+ a=reshape([1,2,3,4],shape(a),order=[2,1])
+ expected = reshape([&
+ 1,2,1,2,&
+ 3,4,3,4,&
+ 1,2,1,2,&
+ 3,4,3,4 ],[4,4],order=[2,1])
+ diff = expected - repmat(a, 2, 2)
+ call check(error, all(abs(diff) == 0), &
+ "all(abs(diff) == 0) failed.")
+ end subroutine test_repmat_rsp
+
+ subroutine test_repmat_rdp(error)
+ !> Error handling
+ type(error_type), allocatable, intent(out) :: error
+ integer, parameter :: n = 2
+ real(dp) :: a(n,n)
+ real(dp),allocatable:: expected(:,:),diff(:,:)
+
+ a=reshape([1,2,3,4],shape(a),order=[2,1])
+ expected = reshape([&
+ 1,2,1,2,&
+ 3,4,3,4,&
+ 1,2,1,2,&
+ 3,4,3,4 ],[4,4],order=[2,1])
+ diff = expected - repmat(a, 2, 2)
+ call check(error, all(abs(diff) == 0), &
+ "all(abs(diff) == 0) failed.")
+ end subroutine test_repmat_rdp
+
+
+ subroutine test_repmat_csp(error)
+ !> Error handling
+ type(error_type), allocatable, intent(out) :: error
+ integer, parameter :: n = 2
+ complex(sp) :: a(n,n)
+ complex(sp),allocatable:: expected(:,:),diff(:,:)
+
+ a=reshape([1,2,3,4],shape(a),order=[2,1])
+ expected = reshape([&
+ 1,2,1,2,&
+ 3,4,3,4,&
+ 1,2,1,2,&
+ 3,4,3,4 ],[4,4],order=[2,1])
+ diff = expected - repmat(a, 2, 2)
+ call check(error, all(abs(diff) == 0), &
+ "all(abs(diff) == 0) failed.")
+ end subroutine test_repmat_csp
+
+ subroutine test_repmat_cdp(error)
+ !> Error handling
+ type(error_type), allocatable, intent(out) :: error
+ integer, parameter :: n = 2
+ complex(dp) :: a(n,n)
+ complex(dp),allocatable:: expected(:,:),diff(:,:)
+
+ a=reshape([1,2,3,4],shape(a),order=[2,1])
+ expected = reshape([&
+ 1,2,1,2,&
+ 3,4,3,4,&
+ 1,2,1,2,&
+ 3,4,3,4 ],[4,4],order=[2,1])
+ diff = expected - repmat(a, 2, 2)
+ call check(error, all(abs(diff) == 0), &
+ "all(abs(diff) == 0) failed.")
+ end subroutine test_repmat_cdp
+
+ subroutine test_repmat_int8(error)
+ !> Error handling
+ type(error_type), allocatable, intent(out) :: error
+ integer, parameter :: n = 2
+ integer(int8) :: a(n,n)
+ integer(int8),allocatable:: expected(:,:),diff(:,:)
+
+ a=reshape([1,2,3,4],shape(a),order=[2,1])
+ expected = reshape([&
+ 1,2,1,2,&
+ 3,4,3,4,&
+ 1,2,1,2,&
+ 3,4,3,4 ],[4,4],order=[2,1])
+ diff = expected - repmat(a, 2, 2)
+ call check(error, all(abs(diff) == 0), &
+ "all(abs(diff) == 0) failed.")
+ end subroutine test_repmat_int8
+
+ subroutine test_repmat_int16(error)
+ !> Error handling
+ type(error_type), allocatable, intent(out) :: error
+ integer, parameter :: n = 2
+ integer(int16) :: a(n,n)
+ integer(int16),allocatable:: expected(:,:),diff(:,:)
+
+ a=reshape([1,2,3,4],shape(a),order=[2,1])
+ expected = reshape([&
+ 1,2,1,2,&
+ 3,4,3,4,&
+ 1,2,1,2,&
+ 3,4,3,4 ],[4,4],order=[2,1])
+ diff = expected - repmat(a, 2, 2)
+ call check(error, all(abs(diff) == 0), &
+ "all(abs(diff) == 0) failed.")
+ end subroutine test_repmat_int16
+
+ subroutine test_repmat_int32(error)
+ !> Error handling
+ type(error_type), allocatable, intent(out) :: error
+ integer, parameter :: n = 2
+ integer(int32) :: a(n,n)
+ integer(int32),allocatable:: expected(:,:),diff(:,:)
+
+ a=reshape([1,2,3,4],shape(a),order=[2,1])
+ expected = reshape([&
+ 1,2,1,2,&
+ 3,4,3,4,&
+ 1,2,1,2,&
+ 3,4,3,4 ],[4,4],order=[2,1])
+ diff = expected - repmat(a, 2, 2)
+ call check(error, all(abs(diff) == 0), &
+ "all(abs(diff) == 0) failed.")
+ end subroutine test_repmat_int32
+
+ subroutine test_repmat_int64(error)
+ !> Error handling
+ type(error_type), allocatable, intent(out) :: error
+ integer, parameter :: n = 2
+ integer(int64) :: a(n,n)
+ integer(int64),allocatable:: expected(:,:),diff(:,:)
+
+ a=reshape([1,2,3,4],shape(a),order=[2,1])
+ expected = reshape([&
+ 1,2,1,2,&
+ 3,4,3,4,&
+ 1,2,1,2,&
+ 3,4,3,4 ],[4,4],order=[2,1])
+ diff = expected - repmat(a, 2, 2)
+ call check(error, all(abs(diff) == 0), &
+ "all(abs(diff) == 0) failed.")
+ end subroutine test_repmat_int64
+
+
+
pure recursive function catalan_number(n) result(value)
integer, intent(in) :: n
integer :: value