summaryrefslogtreecommitdiff
path: root/libm/ldouble/wronkl.c
diff options
context:
space:
mode:
authorEric Andersen <andersen@codepoet.org>2001-05-10 00:40:28 +0000
committerEric Andersen <andersen@codepoet.org>2001-05-10 00:40:28 +0000
commit1077fa4d772832f77a677ce7fb7c2d513b959e3f (patch)
tree579bee13fb0b58d2800206366ec2caecbb15f3fc /libm/ldouble/wronkl.c
parent22358dd7ce7bb49792204b698f01a6f69b9c8e08 (diff)
uClibc now has a math library. muahahahaha!
-Erik
Diffstat (limited to 'libm/ldouble/wronkl.c')
-rw-r--r--libm/ldouble/wronkl.c67
1 files changed, 67 insertions, 0 deletions
diff --git a/libm/ldouble/wronkl.c b/libm/ldouble/wronkl.c
new file mode 100644
index 000000000..bec958f01
--- /dev/null
+++ b/libm/ldouble/wronkl.c
@@ -0,0 +1,67 @@
+/* Wronksian test for Bessel functions. */
+
+long double jnl (), ynl (), floorl ();
+#define PI 3.14159265358979323846L
+
+long double y, Jn, Jnp1, Jmn, Jmnp1, Yn, Ynp1;
+long double w1, w2, err1, max1, err2, max2;
+void wronk ();
+
+main ()
+{
+ long double x, delta;
+ int n, i, j;
+
+ max1 = 0.0L;
+ max2 = 0.0L;
+ delta = 0.6 / PI;
+ for (n = -30; n <= 30; n++)
+ {
+ x = -30.0;
+ while (x < 30.0)
+ {
+ wronk (n, x);
+ x += delta;
+ }
+ delta += .00123456;
+ }
+}
+
+void
+wronk (n, x)
+ int n;
+ long double x;
+{
+
+ Jnp1 = jnl (n + 1, x);
+ Jmn = jnl (-n, x);
+ Jn = jnl (n, x);
+ Jmnp1 = jnl (-(n + 1), x);
+ /* This should be trivially zero. */
+ err1 = Jnp1 * Jmn + Jn * Jmnp1;
+ if (err1 < 0.0)
+ err1 = -err1;
+ if (err1 > max1)
+ {
+ max1 = err1;
+ printf ("1 %3d %.5Le %.3Le\n", n, x, max1);
+ }
+ if (x < 0.0)
+ {
+ x = -x;
+ Jn = jnl (n, x);
+ Jnp1 = jnl (n + 1, x);
+ }
+ Yn = ynl (n, x);
+ Ynp1 = ynl (n + 1, x);
+ /* The Wronksian. */
+ w1 = Jnp1 * Yn - Jn * Ynp1;
+ /* What the Wronksian should be. */
+ w2 = 2.0 / (PI * x);
+ err2 = w1 - w2;
+ if (err2 > max2)
+ {
+ max2 = err2;
+ printf ("2 %3d %.5Le %.3Le\n", n, x, max2);
+ }
+}