--- liboctave/ChangeLog.cvs	Tue Sep  9 11:03:22 2003
+++ liboctave/ChangeLog	Tue Sep  9 11:02:06 2003
@@ -1,3 +1,8 @@
+2003-09-09  David Bateman <address@hidden>
+
+	* lo-specfun.cc (zbesj, zbesy, zbesi, zbesk, zbesh1, zbesh2, 
+	airy, biry): Used scaled versions for overflow/underflow
+
 2003-09-04  John W. Eaton  <address@hidden>
 
 	* lo-specfun.cc (xlgamma): Require nonnegative argument.
--- liboctave/lo-specfun.cc.cvs	Tue Sep  9 10:57:37 2003
+++ liboctave/lo-specfun.cc	Tue Sep  9 11:06:41 2003
@@ -223,7 +223,14 @@
       double zr = z.real ();
       double zi = z.imag ();
 
-      F77_FUNC (zbesj, ZBESJ) (zr, zi, alpha, kode, 1, &yr, &yi, nz, ierr);
+      F77_FUNC (zbesj, ZBESJ) (zr, zi, alpha, 2, 1, &yr, &yi, nz, ierr);
+
+      if (kode != 2)
+	{
+	  double expz = exp(std::abs(zi)); 
+	  yr *= expz;
+	  yi *= expz;
+	}
 
       if (zi == 0.0 && zr >= 0.0)
 	yi = 0.0;
@@ -275,9 +282,16 @@
 	}
       else
 	{
-	  F77_FUNC (zbesy, ZBESY) (zr, zi, alpha, kode, 1, &yr, &yi, nz,
+	  F77_FUNC (zbesy, ZBESY) (zr, zi, alpha, 2, 1, &yr, &yi, nz,
 				  &wr, &wi, ierr);
 
+	  if (kode != 2)
+	    {
+	      double expz = exp(std::abs(zi));
+	      yr *= expz;
+	      yi *= expz;
+	    }
+
 	  if (zi == 0.0 && zr >= 0.0)
 	    yi = 0.0;
 	}
@@ -318,7 +332,14 @@
       double zr = z.real ();
       double zi = z.imag ();
 
-      F77_FUNC (zbesi, ZBESI) (zr, zi, alpha, kode, 1, &yr, &yi, nz, ierr);
+      F77_FUNC (zbesi, ZBESI) (zr, zi, alpha, 2, 1, &yr, &yi, nz, ierr);
+
+      if (kode != 2)
+	{
+	  double expz = exp(std::abs(zr));
+	  yr *= expz;
+	  yi *= expz;
+	}
 
       if (zi == 0.0 && zr >= 0.0)
 	yi = 0.0;
@@ -369,7 +390,17 @@
 	}
       else
 	{
-	  F77_FUNC (zbesk, ZBESK) (zr, zi, alpha, kode, 1, &yr, &yi, nz, ierr);
+	  F77_FUNC (zbesk, ZBESK) (zr, zi, alpha, 2, 1, &yr, &yi, nz, ierr);
+
+	  if (kode != 2)
+	    {
+	      Complex expz = exp(-z);
+	      double rexpz = real(expz);
+	      double iexpz = imag(expz);
+	      double tmp = yr*rexpz - yi*iexpz;
+	      yi = yr*iexpz + yi*rexpz;
+	      yr = tmp;
+	    }
 
 	  if (zi == 0.0 && zr >= 0.0)
 	    yi = 0.0;
@@ -402,7 +433,17 @@
       double zr = z.real ();
       double zi = z.imag ();
 
-      F77_FUNC (zbesh, ZBESH) (zr, zi, alpha, kode, 1, 1, &yr, &yi, nz, ierr);
+      F77_FUNC (zbesh, ZBESH) (zr, zi, alpha, 2, 1, 1, &yr, &yi, nz, ierr);
+
+      if (kode != 2)
+	{
+	  Complex expz = exp( Complex ( 0., 1.) * z);
+	  double rexpz = real(expz);
+	  double iexpz = imag(expz);
+	  double tmp = yr*rexpz - yi*iexpz;
+	  yi = yr*iexpz + yi*rexpz;
+	  yr = tmp;
+	}
 
       retval = bessel_return_value (Complex (yr, yi), ierr);
     }
@@ -435,7 +476,17 @@
       double zr = z.real ();
       double zi = z.imag ();
 
-      F77_FUNC (zbesh, ZBESH) (zr, zi, alpha, kode, 2, 1, &yr, &yi, nz, ierr);
+      F77_FUNC (zbesh, ZBESH) (zr, zi, alpha, 2, 2, 1, &yr, &yi, nz, ierr);
+
+      if (kode != 2)
+	{
+	  Complex expz = exp( - Complex ( 0., 1.) * z);
+	  double rexpz = real(expz);
+	  double iexpz = imag(expz);
+	  double tmp = yr*rexpz - yi*iexpz;
+	  yi = yr*iexpz + yi*rexpz;
+	  yr = tmp;
+	}
 
       retval = bessel_return_value (Complex (yr, yi), ierr);
     }
@@ -620,7 +671,17 @@
 
   int kode = scaled ? 2 : 1;
 
-  F77_FUNC (zairy, ZAIRY) (zr, zi, id, kode, ar, ai, nz, ierr);
+  F77_FUNC (zairy, ZAIRY) (zr, zi, id, 2, ar, ai, nz, ierr);
+
+  if (kode != 2)
+    {
+      Complex expz = exp( - 2. / 3. * z * sqrt(z));
+      double rexpz = real(expz);
+      double iexpz = imag(expz);
+      double tmp = ar*rexpz - ai*iexpz;
+      ai = ar*iexpz + ai*rexpz;
+      ar = tmp;
+    }
 
   if (zi == 0.0 && (! scaled || zr >= 0.0))
     ai = 0.0;
@@ -641,7 +702,17 @@
 
   int kode = scaled ? 2 : 1;
 
-  F77_FUNC (zbiry, ZBIRY) (zr, zi, id, kode, ar, ai, ierr);
+  F77_FUNC (zbiry, ZBIRY) (zr, zi, id, 2, ar, ai, ierr);
+
+  if (kode != 2)
+    {
+      Complex expz = exp( std::abs( real( 2. / 3. * z * sqrt(z))));
+      double rexpz = real(expz);
+      double iexpz = imag(expz);
+      double tmp = ar*rexpz - ai*iexpz;
+      ai = ar*iexpz + ai*rexpz;
+      ar = tmp;
+    }
 
   if (zi == 0.0 && (! scaled || zr >= 0.0))
     ai = 0.0;