diff --git a/doc/specs/stdlib_math.md b/doc/specs/stdlib_math.md
new file mode 100644
index 000000000..3e45b6696
--- /dev/null
+++ b/doc/specs/stdlib_math.md
@@ -0,0 +1,92 @@
+---
+title: math
+---
+
+# The `stdlib_math` module
+
+[TOC]
+
+## Introduction
+
+`stdlib_math` module provides general purpose mathematical functions.
+
+
+## Procedures and Methods provided
+
+
+<!-- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -- -->
+### `clip` function
+
+#### Description
+
+Returns a value which lies in the given interval [`xmin`, `xmax`] (interval is `xmin` and `xmax` inclusive) and is closest to the input value `x`.  
+
+#### Syntax
+
+`res = [[stdlib_math(module):clip(interface)]] (x, xmin, xmax)`
+
+#### Status
+
+Experimental
+
+#### Class
+
+Elemental function.
+
+#### Argument(s)
+
+`x`: scalar of either `integer` or `real` type. This argument is `intent(in)`.  
+`xmin`: scalar of either `integer` or `real` type. This argument is `intent(in)`.  
+`xmax`: scalar of either `integer` or `real` type, which must be greater than or equal to `xmin`. This argument is `intent(in)`.
+
+Note: All arguments must have same `type` and same `kind`.
+
+#### Output value or Result value
+
+The output is a scalar of `type` and `kind` same as to that of the arguments.
+
+#### Examples
+
+##### Example 1:
+
+Here inputs are of type `integer` and kind `int32`
+```fortran
+program demo_clip_integer
+  use stdlib_math, only: clip
+  use stdlib_kinds, only: int32
+  implicit none
+  integer(int32) :: x
+  integer(int32) :: xmin
+  integer(int32) :: xmax
+  integer(int32) :: clipped_value
+
+  xmin = -5_int32
+  xmax = 5_int32
+  x = 12_int32
+
+  clipped_value = clip(x, xmin, xmax)
+  ! clipped_value <- 5
+end program demo_clip_integer
+```
+
+##### Example 2:
+
+Here inputs are of type `real` and kind `sp`
+```fortran
+program demo_clip_real
+  use stdlib_math, only: clip
+  use stdlib_kinds, only: sp
+  implicit none
+  real(sp) :: x
+  real(sp) :: xmin
+  real(sp) :: xmax
+  real(sp) :: clipped_value
+
+  xmin = -5.769_sp
+  xmax = 3.025_sp
+  x = 3.025_sp
+
+  clipped_value = clip(x, xmin, xmax)
+  ! clipped_value <- 3.02500010
+end program demo_clip_real
+```
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index 91be20685..a1df4d52c 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -23,6 +23,7 @@ set(fppFiles
     stdlib_quadrature_trapz.fypp
     stdlib_quadrature_simps.fypp
     stdlib_stats_distribution_PRNG.fypp
+    stdlib_math.fypp
     stdlib_string_type.fypp
 )
 
diff --git a/src/Makefile.manual b/src/Makefile.manual
index 37d54dc3c..8ad3fe61a 100644
--- a/src/Makefile.manual
+++ b/src/Makefile.manual
@@ -19,6 +19,7 @@ SRCFYPP =\
         stdlib_stats_moment_mask.fypp \
         stdlib_stats_moment_scalar.fypp \
         stdlib_stats_var.fypp \
+        stdlib_math.fypp \
         stdlib_stats_distribution_PRNG.fypp \
         stdlib_string_type.fypp
 
