Mercurial > octave-nkf

changeset 4506:3c82fc8f822c

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 2003-09-10 13:56:57 by jwe]
author jwe
date Wed, 10 Sep 2003 13:56:57 +0000
parents e944fbe3fff2
children 65f47f8a92a2
files liboctave/ChangeLog liboctave/lo-specfun.cc src/ChangeLog
diffstat 3 files changed, 105 insertions(+), 12 deletions(-) [+]
line wrap: on
line diff
--- a/liboctave/ChangeLog	Tue Sep 09 19:14:06 2003 +0000
+++ b/liboctave/ChangeLog	Wed Sep 10 13:56:57 2003 +0000
@@ -1,3 +1,14 @@
+2003-09-09  David Bateman <dbateman@free.fr>
+
+	* lo-specfun.cc (zbesj, zbesy, zbesi, zbesk, zbesh1, zbesh2, airy,
+	biry): Always request scaled results from AMOS functions and
+	perform reverse scaling on results if scaled result not requested
+	by user.
+
+2003-09-04  John W. Eaton  <jwe@bevo.che.wisc.edu>
+ 
+ 	* lo-specfun.cc (xlgamma): Require nonnegative argument.
+
 2003-09-09  John W. Eaton  <jwe@bevo.che.wisc.edu>
 
 	* Array-d.cc: Instantiate assign functions.
--- a/liboctave/lo-specfun.cc	Tue Sep 09 19:14:06 2003 +0000
+++ b/liboctave/lo-specfun.cc	Wed Sep 10 13:56:57 2003 +0000
@@ -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,8 +282,15 @@
 	}
       else
 	{
-	  F77_FUNC (zbesy, ZBESY) (zr, zi, alpha, kode, 1, &yr, &yi, nz,
-				  &wr, &wi, ierr);
+	  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,20 @@
 	}
       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 +436,20 @@
       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.0, 1.0) * 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 +482,20 @@
       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.0, 1.0) * 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);
     }
@@ -618,9 +678,20 @@
 
   int id = deriv ? 1 : 0;
 
-  int kode = scaled ? 2 : 1;
+  F77_FUNC (zairy, ZAIRY) (zr, zi, id, 2, ar, ai, nz, ierr);
+
+  if (! scaled)
+    {
+      Complex expz = exp (- 2.0 / 3.0 * z * sqrt(z));
 
-  F77_FUNC (zairy, ZAIRY) (zr, zi, id, kode, ar, ai, nz, ierr);
+      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;
@@ -639,9 +710,20 @@
 
   int id = deriv ? 1 : 0;
 
-  int kode = scaled ? 2 : 1;
+  F77_FUNC (zbiry, ZBIRY) (zr, zi, id, 2, ar, ai, ierr);
+
+  if (! scaled)
+    {
+      Complex expz = exp (std::abs (real (2.0 / 3.0 * z * sqrt (z))));
 
-  F77_FUNC (zbiry, ZBIRY) (zr, zi, id, kode, ar, ai, ierr);
+      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;
--- a/src/ChangeLog	Tue Sep 09 19:14:06 2003 +0000
+++ b/src/ChangeLog	Wed Sep 10 13:56:57 2003 +0000
@@ -29,7 +29,7 @@
 	New files.
 	* Makefile.in (OP_XSRC): Add them to the list.
 
-2003-09-08  D.  <dbateman@free.fr>
+2003-09-09  David Bateman <dbateman@free.fr>
 
 	* OPERATORS/op-cs-s.cc (DEFBINOP): First arg is complex, second is
 	double.