mirror of
https://git.savannah.gnu.org/git/guile.git
synced 2025-05-29 16:30:19 +02:00
(scm_divide): Adapt code from libstdc++/f2c to void potential overflow
problems. Thanks to John W Eaton!
This commit is contained in:
parent
1b74237d4f
commit
5c90824dd9
1 changed files with 58 additions and 7 deletions
|
@ -1,4 +1,8 @@
|
|||
/* Copyright (C) 1995,1996,1997,1998,1999,2000,2001 Free Software Foundation, Inc.
|
||||
*
|
||||
* Portions Copyright 1990, 1991, 1992, 1993 by AT&T Bell Laboratories
|
||||
* and Bellcore. See scm_divide.
|
||||
*
|
||||
*
|
||||
* This program is free software; you can redistribute it and/or modify
|
||||
* it under the terms of the GNU General Public License as published by
|
||||
|
@ -3645,6 +3649,33 @@ scm_num2dbl (SCM a, const char *why)
|
|||
#undef FUNC_NAME
|
||||
|
||||
|
||||
/* The code below for complex division is adapted from the GNU
|
||||
libstdc++, which adapted it from f2c's libF77, and is subject to
|
||||
this copyright: */
|
||||
|
||||
/****************************************************************
|
||||
Copyright 1990, 1991, 1992, 1993 by AT&T Bell Laboratories and Bellcore.
|
||||
|
||||
Permission to use, copy, modify, and distribute this software
|
||||
and its documentation for any purpose and without fee is hereby
|
||||
granted, provided that the above copyright notice appear in all
|
||||
copies and that both that the copyright notice and this
|
||||
permission notice and warranty disclaimer appear in supporting
|
||||
documentation, and that the names of AT&T Bell Laboratories or
|
||||
Bellcore or any of their entities not be used in advertising or
|
||||
publicity pertaining to distribution of the software without
|
||||
specific, written prior permission.
|
||||
|
||||
AT&T and Bellcore disclaim all warranties with regard to this
|
||||
software, including all implied warranties of merchantability
|
||||
and fitness. In no event shall AT&T or Bellcore be liable for
|
||||
any special, indirect or consequential damages or any damages
|
||||
whatsoever resulting from loss of use, data or profits, whether
|
||||
in an action of contract, negligence or other tortious action,
|
||||
arising out of or in connection with the use or performance of
|
||||
this software.
|
||||
****************************************************************/
|
||||
|
||||
SCM_GPROC1 (s_divide, "/", scm_tc7_asubr, scm_divide, g_divide);
|
||||
/* Divide the first argument by the product of the remaining
|
||||
arguments. If called with one argument @var{z1}, 1/@var{z1} is
|
||||
|
@ -3671,8 +3702,15 @@ scm_divide (SCM x, SCM y)
|
|||
} else if (SCM_COMPLEXP (x)) {
|
||||
double r = SCM_COMPLEX_REAL (x);
|
||||
double i = SCM_COMPLEX_IMAG (x);
|
||||
double d = r * r + i * i;
|
||||
return scm_make_complex (r / d, -i / d);
|
||||
if (r <= i) {
|
||||
double t = r / i;
|
||||
double d = i * (1.0 + t * t);
|
||||
return scm_make_complex (t / d, -1.0 / d);
|
||||
} else {
|
||||
double t = i / r;
|
||||
double d = r * (1.0 + t * t);
|
||||
return scm_make_complex (1.0 / d, -t / d);
|
||||
}
|
||||
} else {
|
||||
SCM_WTA_DISPATCH_1 (g_divide, x, SCM_ARG1, s_divide);
|
||||
}
|
||||
|
@ -3708,8 +3746,15 @@ scm_divide (SCM x, SCM y)
|
|||
{
|
||||
double r = SCM_COMPLEX_REAL (y);
|
||||
double i = SCM_COMPLEX_IMAG (y);
|
||||
double d = r * r + i * i;
|
||||
return scm_make_complex ((a * r) / d, (-a * i) / d);
|
||||
if (r <= i) {
|
||||
double t = r / i;
|
||||
double d = i * (1.0 + t * t);
|
||||
return scm_make_complex ((a * t) / d, -a / d);
|
||||
} else {
|
||||
double t = i / r;
|
||||
double d = r * (1.0 + t * t);
|
||||
return scm_make_complex (a / d, -(a * t) / d);
|
||||
}
|
||||
}
|
||||
} else {
|
||||
SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
|
||||
|
@ -3792,9 +3837,15 @@ scm_divide (SCM x, SCM y)
|
|||
} else if (SCM_COMPLEXP (y)) {
|
||||
double ry = SCM_COMPLEX_REAL (y);
|
||||
double iy = SCM_COMPLEX_IMAG (y);
|
||||
double d = ry * ry + iy * iy;
|
||||
return scm_make_complex ((rx * ry + ix * iy) / d,
|
||||
(ix * ry - rx * iy) / d);
|
||||
if (ry <= iy) {
|
||||
double t = ry / iy;
|
||||
double d = iy * (1.0 + t * t);
|
||||
return scm_make_complex ((rx * t + ix) / d, (ix * t - rx) / d);
|
||||
} else {
|
||||
double t = iy / ry;
|
||||
double d = ry * (1.0 + t * t);
|
||||
return scm_make_complex ((rx + ix * t) / d, (ix - rx * t) / d);
|
||||
}
|
||||
} else {
|
||||
SCM_WTA_DISPATCH_2 (g_divide, x, y, SCM_ARGn, s_divide);
|
||||
}
|
||||
|
|
Loading…
Add table
Add a link
Reference in a new issue