diff --git a/src/stdlib_math.fypp b/src/stdlib_math.fypp
new file mode 100644
index 000000000..fc00d94f8
--- /dev/null
+++ b/src/stdlib_math.fypp
@@ -0,0 +1,30 @@
+#:include "common.fypp"
+#:set IR_KINDS_TYPES = INT_KINDS_TYPES + REAL_KINDS_TYPES
+
+module stdlib_math
+    use stdlib_kinds, only: int8, int16, int32, int64, sp, dp, qp
+
+    implicit none
+    private
+    public :: clip
+
+    interface clip
+    #:for k1, t1 in IR_KINDS_TYPES
+        module procedure clip_${k1}$
+    #:endfor
+    end interface clip
+
+contains
+    
+    #:for k1, t1 in IR_KINDS_TYPES
+    elemental function clip_${k1}$(x, xmin, xmax) result(res)
+        ${t1}$, intent(in) :: x
+        ${t1}$, intent(in) :: xmin
+        ${t1}$, intent(in) :: xmax
+        ${t1}$ :: res
+
+        res = max(min(x, xmax), xmin)
+    end function clip_${k1}$
+    
+    #:endfor
+end module stdlib_math
diff --git a/src/tests/CMakeLists.txt b/src/tests/CMakeLists.txt
index 288445de9..71275e021 100644
--- a/src/tests/CMakeLists.txt
+++ b/src/tests/CMakeLists.txt
@@ -16,6 +16,7 @@ add_subdirectory(stats)
 add_subdirectory(string)
 add_subdirectory(system)
 add_subdirectory(quadrature)
+add_subdirectory(math)
 
 ADDTEST(always_skip)
 set_tests_properties(always_skip PROPERTIES SKIP_RETURN_CODE 77)
diff --git a/src/tests/Makefile.manual b/src/tests/Makefile.manual
index 773066b1c..e8fdbc139 100644
--- a/src/tests/Makefile.manual
+++ b/src/tests/Makefile.manual
@@ -9,3 +9,4 @@ all test clean:
 	$(MAKE) -f Makefile.manual --directory=quadrature $@
 	$(MAKE) -f Makefile.manual --directory=stats $@
 	$(MAKE) -f Makefile.manual --directory=string $@
