Test_extended_meson_field: add --Ni and --Nj command line options

Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>

Grid/lattice/Lattice_slicesum_core.h

@@ -59,7 +59,7 @@ inline void sliceSumReduction_cub_small(const vobj *Data,
 #if defined(__CUDACC__) && (__CUDACC_VER_MAJOR__ >= 13)
   #define GRID_CUB_SUM_OP ::cuda::std::plus<>{}
 #else
-  #define GRID_CUB_SUM_OP ::cub::Sum()
+  #define GRID_CUB_SUM_OP ::gpucub::Sum()
 #endif
   
   gpuError_t gpuErr = gpucub::DeviceSegmentedReduce::Reduce(temp_storage_array, temp_storage_bytes, rb_p,d_out, rd, d_offsets, d_offsets+1, GRID_CUB_SUM_OP, zero_init, computeStream);

Grid/simd/Simd.h

@@ -113,6 +113,14 @@ accelerator_inline RealD    adj(const RealD  & r){ return r; }
 accelerator_inline ComplexD adj(const ComplexD& r){ return(conjugate(r)); }
 accelerator_inline ComplexF adj(const ComplexF& r ){ return(conjugate(r)); }
 
+#if defined(GRID_CUDA) || defined(GRID_HIP)
+//Provide for convenience
+accelerator_inline std::complex<double> conjugate(const std::complex<double>& r){ return(conj(r)); }
+accelerator_inline std::complex<float>  conjugate(const std::complex<float>& r) { return(conj(r)); }
+accelerator_inline std::complex<double> adj(const std::complex<double>& r)      { return(conj(r)); }
+accelerator_inline std::complex<float>  adj(const std::complex<float>& r)       { return(conj(r)); }
+#endif
+
 accelerator_inline RealF real(const RealF  & r){ return r; }
 accelerator_inline RealD real(const RealD  & r){ return r; }
 accelerator_inline RealF real(const ComplexF  & r){ return r.real(); }

Grid/tensors/Tensor_inner.h

@@ -32,6 +32,33 @@ Author: Christoph Lehner <christoph@lhnr.de>
 
 NAMESPACE_BEGIN(Grid);
 
+//////////////////////////////////////
+// innerProductD scalar overloads must be visible before norm2 is defined
+//////////////////////////////////////
+accelerator_inline ComplexD innerProductD(const ComplexF &l,const ComplexF &r){  return innerProduct(l,r); }
+accelerator_inline ComplexD innerProductD(const ComplexD &l,const ComplexD &r){  return innerProduct(l,r); }
+accelerator_inline RealD    innerProductD(const RealD    &l,const RealD    &r){  return innerProduct(l,r); }
+accelerator_inline RealD    innerProductD(const RealF    &l,const RealF    &r){  return innerProduct(l,r); }
+
+accelerator_inline vComplexD innerProductD(const vComplexD &l,const vComplexD &r){  return innerProduct(l,r); }
+accelerator_inline vRealD    innerProductD(const vRealD    &l,const vRealD    &r){  return innerProduct(l,r); }
+accelerator_inline vComplexD innerProductD(const vComplexF &l,const vComplexF &r)
+{
+  vComplexD la,lb;
+  vComplexD ra,rb;
+  Optimization::PrecisionChange::StoD(l.v,la.v,lb.v);
+  Optimization::PrecisionChange::StoD(r.v,ra.v,rb.v);
+  return innerProduct(la,ra) + innerProduct(lb,rb);
+}
+accelerator_inline vRealD innerProductD(const vRealF &l,const vRealF &r)
+{
+  vRealD la,lb;
+  vRealD ra,rb;
+  Optimization::PrecisionChange::StoD(l.v,la.v,lb.v);
+  Optimization::PrecisionChange::StoD(r.v,ra.v,rb.v);
+  return innerProduct(la,ra) + innerProduct(lb,rb);
+}
+
 ///////////////////////////////////////////////////////////////////////////////////////
 // innerProduct Scalar x Scalar -> Scalar
 // innerProduct Vector x Vector -> Scalar
@@ -138,30 +165,8 @@ auto Reduce (const iScalar<l>& lhs) -> iScalar<decltype(Reduce(lhs._internal))>
 
 //////////////////////////////////////
 // innerProductD : if single promote to double and evaluate with sum 2x
+// (scalar/vector overloads are declared above norm2 for ADL visibility)
 //////////////////////////////////////
-accelerator_inline ComplexD innerProductD(const ComplexF &l,const ComplexF &r){  return innerProduct(l,r); }
-accelerator_inline ComplexD innerProductD(const ComplexD &l,const ComplexD &r){  return innerProduct(l,r); }
-accelerator_inline RealD    innerProductD(const RealD    &l,const RealD    &r){  return innerProduct(l,r); }
-accelerator_inline RealD    innerProductD(const RealF    &l,const RealF    &r){  return innerProduct(l,r); }
-
-accelerator_inline vComplexD innerProductD(const vComplexD &l,const vComplexD &r){  return innerProduct(l,r); }
-accelerator_inline vRealD    innerProductD(const vRealD    &l,const vRealD    &r){  return innerProduct(l,r); }
-accelerator_inline vComplexD innerProductD(const vComplexF &l,const vComplexF &r)
-{  
-  vComplexD la,lb;
-  vComplexD ra,rb;
-  Optimization::PrecisionChange::StoD(l.v,la.v,lb.v);
-  Optimization::PrecisionChange::StoD(r.v,ra.v,rb.v);
-  return innerProduct(la,ra) + innerProduct(lb,rb); 
-}
-accelerator_inline vRealD innerProductD(const vRealF &l,const vRealF &r)
-{  
-  vRealD la,lb;
-  vRealD ra,rb;
-  Optimization::PrecisionChange::StoD(l.v,la.v,lb.v);
-  Optimization::PrecisionChange::StoD(r.v,ra.v,rb.v);
-  return innerProduct(la,ra) + innerProduct(lb,rb); 
-}
 
 // Now do it for vector, matrix, scalar
 template<class l,class r,int N> accelerator_inline

systems/Frontier/config-command

@@ -1,4 +1,3 @@
-CLIME=`spack find --paths c-lime@2-3-9 | grep c-lime| cut -c 15-`
 ../../configure --enable-comms=mpi-auto \
 --with-lime=$CLIME \
 --enable-unified=no \
@@ -9,12 +8,13 @@ CLIME=`spack find --paths c-lime@2-3-9 | grep c-lime| cut -c 15-`
 --disable-gparity \
 --disable-fermion-reps \
 --enable-simd=GPU \
---with-gmp=$OLCF_GMP_ROOT \
---with-mpfr=/opt/cray/pe/gcc/mpfr/3.1.4/ \
+--with-gmp=$GMP \
+--with-mpfr=$MPFR \
+--with-openssl=$OPENSSL \
 --disable-fermion-reps \
 CXX=hipcc MPICXX=mpicxx \
-CXXFLAGS="-fPIC -I${ROCM_PATH}/include/ -I${MPICH_DIR}/include -L/lib64 " \
- LDFLAGS="-L/lib64 -L${ROCM_PATH}/lib -L${MPICH_DIR}/lib -lmpi -L${CRAY_MPICH_ROOTDIR}/gtl/lib -lmpi_gtl_hsa -lhipblas -lrocblas -lhipfft"
+CXXFLAGS="-fPIC -I${ROCM_PATH}/include/ -I${MPICH_DIR}/include " \
+ LDFLAGS="-L${ROCM_PATH}/lib -L${MPICH_DIR}/lib -lmpi -lmpi_gtl_hsa -lhipblas -lrocblas -lhipfft -lamdhip64"
 
 
 

systems/Frontier/sourceme.sh

@@ -1,16 +1,14 @@
 
 echo spack
-. /autofs/nccs-svm1_home1/paboyle/Crusher/Grid/spack/share/spack/setup-env.sh
+. /autofs/nccs-svm1_home1/paboyle/spack/share/spack/setup-env.sh
 
-module load amd/7.0.2
-module load cray-fftw
-module load craype-accel-amd-gfx90a
-mkdir $HOME/LD_PATH
-ln -s /opt/rocm-6.4.2/lib/libamdhip* $HOME/LD_PATH
+export CLIME=`spack find --paths c-lime | grep ^c-lime | awk '{print $2}' `
+export MPFR=`spack find --paths mpfr    | grep ^mpfr  | awk '{print $2}' `
+export OPENSSL=`spack find --paths openssl | grep openssl | awk  '{print $2}' `
+export GMP=`spack find --paths gmp      | grep ^gmp | awk '{print $2}' `
 
-#Ugly hacks to get down level software working on current system
-export LD_LIBRARY_PATH=/opt/cray/libfabric/1.20.1/lib64/:$LD_LIBRARY_PATH
-export LD_LIBRARY_PATH=/opt/gcc/mpfr/3.1.4/lib:$LD_LIBRARY_PATH
-#export LD_LIBRARY_PATH=`pwd`/:$LD_LIBRARY_PATH
-export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:$HOME/LD_PATH/
-export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/opt/rocm-7.0.2/lib
+module load cce/21.0.0
+module load cpe/26.03
+module load rocm/7.0.2
+export LD_LIBRARY_PATH=$CRAY_LD_LIBRARY_PATH:$LD_LIBRARY_PATH
+export LD_LIBRARY_PATH=/opt/rocm-7.0.2/lib/llvm/lib/:$LD_LIBRARY_PATH

systems/Perlmutter/config-command

@@ -1,12 +1,12 @@
 DIR=`pwd`
-PREFIX=$HOME/DDHMC/Grid/systems/Prerequisites/install/
+
 ../../configure \
     --enable-comms=mpi \
     --enable-simd=GPU \
     --enable-shm=nvlink \
     --enable-gen-simd-width=64 \
-    --with-gmp=$PREFIX \
-    --with-mpfr=$PREFIX \
+    --with-gmp=$GMP \
+    --with-mpfr=$MPFR \
     --enable-accelerator=cuda \
     --disable-fermion-reps \
     --disable-unified \

systems/Perlmutter/sourceme.sh

@@ -1,4 +1,6 @@
-
 export CRAY_ACCEL_TARGET=nvidia80
+source /global/homes/p/pboyle/spack/share/spack/setup-env.sh
+export MPFR=`spack find --paths mpfr | grep mpfr | cut -c 13-`
+export GMP=`spack find --paths gmp   | grep gmp  | cut -c 12-`
 
-module load PrgEnv-gnu cpe-cuda cudatoolkit/11.4
+module load PrgEnv-gnu cpe-cuda cudatoolkit/12.0

systems/mac-arm/sourceme.sh

@@ -1,11 +1,7 @@
 source /Users/peterboyle/QCD//Spack/spack//share/spack/setup-env.sh
 export FFTW=`spack find --paths fftw    | grep ^fftw   | awk '{print $2}' `
-#export HDF5=`spack find --paths hdf5+cxx   | grep ^hdf5   | awk '{print $2}' `
 export CLIME=`spack find --paths c-lime | grep ^c-lime | awk '{print $2}' `
 export MPFR=`spack find --paths mpfr    | grep ^mpfr  | awk '{print $2}' `
 export OPENSSL=`spack find --paths openssl | grep openssl | awk  '{print $2}' `
 export GMP=`spack find --paths gmp      | grep ^gmp | awk '{print $2}' `
 
-export LD_LIBRARY_PATH=$MPFR/lib:$LD_LIBRARY_PATH
-export LD_LIBRARY_PATH=$GMP/lib:$LD_LIBRARY_PATH
-

tests/Test_extended_meson_field.cc

@@ -481,14 +481,16 @@ void A2ALoopLeftContractionType3(PropagatorField &tloop, const PropagatorField &
   uint64_t Osites = loop.Grid()->oSites();
   int Nsimd = SpinColourMatrix_v::Nsimd();
   accelerator_for(ss, Osites, Nsimd, {
+    typedef decltype(coalescedRead(loopv[0])) calcSCMatrix;
+    typedef iSpinMatrix<typename calcSCMatrix::vector_type> calcSpinMatrix;
     auto l = loopv(ss);
-    SpinMatrix_v spinLoop; spinLoop = Zero();
+    calcSpinMatrix spinLoop; spinLoop = Zero();
     for (int s1 = 0; s1 < Ns; ++s1)
     for (int s2 = 0; s2 < Ns; ++s2)
       spinLoop()(s1,s2)() = l()(s1,s2)(0,0) + l()(s1,s2)(1,1) + l()(s1,s2)(2,2);
     auto tmp = l; tmp = Zero();
     for (int mu = 0; mu < ng; ++mu) {
-      SpinMatrix_v tmp2 = Gamma(g1[mu]) * spinLoop * Gamma(g2[mu]);
+      calcSpinMatrix tmp2 = Gamma(g1[mu]) * spinLoop * Gamma(g2[mu]);
       for (int s1 = 0; s1 < Ns; ++s1)
       for (int s2 = 0; s2 < Ns; ++s2)
         tmp()(s1,s2)(0,0) = tmp()(s1,s2)(0,0) + tmp2()(s1,s2)();
@@ -514,13 +516,16 @@ void A2ALoopRightContractionType0(FermionField &loopRight,
   uint64_t Osites = right.Grid()->oSites();
   int Nsimd = SpinColourVector_v::Nsimd();
   accelerator_for(ss, Osites, Nsimd, {
+    typedef decltype(coalescedRead(rv[0]))  calcSCVector;
+    typedef decltype(coalescedRead(tlv[0])) calcSCMatrix;
+    typedef iSpinMatrix<typename calcSCMatrix::vector_type> calcSpinMatrix;
     auto loopm  = tlv(ss);
     auto rightv = rv(ss);
-    SpinMatrix_v spinLoop; spinLoop = Zero();
+    calcSpinMatrix spinLoop; spinLoop = Zero();
     for (int s1 = 0; s1 < Ns; ++s1)
     for (int s2 = 0; s2 < Ns; ++s2)
       spinLoop()(s1,s2)() = loopm()(s1,s2)(0,0);
-    SpinColourVector_v lR; lR = Zero();
+    calcSCVector lR; lR = Zero();
     for (int mu = 0; mu < ng; ++mu) {
       auto GLoop   = Gamma(g2[mu]) * spinLoop;
       auto trGLoop = GLoop()(0,0)() + GLoop()(1,1)() + GLoop()(2,2)() + GLoop()(3,3)();
@@ -562,9 +567,10 @@ void A2ALoopRightContractionType2(FermionField &loopRight,
   uint64_t Osites = right.Grid()->oSites();
   int Nsimd = SpinColourVector_v::Nsimd();
   accelerator_for(ss, Osites, Nsimd, {
+    typedef decltype(coalescedRead(rv[0])) calcSCVector;
     auto loopm  = tlv(ss);
     auto rightv = rv(ss);
-    SpinColourVector_v lR; lR = Zero();
+    calcSCVector lR; lR = Zero();
     for (int mu = 0; mu < ng; ++mu) {
       int s1 = mu / Ns;
       int s2 = mu % Ns;
@@ -591,9 +597,10 @@ void A2ALoopRightContractionType3(FermionField &loopRight,
   uint64_t Osites = right.Grid()->oSites();
   int Nsimd = SpinColourVector_v::Nsimd();
   accelerator_for(ss, Osites, Nsimd, {
+    typedef decltype(coalescedRead(rv[0])) calcSCVector;
     auto loopm  = tlv(ss);
     auto rightv = rv(ss);
-    SpinColourVector_v lR; lR = Zero();
+    calcSCVector lR; lR = Zero();
     for (int s = 0; s < Ns; ++s)
     for (int c = 0; c < Nc; ++c)
       lR()(s)(c) = loopm()(s,0)(0,0) * rightv()(0)(c)
@@ -703,6 +710,11 @@ int main(int argc, char *argv[])
   int N_j   = 8;
   int Nloop = 4;
 
+  if (GridCmdOptionExists(argv, argv+argc, "--Ni"))
+    N_i = std::stoi(GridCmdOptionPayload(argv, argv+argc, "--Ni"));
+  if (GridCmdOptionExists(argv, argv+argc, "--Nj"))
+    N_j = std::stoi(GridCmdOptionPayload(argv, argv+argc, "--Nj"));
+
   GridParallelRNG pRNG(&grid);
   pRNG.SeedFixedIntegers({1, 2, 3, 4});
 
@@ -753,13 +765,13 @@ int main(int argc, char *argv[])
     for (int t  = 0; t  < Nt;  t++)
     for (int ii = 0; ii < N_i; ii++)
     for (int jj = 0; jj < N_j; jj++) {
-      norm2_ref  += std::norm(result_ref(t, ii, jj));
-      norm2_blas += std::norm(result_blas(t, ii, jj));
-      norm2_gpu  += std::norm(result_gpu(t, ii, jj));
+      norm2_ref  += norm2(result_ref(t, ii, jj));
+      norm2_blas += norm2(result_blas(t, ii, jj));
+      norm2_gpu  += norm2(result_gpu(t, ii, jj));
       ComplexD diff_blas = result_ref(t, ii, jj) - result_blas(t, ii, jj);
       ComplexD diff_gpu  = result_ref(t, ii, jj) - result_gpu(t, ii, jj);
-      norm2_diff_blas += std::norm(diff_blas);
-      norm2_diff_gpu  += std::norm(diff_gpu);
+      norm2_diff_blas += norm2(diff_blas);
+      norm2_diff_gpu  += norm2(diff_gpu);
     }
 
     double rel_blas = (norm2_ref > 0) ? std::sqrt(norm2_diff_blas / norm2_ref) : 0.0;

tests/debug/Test_hipfft_bug_fail.cc

@@ -0,0 +1,76 @@
+/*
+ * Isolating the hipfft HIPFFT_PARSE_ERROR on ROCm 7 / hipFFT 1.0.20.
+ *
+ * Tests three orderings with an empty rocFFT cache to find which GPU
+ * operation before plan creation triggers the failure:
+ *   A) hipMalloc only            — hypothesis: passes (no async GPU work)
+ *   B) hipMalloc + hipMemset     — hypothesis: fails  (async GPU work in flight)
+ *   C) hipMalloc + hipMemset     — hypothesis: passes (work completed before plan)
+ *      + hipDeviceSynchronize
+ *
+ * Compile:
+ *   hipcc -o Test_hipfft_bug_fail Test_hipfft_bug_fail.cc -lhipfft
+ *
+ * Run with empty cache:
+ *   rm -rf ~/.cache/
+ *   ./Test_hipfft_bug_fail
+ */
+
+#include <cstdio>
+#include <hipfft/hipfft.h>
+#include <hip/hip_runtime.h>
+
+static const char *res(hipfftResult rv) {
+  return rv == HIPFFT_SUCCESS ? "SUCCESS" : "PARSE_ERROR";
+}
+
+static hipfftResult makePlan(int G, int howmany) {
+  int n[] = {G};
+  hipfftHandle p;
+  size_t workSize = 0;
+  hipfftCreate(&p);
+  hipfftResult rv = hipfftMakePlanMany(p, 1, n,
+                                       nullptr, 1, G, nullptr, 1, G,
+                                       HIPFFT_Z2Z, howmany, &workSize);
+  hipfftDestroy(p);
+  return rv;
+}
+
+int main(void) {
+  hipDeviceProp_t prop;
+  hipGetDeviceProperties(&prop, 0);
+  printf("Device: %s\n", prop.name);
+#ifdef hipfftVersionMinor
+  printf("hipFFT version: %d.%d.%d\n\n",
+         hipfftVersionMajor, hipfftVersionMinor, hipfftVersionPatch);
+#endif
+
+  for (int G : {4, 8, 16, 32}) {
+    int howmany = 512;
+    long nelems = (long)G * howmany;
+    hipfftDoubleComplex *buf = nullptr;
+    hipMalloc(&buf, nelems * sizeof(hipfftDoubleComplex));
+
+    // Tests ordered so each runs before a prior success can populate the cache.
+
+    // B first: hipMalloc + hipMemset (async GPU work in flight)
+    // If this fails, A (no hipMemset) will pass, confirming hipMemset is the trigger.
+    hipMemset(buf, 0, nelems * sizeof(hipfftDoubleComplex));
+    hipfftResult rvB = makePlan(G, howmany);
+    printf("G=%-4d  B) hipMalloc + hipMemset       : %s\n", G, res(rvB));
+
+    // C: hipMalloc + hipMemset + sync — does syncing before plan creation fix it?
+    hipMemset(buf, 0, nelems * sizeof(hipfftDoubleComplex));
+    hipDeviceSynchronize();
+    hipfftResult rvC = makePlan(G, howmany);
+    printf("G=%-4d  C) hipMalloc + hipMemset + sync: %s\n", G, res(rvC));
+
+    // A last: hipMalloc only, no async GPU work — should always pass
+    hipfftResult rvA = makePlan(G, howmany);
+    printf("G=%-4d  A) hipMalloc only             : %s\n\n", G, res(rvA));
+
+    hipFree(buf);
+  }
+
+  return 0;
+}

tests/debug/Test_hipfft_bug_pass.cc

@@ -0,0 +1,61 @@
+/*
+ * Minimal program demonstrating the workaround for the hipfft ROCm 7 bug.
+ *
+ * Workaround: create the hipfft plan BEFORE any hipMalloc.  Plan creation
+ * for G < 32 then succeeds even with an empty rocFFT cache.
+ *
+ * Compile:
+ *   hipcc -o Test_hipfft_bug_pass Test_hipfft_bug_pass.cc -lhipfft
+ *
+ * Run:
+ *   rm -rf ~/.cache/rocfft
+ *   ./Test_hipfft_bug_pass
+ *
+ * Expected: all G values succeed.
+ * Compare with Test_hipfft_bug_fail.cc which uses the opposite ordering.
+ */
+
+#include <cstdio>
+#include <hipfft/hipfft.h>
+#include <hip/hip_runtime.h>
+
+int main(void) {
+  hipDeviceProp_t prop;
+  hipGetDeviceProperties(&prop, 0);
+  printf("Device: %s\n", prop.name);
+#ifdef hipfftVersionMinor
+  printf("hipFFT version: %d.%d.%d\n\n",
+         hipfftVersionMajor, hipfftVersionMinor, hipfftVersionPatch);
+#endif
+
+  for (int G : {8, 16, 32}) {
+    int howmany = 512;
+    int n[] = {G};
+    long nelems = (long)G * howmany;
+
+    // Plan created BEFORE hipMalloc — succeeds for all G
+    hipfftHandle p;
+    size_t workSize = 0;
+    hipfftCreate(&p);
+    hipfftResult rv = hipfftMakePlanMany(p, 1, n,
+                                         nullptr, 1, G, nullptr, 1, G,
+                                         HIPFFT_Z2Z, howmany, &workSize);
+    printf("G=%-4d  plan-then-hipMalloc: %d (%s)\n",
+           G, (int)rv, rv == HIPFFT_SUCCESS ? "HIPFFT_SUCCESS" : "HIPFFT_PARSE_ERROR");
+
+    if (rv == HIPFFT_SUCCESS) {
+      hipfftDoubleComplex *buf = nullptr;
+      hipMalloc(&buf, nelems * sizeof(hipfftDoubleComplex));
+      hipMemset(buf, 0, nelems * sizeof(hipfftDoubleComplex));
+      rv = hipfftExecZ2Z(p, buf, buf, HIPFFT_FORWARD);
+      hipDeviceSynchronize();
+      printf("G=%-4d  execFwd:             %d (%s)\n",
+             G, (int)rv, rv == HIPFFT_SUCCESS ? "HIPFFT_SUCCESS" : "FAILED");
+      hipFree(buf);
+    }
+
+    hipfftDestroy(p);
+  }
+
+  return 0;
+}

tests/debug/Test_hipfft_minimal.cc

@@ -0,0 +1,161 @@
+/*
+ * Minimal reproducer for hipfftMakePlanMany / hipfftPlanMany failures.
+ *
+ * Compile on Frontier (no Grid headers needed):
+ *   hipcc -o Test_hipfft_minimal Test_hipfft_minimal.cc -lhipfft
+ *
+ * Run:
+ *   ./Test_hipfft_minimal
+ */
+
+#include <cstdio>
+#include <cstdlib>
+#include <hipfft/hipfft.h>
+#include <hip/hip_runtime.h>
+
+static const char *hipfftResultString(hipfftResult r) {
+  switch (r) {
+  case HIPFFT_SUCCESS:                  return "HIPFFT_SUCCESS";
+  case HIPFFT_INVALID_PLAN:             return "HIPFFT_INVALID_PLAN";
+  case HIPFFT_ALLOC_FAILED:             return "HIPFFT_ALLOC_FAILED";
+  case HIPFFT_INVALID_TYPE:             return "HIPFFT_INVALID_TYPE";
+  case HIPFFT_INVALID_VALUE:            return "HIPFFT_INVALID_VALUE";
+  case HIPFFT_INTERNAL_ERROR:           return "HIPFFT_INTERNAL_ERROR";
+  case HIPFFT_EXEC_FAILED:              return "HIPFFT_EXEC_FAILED";
+  case HIPFFT_SETUP_FAILED:             return "HIPFFT_SETUP_FAILED";
+  case HIPFFT_INVALID_SIZE:             return "HIPFFT_INVALID_SIZE";
+  case HIPFFT_UNALIGNED_DATA:           return "HIPFFT_UNALIGNED_DATA";
+  case HIPFFT_INCOMPLETE_PARAMETER_LIST:return "HIPFFT_INCOMPLETE_PARAMETER_LIST";
+  case HIPFFT_INVALID_DEVICE:           return "HIPFFT_INVALID_DEVICE";
+  case HIPFFT_PARSE_ERROR:              return "HIPFFT_PARSE_ERROR";
+  case HIPFFT_NO_WORKSPACE:             return "HIPFFT_NO_WORKSPACE";
+  case HIPFFT_NOT_IMPLEMENTED:          return "HIPFFT_NOT_IMPLEMENTED";
+  case HIPFFT_NOT_SUPPORTED:            return "HIPFFT_NOT_SUPPORTED";
+  default:                              return "UNKNOWN";
+  }
+}
+
+// Plan creation + execution for (G, howmany).
+// Tests two orderings to isolate whether a prior hipMalloc poisons hipfft
+// plan creation for small G on ROCm 7:
+//   A) plan BEFORE hipMalloc  — hypothesis: succeeds
+//   B) hipMalloc BEFORE plan  — hypothesis: fails for G < 32
+static void tryPlanAndExec(int G, long howmany) {
+  int n[] = {G};
+  long nelems = (long)G * howmany;
+
+  printf("--- G=%-4d  howmany=%-10ld  total_elems=%-12ld ---\n",
+         G, howmany, nelems);
+
+  // Allocate device buffer (hipfftDoubleComplex = 16 bytes each)
+  hipfftDoubleComplex *dbuf = nullptr;
+  hipError_t herr = hipMalloc(&dbuf, nelems * sizeof(hipfftDoubleComplex));
+  if (herr != hipSuccess) {
+    printf("  hipMalloc failed (%d) for %ld elems — skipping\n\n", (int)herr, nelems);
+    return;
+  }
+  hipMemset(dbuf, 0, nelems * sizeof(hipfftDoubleComplex));
+
+  // 1. hipfftPlanMany (one-step, nullptr embed) — current Grid path
+  {
+    hipfftHandle p;
+    hipfftResult rv = hipfftPlanMany(&p, 1, n,
+                                     nullptr, 1, G,
+                                     nullptr, 1, G,
+                                     HIPFFT_Z2Z, (int)howmany);
+    printf("  hipfftPlanMany   create : %d (%s)\n", (int)rv, hipfftResultString(rv));
+    if (rv == HIPFFT_SUCCESS) {
+      rv = hipfftExecZ2Z(p, dbuf, dbuf, HIPFFT_FORWARD);
+      hipDeviceSynchronize();
+      printf("  hipfftPlanMany   execFwd: %d (%s)\n", (int)rv, hipfftResultString(rv));
+      hipfftDestroy(p);
+    }
+  }
+
+  // 2. hipfftCreate + hipfftMakePlanMany (two-step) — also current Grid path
+  {
+    hipfftHandle p;
+    size_t workSize = 0;
+    hipfftResult rc = hipfftCreate(&p);
+    if (rc == HIPFFT_SUCCESS) {
+      hipfftResult rv = hipfftMakePlanMany(p, 1, n,
+                                           nullptr, 1, G,
+                                           nullptr, 1, G,
+                                           HIPFFT_Z2Z, (int)howmany, &workSize);
+      printf("  hipfftMakePlanMany      : %d (%s)  workSize=%zu\n",
+             (int)rv, hipfftResultString(rv), workSize);
+      if (rv == HIPFFT_SUCCESS) {
+        rv = hipfftExecZ2Z(p, dbuf, dbuf, HIPFFT_FORWARD);
+        hipDeviceSynchronize();
+        printf("  hipfftMakePlanMany exec : %d (%s)\n", (int)rv, hipfftResultString(rv));
+      }
+      hipfftDestroy(p);
+    } else {
+      printf("  hipfftCreate            : %d (%s)\n", (int)rc, hipfftResultString(rc));
+    }
+  }
+
+  // 3. hipfftPlan1d (simplest API, batch = howmany)
+  {
+    hipfftHandle p;
+    hipfftResult rv = hipfftPlan1d(&p, G, HIPFFT_Z2Z, (int)howmany);
+    printf("  hipfftPlan1d     create : %d (%s)\n", (int)rv, hipfftResultString(rv));
+    if (rv == HIPFFT_SUCCESS) {
+      rv = hipfftExecZ2Z(p, dbuf, dbuf, HIPFFT_FORWARD);
+      hipDeviceSynchronize();
+      printf("  hipfftPlan1d     execFwd: %d (%s)\n", (int)rv, hipfftResultString(rv));
+      hipfftDestroy(p);
+    }
+  }
+
+  hipFree(dbuf);
+  printf("\n");
+}
+
+int main(void) {
+  // Print HIP device info
+  int device = 0;
+  hipGetDevice(&device);
+  hipDeviceProp_t prop;
+  hipGetDeviceProperties(&prop, device);
+  printf("Device %d: %s  warpSize=%d\n\n", device, prop.name, prop.warpSize);
+
+#ifdef hipfftVersionMinor
+  printf("hipFFT version: %d.%d.%d\n\n",
+         hipfftVersionMajor, hipfftVersionMinor, hipfftVersionPatch);
+#endif
+
+  // Original sweep with small howmany (these passed first time)
+  printf("=== Small howmany (original sweep) ===\n\n");
+  for (int G : {4, 8, 12, 16, 24, 32, 48, 64})
+    tryPlanAndExec(G, 512);
+
+  // Grid-realistic howmany values derived from actual lattice geometries.
+  // howmany = Ncomp * product(ldimensions[d] for d != dim)
+  // For LatticeComplexD: Ncomp=1.
+  printf("=== Grid-realistic parameters ===\n\n");
+
+  // --grid 16.16.16.16  4D FFT  (KNOWN TO FAIL in Grid)
+  // Each dim: G=16, Nperp=16^3=4096
+  tryPlanAndExec(16, 4096);
+
+  // --grid 32.32.32.32  4D FFT  (KNOWN TO SUCCEED in Grid)
+  // Each dim: G=32, Nperp=32^3=32768
+  tryPlanAndExec(32, 32768);
+
+  // --grid 32.32.32.32 Ls=8  5D DWF FFT  (KNOWN TO FAIL on dim 0 in Grid)
+  // dim 0: G=8,  Nperp=32^4=1048576
+  tryPlanAndExec(8, 1048576);
+  // dim 1-4: G=32, Nperp=8*32^3=262144
+  tryPlanAndExec(32, 262144);
+
+  // Extra intermediate cases to bracket the failure
+  tryPlanAndExec(16, 1024);
+  tryPlanAndExec(16, 2048);
+  tryPlanAndExec(16, 8192);
+  tryPlanAndExec(8,  4096);
+  tryPlanAndExec(8,  65536);
+  tryPlanAndExec(8,  262144);
+
+  return 0;
+}

tests/debug/Test_hipfft_repro.cc

@@ -0,0 +1,168 @@
+/*
+ * Reproducer for HIPFFT_PARSE_ERROR (error 12) from hipfftMakePlanMany on
+ * ROCm 7 / hipFFT 1.0.20 (Frontier, MI210 login and MI250X compute nodes).
+ *
+ * Observed failure: G < 32 returns HIPFFT_PARSE_ERROR from all three plan
+ * creation APIs (hipfftPlanMany, hipfftMakePlanMany, hipfftPlan1d) when a
+ * device buffer is allocated and zeroed with hipMalloc+hipMemset before the
+ * plan creation call.  G >= 32 succeeds.
+ *
+ * Contrast with Test_hipfft_minimal.cc (plan-first ordering) which passes
+ * for all G even with an empty rocFFT cache.
+ *
+ * Compile on Frontier (no Grid headers needed):
+ *   hipcc -o Test_hipfft_repro Test_hipfft_repro.cc -lhipfft
+ *
+ * Run with empty cache to reproduce the failure:
+ *   rm -rf ~/.cache/rocfft
+ *   ./Test_hipfft_repro
+ */
+
+#include <cstdio>
+#include <cstdlib>
+#include <hipfft/hipfft.h>
+#include <hip/hip_runtime.h>
+
+static const char *hipfftResultString(hipfftResult r) {
+  switch (r) {
+  case HIPFFT_SUCCESS:                  return "HIPFFT_SUCCESS";
+  case HIPFFT_INVALID_PLAN:             return "HIPFFT_INVALID_PLAN";
+  case HIPFFT_ALLOC_FAILED:             return "HIPFFT_ALLOC_FAILED";
+  case HIPFFT_INVALID_TYPE:             return "HIPFFT_INVALID_TYPE";
+  case HIPFFT_INVALID_VALUE:            return "HIPFFT_INVALID_VALUE";
+  case HIPFFT_INTERNAL_ERROR:           return "HIPFFT_INTERNAL_ERROR";
+  case HIPFFT_EXEC_FAILED:              return "HIPFFT_EXEC_FAILED";
+  case HIPFFT_SETUP_FAILED:             return "HIPFFT_SETUP_FAILED";
+  case HIPFFT_INVALID_SIZE:             return "HIPFFT_INVALID_SIZE";
+  case HIPFFT_UNALIGNED_DATA:           return "HIPFFT_UNALIGNED_DATA";
+  case HIPFFT_INCOMPLETE_PARAMETER_LIST:return "HIPFFT_INCOMPLETE_PARAMETER_LIST";
+  case HIPFFT_INVALID_DEVICE:           return "HIPFFT_INVALID_DEVICE";
+  case HIPFFT_PARSE_ERROR:              return "HIPFFT_PARSE_ERROR";
+  case HIPFFT_NO_WORKSPACE:             return "HIPFFT_NO_WORKSPACE";
+  case HIPFFT_NOT_IMPLEMENTED:          return "HIPFFT_NOT_IMPLEMENTED";
+  case HIPFFT_NOT_SUPPORTED:            return "HIPFFT_NOT_SUPPORTED";
+  default:                              return "UNKNOWN";
+  }
+}
+
+// Plan creation + execution for (G, howmany) using hipfftCreate+hipfftMakePlanMany.
+// This is the path Grid's FFT.h now uses.
+static void tryPlanAndExec(int G, long howmany) {
+  int n[] = {G};
+  long nelems = (long)G * howmany;
+
+  printf("--- G=%-4d  howmany=%-10ld  total_elems=%-12ld ---\n",
+         G, howmany, nelems);
+
+  // Allocate device buffer (hipfftDoubleComplex = 16 bytes each)
+  hipfftDoubleComplex *dbuf = nullptr;
+  hipError_t herr = hipMalloc(&dbuf, nelems * sizeof(hipfftDoubleComplex));
+  if (herr != hipSuccess) {
+    printf("  hipMalloc failed (%d) for %ld elems — skipping\n\n", (int)herr, nelems);
+    return;
+  }
+  hipMemset(dbuf, 0, nelems * sizeof(hipfftDoubleComplex));
+
+  // 1. hipfftPlanMany (one-step, nullptr embed) — current Grid path
+  {
+    hipfftHandle p;
+    hipfftResult rv = hipfftPlanMany(&p, 1, n,
+                                     nullptr, 1, G,
+                                     nullptr, 1, G,
+                                     HIPFFT_Z2Z, (int)howmany);
+    printf("  hipfftPlanMany   create : %d (%s)\n", (int)rv, hipfftResultString(rv));
+    if (rv == HIPFFT_SUCCESS) {
+      rv = hipfftExecZ2Z(p, dbuf, dbuf, HIPFFT_FORWARD);
+      hipDeviceSynchronize();
+      printf("  hipfftPlanMany   execFwd: %d (%s)\n", (int)rv, hipfftResultString(rv));
+      hipfftDestroy(p);
+    }
+  }
+
+  // 2. hipfftCreate + hipfftMakePlanMany (two-step) — also current Grid path
+  {
+    hipfftHandle p;
+    size_t workSize = 0;
+    hipfftResult rc = hipfftCreate(&p);
+    if (rc == HIPFFT_SUCCESS) {
+      hipfftResult rv = hipfftMakePlanMany(p, 1, n,
+                                           nullptr, 1, G,
+                                           nullptr, 1, G,
+                                           HIPFFT_Z2Z, (int)howmany, &workSize);
+      printf("  hipfftMakePlanMany      : %d (%s)  workSize=%zu\n",
+             (int)rv, hipfftResultString(rv), workSize);
+      if (rv == HIPFFT_SUCCESS) {
+        rv = hipfftExecZ2Z(p, dbuf, dbuf, HIPFFT_FORWARD);
+        hipDeviceSynchronize();
+        printf("  hipfftMakePlanMany exec : %d (%s)\n", (int)rv, hipfftResultString(rv));
+      }
+      hipfftDestroy(p);
+    } else {
+      printf("  hipfftCreate            : %d (%s)\n", (int)rc, hipfftResultString(rc));
+    }
+  }
+
+  // 3. hipfftPlan1d (simplest API, batch = howmany)
+  {
+    hipfftHandle p;
+    hipfftResult rv = hipfftPlan1d(&p, G, HIPFFT_Z2Z, (int)howmany);
+    printf("  hipfftPlan1d     create : %d (%s)\n", (int)rv, hipfftResultString(rv));
+    if (rv == HIPFFT_SUCCESS) {
+      rv = hipfftExecZ2Z(p, dbuf, dbuf, HIPFFT_FORWARD);
+      hipDeviceSynchronize();
+      printf("  hipfftPlan1d     execFwd: %d (%s)\n", (int)rv, hipfftResultString(rv));
+      hipfftDestroy(p);
+    }
+  }
+
+  hipFree(dbuf);
+  printf("\n");
+}
+
+int main(void) {
+  // Print HIP device info
+  int device = 0;
+  hipGetDevice(&device);
+  hipDeviceProp_t prop;
+  hipGetDeviceProperties(&prop, device);
+  printf("Device %d: %s  warpSize=%d\n\n", device, prop.name, prop.warpSize);
+
+#ifdef hipfftVersionMinor
+  printf("hipFFT version: %d.%d.%d\n\n",
+         hipfftVersionMajor, hipfftVersionMinor, hipfftVersionPatch);
+#endif
+
+  // Original sweep with small howmany (these passed first time)
+  printf("=== Small howmany (original sweep) ===\n\n");
+  for (int G : {4, 8, 12, 16, 24, 32, 48, 64})
+    tryPlanAndExec(G, 512);
+
+  // Grid-realistic howmany values derived from actual lattice geometries.
+  // howmany = Ncomp * product(ldimensions[d] for d != dim)
+  // For LatticeComplexD: Ncomp=1.
+  printf("=== Grid-realistic parameters ===\n\n");
+
+  // --grid 16.16.16.16  4D FFT  (KNOWN TO FAIL in Grid)
+  // Each dim: G=16, Nperp=16^3=4096
+  tryPlanAndExec(16, 4096);
+
+  // --grid 32.32.32.32  4D FFT  (KNOWN TO SUCCEED in Grid)
+  // Each dim: G=32, Nperp=32^3=32768
+  tryPlanAndExec(32, 32768);
+
+  // --grid 32.32.32.32 Ls=8  5D DWF FFT  (KNOWN TO FAIL on dim 0 in Grid)
+  // dim 0: G=8,  Nperp=32^4=1048576
+  tryPlanAndExec(8, 1048576);
+  // dim 1-4: G=32, Nperp=8*32^3=262144
+  tryPlanAndExec(32, 262144);
+
+  // Extra intermediate cases to bracket the failure
+  tryPlanAndExec(16, 1024);
+  tryPlanAndExec(16, 2048);
+  tryPlanAndExec(16, 8192);
+  tryPlanAndExec(8,  4096);
+  tryPlanAndExec(8,  65536);
+  tryPlanAndExec(8,  262144);
+
+  return 0;
+}