diff --git a/.clang-tidy b/.clang-tidy
index f98062b30ed..f41c44efa2a 100644
--- a/.clang-tidy
+++ b/.clang-tidy
@@ -43,7 +43,6 @@ Checks: >
        -modernize-use-integer-sign-comparison,
        -modernize-use-ranges,
        -modernize-use-starts-ends-with,
-       -modernize-use-std-numbers,
        -modernize-use-trailing-return-type,
      performance-*,
        -performance-enum-size,
@@ -81,7 +80,6 @@ HeaderFileExtensions: ['', "H", 'h', 'hh', 'hpp', 'hxx']
 # We will try to modernize this after switching to C++20
 # -modernize-use-constraints
 # -modernize-use-designated-initializers
-# -modernize-use-std-numbers
 # -modernize-use-ranges
 # -modernize-use-integer-sign-comparison
 # -modernize-use-starts-ends-with
diff --git a/Src/Base/AMReX_COORDSYS_1D_C.H b/Src/Base/AMReX_COORDSYS_1D_C.H
index 5274725fde5..63a03dc677a 100644
--- a/Src/Base/AMReX_COORDSYS_1D_C.H
+++ b/Src/Base/AMReX_COORDSYS_1D_C.H
@@ -6,6 +6,7 @@
 #include <AMReX_Array.H>
 #include <AMReX_FArrayBox.H>
 #include <cmath>
