Allow for testing using long double instead of pari

tests/trigtest.c

@@ -1,12 +1,14 @@
 /*
- * program to test the accuracy of trig routines. Requires libpari.
+ * program to test the accuracy of trig routines. Requires libpari or
+ * extended precision.
  * This program is not part of fftw.
  */
 
 #include <math.h>
 #include <stdio.h>
 #include <stdlib.h>
-#include <pari/pari.h>
+
+#define USE_PARI 1
 
 typedef double trigreal;
 #  define COS cos
@@ -57,6 +59,9 @@ trigreal sin2pi(int m, uint n)
      return sin2pi0((trigreal)m, (trigreal)n);
 }
 
+#if USE_PARI
+#include <pari/pari.h>
+
 long prec = 25;
 
 double ck(long m, long n, double (*cf) (int, uint), GEN(*gf) (GEN, long))
@@ -77,13 +82,28 @@ double ck(long m, long n, double (*cf) (int, uint), GEN(*gf) (GEN, long))
      avma = ltop;
      return cerr;
 }
+#else
+
+double ck(long m, long n, double (*cf) (int, uint), 
+	  long double(*gf) (long double))
+{
+     long double l2pi = 6.2831853071795864769252867665590057683943388L;
+     return cf(m, n) - gf(l2pi * (long double)m / (long double)n);
+}
+
+#define gsin sinl
+#define gcos cosl
+#endif
 
 int main(int argc, char *argv[])
 {
      long nmin, nmax;
      long n, m;
+
+#if USE_PARI
      pari_init(500000, 2);
      mppi(prec);
+#endif
 
      if (argc > 1)
 	  nmin = atoi(argv[1]);
@@ -98,20 +118,12 @@ int main(int argc, char *argv[])
 	  double maxe = 0.0, nmaxe = 0.0;;
 	  for (m = 0; m < n; ++m) {
 	       double e;
-	       e = ck(m, n, sin2pi, gsin);
-	       if (e > maxe)
-		    maxe = e;
-	       e = ck(m, n, cos2pi, gcos);
-	       if (e > maxe)
-		    maxe = e;
-	       e = ck(m, n, naive_sin2pi, gsin);
-	       if (e > nmaxe)
-		    nmaxe = e;
-	       e = ck(m, n, naive_cos2pi, gcos);
-	       if (e > nmaxe)
-		    nmaxe = e;
-
+	       e = ck(m, n, sin2pi, gsin); if (e > maxe) maxe = e;
+	       e = ck(m, n, cos2pi, gcos); if (e > maxe) maxe = e;
+	       e = ck(m, n, naive_sin2pi, gsin); if (e > nmaxe) nmaxe = e;
+	       e = ck(m, n, naive_cos2pi, gcos); if (e > nmaxe) nmaxe = e;
 	  }
 	  printf("%ld %g %g\n", n, maxe, nmaxe);
      }
+     return 0;
 }