Parallelize AEILocalInterp with OpenMP

This leads to a slight change in behaviour. Currently, AEILocalInterp
traverses the list of points sequentially, and aborts when the first
error is encountered. After parallelisation, there is no fixed order
in which the points are traversed, and if several errors are
encountered, any one of the errors may be returned, not necessarily
the first. I am not aware of any thorn that would or should rely on
such an ordering.

This patch also adds "restrict" and "const" statements that may
improve performance as it gives the compiler more information about
dependencies between pointers.


git-svn-id: http://svn.aei.mpg.de/numrel/AEIThorns/AEILocalInterp/trunk@55 0f49ee68-0e4f-0410-9b9c-b2c123ded7ef

6 files changed, 156 insertions, 170 deletions

diff --git a/src/common/evaluate.c b/src/common/evaluate.c
index 60156b3..3c814e4 100644
--- a/src/common/evaluate.c
+++ b/src/common/evaluate.c
@@ -13,43 +13,43 @@
  * 1-D routines
  */
 
-fp AEILocalInterp_eval_1dcube2(const struct coeffs_struct_1d_cube_size2 *coeffs,
-			       const struct data_struct_1d_cube_size2 *data)
+fp AEILocalInterp_eval_1dcube2(const struct coeffs_struct_1d_cube_size2 *restrict const coeffs,
+			       const struct data_struct_1d_cube_size2 *restrict const data)
 {
 return
   #include "1d.cube.size2/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_1dcube3(const struct coeffs_struct_1d_cube_size3 *coeffs,
-			       const struct data_struct_1d_cube_size3 *data)
+fp AEILocalInterp_eval_1dcube3(const struct coeffs_struct_1d_cube_size3 *restrict const coeffs,
+			       const struct data_struct_1d_cube_size3 *restrict const data)
 {
 return
   #include "1d.cube.size3/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_1dcube4(const struct coeffs_struct_1d_cube_size4 *coeffs,
-			       const struct data_struct_1d_cube_size4 *data)
+fp AEILocalInterp_eval_1dcube4(const struct coeffs_struct_1d_cube_size4 *restrict const coeffs,
+			       const struct data_struct_1d_cube_size4 *restrict const data)
 {
 return
   #include "1d.cube.size4/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_1dcube5(const struct coeffs_struct_1d_cube_size5 *coeffs,
-			       const struct data_struct_1d_cube_size5 *data)
+fp AEILocalInterp_eval_1dcube5(const struct coeffs_struct_1d_cube_size5 *restrict const coeffs,
+			       const struct data_struct_1d_cube_size5 *restrict const data)
 {
 return
   #include "1d.cube.size5/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_1dcube6(const struct coeffs_struct_1d_cube_size6 *coeffs,
-			       const struct data_struct_1d_cube_size6 *data)
+fp AEILocalInterp_eval_1dcube6(const struct coeffs_struct_1d_cube_size6 *restrict const coeffs,
+			       const struct data_struct_1d_cube_size6 *restrict const data)
 {
 return
   #include "1d.cube.size6/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_1dcube7(const struct coeffs_struct_1d_cube_size7 *coeffs,
-			       const struct data_struct_1d_cube_size7 *data)
+fp AEILocalInterp_eval_1dcube7(const struct coeffs_struct_1d_cube_size7 *restrict const coeffs,
+			       const struct data_struct_1d_cube_size7 *restrict const data)
 {
 return
   #include "1d.cube.size7/evaluate-molecule.c"
@@ -61,36 +61,36 @@ return
  * 2-D routines
  */
 
-fp AEILocalInterp_eval_2dcube2(const struct coeffs_struct_2d_cube_size2 *coeffs,
-			       const struct data_struct_2d_cube_size2 *data)
+fp AEILocalInterp_eval_2dcube2(const struct coeffs_struct_2d_cube_size2 *restrict const coeffs,
+			       const struct data_struct_2d_cube_size2 *restrict const data)
 {
 return
   #include "2d.cube.size2/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_2dcube3(const struct coeffs_struct_2d_cube_size3 *coeffs,
-			       const struct data_struct_2d_cube_size3 *data)
+fp AEILocalInterp_eval_2dcube3(const struct coeffs_struct_2d_cube_size3 *restrict const coeffs,
+			       const struct data_struct_2d_cube_size3 *restrict const data)
 {
 return
   #include "2d.cube.size3/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_2dcube4(const struct coeffs_struct_2d_cube_size4 *coeffs,
-			       const struct data_struct_2d_cube_size4 *data)
+fp AEILocalInterp_eval_2dcube4(const struct coeffs_struct_2d_cube_size4 *restrict const coeffs,
+			       const struct data_struct_2d_cube_size4 *restrict const data)
 {
 return
   #include "2d.cube.size4/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_2dcube5(const struct coeffs_struct_2d_cube_size5 *coeffs,
-			       const struct data_struct_2d_cube_size5 *data)
+fp AEILocalInterp_eval_2dcube5(const struct coeffs_struct_2d_cube_size5 *restrict const coeffs,
+			       const struct data_struct_2d_cube_size5 *restrict const data)
 {
 return
   #include "2d.cube.size5/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_2dcube6(const struct coeffs_struct_2d_cube_size6 *coeffs,
-			       const struct data_struct_2d_cube_size6 *data)
+fp AEILocalInterp_eval_2dcube6(const struct coeffs_struct_2d_cube_size6 *restrict const coeffs,
+			       const struct data_struct_2d_cube_size6 *restrict const data)
 {
 return
   #include "2d.cube.size6/evaluate-molecule.c"
@@ -102,36 +102,36 @@ return
  * 3-D routines
  */
 
-fp AEILocalInterp_eval_3dcube2(const struct coeffs_struct_3d_cube_size2 *coeffs,
-			       const struct data_struct_3d_cube_size2 *data)
+fp AEILocalInterp_eval_3dcube2(const struct coeffs_struct_3d_cube_size2 *restrict const coeffs,
+			       const struct data_struct_3d_cube_size2 *restrict const data)
 {
 return
   #include "3d.cube.size2/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_3dcube3(const struct coeffs_struct_3d_cube_size3 *coeffs,
-			       const struct data_struct_3d_cube_size3 *data)
+fp AEILocalInterp_eval_3dcube3(const struct coeffs_struct_3d_cube_size3 *restrict const coeffs,
+			       const struct data_struct_3d_cube_size3 *restrict const data)
 {
 return
   #include "3d.cube.size3/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_3dcube4(const struct coeffs_struct_3d_cube_size4 *coeffs,
-			       const struct data_struct_3d_cube_size4 *data)
+fp AEILocalInterp_eval_3dcube4(const struct coeffs_struct_3d_cube_size4 *restrict const coeffs,
+			       const struct data_struct_3d_cube_size4 *restrict const data)
 {
 return
   #include "3d.cube.size4/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_3dcube5(const struct coeffs_struct_3d_cube_size5 *coeffs,
-			       const struct data_struct_3d_cube_size5 *data)
+fp AEILocalInterp_eval_3dcube5(const struct coeffs_struct_3d_cube_size5 *restrict const coeffs,
+			       const struct data_struct_3d_cube_size5 *restrict const data)
 {
 return
   #include "3d.cube.size5/evaluate-molecule.c"
 }
 
-fp AEILocalInterp_eval_3dcube6(const struct coeffs_struct_3d_cube_size6 *coeffs,
-			       const struct data_struct_3d_cube_size6 *data)
+fp AEILocalInterp_eval_3dcube6(const struct coeffs_struct_3d_cube_size6 *restrict const coeffs,
+			       const struct data_struct_3d_cube_size6 *restrict const data)
 {
 return
   #include "3d.cube.size6/evaluate-molecule.c"
diff --git a/src/common/evaluate.h b/src/common/evaluate.h
index 4e2046a..a7c6e53 100644
--- a/src/common/evaluate.h
+++ b/src/common/evaluate.h
@@ -8,37 +8,37 @@
  *	"structs.h"
  */
 
-fp AEILocalInterp_eval_1dcube2(const struct coeffs_struct_1d_cube_size2 *coeffs,
-			       const struct data_struct_1d_cube_size2 *data);
-fp AEILocalInterp_eval_1dcube3(const struct coeffs_struct_1d_cube_size3 *coeffs,
-			       const struct data_struct_1d_cube_size3 *data);
-fp AEILocalInterp_eval_1dcube4(const struct coeffs_struct_1d_cube_size4 *coeffs,
-			       const struct data_struct_1d_cube_size4 *data);
-fp AEILocalInterp_eval_1dcube5(const struct coeffs_struct_1d_cube_size5 *coeffs,
-			       const struct data_struct_1d_cube_size5 *data);
-fp AEILocalInterp_eval_1dcube6(const struct coeffs_struct_1d_cube_size6 *coeffs,
-			       const struct data_struct_1d_cube_size6 *data);
-fp AEILocalInterp_eval_1dcube7(const struct coeffs_struct_1d_cube_size7 *coeffs,
-			       const struct data_struct_1d_cube_size7 *data);
+fp AEILocalInterp_eval_1dcube2(const struct coeffs_struct_1d_cube_size2 *restrict coeffs,
+			       const struct data_struct_1d_cube_size2 *restrict data);
+fp AEILocalInterp_eval_1dcube3(const struct coeffs_struct_1d_cube_size3 *restrict coeffs,
+			       const struct data_struct_1d_cube_size3 *restrict data);
+fp AEILocalInterp_eval_1dcube4(const struct coeffs_struct_1d_cube_size4 *restrict coeffs,
+			       const struct data_struct_1d_cube_size4 *restrict data);
+fp AEILocalInterp_eval_1dcube5(const struct coeffs_struct_1d_cube_size5 *restrict coeffs,
+			       const struct data_struct_1d_cube_size5 *restrict data);
+fp AEILocalInterp_eval_1dcube6(const struct coeffs_struct_1d_cube_size6 *restrict coeffs,
+			       const struct data_struct_1d_cube_size6 *restrict data);
+fp AEILocalInterp_eval_1dcube7(const struct coeffs_struct_1d_cube_size7 *restrict coeffs,
+			       const struct data_struct_1d_cube_size7 *restrict data);
 
-fp AEILocalInterp_eval_2dcube2(const struct coeffs_struct_2d_cube_size2 *coeffs,
-			       const struct data_struct_2d_cube_size2 *data);
-fp AEILocalInterp_eval_2dcube3(const struct coeffs_struct_2d_cube_size3 *coeffs,
-			       const struct data_struct_2d_cube_size3 *data);
-fp AEILocalInterp_eval_2dcube4(const struct coeffs_struct_2d_cube_size4 *coeffs,
-			       const struct data_struct_2d_cube_size4 *data);
-fp AEILocalInterp_eval_2dcube5(const struct coeffs_struct_2d_cube_size5 *coeffs,
-			       const struct data_struct_2d_cube_size5 *data);
-fp AEILocalInterp_eval_2dcube6(const struct coeffs_struct_2d_cube_size6 *coeffs,
-			       const struct data_struct_2d_cube_size6 *data);
+fp AEILocalInterp_eval_2dcube2(const struct coeffs_struct_2d_cube_size2 *restrict coeffs,
+			       const struct data_struct_2d_cube_size2 *restrict data);
+fp AEILocalInterp_eval_2dcube3(const struct coeffs_struct_2d_cube_size3 *restrict coeffs,
+			       const struct data_struct_2d_cube_size3 *restrict data);
+fp AEILocalInterp_eval_2dcube4(const struct coeffs_struct_2d_cube_size4 *restrict coeffs,
+			       const struct data_struct_2d_cube_size4 *restrict data);
+fp AEILocalInterp_eval_2dcube5(const struct coeffs_struct_2d_cube_size5 *restrict coeffs,
+			       const struct data_struct_2d_cube_size5 *restrict data);
+fp AEILocalInterp_eval_2dcube6(const struct coeffs_struct_2d_cube_size6 *restrict coeffs,
+			       const struct data_struct_2d_cube_size6 *restrict data);
 
-fp AEILocalInterp_eval_3dcube2(const struct coeffs_struct_3d_cube_size2 *coeffs,
-			       const struct data_struct_3d_cube_size2 *data);
-fp AEILocalInterp_eval_3dcube3(const struct coeffs_struct_3d_cube_size3 *coeffs,
-			       const struct data_struct_3d_cube_size3 *data);
-fp AEILocalInterp_eval_3dcube4(const struct coeffs_struct_3d_cube_size4 *coeffs,
-			       const struct data_struct_3d_cube_size4 *data);
-fp AEILocalInterp_eval_3dcube5(const struct coeffs_struct_3d_cube_size5 *coeffs,
-			       const struct data_struct_3d_cube_size5 *data);
-fp AEILocalInterp_eval_3dcube6(const struct coeffs_struct_3d_cube_size6 *coeffs,
-			       const struct data_struct_3d_cube_size6 *data);
+fp AEILocalInterp_eval_3dcube2(const struct coeffs_struct_3d_cube_size2 *restrict coeffs,
+			       const struct data_struct_3d_cube_size2 *restrict data);
+fp AEILocalInterp_eval_3dcube3(const struct coeffs_struct_3d_cube_size3 *restrict coeffs,
+			       const struct data_struct_3d_cube_size3 *restrict data);
+fp AEILocalInterp_eval_3dcube4(const struct coeffs_struct_3d_cube_size4 *restrict coeffs,
+			       const struct data_struct_3d_cube_size4 *restrict data);
+fp AEILocalInterp_eval_3dcube5(const struct coeffs_struct_3d_cube_size5 *restrict coeffs,
+			       const struct data_struct_3d_cube_size5 *restrict data);
+fp AEILocalInterp_eval_3dcube6(const struct coeffs_struct_3d_cube_size6 *restrict coeffs,
+			       const struct data_struct_3d_cube_size6 *restrict data);
diff --git a/src/common/load-template.c b/src/common/load-template.c
index 5caa23b..5e197f7 100644
--- a/src/common/load-template.c
+++ b/src/common/load-template.c
@@ -16,35 +16,28 @@
 #undef  DATA
 #define DATA	DATA_REAL
 
-void LOAD_FUNCTION_NAME(r)(const CCTK_REAL *ptr,
-			   INT_STRIDE_IJK,
-			  struct DATA_STRUCT *data)
-{
-#include LOAD_DATA_FILE_NAME
-}
-
 #ifdef HAVE_CCTK_REAL4
-  void LOAD_FUNCTION_NAME(r4)(const CCTK_REAL4 *ptr,
+  void LOAD_FUNCTION_NAME(r4)(const CCTK_REAL4 *restrict const ptr,
 			      INT_STRIDE_IJK,
-			      struct DATA_STRUCT *data)
+			      struct DATA_STRUCT *restrict const data)
   {
   #include LOAD_DATA_FILE_NAME
   }
 #endif
 
 #ifdef HAVE_CCTK_REAL8
-  void LOAD_FUNCTION_NAME(r8)(const CCTK_REAL8 *ptr,
+  void LOAD_FUNCTION_NAME(r8)(const CCTK_REAL8 *restrict const ptr,
 			      INT_STRIDE_IJK,
-			      struct DATA_STRUCT *data)
+			      struct DATA_STRUCT *restrict const data)
   {
   #include LOAD_DATA_FILE_NAME
   }
 #endif
 
 #ifdef HAVE_CCTK_REAL16
-  void LOAD_FUNCTION_NAME(r16)(const CCTK_REAL16 *ptr,
+  void LOAD_FUNCTION_NAME(r16)(const CCTK_REAL16 *restrict const ptr,
 				INT_STRIDE_IJK,
-				struct DATA_STRUCT *data)
+				struct DATA_STRUCT *restrict const data)
   {
   #include LOAD_DATA_FILE_NAME
   }
@@ -59,35 +52,28 @@ void LOAD_FUNCTION_NAME(r)(const CCTK_REAL *ptr,
 #undef  DATA
 #define DATA	DATA_COMPLEX
 
-void LOAD_FUNCTION_NAME(c)(const CCTK_REAL (*ptr)[COMPLEX_N_PARTS],
-			   INT_STRIDE_IJK, int part,
-			   struct DATA_STRUCT *data)
-{
-#include LOAD_DATA_FILE_NAME
-}
-
 #ifdef HAVE_CCTK_COMPLEX8
-  void LOAD_FUNCTION_NAME(c8)(const CCTK_REAL4 (*ptr)[COMPLEX_N_PARTS],
+  void LOAD_FUNCTION_NAME(c8)(const CCTK_REAL4 (*restrict const ptr)[COMPLEX_N_PARTS],
 			      INT_STRIDE_IJK, int part,
-			      struct DATA_STRUCT *data)
+			      struct DATA_STRUCT *restrict const data)
   {
   #include LOAD_DATA_FILE_NAME
   }
 #endif
 
 #ifdef HAVE_CCTK_COMPLEX16
-  void LOAD_FUNCTION_NAME(c16)(const CCTK_REAL8 (*ptr)[COMPLEX_N_PARTS],
+  void LOAD_FUNCTION_NAME(c16)(const CCTK_REAL8 (*restrict const ptr)[COMPLEX_N_PARTS],
 				INT_STRIDE_IJK, int part,
-				struct DATA_STRUCT *data)
+				struct DATA_STRUCT *restrict const data)
   {
   #include LOAD_DATA_FILE_NAME
   }
 #endif
 
 #ifdef HAVE_CCTK_COMPLEX32
-  void LOAD_FUNCTION_NAME(c32)(const CCTK_REAL16 (*ptr)[COMPLEX_N_PARTS],
+  void LOAD_FUNCTION_NAME(c32)(const CCTK_REAL16 (*restrict const ptr)[COMPLEX_N_PARTS],
 				INT_STRIDE_IJK, int part,
-				struct DATA_STRUCT *data)
+				struct DATA_STRUCT *restrict const data)
   {
   #include LOAD_DATA_FILE_NAME
   }
diff --git a/src/common/load-template.h b/src/common/load-template.h
index 08cdfd4..f4b7bf3 100644
--- a/src/common/load-template.h
+++ b/src/common/load-template.h
@@ -18,26 +18,26 @@
  * load-data routines for real datatypes
  */
 
-void LOAD_FUNCTION_NAME(r)(const CCTK_REAL *ptr,
+void LOAD_FUNCTION_NAME(r)(const CCTK_REAL *restrict ptr,
 			   INT_STRIDE_IJK,
-			   struct DATA_STRUCT *data);
+			   struct DATA_STRUCT *restrict data);
 
 #ifdef HAVE_CCTK_REAL4
-  void LOAD_FUNCTION_NAME(r4)(const CCTK_REAL4 *ptr,
+  void LOAD_FUNCTION_NAME(r4)(const CCTK_REAL4 *restrict ptr,
 			      INT_STRIDE_IJK,
-			      struct DATA_STRUCT *data);
+			      struct DATA_STRUCT *restrict data);
 #endif
 
 #ifdef HAVE_CCTK_REAL8
-  void LOAD_FUNCTION_NAME(r8)(const CCTK_REAL8 *ptr,
+  void LOAD_FUNCTION_NAME(r8)(const CCTK_REAL8 *restrict ptr,
 			      INT_STRIDE_IJK,
-			      struct DATA_STRUCT *data);
+			      struct DATA_STRUCT *restrict data);
 #endif
 
 #ifdef HAVE_CCTK_REAL16
-  void LOAD_FUNCTION_NAME(r16)(const CCTK_REAL16 *ptr,
+  void LOAD_FUNCTION_NAME(r16)(const CCTK_REAL16 *restrict ptr,
 				INT_STRIDE_IJK,
-				struct DATA_STRUCT *data);
+				struct DATA_STRUCT *restrict data);
 #endif
 
 /******************************************************************************/
@@ -46,24 +46,24 @@ void LOAD_FUNCTION_NAME(r)(const CCTK_REAL *ptr,
  * load-data routines for complex datatypes
  */
 
-void LOAD_FUNCTION_NAME(c)(const CCTK_REAL (*ptr)[COMPLEX_N_PARTS],
+void LOAD_FUNCTION_NAME(c)(const CCTK_REAL (*restrict ptr)[COMPLEX_N_PARTS],
 			   INT_STRIDE_IJK, int part,
-			   struct DATA_STRUCT *data);
+			   struct DATA_STRUCT *restrict data);
 
 #ifdef HAVE_CCTK_COMPLEX8
-  void LOAD_FUNCTION_NAME(c8)(const CCTK_REAL4 (*ptr)[COMPLEX_N_PARTS],
+  void LOAD_FUNCTION_NAME(c8)(const CCTK_REAL4 (*restrict ptr)[COMPLEX_N_PARTS],
 			      INT_STRIDE_IJK, int part,
-			      struct DATA_STRUCT *data);
+			      struct DATA_STRUCT *restrict data);
 #endif
 
 #ifdef HAVE_CCTK_COMPLEX16
-  void LOAD_FUNCTION_NAME(c16)(const CCTK_REAL8 (*ptr)[COMPLEX_N_PARTS],
+  void LOAD_FUNCTION_NAME(c16)(const CCTK_REAL8 (*restrict ptr)[COMPLEX_N_PARTS],
 				INT_STRIDE_IJK, int part,
-				struct DATA_STRUCT *data);
+				struct DATA_STRUCT *restrict data);
 #endif
 
 #ifdef HAVE_CCTK_COMPLEX32
-  void LOAD_FUNCTION_NAME(c32)(const CCTK_REAL16 (*ptr)[COMPLEX_N_PARTS],
+  void LOAD_FUNCTION_NAME(c32)(const CCTK_REAL16 (*restrict ptr)[COMPLEX_N_PARTS],
 				INT_STRIDE_IJK, int part,
-				struct DATA_STRUCT *data);
+				struct DATA_STRUCT *restrict data);
 #endif
diff --git a/src/common/store.c b/src/common/store.c
index 3ed6003..1d414fe 100644
--- a/src/common/store.c
+++ b/src/common/store.c
@@ -17,48 +17,48 @@
 #define COEFF(mi)	Jacobian_ptr[Jacobian_mi_stride*mi]
 
 void AEILocalInterp_store_1dcube2
-      (fp factor, const struct coeffs_struct_1d_cube_size2 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size2 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride)
 {
 #include "1d.cube.size2/store-coeffs.c"
 }
 
 void AEILocalInterp_store_1dcube3
-      (fp factor, const struct coeffs_struct_1d_cube_size3 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size3 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride)
 {
 #include "1d.cube.size3/store-coeffs.c"
 }
 
 void AEILocalInterp_store_1dcube4
-      (fp factor, const struct coeffs_struct_1d_cube_size4 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size4 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride)
 {
 #include "1d.cube.size4/store-coeffs.c"
 }
 
 void AEILocalInterp_store_1dcube5
-      (fp factor, const struct coeffs_struct_1d_cube_size5 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size5 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride)
 {
 #include "1d.cube.size5/store-coeffs.c"
 }
 
 void AEILocalInterp_store_1dcube6
-      (fp factor, const struct coeffs_struct_1d_cube_size6 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size6 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride)
 {
 #include "1d.cube.size6/store-coeffs.c"
 }
 
 void AEILocalInterp_store_1dcube7
-      (fp factor, const struct coeffs_struct_1d_cube_size7 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size7 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride)
 {
 #include "1d.cube.size7/store-coeffs.c"
@@ -75,40 +75,40 @@ void AEILocalInterp_store_1dcube7
 				      + Jacobian_mj_stride*mj ]
 
 void AEILocalInterp_store_2dcube2
-      (fp factor, const struct coeffs_struct_2d_cube_size2 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_2d_cube_size2 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride)
 {
 #include "2d.cube.size2/store-coeffs.c"
 }
 
 void AEILocalInterp_store_2dcube3
-      (fp factor, const struct coeffs_struct_2d_cube_size3 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_2d_cube_size3 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride)
 {
 #include "2d.cube.size3/store-coeffs.c"
 }
 
 void AEILocalInterp_store_2dcube4
-      (fp factor, const struct coeffs_struct_2d_cube_size4 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_2d_cube_size4 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride)
 {
 #include "2d.cube.size4/store-coeffs.c"
 }
 
 void AEILocalInterp_store_2dcube5
-      (fp factor, const struct coeffs_struct_2d_cube_size5 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_2d_cube_size5 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride)
 {
 #include "2d.cube.size5/store-coeffs.c"
 }
 
 void AEILocalInterp_store_2dcube6
-      (fp factor, const struct coeffs_struct_2d_cube_size6 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_2d_cube_size6 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride)
 {
 #include "2d.cube.size6/store-coeffs.c"
@@ -126,40 +126,40 @@ void AEILocalInterp_store_2dcube6
 				      + Jacobian_mk_stride*mk ]
 
 void AEILocalInterp_store_3dcube2
-      (fp factor, const struct coeffs_struct_3d_cube_size2 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_3d_cube_size2 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride, int Jacobian_mk_stride)
 {
 #include "3d.cube.size2/store-coeffs.c"
 }
 
 void AEILocalInterp_store_3dcube3
-      (fp factor, const struct coeffs_struct_3d_cube_size3 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_3d_cube_size3 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride, int Jacobian_mk_stride)
 {
 #include "3d.cube.size3/store-coeffs.c"
 }
 
 void AEILocalInterp_store_3dcube4
-      (fp factor, const struct coeffs_struct_3d_cube_size4 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_3d_cube_size4 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride, int Jacobian_mk_stride)
 {
 #include "3d.cube.size4/store-coeffs.c"
 }
 
 void AEILocalInterp_store_3dcube5
-      (fp factor, const struct coeffs_struct_3d_cube_size5 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_3d_cube_size5 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride, int Jacobian_mk_stride)
 {
 #include "3d.cube.size5/store-coeffs.c"
 }
 
 void AEILocalInterp_store_3dcube6
-      (fp factor, const struct coeffs_struct_3d_cube_size6 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_3d_cube_size6 *restrict const coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride, int Jacobian_mk_stride)
 {
 #include "3d.cube.size6/store-coeffs.c"
diff --git a/src/common/store.h b/src/common/store.h
index e042552..678d663 100644
--- a/src/common/store.h
+++ b/src/common/store.h
@@ -9,68 +9,68 @@
  */
 
 void AEILocalInterp_store_1dcube2
-      (fp factor, const struct coeffs_struct_1d_cube_size2 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size2 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride);
 void AEILocalInterp_store_1dcube3
-      (fp factor, const struct coeffs_struct_1d_cube_size3 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size3 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride);
 void AEILocalInterp_store_1dcube4
-      (fp factor, const struct coeffs_struct_1d_cube_size4 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size4 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride);
 void AEILocalInterp_store_1dcube5
-      (fp factor, const struct coeffs_struct_1d_cube_size5 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size5 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride);
 void AEILocalInterp_store_1dcube6
-      (fp factor, const struct coeffs_struct_1d_cube_size6 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size6 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride);
 void AEILocalInterp_store_1dcube7
-      (fp factor, const struct coeffs_struct_1d_cube_size7 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_1d_cube_size7 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride);
 
 void AEILocalInterp_store_2dcube2
-      (fp factor, const struct coeffs_struct_2d_cube_size2 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_2d_cube_size2 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride);
 void AEILocalInterp_store_2dcube3
-      (fp factor, const struct coeffs_struct_2d_cube_size3 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_2d_cube_size3 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride);
 void AEILocalInterp_store_2dcube4
-      (fp factor, const struct coeffs_struct_2d_cube_size4 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_2d_cube_size4 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride);
 void AEILocalInterp_store_2dcube5
-      (fp factor, const struct coeffs_struct_2d_cube_size5 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_2d_cube_size5 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride);
 void AEILocalInterp_store_2dcube6
-      (fp factor, const struct coeffs_struct_2d_cube_size6 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_2d_cube_size6 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride);
 
 void AEILocalInterp_store_3dcube2
-      (fp factor, const struct coeffs_struct_3d_cube_size2 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_3d_cube_size2 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride, int Jacobian_mk_stride);
 void AEILocalInterp_store_3dcube3
-      (fp factor, const struct coeffs_struct_3d_cube_size3 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_3d_cube_size3 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride, int Jacobian_mk_stride);
 void AEILocalInterp_store_3dcube4
-      (fp factor, const struct coeffs_struct_3d_cube_size4 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_3d_cube_size4 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride, int Jacobian_mk_stride);
 void AEILocalInterp_store_3dcube5
-      (fp factor, const struct coeffs_struct_3d_cube_size5 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_3d_cube_size5 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride, int Jacobian_mk_stride);
 void AEILocalInterp_store_3dcube6
-      (fp factor, const struct coeffs_struct_3d_cube_size6 *coeffs,
-       fp Jacobian_ptr[],
+      (fp factor, const struct coeffs_struct_3d_cube_size6 *restrict coeffs,
+       fp *restrict Jacobian_ptr,
        int Jacobian_mi_stride, int Jacobian_mj_stride, int Jacobian_mk_stride);