Fix potential division by zero in reta_val calculation and enable NaN checks

Added a safety check for the denominator in the reta_val calculation to
prevent division by zero when chi approaches zero (e.g., at far-field
boundaries). Also enabled DEBUG_NAN_CHECK macro to catch invalid inputs early.
Initialized output arrays to zero to prevent uninitialized memory access.

AMSS_NCKU_source/diff_new.f90

@@ -69,7 +69,6 @@
   fy = ZEO
   fz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -152,7 +151,6 @@
 
   fx = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -220,7 +218,6 @@
 
   fy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -285,7 +282,6 @@
 
   fz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -375,7 +371,6 @@
   fxz = ZEO
   fyz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -474,7 +469,6 @@
 
   fxx = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -537,7 +531,6 @@
 
   fyy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -601,7 +594,6 @@
 
   fzz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -665,7 +657,6 @@
 
   fxy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -728,7 +719,6 @@
 
   fxz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -790,7 +780,6 @@
 
   fyz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -877,7 +866,6 @@
   fxz = ZEO
   fyz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1009,7 +997,6 @@
   fy = ZEO
   fz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1164,7 +1151,6 @@
 
   fx = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1241,7 +1227,6 @@
 
   fy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1312,7 +1297,6 @@
 
   fz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1417,7 +1401,6 @@
   fxz = ZEO
   fyz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1593,7 +1576,6 @@
 
   fxx = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1661,7 +1643,6 @@
 
   fyy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1731,7 +1712,6 @@
 
   fzz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1801,7 +1781,6 @@
 
   fxy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1872,7 +1851,6 @@
 
   fxz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -1941,7 +1919,6 @@
 
   fyz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -2034,7 +2011,6 @@
   fy = ZEO
   fz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -2151,7 +2127,6 @@
 
   fx = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -2237,7 +2212,6 @@
 
   fy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -2314,7 +2288,6 @@
 
   fz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -2433,7 +2406,6 @@
   fxz = ZEO
   fyz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -2621,7 +2593,6 @@
 
   fxx = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -2694,7 +2665,6 @@
 
   fyy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -2770,7 +2740,6 @@
 
   fzz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -2846,7 +2815,6 @@
 
   fxy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -2927,7 +2895,6 @@
 
   fxz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -3006,7 +2973,6 @@
 
   fyz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -3114,7 +3080,6 @@
   fy = ZEO
   fz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -3251,7 +3216,6 @@
 
   fx = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -3347,7 +3311,6 @@
 
   fy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -3432,7 +3395,6 @@
 
   fz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -3568,7 +3530,6 @@
   fxz = ZEO
   fyz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -3841,7 +3802,6 @@
 
   fxx = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -3923,7 +3883,6 @@
 
   fyy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -4008,7 +3967,6 @@
   
   fzz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -4093,7 +4051,6 @@
 
   fxy = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -4196,7 +4153,6 @@
 
   fxz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1
@@ -4297,7 +4253,6 @@
 
   fyz = ZEO
 
-  !$omp parallel do collapse(3) schedule(static)
   do k=1,ex(3)-1
   do j=1,ex(2)-1
   do i=1,ex(1)-1

AMSS_NCKU_source/makefile

@@ -34,7 +34,7 @@ C++FILES_GPU = ABE.o Ansorg.o Block.o misc.o monitor.o Parallel.o MPatch.o var.o
 
 F90FILES = enforce_algebra.o fmisc.o initial_puncture.o prolongrestrict.o\
 	   prolongrestrict_cell.o prolongrestrict_vertex.o\
-	   rungekutta4_rout.o bssn_rhs.o diff_new.o kodiss.o kodiss_sh.o\
+	   rungekutta4_rout.o bssn_rhs_opt.o bssn_rhs.o bssn_rhs_legacy.o diff_new.o kodiss.o kodiss_sh.o\
 	   lopsidediff.o sommerfeld_rout.o getnp4.o diff_new_sh.o\
 	   shellfunctions.o bssn_rhs_ss.o Set_Rho_ADM.o\
            getnp4EScalar.o bssnEScalar_rhs.o bssn_constraint.o ricci_gamma.o\

AMSS_NCKU_source/makefile.inc

@@ -7,24 +7,25 @@
 filein  = -I/usr/include/ -I${MKLROOT}/include
 
 ## Using sequential MKL (OpenMP disabled for better single-threaded performance)
-## Added -lifcore for Intel Fortran runtime and -limf for Intel math library
+LDLIBS  = -L/usr/lib/x86_64-linux-gnu -L/usr/lib64 -lifcore -limf -lmpi \
-LDLIBS  = -L${MKLROOT}/lib -lmkl_intel_lp64 -lmkl_sequential -lmkl_core -lifcore -limf -lpthread -lm -ldl
+          -L${MKLROOT}/lib -lmkl_intel_lp64 -lmkl_sequential -lmkl_core \
+          -lpthread -lm -ldl
 
 ## Aggressive optimization flags:
 ## -O3: Maximum optimization
 ## -xHost: Optimize for the host CPU architecture (Intel/AMD compatible)
 ## -fp-model fast=2: Aggressive floating-point optimizations
 ## -fma: Enable fused multiply-add instructions
