Mercurial > octave-nkf

changeset 8279:b3734f1cb592

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
lo-specfun.cc: fix prototypes and calls to cbes{h,i,j,k,y} subroutines
author John W. Eaton <jwe@octave.org>
date Tue, 28 Oct 2008 11:45:57 -0400
parents ab0674a8b345
children 5ee11a81688e
files liboctave/ChangeLog liboctave/lo-specfun.cc
diffstat 2 files changed, 77 insertions(+), 96 deletions(-) [+]
line wrap: on
line diff
--- a/liboctave/ChangeLog	Tue Oct 28 10:47:26 2008 -0400
+++ b/liboctave/ChangeLog	Tue Oct 28 11:45:57 2008 -0400
@@ -1,3 +1,10 @@
+2008-10-28  John W. Eaton  <jwe@octave.org>
+
+	* lo-specfun.cc: Fix prototypes for the Fortran subroutines cbesh,
+	cbesi, cbesj, cbesk, and cbesy.
+	(cbesh, cbesi, cbesj, cbesk, cbesy): Fix calls to Fortran
+	subroutines.
+
 2008-10-28  Brian Gough  <bjg@gnu.org>
 
 	* lo-specfun.cc (zbesi): Fix scaling factor for negative alpha.
--- a/liboctave/lo-specfun.cc	Tue Oct 28 10:47:26 2008 -0400
+++ b/liboctave/lo-specfun.cc	Tue Oct 28 11:45:57 2008 -0400
@@ -77,29 +77,31 @@
 			   double*, octave_idx_type&, octave_idx_type&);
 
   F77_RET_T
-  F77_FUNC (cbesj, cBESJ) (const float&, const float&, const float&,
-			   const octave_idx_type&, const octave_idx_type&, float*, float*,
-			   octave_idx_type&, octave_idx_type&);
+  F77_FUNC (cbesj, cBESJ) (const FloatComplex&, const float&,
+			   const octave_idx_type&, const octave_idx_type&,
+			   FloatComplex*, octave_idx_type&, octave_idx_type&);
 
   F77_RET_T
-  F77_FUNC (cbesy, CBESY) (const float&, const float&, const float&,
-			   const octave_idx_type&, const octave_idx_type&, float*, float*,
-			   octave_idx_type&, float*, float*, octave_idx_type&);
+  F77_FUNC (cbesy, CBESY) (const FloatComplex&, const float&,
+			   const octave_idx_type&, const octave_idx_type&,
+			   FloatComplex*, octave_idx_type&,
+			   FloatComplex*, octave_idx_type&);
 
   F77_RET_T
-  F77_FUNC (cbesi, CBESI) (const float&, const float&, const float&,
-			   const octave_idx_type&, const octave_idx_type&, float*, float*,
-			   octave_idx_type&, octave_idx_type&);
+  F77_FUNC (cbesi, CBESI) (const FloatComplex&, const float&,
+			   const octave_idx_type&, const octave_idx_type&,
+			   FloatComplex*, octave_idx_type&, octave_idx_type&);
 
   F77_RET_T
-  F77_FUNC (cbesk, CBESK) (const float&, const float&, const float&,
-			   const octave_idx_type&, const octave_idx_type&, float*, float*,
-			   octave_idx_type&, octave_idx_type&);
+  F77_FUNC (cbesk, CBESK) (const FloatComplex&, const float&,
+			   const octave_idx_type&, const octave_idx_type&,
+			   FloatComplex*, octave_idx_type&, octave_idx_type&);
 
   F77_RET_T