+#include <numbers>
 
 namespace amrex {
 
@@ -29,7 +30,7 @@ amrex_setvol (Box const& bx, Array4<Real> const& vol,
     }
     else if (coord == 1)
     {
-        const Real pi = Real(3.1415926535897932);
+        const Real pi = std::numbers::pi_v<Real>;
         AMREX_PRAGMA_SIMD
         for (int i = lo.x; i <= hi.x; ++i) {
             Real ri = offset[0] + dx[0]*i;
@@ -40,7 +41,7 @@ amrex_setvol (Box const& bx, Array4<Real> const& vol,
     }
     else
     {
-        const Real pi = Real(3.1415926535897932);
+        const Real pi = std::numbers::pi_v<Real>;
         AMREX_PRAGMA_SIMD
         for (int i = lo.x; i <= hi.x; ++i) {
             Real ri = offset[0] + dx[0]*i;
@@ -70,7 +71,7 @@ amrex_setarea (Box const& bx, Array4<Real> const& area,
     }
     else if (coord == 1)
     {
-        const Real pi = Real(3.1415926535897932);
+        const Real pi = std::numbers::pi_v<Real>;
         AMREX_PRAGMA_SIMD
         for (int i = lo.x; i <= hi.x; ++i) {
             Real ri = offset[0] + dx[0]*i;
@@ -80,7 +81,7 @@ amrex_setarea (Box const& bx, Array4<Real> const& area,
     }
     else
     {
-        const Real pi = Real(3.1415926535897932);
+        const Real pi = std::numbers::pi_v<Real>;
         AMREX_PRAGMA_SIMD
         for (int i = lo.x; i <= hi.x; ++i) {
             Real ri = offset[0] + dx[0]*i;
diff --git a/Src/Base/AMReX_COORDSYS_2D_C.H b/Src/Base/AMReX_COORDSYS_2D_C.H
index 45f3f25052b..e7a790c4d59 100644
--- a/Src/Base/AMReX_COORDSYS_2D_C.H
+++ b/Src/Base/AMReX_COORDSYS_2D_C.H
@@ -6,6 +6,7 @@
 #include <AMReX_Array.H>
 #include <AMReX_FArrayBox.H>
 #include <cmath>
+#include <numbers>
 
 namespace amrex {
 
@@ -31,7 +32,7 @@ amrex_setvol (Box const& bx, Array4<Real> const& vol,
     }
     else if (coord == 1) // r-z
     {
-        const Real pi = Real(3.1415926535897932);
+        const Real pi = std::numbers::pi_v<Real>;
         for     (int j = lo.y; j <= hi.y; ++j) {
             AMREX_PRAGMA_SIMD
             for (int i = lo.x; i <= hi.x; ++i) {
@@ -44,7 +45,7 @@ amrex_setvol (Box const& bx, Array4<Real> const& vol,
     }
     else // r-theta
     {
-        const Real pi = Real(3.1415926535897932);
+        const Real pi = std::numbers::pi_v<Real>;
         for     (int j = lo.y; j <= hi.y; ++j) {
             Real ti = offset[1] + dx[1]*(j);
             Real to = ti + dx[1];
@@ -82,7 +83,7 @@ amrex_setarea (Box const& bx, Array4<Real> const& area,
     }
     else if (coord == 1)
     {
-        const Real pi = Real(3.1415926535897932);
+        const Real pi = std::numbers::pi_v<Real>;
         if (dir == 0)
         {
             for     (int j = lo.y; j <= hi.y; ++j) {
@@ -108,7 +109,7 @@ amrex_setarea (Box const& bx, Array4<Real> const& area,
     }
     else
     {
-        const Real pi = Real(3.1415926535897932);
+        const Real pi = std::numbers::pi_v<Real>;
         if (dir == 0)
         {
             for     (int j = lo.y; j <= hi.y; ++j) {
diff --git a/Src/Base/AMReX_CoordSys.cpp b/Src/Base/AMReX_CoordSys.cpp
index a958d33a5b6..2452f0df64e 100644
--- a/Src/Base/AMReX_CoordSys.cpp
+++ b/Src/Base/AMReX_CoordSys.cpp
@@ -5,12 +5,13 @@
 #include <AMReX_ParallelDescriptor.H>
 
 #include <iostream>
+#include <numbers>
 
 namespace {
 #if (AMREX_SPACEDIM == 2)
-    constexpr double  TWOPI = 2.*3.14159265358979323846264338327950288;
+    constexpr double TWOPI = 2.0 * std::numbers::pi_v<double>;
 #elif (AMREX_SPACEDIM == 1)
-    constexpr double FOURPI = 4.*3.14159265358979323846264338327950288;
+    constexpr double FOURPI = 4.0 * std::numbers::pi_v<double>;
 #endif
 }
 
diff --git a/Src/Base/AMReX_Geometry.H b/Src/Base/AMReX_Geometry.H
index a1a8cc6fb6d..721ec9e06f4 100644
--- a/Src/Base/AMReX_Geometry.H
+++ b/Src/Base/AMReX_Geometry.H
@@ -14,6 +14,7 @@
 
 #include <iosfwd>
 #include <map>
+#include <numbers>
 
 namespace amrex {
 
@@ -252,7 +253,7 @@ public:
             Real rl = geomdata.ProbLo()[0] + static_cast<Real>(point[0]) * dx[0];
             Real rr = rl + dx[0];
 
-            constexpr Real pi = Real(3.1415926535897932);
+            constexpr Real pi = std::numbers::pi_v<Real>;
             vol = (4.0_rt / 3.0_rt) * pi * dx[0] * (rl * rl + rl * rr + rr * rr);
         }
 
@@ -271,7 +272,7 @@ public:
             Real r_l = geomdata.ProbLo()[0] + static_cast<Real>(point[0]) * dx[0];
             Real r_r = geomdata.ProbLo()[0] + static_cast<Real>(point[0]+1) * dx[0];
 
-            constexpr Real pi = Real(3.1415926535897932);
+            constexpr Real pi = std::numbers::pi_v<Real>;
             vol = pi * (r_l + r_r) * dx[0] * dx[1];
         }
         else {
@@ -283,7 +284,7 @@ public:
             Real theta_l = geomdata.ProbLo()[1] + static_cast<Real>(point[1]) * dx[1];
             Real theta_r = geomdata.ProbLo()[1] + static_cast<Real>(point[1]+1) * dx[1];
 
-            constexpr Real twoThirdsPi = static_cast<Real>(2.0 * 3.1415926535897932 / 3.0);
+            constexpr Real twoThirdsPi = static_cast<Real>((2.0 / 3.0) * std::numbers::pi_v<Real>);
             vol = twoThirdsPi * (std::cos(theta_l) - std::cos(theta_r)) * dx[0] *
                 (r_r*r_r + r_r*r_l + r_l*r_l);
         }
diff --git a/Src/Base/AMReX_Math.H b/Src/Base/AMReX_Math.H
index bb943a04cbf..a4124a7ee40 100644
--- a/Src/Base/AMReX_Math.H
+++ b/Src/Base/AMReX_Math.H
@@ -10,6 +10,7 @@
 #include <cmath>
 #include <cstdlib>
 #include <limits>
+#include <numbers>
 #include <type_traits>
 #include <utility>
 
@@ -67,7 +68,7 @@ using std::isinf;
 template <typename T>
 constexpr std::enable_if_t<std::is_floating_point_v<T>,T> pi ()
 {
-    return T(3.1415926535897932384626433832795029L);
+    return std::numbers::pi_v<T>;
 }
 
 //! Return cos(x*pi) given x
diff --git a/Src/Base/AMReX_MultiFabUtil.cpp b/Src/Base/AMReX_MultiFabUtil.cpp
index 1a64ae5f604..c51facc7a6f 100644
--- a/Src/Base/AMReX_MultiFabUtil.cpp
+++ b/Src/Base/AMReX_MultiFabUtil.cpp
@@ -1,6 +1,7 @@
 
 #include <AMReX_MultiFabUtil.H>
 #include <AMReX_Random.H>
+#include <numbers>
 #include <sstream>
 #include <iostream>
 
@@ -1052,7 +1053,7 @@ namespace amrex
                         if (m[box_no](i,j,k)) {
                             return Real(0.);
                         } else {
-                            constexpr Real pi = Real(3.1415926535897932);
+                            constexpr Real pi = std::numbers::pi_v<Real>;
                             Real ri = rlo + dx[0]*i;
                             Real ro = ri + dx[0];
                             return Real(4./3.)*pi*(ro-ri)*(ro*ro+ro*ri+ri*ri)
@@ -1072,7 +1073,7 @@ namespace amrex
                         } else {
                             Real ri = rlo + dx[0]*i;
                             Real ro = ri + dx[0];
-                            constexpr Real pi = Real(3.1415926535897932);
+                            constexpr Real pi = std::numbers::pi_v<Real>;
                             return pi*dx[1]*dx[0]*(ro+ri)
                                 * a[box_no](i,j,k,icomp);
                         }
@@ -1117,7 +1118,7 @@ namespace amrex
                 [=] AMREX_GPU_DEVICE (int box_no, int i, int j, int k)
                                noexcept -> Real
                 {
-                    constexpr Real pi = Real(3.1415926535897932);
+                    constexpr Real pi = std::numbers::pi_v<Real>;
                     Real ri = rlo + dx[0]*i;
                     Real ro = ri + dx[0];
                     return Real(4./3.)*pi*(ro-ri)*(ro*ro+ro*ri+ri*ri)
@@ -1133,7 +1134,7 @@ namespace amrex
                 {
                     Real ri = rlo + dx[0]*i;
                     Real ro = ri + dx[0];
-                    constexpr Real pi = Real(3.1415926535897932);
+                    constexpr Real pi = std::numbers::pi_v<Real>;
                     return pi*dx[1]*dx[0]*(ro+ri)
                         * a[box_no](i,j,k,icomp);
                 });
diff --git a/Src/LinearSolvers/OpenBC/AMReX_OpenBC_K.H b/Src/LinearSolvers/OpenBC/AMReX_OpenBC_K.H
index b9a08f79ac8..1f0aab7ecbb 100644
--- a/Src/LinearSolvers/OpenBC/AMReX_OpenBC_K.H
+++ b/Src/LinearSolvers/OpenBC/AMReX_OpenBC_K.H
@@ -3,6 +3,7 @@
 
 #include <AMReX_OpenBC.H>
 #include <AMReX_LOUtil_K.H>
+#include <numbers>
 
 namespace amrex::openbc {
 
@@ -46,7 +47,7 @@ void scale_moments (openbc::Moments::array_type& mom)
 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
 Real block_potential (openbc::Moments const& mom, Real xb, Real yb, Real zb)
 {
-    constexpr Real oneover4pi = Real(1.)/Real(4.*3.1415926535897932);
+    constexpr Real oneover4pi = Real(1.) / (Real(4.) * std::numbers::pi_v<Real>);
 
     xb -= mom.x;
     yb -= mom.y;
diff --git a/Tests/Amr/Advection_AmrCore/Source/face_velocity.H b/Tests/Amr/Advection_AmrCore/Source/face_velocity.H
index a06392c2a69..1241d999f71 100644
--- a/Tests/Amr/Advection_AmrCore/Source/face_velocity.H
+++ b/Tests/Amr/Advection_AmrCore/Source/face_velocity.H
@@ -4,6 +4,7 @@
 #include <AMReX_Geometry.H>
 #include <AMReX_FArrayBox.H>
 #include <AMReX_REAL.H>
+#include <numbers>
 
 AMREX_GPU_DEVICE
 AMREX_FORCE_INLINE
@@ -13,7 +14,7 @@ void get_face_velocity_psi(int i, int j,
                            amrex::GeometryData const& geomdata)
 {
     using namespace amrex;
-    constexpr Real PI = 3.1415926535897932384626;
+    constexpr Real PI = std::numbers::pi_v<Real>;
 
     const Real* AMREX_RESTRICT prob_lo = geomdata.ProbLo();
     const Real* AMREX_RESTRICT dx      = geomdata.CellSize();
diff --git a/Tests/Amr/Advection_AmrLevel/Exec/SingleVortex/face_velocity_2d_K.H b/Tests/Amr/Advection_AmrLevel/Exec/SingleVortex/face_velocity_2d_K.H
index a6c3858e312..44013af8cec 100644
--- a/Tests/Amr/Advection_AmrLevel/Exec/SingleVortex/face_velocity_2d_K.H
+++ b/Tests/Amr/Advection_AmrLevel/Exec/SingleVortex/face_velocity_2d_K.H
@@ -3,6 +3,7 @@
 
 #include <AMReX_FArrayBox.H>
 #include <AMReX_Array.H>
+#include <numbers>
 
 AMREX_GPU_DEVICE
 AMREX_FORCE_INLINE
@@ -15,7 +16,7 @@ void get_face_velocity_psi(int i, int j,
                            amrex::Real prob_lo_y)
 {
     using namespace amrex;
-    constexpr Real PI = 3.141592653589793238462643383279502884197;
+    constexpr Real PI = std::numbers::pi_v<Real>;
 
     // Compute streamfunction psi
     Real y = ( (Real)j + 0.5 )*dy + prob_lo_y;
diff --git a/Tests/Amr/Advection_AmrLevel/Exec/SingleVortex/face_velocity_3d_K.H b/Tests/Amr/Advection_AmrLevel/Exec/SingleVortex/face_velocity_3d_K.H
index 73ff18c48c1..a49bab2bc50 100644
--- a/Tests/Amr/Advection_AmrLevel/Exec/SingleVortex/face_velocity_3d_K.H
+++ b/Tests/Amr/Advection_AmrLevel/Exec/SingleVortex/face_velocity_3d_K.H
@@ -3,6 +3,7 @@
 
 #include <AMReX_FArrayBox.H>
 #include <AMReX_Array.H>
+#include <numbers>
 
 AMREX_GPU_DEVICE
 AMREX_FORCE_INLINE
@@ -17,7 +18,7 @@ void get_face_velocity_psi(int i, int j,
                            amrex::Real /*prob_lo_z*/)
 {
     using namespace amrex;
-    constexpr Real PI = 3.141592653589793238462643383279502884197;
+    constexpr Real PI = std::numbers::pi_v<Real>;
 
     // Compute streamfunction psi
     Real y = ( (Real)j + 0.5 )*dy + prob_lo_y;
diff --git a/Tests/EB_CNS/Exec/Pulse/cns_prob.H b/Tests/EB_CNS/Exec/Pulse/cns_prob.H
index 10f2ec40672..82bf4298a2b 100644
--- a/Tests/EB_CNS/Exec/Pulse/cns_prob.H
+++ b/Tests/EB_CNS/Exec/Pulse/cns_prob.H
@@ -3,6 +3,7 @@
 
 #include <AMReX_Geometry.H>
 #include <AMReX_FArrayBox.H>
+#include <numbers>
 #include "CNS_index_macros.H"
 #include "CNS_parm.H"
 #include "cns_prob_parm.H"
@@ -15,7 +16,7 @@ cns_initdata (int i, int j, int k, amrex::Array4<amrex::Real> const& state,
 {
     using amrex::Real;
 
-    constexpr Real pi = 3.1415926535897932;
+    constexpr Real pi = std::numbers::pi_v<Real>;
 
     const Real* prob_lo = geomdata.ProbLo();
     const Real* prob_hi = geomdata.ProbHi();
diff --git a/Tests/GPU/CNS/Exec/RT/cns_prob.H b/Tests/GPU/CNS/Exec/RT/cns_prob.H
index a5590ba95d6..9c947b68e9d 100644
--- a/Tests/GPU/CNS/Exec/RT/cns_prob.H
+++ b/Tests/GPU/CNS/Exec/RT/cns_prob.H
@@ -3,6 +3,7 @@
 
 #include <AMReX_Geometry.H>
 #include <AMReX_FArrayBox.H>
+#include <numbers>
 #include "CNS_index_macros.H"
 #include "CNS_parm.H"
 #include "cns_prob_parm.H"
@@ -19,7 +20,7 @@ cns_initdata (int i, int j, int k, amrex::Array4<amrex::Real> const& state,
     const Real* prob_hi = geomdata.ProbHi();
     const Real* dx      = geomdata.CellSize();
 
-    constexpr Real pi = Real(3.14159265358979323846264338327950288);
+    constexpr Real pi = std::numbers::pi_v<Real>;
     const Real splitx = Real(0.5)*(prob_lo[0]+prob_hi[0]);
     const Real splity = Real(0.5)*(prob_lo[1]+prob_hi[1]);
     const Real splitz = Real(0.5)*(prob_lo[2]+prob_hi[2]);
diff --git a/Tests/LinearSolvers/ABecLap_SP/initProb_K.H b/Tests/LinearSolvers/ABecLap_SP/initProb_K.H
index d593cb2cbce..c3cfa056b21 100644
--- a/Tests/LinearSolvers/ABecLap_SP/initProb_K.H
+++ b/Tests/LinearSolvers/ABecLap_SP/initProb_K.H
@@ -2,6 +2,7 @@
 #define INIT_PROB_K_H_
 
 #include <AMReX_FArrayBox.H>
+#include <numbers>
 
 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
 void actual_init_poisson (int i, int j, int k,
@@ -10,8 +11,9 @@ void actual_init_poisson (int i, int j, int k,
                           amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& prob_lo,
                           amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dx)
 {
-    constexpr amrex::Real tpi = 2.*3.1415926535897932;
-    constexpr amrex::Real fpi = 4.*3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
+    constexpr amrex::Real tpi = 2.*pi;
+    constexpr amrex::Real fpi = 4.*pi;
     constexpr amrex::Real fac = tpi*tpi*static_cast<amrex::Real>(AMREX_SPACEDIM);
     amrex::Real x = prob_lo[0] + dx[0] * (i + 0.5);
     amrex::Real y = prob_lo[1] + dx[1] * (j + 0.5);
@@ -50,7 +52,7 @@ void actual_init_bcoef (int i, int j, int k,
     constexpr amrex::Real sigma = 10.;
     const amrex::Real theta = 0.5*std::log(3.) / (w + 1.e-50);
 
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
 
@@ -99,7 +101,7 @@ void actual_init_abeclap (int i, int j, int k,
     constexpr amrex::Real sigma = 10.;
     const amrex::Real theta = 0.5*std::log(3.) / (w + 1.e-50);
 
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
     constexpr amrex::Real fac = static_cast<amrex::Real>(AMREX_SPACEDIM*4)*pi*pi;
diff --git a/Tests/LinearSolvers/ABecLaplacian_C/initProb_K.H b/Tests/LinearSolvers/ABecLaplacian_C/initProb_K.H
index 3dba4a77124..3989de179f2 100644
--- a/Tests/LinearSolvers/ABecLaplacian_C/initProb_K.H
+++ b/Tests/LinearSolvers/ABecLaplacian_C/initProb_K.H
@@ -3,6 +3,7 @@
 
 #include <AMReX_FArrayBox.H>
 #include <algorithm>
+#include <numbers>
 
 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
 void actual_init_poisson (int i, int j, int k,
@@ -11,8 +12,9 @@ void actual_init_poisson (int i, int j, int k,
                           amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& prob_lo,
                           amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dx)
 {
-    constexpr amrex::Real tpi = 2.*3.1415926535897932;
-    constexpr amrex::Real fpi = 4.*3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
+    constexpr amrex::Real tpi = 2.*pi;
+    constexpr amrex::Real fpi = 4.*pi;
     constexpr amrex::Real fac = tpi*tpi*static_cast<amrex::Real>(AMREX_SPACEDIM);
     amrex::Real x = prob_lo[0] + dx[0] * (i + 0.5);
     amrex::Real y = prob_lo[1] + dx[1] * (j + 0.5);
@@ -86,7 +88,7 @@ void actual_init_abeclap (int i, int j, int k,
     constexpr amrex::Real sigma = 10.;
     const amrex::Real theta = 0.5*std::log(3.) / (w + 1.e-50);
 
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
     constexpr amrex::Real fac = static_cast<amrex::Real>(AMREX_SPACEDIM*4)*pi*pi;
@@ -162,7 +164,7 @@ void actual_init_abeclap_in (int i, int j, int k,
     constexpr amrex::Real sigma = 10.;
     const amrex::Real theta = 0.5*std::log(3.) / (w + 1.e-50);
 
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
     constexpr amrex::Real fac = static_cast<amrex::Real>(AMREX_SPACEDIM*4)*pi*pi;
@@ -210,7 +212,7 @@ void actual_init_dphi_dx_lo (int i, int j, int k, amrex::Array4<amrex::Real> con
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& prob_lo,
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dx)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
 
@@ -231,7 +233,7 @@ void actual_init_dphi_dx_hi (int i, int j, int k, amrex::Array4<amrex::Real> con
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& prob_lo,
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dx)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
 
@@ -252,7 +254,7 @@ void actual_init_dphi_dy_lo (int i, int j, int k, amrex::Array4<amrex::Real> con
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& prob_lo,
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dx)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
 
@@ -273,7 +275,7 @@ void actual_init_dphi_dy_hi (int i, int j, int k, amrex::Array4<amrex::Real> con
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& prob_lo,
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dx)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
 
@@ -294,7 +296,7 @@ void actual_init_dphi_dz_lo (int i, int j, int k, amrex::Array4<amrex::Real> con
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& prob_lo,
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dx)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
 
@@ -315,7 +317,7 @@ void actual_init_dphi_dz_hi (int i, int j, int k, amrex::Array4<amrex::Real> con
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& prob_lo,
                              amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dx)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
 
@@ -348,7 +350,7 @@ void actual_init_nodeabeclap(int i, int j, int k,
     constexpr amrex::Real sigma = 10.;
     const amrex::Real theta = 0.5*std::log(3.) / (w + 1.e-50);
 
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real tpi =  2.*pi;
     constexpr amrex::Real fpi =  4.*pi;
     constexpr amrex::Real fac = static_cast<amrex::Real>(AMREX_SPACEDIM*4)*pi*pi;
diff --git a/Tests/LinearSolvers/CellEB2/MyTest_K.H b/Tests/LinearSolvers/CellEB2/MyTest_K.H
index deafabf9d7d..d9e98670e99 100644
--- a/Tests/LinearSolvers/CellEB2/MyTest_K.H
+++ b/Tests/LinearSolvers/CellEB2/MyTest_K.H
@@ -3,6 +3,7 @@
 
 #include <AMReX_FArrayBox.H>
 #include <AMReX_EBCellFlag.H>
+#include <numbers>
 
 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
 void mytest_set_phi_reg (int i, int j, int k, amrex::Array4<amrex::Real> const& phie,
@@ -13,7 +14,7 @@ void mytest_set_phi_reg (int i, int j, int k, amrex::Array4<amrex::Real> const&
                          amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> const& dx,
                          int prob_type, amrex::Box const& vbx)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     amrex::IntVect iv(AMREX_D_DECL(i,j,k));
 
     if (vbx.contains(iv)) {
@@ -70,7 +71,7 @@ void mytest_set_phi_eb (int i, int j, int k, amrex::Array4<amrex::Real> const& p
                         amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> const& dx,
                         int prob_type, amrex::Box const& vbx)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     amrex::IntVect iv(AMREX_D_DECL(i,j,k));
 
     if (vbx.contains(iv)) {
@@ -141,7 +142,7 @@ void mytest_set_phi_boundary (int i, int j, int k, amrex::Array4<amrex::Real> co
                               amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> const& dx,
                               amrex::Box const& dbx)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
 
     if (!dbx.contains(amrex::IntVect(AMREX_D_DECL(i,j,k)))) {
         amrex::Real x, y;
diff --git a/Tests/LinearSolvers/CellOverset/MyTest.cpp b/Tests/LinearSolvers/CellOverset/MyTest.cpp
index 4d45b8d7d85..d84e39b617b 100644
--- a/Tests/LinearSolvers/CellOverset/MyTest.cpp
+++ b/Tests/LinearSolvers/CellOverset/MyTest.cpp
@@ -5,6 +5,8 @@
 #include <AMReX_MultiFabUtil.H>
 #include <AMReX_PlotFileUtil.H>
 
+#include <numbers>
+
 using namespace amrex;
 
 MyTest::MyTest ()
@@ -173,7 +175,7 @@ MyTest::initData ()
             constexpr amrex::Real sigma = 10.;
             const amrex::Real theta = 0.5*std::log(3.) / (w + 1.e-50);
 
-            constexpr amrex::Real pi = 3.1415926535897932;
+            constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
             constexpr amrex::Real tpi =  2.*pi;
             constexpr amrex::Real fpi =  4.*pi;
             constexpr amrex::Real fac = static_cast<amrex::Real>(AMREX_SPACEDIM)*4.*pi*pi;
diff --git a/Tests/LinearSolvers/EBTensor/MyTest_2D_K.H b/Tests/LinearSolvers/EBTensor/MyTest_2D_K.H
index 5087c5ed4d1..2feb0ce1a71 100644
--- a/Tests/LinearSolvers/EBTensor/MyTest_2D_K.H
+++ b/Tests/LinearSolvers/EBTensor/MyTest_2D_K.H
@@ -2,6 +2,7 @@
 #define MY_TEST_K_H_
 
 #include <AMReX_FArrayBox.H>
+#include <numbers>
 
 AMREX_GPU_DEVICE
 inline void init (amrex::Real x, amrex::Real y, amrex::Real R2,
@@ -9,7 +10,7 @@ inline void init (amrex::Real x, amrex::Real y, amrex::Real R2,
                   amrex::Real& urhs, amrex::Real& vrhs,
                   amrex::Real& eta)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     amrex::Real x2 = x*x;
     amrex::Real y2 = y*y;
     amrex::Real r2 = x2+y2;
diff --git a/Tests/LinearSolvers/EBTensor/MyTest_3D_K.H b/Tests/LinearSolvers/EBTensor/MyTest_3D_K.H
index db3b726d8fe..3b69b0eddc5 100644
--- a/Tests/LinearSolvers/EBTensor/MyTest_3D_K.H
+++ b/Tests/LinearSolvers/EBTensor/MyTest_3D_K.H
@@ -2,6 +2,7 @@
 #define MY_TEST_K_H_
 
 #include <AMReX_FArrayBox.H>
+#include <numbers>
 
 AMREX_GPU_DEVICE
 inline void init (amrex::Real x, amrex::Real y, amrex::Real z, amrex::Real R2,
@@ -9,7 +10,7 @@ inline void init (amrex::Real x, amrex::Real y, amrex::Real z, amrex::Real R2,
                   amrex::Real& urhs, amrex::Real& vrhs, amrex::Real& wrhs,
                   amrex::Real& eta)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     amrex::Real x2 = x*x;
     amrex::Real y2 = y*y;
     amrex::Real r2 = x2+y2;
diff --git a/Tests/LinearSolvers/LeastSquares/initEB.cpp b/Tests/LinearSolvers/LeastSquares/initEB.cpp
index a5fc5617856..93da9d0812c 100644
--- a/Tests/LinearSolvers/LeastSquares/initEB.cpp
+++ b/Tests/LinearSolvers/LeastSquares/initEB.cpp
@@ -6,6 +6,7 @@
 
 #include <cmath>
 #include <algorithm>
+#include <numbers>
 
 #include "MyTest.H"
 #include "MyEB.H"
@@ -15,7 +16,7 @@ using namespace amrex;
 void
 MyTest::initializeEB ()
 {
-    constexpr Real pi = 3.1415926535897932;
+    constexpr Real pi = std::numbers::pi_v<Real>;
 
     ParmParse pp("eb2");
     std::string geom_type;
diff --git a/Tests/LinearSolvers/LeastSquares/initPoiseuilleDataFor2D.cpp b/Tests/LinearSolvers/LeastSquares/initPoiseuilleDataFor2D.cpp
index 5bc85363840..6a2a6044b24 100644
--- a/Tests/LinearSolvers/LeastSquares/initPoiseuilleDataFor2D.cpp
+++ b/Tests/LinearSolvers/LeastSquares/initPoiseuilleDataFor2D.cpp
@@ -3,6 +3,8 @@
 #include <AMReX_EBMultiFabUtil.H>
 #include <AMReX_MultiFabUtil.H>
 
+#include <numbers>
+
 using namespace amrex;
 
 #if (AMREX_SPACEDIM == 2)
@@ -40,7 +42,7 @@ void MyTest::initializePoiseuilleDataFor2D(int ilev) {
 
     amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
       Real H = poiseuille_height;
-      constexpr Real pi = 3.1415926535897932;
+      constexpr Real pi = std::numbers::pi_v<Real>;
       Real t = (poiseuille_rotation / 180.) * pi;
 
       Real a = std::tan(t);
diff --git a/Tests/LinearSolvers/LeastSquares/initPoiseuilleDataFor3D.cpp b/Tests/LinearSolvers/LeastSquares/initPoiseuilleDataFor3D.cpp
index 98c8c4d051d..a14c4cd11b9 100644
--- a/Tests/LinearSolvers/LeastSquares/initPoiseuilleDataFor3D.cpp
+++ b/Tests/LinearSolvers/LeastSquares/initPoiseuilleDataFor3D.cpp
@@ -3,6 +3,8 @@
 #include <AMReX_EBMultiFabUtil.H>
 #include <AMReX_MultiFabUtil.H>
 
+#include <numbers>
+
 using namespace amrex;
 
 #if (AMREX_SPACEDIM == 3)
@@ -49,7 +51,7 @@ void MyTest::initializePoiseuilleDataFor3D(int ilev) {
                                                   int k) noexcept {
         Real H = poiseuille_height;
         int nfdir = poiseuille_no_flow_dir;
-        constexpr Real pi = 3.1415926535897932;
+        constexpr Real pi = std::numbers::pi_v<Real>;
         Real alpha = (poiseuille_askew_rotation[0] / 180.) * pi;
         Real gamma = (poiseuille_askew_rotation[1] / 180.) * pi;
 
diff --git a/Tests/LinearSolvers/NodalOverset/MyTest.cpp b/Tests/LinearSolvers/NodalOverset/MyTest.cpp
index 472e75cf7db..6a539e6dddd 100644
--- a/Tests/LinearSolvers/NodalOverset/MyTest.cpp
+++ b/Tests/LinearSolvers/NodalOverset/MyTest.cpp
@@ -5,6 +5,8 @@
 #include <AMReX_MultiFabUtil.H>
 #include <AMReX_PlotFileUtil.H>
 
+#include <numbers>
+
 using namespace amrex;
 
 MyTest::MyTest ()
@@ -132,7 +134,7 @@ MyTest::initData ()
         {
             Real x = i*dx[0] - 0.5;
             Real y = j*dx[1] - 0.5;
-            Real theta = std::atan2(x,y) + 0.5*3.1415926535897932;
+            Real theta = std::atan2(x,y) + Real(0.5)*std::numbers::pi_v<Real>;
             Real r2 = x*x + y*y;
             pa(i,j,k) =  r2*r2*std::cos(3.0*theta);
             ra(i,j,k) = 7.0*r2*std::cos(3.0*theta);
diff --git a/Tests/LinearSolvers/NodalPoisson/MyTest.cpp b/Tests/LinearSolvers/NodalPoisson/MyTest.cpp
index a8b6f825862..66de157c77c 100644
--- a/Tests/LinearSolvers/NodalPoisson/MyTest.cpp
+++ b/Tests/LinearSolvers/NodalPoisson/MyTest.cpp
@@ -6,6 +6,8 @@
 #include <AMReX_FillPatchUtil.H>
 #include <AMReX_PlotFileUtil.H>
 
+#include <numbers>
+
 using namespace amrex;
 
 MyTest::MyTest ()
@@ -251,7 +253,7 @@ MyTest::initData ()
             Array4<Real> const rh  = rhs[ilev].array(mfi);
             amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
             {
-                constexpr Real pi = 3.1415926535897932;
+                constexpr Real pi = std::numbers::pi_v<Real>;
                 constexpr Real tpi = 2.*pi;
                 constexpr Real fpi = 4.*pi;
                 constexpr Real fac = tpi*tpi*AMREX_SPACEDIM;
diff --git a/Tests/LinearSolvers/NodeTensorLap/MyTest.cpp b/Tests/LinearSolvers/NodeTensorLap/MyTest.cpp
index 46d52a4c109..0385b1b4522 100644
--- a/Tests/LinearSolvers/NodeTensorLap/MyTest.cpp
+++ b/Tests/LinearSolvers/NodeTensorLap/MyTest.cpp
@@ -5,6 +5,8 @@
 #include <AMReX_FillPatchUtil.H>
 #include <AMReX_PlotFileUtil.H>
 
+#include <numbers>
+
 using namespace amrex;
 
 MyTest::MyTest ()
@@ -161,7 +163,7 @@ MyTest::initData ()
             Array4<Real> const rh  = rhs[ilev].array(mfi);
             amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
             {
-                constexpr Real pi = 3.1415926535897932;
+                constexpr Real pi = std::numbers::pi_v<Real>;
                 constexpr Real tpi = 2.*pi;
                 constexpr Real fpi = 4.*pi;
                 constexpr Real fac = 4.*pi*pi;
diff --git a/Tests/LinearSolvers/TensorOverset/MyTest_K.H b/Tests/LinearSolvers/TensorOverset/MyTest_K.H
index 12fe94bb057..87470c115a8 100644
--- a/Tests/LinearSolvers/TensorOverset/MyTest_K.H
+++ b/Tests/LinearSolvers/TensorOverset/MyTest_K.H
@@ -2,6 +2,7 @@
 #define MY_TEST_K_H_
 
 #include <AMReX_FArrayBox.H>
+#include <numbers>
 
 AMREX_GPU_DEVICE AMREX_FORCE_INLINE
 void init (amrex::Real x, amrex::Real y, amrex::Real z, amrex::Real R2,
@@ -9,7 +10,7 @@ void init (amrex::Real x, amrex::Real y, amrex::Real z, amrex::Real R2,
            amrex::Real& urhs, amrex::Real& vrhs, amrex::Real& wrhs,
            amrex::Real& eta)
 {
-    constexpr amrex::Real pi = 3.1415926535897932;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     amrex::Real x2 = x*x;
     amrex::Real y2 = y*y;
     amrex::Real r2 = x2+y2;
diff --git a/Tests/ParmParse/main.cpp b/Tests/ParmParse/main.cpp
index 8bd54e37c44..141debf3ebc 100644
--- a/Tests/ParmParse/main.cpp
+++ b/Tests/ParmParse/main.cpp
@@ -1,4 +1,5 @@
 #include <AMReX.H>
+#include <numbers>
 #include <AMReX_Box.H>
 #include <AMReX_Utility.H>
 #include <AMReX_Print.H>
@@ -21,7 +22,7 @@ int main(int argc, char* argv[])
         ParmParse::SetParserPrefix("physical_constants");
         ParmParse pp("physical_constants");
         pp.add("c", 299792458.);
-        pp.add("pi", 3.14159265358979323846);
+        pp.add("pi", std::numbers::pi_v<double>);
     }
     {
         ParmParse pp;
@@ -106,9 +107,9 @@ int main(int argc, char* argv[])
         ParmParse pp;
         auto parser = pp.makeParser("pi*x+c*y", {"x","y"});
         auto exe = parser.compile<2>();
-        AMREX_ALWAYS_ASSERT(amrex::almostEqual(3.14159265358979323846+299792458.,
+        AMREX_ALWAYS_ASSERT(amrex::almostEqual(std::numbers::pi_v<double>+299792458.,
                                                exe(1.0,1.0)) &&
-                            amrex::almostEqual(3.14159265358979323846, exe(1.0,0.0)) &&
+                            amrex::almostEqual(std::numbers::pi_v<double>, exe(1.0,0.0)) &&
                             amrex::almostEqual(299792458., exe(0.0, 1.0)));
     }
     {
diff --git a/Tests/Parser/main.cpp b/Tests/Parser/main.cpp
index 188ac9ed17b..9ca46e67932 100644
--- a/Tests/Parser/main.cpp
+++ b/Tests/Parser/main.cpp
@@ -2,6 +2,7 @@
 #include <AMReX_Parser.H>
 #include <AMReX_IParser.H>
 #include <map>
+#include <numbers>
 
 using namespace amrex;
 
@@ -311,7 +312,7 @@ int main (int argc, char* argv[])
                         {"x","y","z"},
                         [=] (double x, double, double z) -> double {
                             double pi = 3.14;
-                            return 2.*std::sqrt(2.)+std::sqrt(-std::log(x))*std::cos(2*pi*z);
+                            return 2.*std::numbers::sqrt2+std::sqrt(-std::log(x))*std::cos(2*pi*z);
                         },
                         {0.5, 0.8, 0.3}, {16, 16, 16}, 100,
                         1.e-12, 1.e-15);
diff --git a/Tools/Plotfile/fvolumesum.cpp b/Tools/Plotfile/fvolumesum.cpp
index 398d5831c85..71129371a89 100644
--- a/Tools/Plotfile/fvolumesum.cpp
+++ b/Tools/Plotfile/fvolumesum.cpp
@@ -1,4 +1,5 @@
 #include <AMReX.H>
+#include <numbers>
 #include <AMReX_Print.H>
 #include <AMReX_PlotFileUtil.H>
 #include <AMReX_MultiFabUtil.H>
@@ -89,7 +90,7 @@ void main_main()
 
     int coord = pf.coordSys();
 
-    constexpr Real pi = 3.1415926535897932;
+    constexpr Real pi = std::numbers::pi_v<Real>;
 
     for (int ilev = 0; ilev <= fine_level; ++ilev) {