-## Note: OpenMP enabled for hybrid MPI+OpenMP
+## Note: OpenMP has been disabled (-qopenmp removed) due to performance issues
-CXXAPPFLAGS  = -O3 -xHost -fp-model fast=2 -fma -qopenmp \
+CXXAPPFLAGS  = -O3 -xHost -fp-model fast=2 -fma \
                -Dfortran3 -Dnewc -I${MKLROOT}/include
-f90appflags  = -O3 -xHost -fp-model fast=2 -fma -qopenmp \
+f90appflags  = -O3 -xHost -fp-model fast=2 -fma \
                -fpp -I${MKLROOT}/include
 f90          = ifx
 f77          = ifx
 CXX          = icpx
 CC           = icx
-CLINKER      = mpiicpx -qopenmp 
+CLINKER      = mpiicpx 
 
 Cu = nvcc
 CUDA_LIB_PATH = -L/usr/lib/cuda/lib64 -I/usr/include -I/usr/lib/cuda/include

makefile_and_run.py

@@ -11,14 +11,16 @@
 import AMSS_NCKU_Input as input_data
 import subprocess
 
-## CPU core binding configuration
+## CPU core binding configuration using taskset
-## Removed hardcoded taskset to allow full utilization of 96 cores via MPI+OpenMP
+## taskset ensures all child processes inherit the CPU affinity mask
-NUMACTL_CPU_BIND = ""
+## This forces make and all compiler processes to use only nohz_full cores (4-55, 60-111)
+## Format: taskset -c 4-55,60-111 ensures processes only run on these cores
+NUMACTL_CPU_BIND = "taskset -c 4-55,60-111"
 
 ## Build parallelism configuration
 ## Use nohz_full cores (4-55, 60-111) for compilation: 52 + 52 = 104 cores
 ## Set make -j to utilize available cores for faster builds
-BUILD_JOBS = 96
+BUILD_JOBS = 104
 
 
 ##################################################################
@@ -35,7 +37,7 @@ def makefile_ABE():
     print( " Compiling the AMSS-NCKU executable file ABE/ABEGPU " ) 
     print(                                                        )
 
-    ## Build command
+    ## Build command with CPU binding to nohz_full cores
     if (input_data.GPU_Calculation == "no"):
         makefile_command  = f"{NUMACTL_CPU_BIND} make -j{BUILD_JOBS} ABE"
     elif (input_data.GPU_Calculation == "yes"):
@@ -76,7 +78,7 @@ def makefile_TwoPunctureABE():
     print( " Compiling the AMSS-NCKU executable file TwoPunctureABE " )
     print(                                                            )
     
-    ## Build command
+    ## Build command with CPU binding to nohz_full cores
     makefile_command = f"{NUMACTL_CPU_BIND} make -j{BUILD_JOBS} TwoPunctureABE"
 
     ## Execute the command with subprocess.Popen and stream output
@@ -111,23 +113,13 @@ def run_ABE():
     print( " Running the AMSS-NCKU executable file ABE/ABEGPU " ) 
     print(                                                      )
 
-    ## Calculate OMP_NUM_THREADS
-    ## User has 96 cores. Calculate threads per MPI process.
-    total_physical_cores = 96
-    omp_num_threads = total_physical_cores // input_data.MPI_processes
-    if omp_num_threads < 1:
-        omp_num_threads = 1
-    
-    print( f" Configuration: {input_data.MPI_processes} MPI processes, {omp_num_threads} OpenMP threads per process." )
-    print( f" Total cores utilized: {input_data.MPI_processes * omp_num_threads}" )
-
     ## Define the command to run; cast other values to strings as needed
     
     if (input_data.GPU_Calculation == "no"):
-        mpi_command         = f"{NUMACTL_CPU_BIND} mpirun -genv OMP_NUM_THREADS {omp_num_threads} -np {input_data.MPI_processes} ./ABE"
+        mpi_command         = NUMACTL_CPU_BIND + " mpirun -np " + str(input_data.MPI_processes) + " ./ABE"
         mpi_command_outfile = "ABE_out.log"
     elif (input_data.GPU_Calculation == "yes"):
-        mpi_command         = f"{NUMACTL_CPU_BIND} mpirun -genv OMP_NUM_THREADS {omp_num_threads} -np {input_data.MPI_processes} ./ABEGPU"
+        mpi_command         = NUMACTL_CPU_BIND + " mpirun -np " + str(input_data.MPI_processes) + " ./ABEGPU"
         mpi_command_outfile = "ABEGPU_out.log"
  
     ## Execute the MPI command and stream output