-  F77_FUNC (cbesh, CBESH) (const float&, const float&, const float&,
-			   const octave_idx_type&, const octave_idx_type&, const octave_idx_type&, float*,
-			   float*, octave_idx_type&, octave_idx_type&);
+  F77_FUNC (cbesh, CBESH) (const FloatComplex&, const float&,
+			   const octave_idx_type&, const octave_idx_type&,
+			   const octave_idx_type&, FloatComplex*,
+			   octave_idx_type&, octave_idx_type&);
 
   F77_RET_T
   F77_FUNC (zairy, ZAIRY) (const double&, const double&, const octave_idx_type&,
@@ -1289,27 +1291,23 @@
 
   if (alpha >= 0.0)
     {
-      float yr = 0.0;
-      float yi = 0.0;
+      FloatComplex y = 0.0;
 
       octave_idx_type nz;
 
-      float zr = z.real ();
-      float zi = z.imag ();
-
-      F77_FUNC (cbesj, CBESJ) (zr, zi, alpha, 2, 1, &yr, &yi, nz, ierr);
+      F77_FUNC (cbesj, CBESJ) (z, alpha, 2, 1, &y, nz, ierr);
 
       if (kode != 2)
 	{
-	  float expz = exp (std::abs (zi)); 
-	  yr *= expz;
-	  yi *= expz;
+	  float expz = exp (std::abs (imag (z))); 
+	  y.real () *= expz;
+	  y.imag () *= expz;
 	}
 
-      if (zi == 0.0 && zr >= 0.0)
-	yi = 0.0;
-
-      retval = bessel_return_value (FloatComplex (yr, yi), ierr);
+      if (imag (z) == 0.0 && real (z) >= 0.0)
+	y.imag () = 0.0;
+
+      retval = bessel_return_value (y, ierr);
     }
   else if (is_integer_value (alpha))
     {
@@ -1346,40 +1344,34 @@
 
   if (alpha >= 0.0)
     {
-      float yr = 0.0;
-      float yi = 0.0;
+      FloatComplex y = 0.0;
 
       octave_idx_type nz;
 
-      float wr, wi;
-
-      float zr = z.real ();
-      float zi = z.imag ();
+      FloatComplex w;
 
       ierr = 0;
 
-      if (zr == 0.0 && zi == 0.0)
+      if (real (z) == 0.0 && imag (z) == 0.0)
 	{
-	  yr = -octave_Float_Inf;
-	  yi = 0.0;
+	  y = FloatComplex (-octave_Float_Inf, 0.0);
 	}
       else
 	{
-	  F77_FUNC (cbesy, CBESY) (zr, zi, alpha, 2, 1, &yr, &yi, nz,
-				   &wr, &wi, ierr);
+	  F77_FUNC (cbesy, CBESY) (z, alpha, 2, 1, &y, nz, &w, ierr);
 
 	  if (kode != 2)
 	    {
-	      float expz = exp (std::abs (zi));
-	      yr *= expz;
-	      yi *= expz;
+	      float expz = exp (std::abs (imag (z)));
+	      y.real () *= expz;
+	      y.imag () *= expz;
 	    }
 
-	  if (zi == 0.0 && zr >= 0.0)
-	    yi = 0.0;
+	  if (imag (z) == 0.0 && real (z) >= 0.0)
+	    y.imag () = 0.0;
 	}
 
-      return bessel_return_value (FloatComplex (yr, yi), ierr);
+      return bessel_return_value (y, ierr);
     }
   else if (is_integer_value (alpha - 0.5))
     {
@@ -1416,27 +1408,22 @@
 
   if (alpha >= 0.0)
     {
-      float yr = 0.0;
-      float yi = 0.0;
+      FloatComplex y = 0.0;
 
       octave_idx_type nz;
 
-      float zr = z.real ();
-      float zi = z.imag ();
-
-      F77_FUNC (cbesi, CBESI) (zr, zi, alpha, 2, 1, &yr, &yi, nz, ierr);
+      F77_FUNC (cbesi, CBESI) (z, alpha, 2, 1, &y, nz, ierr);
 
       if (kode != 2)
 	{
-	  float expz = exp (std::abs (zr));
-	  yr *= expz;
-	  yi *= expz;
+	  float expz = exp (std::abs (real (z)));
+	  y *= expz;
 	}
 
-      if (zi == 0.0 && zr >= 0.0)
-	yi = 0.0;
-
-      retval = bessel_return_value (FloatComplex (yr, yi), ierr);
+      if (imag (z) == 0.0 && real (z) >= 0.0)
+	y.imag () = 0.0;
+
+      retval = bessel_return_value (y, ierr);
     }
   else
     {
@@ -1473,24 +1460,19 @@
 
   if (alpha >= 0.0)
     {
-      float yr = 0.0;
-      float yi = 0.0;
+      FloatComplex y = 0.0;
 
       octave_idx_type nz;
 
-      float zr = z.real ();
-      float zi = z.imag ();
-
       ierr = 0;
 
-      if (zr == 0.0 && zi == 0.0)
+      if (real (z) == 0.0 && imag (z) == 0.0)
 	{
-	  yr = octave_Float_Inf;
-	  yi = 0.0;
+	  y = FloatComplex (octave_Float_Inf, 0.0);
 	}
       else
 	{
-	  F77_FUNC (cbesk, CBESK) (zr, zi, alpha, 2, 1, &yr, &yi, nz, ierr);
+	  F77_FUNC (cbesk, CBESK) (z, alpha, 2, 1, &y, nz, ierr);
 
 	  if (kode != 2)
 	    {
@@ -1499,17 +1481,17 @@
 	      float rexpz = real (expz);
 	      float iexpz = imag (expz);
 
-	      float tmp = yr*rexpz - yi*iexpz;
-
-	      yi = yr*iexpz + yi*rexpz;
-	      yr = tmp;
+	      float tmp = real (y) * rexpz - imag (y) * iexpz;
+
+	      y.imag () = real (y) * iexpz + imag (y) * rexpz;
+	      y.real () = tmp;
 	    }
 
-	  if (zi == 0.0 && zr >= 0.0)
-	    yi = 0.0;
+	  if (imag (z) == 0.0 && real (z) >= 0.0)
+	    y.imag () = 0.0;
 	}
 
-      retval = bessel_return_value (FloatComplex (yr, yi), ierr);
+      retval = bessel_return_value (y, ierr);
     }
   else
     {
@@ -1528,15 +1510,11 @@
 
   if (alpha >= 0.0)
     {
-      float yr = 0.0;
-      float yi = 0.0;
+      FloatComplex y = 0.0;
 
       octave_idx_type nz;
 
-      float zr = z.real ();
-      float zi = z.imag ();
-
-      F77_FUNC (cbesh, CBESH) (zr, zi, alpha, 2, 1, 1, &yr, &yi, nz, ierr);
+      F77_FUNC (cbesh, CBESH) (z, alpha, 2, 1, 1, &y, nz, ierr);
 
       if (kode != 2)
 	{
@@ -1545,13 +1523,13 @@
 	  float rexpz = real (expz);
 	  float iexpz = imag (expz);
 
-	  float tmp = yr*rexpz - yi*iexpz;
-
-	  yi = yr*iexpz + yi*rexpz;
-	  yr = tmp;
+	  float tmp = real (y) * rexpz - imag (y) * iexpz;
+
+	  y.imag () = real (y) * iexpz + imag (y) * rexpz;
+	  y.real () = tmp;
 	}
 
-      retval = bessel_return_value (FloatComplex (yr, yi), ierr);
+      retval = bessel_return_value (y, ierr);
     }
   else
     {
@@ -1574,15 +1552,11 @@
 
   if (alpha >= 0.0)
     {
-      float yr = 0.0;
-      float yi = 0.0;
+      FloatComplex y = 0.0;
 
       octave_idx_type nz;
 
-      float zr = z.real ();
-      float zi = z.imag ();
-
-      F77_FUNC (cbesh, CBESH) (zr, zi, alpha, 2, 2, 1, &yr, &yi, nz, ierr);
+      F77_FUNC (cbesh, CBESH) (z, alpha, 2, 2, 1, &y, nz, ierr);
 
       if (kode != 2)
 	{
@@ -1591,13 +1565,13 @@
 	  float rexpz = real (expz);
 	  float iexpz = imag (expz);
 
-	  float tmp = yr*rexpz - yi*iexpz;
-
-	  yi = yr*iexpz + yi*rexpz;
-	  yr = tmp;
+	  float tmp = real (y) * rexpz - imag (y) * iexpz;
+
+	  y.imag () = real (y) * iexpz + imag (y) * rexpz;
+	  y.real () = tmp;
 	}
 
-      retval = bessel_return_value (FloatComplex (yr, yi), ierr);
+      retval = bessel_return_value (y, ierr);
     }
   else
     {