+	$(MAKE) -f Makefile.manual --directory=math $@
diff --git a/src/tests/math/CMakeLists.txt b/src/tests/math/CMakeLists.txt
new file mode 100644
index 000000000..ed5f32894
--- /dev/null
+++ b/src/tests/math/CMakeLists.txt
@@ -0,0 +1 @@
+ADDTEST(stdlib_math)
diff --git a/src/tests/math/Makefile.manual b/src/tests/math/Makefile.manual
new file mode 100644
index 000000000..de5f87d26
--- /dev/null
+++ b/src/tests/math/Makefile.manual
@@ -0,0 +1,4 @@
+PROGS_SRC = test_stdlib_math.f90
+
+
+include ../Makefile.manual.test.mk
diff --git a/src/tests/math/test_stdlib_math.f90 b/src/tests/math/test_stdlib_math.f90
new file mode 100644
index 000000000..7fafc6bfe
--- /dev/null
+++ b/src/tests/math/test_stdlib_math.f90
@@ -0,0 +1,98 @@
+! SPDX-Identifier: MIT
+
+program test_stdlib_math
+    use stdlib_math, only: clip
+    use stdlib_error, only: check
+    use stdlib_kinds, only: int8, int16, int32, int64, sp, dp, qp
+    implicit none
+
+    ! clip function
+    ! testing format: check(clip(x, xmin, xmax) == correct answer)
+    ! valid case: xmin is not greater than xmax
+    ! invalid case: xmin is greater than xmax
+
+    ! type: integer(int8), kind: int8
+    ! valid test case
+    call check(clip(2_int8, -2_int8, 5_int8) == 2_int8, &
+               'clip_int8 failed for valid case', warn=.true.)
+    call check(clip(127_int8, -127_int8, 0_int8) == 0_int8, &
+               'clip_int8 failed for valid case', warn=.true.)
+    ! invalid test case
+    call check(clip(2_int8, 5_int8, -2_int8) == 5_int8, &
+               'clip_int8 failed for invalid case', warn=.true.)
+    call check(clip(127_int8, 0_int8, -127_int8) == 0_int8, &
+               'clip_int8 failed for invalid case', warn=.true.)
+
+    ! type: integer(int16), kind: int16
+    ! valid test case
+    call check(clip(2_int16, -2_int16, 5_int16) == 2_int16, &
+               'clip_int16 failed for valid case', warn=.true.)
+    call check(clip(32767_int16, -32767_int16, 0_int16) == 0_int16, &
+               'clip_int16 failed for valid case', warn=.true.)
+    ! invalid test case
+    call check(clip(2_int16, 5_int16, -2_int16) == 5_int16, &
+               'clip_int16 failed for invalid case', warn=.true.)
+    call check(clip(32767_int16, 0_int16, -32767_int16) == 0_int16, &
+               'clip_int16 failed for invalid case', warn=.true.)
+
+    ! type: integer(int32), kind: int32
+    ! valid test case
+    call check(clip(2_int32, -2_int32, 5_int32) == 2_int32, &
+               'clip_int32 failed for valid case', warn=.true.)
+    call check(clip(-2147483647_int32, 0_int32, 2147483647_int32) == 0_int32, &
+               'clip_int32 failed for valid case', warn=.true.)
+    ! invalid test case
+    call check(clip(2_int32, 5_int32, -2_int32) == 5_int32, &
+               'clip_int32 failed for invalid case', warn=.true.)
+    call check(clip(-2147483647_int32, 2147483647_int32, 0_int32) == 2147483647_int32, &
+               'clip_int32 failed for invalid case', warn=.true.)
+
+    ! type: integer(int64), kind: int64
+    ! valid test case
+    call check(clip(2_int64, -2_int64, 5_int64) == 2_int64, &
+               'clip_int64 failed for valid case', warn=.true.)
+    call check(clip(-922337203_int64, -10_int64, 25_int64) == -10_int64, &
+               'clip_int64 failed for valid case', warn=.true.)
+    ! invalid test case
+    call check(clip(2_int64, 5_int64, -2_int64) == 5_int64, &
+               'clip_int64 failed for invalid case', warn=.true.)
+    call check(clip(-922337203_int64, 25_int64, -10_int64) == 25_int64, &
+               'clip_int64 failed for invalid case', warn=.true.)
+
+    ! type: real(sp), kind: sp
+    ! valid test case
+    call check(clip(3.025_sp, -5.77_sp, 3.025_sp) == 3.025_sp, &
+               'clip_sp failed for valid case', warn=.true.)
+    call check(clip(0.0_sp, -1578.025_sp, -59.68_sp) == -59.68_sp, &
+               'clip_sp failed for valid case', warn=.true.)
+    ! invalid test case
+    call check(clip(3.025_sp, 3.025_sp, -5.77_sp) == 3.025_sp, &
+               'clip_sp failed for invalid case', warn=.true.)
+    call check(clip(0.0_sp, -59.68_sp, -1578.025_sp) == -59.68_sp, &
+               'clip_sp failed for invalid case', warn=.true.)
+
+    ! type: real(dp), kind: dp
+    ! valid test case
+    call check(clip(3.025_dp, -5.77_dp, 3.025_dp) == 3.025_dp, &
+               'clip_dp failed for valid case', warn=.true.)
+    call check(clip(-7.0_dp, 0.059668_dp, 1.00268_dp) == 0.059668_dp, &
+               'clip_dp failed for valid case', warn=.true.)
+    ! invalid test case
+    call check(clip(3.025_dp, 3.025_dp, -5.77_dp) == 3.025_dp, &
+               'clip_dp failed for invalid case', warn=.true.)
+    call check(clip(-7.0_dp, 1.00268_dp, 0.059668_dp) == 1.00268_dp, &
+               'clip_dp failed for invalid case', warn=.true.)
+
+    ! type: real(qp), kind: qp
+    ! valid test case
+    call check(clip(3.025_qp, -5.77_qp, 3.025_qp) == 3.025_qp, &
+               'clip_qp failed for valid case', warn=.true.)
+    call check(clip(-55891546.2_qp, -8958133457.23_qp, -689712245.23_qp) == -689712245.23_qp, &
+               'clip_qp failed for valid case', warn=.true.)
+    ! invalid test case
+    call check(clip(3.025_qp, 3.025_qp, -5.77_qp) == 3.025_qp, &
+               'clip_qp failed for invalid case', warn=.true.)
+    call check(clip(-55891546.2_qp, -689712245.23_qp, -8958133457.23_qp) == -689712245.23_qp, &
+               'clip_qp failed for invalid case', warn=.true.)
+
+end program test_stdlib_math