Mercurial > hg > octave-lyh
changeset 4506:3c82fc8f822c
[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 +++ b/liboctave/ChangeLog @@ -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 +++ b/liboctave/lo-specfun.cc @@ -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;