Moving source code

1 files changed, 0 insertions, 241 deletions

diff --git a/src/elementaryFunctions/atan/zatans.c b/src/elementaryFunctions/atan/zatans.c
deleted file mode 100644
index c511d790..00000000
--- a/src/elementaryFunctions/atan/zatans.c
+++ /dev/null
@@ -1,241 +0,0 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2006-2008 - INRIA - Bruno JOFRET
- *
- *  This file must be used under the terms of the CeCILL.
- *  This source file is licensed as described in the file COPYING, which
- *  you should have received as part of this distribution.  The terms
- *  are also available at
- *  http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
- *
- */
-
-/*
-  PURPOSE
-  watan compute the arctangent of a complex number
-  y = yr + i yi = atan(x), x = xr + i xi
-
-  CALLING LIST / PARAMETERS
-  subroutine watan(xr,xi,yr,yi)
-  double precision xr,xi,yr,yi
-
-  xr,xi: real and imaginary parts of the complex number
-  yr,yi: real and imaginary parts of the result
-  yr,yi may have the same memory cases than xr et xi
-
-  COPYRIGHT (C) 2001 Bruno Pincon and Lydia van Dijk
-  Written by Bruno Pincon <Bruno.Pincon@iecn.u-nancy.fr> so
-  as to get more precision.  Also to fix the
-  behavior at the singular points and at the branch cuts.
-  Polished by Lydia van Dijk
-  <lvandijk@hammersmith-consulting.com>
-
-  CHANGES : - (Bruno on 2001 May 22) for ysptrk use a
-  minimax polynome to enlarge the special
-  evaluation zone |s| < SLIM. Also rename
-  this function as lnp1m1.
-  - (Bruno on 2001 June 7) better handling
-  of spurious over/underflow ; remove
-  the call to pythag ; better accuracy
-  in the real part for z near +-i
-
-  EXTERNALS FUNCTIONS
-  dlamch
-  lnp1m1  (at the end of this file)
-
-  ALGORITHM : noting z = a + i*b, we have:
-  Z = yr + yi*b = arctan(z) = (i/2) * log( (i+z)/(i-z) )
-
-  This function has two branch points at +i and -i and the
-  chosen  branch cuts are the two half-straight lines
-  D1 = [i, i*oo) and D2 = (-i*oo, i].  The function is then
-  analytic in C \ (D1 U D2)).
-
-  From the definition it follows that:
-
-  yr = 0.5 Arg ( (i+z)/(i-z) )                   (1)
-  yi = 0.5 log (|(i+z)/(i-z)|)                   (2)
-
-  so lim (z -> +- i) yr = undefined (and Nan is logical)
-  lim (z -> +i)   yi = +oo
-  lim (z -> -i)   yi = -oo
-
-  The real part of arctan(z) is discontinuous across D1 and D2
-  and we impose the following definitions:
-  if imag(z) > 1 then
-  Arg(arctan(z)) =  pi/2 (=lim real(z) -> 0+)
-  if imag(z) < 1 then
-  Arg(arctan(z)) = -pi/2 (=lim real(z) -> 0-)
-
-
-  Basic evaluation: if we write (i+z)/(i-z) using
-  z = a + i*b, we get:
-
-  i+z    1-(a**2+b**2) + i*(2a)
-  --- =  ----------------------
-  i-z       a**2 + (1-b)**2
-
-  then, with r2 = |z|^2 = a**2 + b**2 :
-
-  yr = 0.5 * Arg(1-r2 + (2*a)*i)
-  = 0.5 * atan2(2a, (1-r2))                      (3)
-
-  This formula is changed when r2 > RMAX (max pos float)
-  and also when |1-r2| and |a| are near 0 (see comments
-  in the code).
-
-  After some math:
-
-  yi = 0.25 * log( (a**2 + (b + 1)**2) /
-  (a**2 + (b - 1)**2) )            (4)
-
-  Evaluation for "big" |z|
-  ------------------------
-
-  If |z| is "big", the direct evaluation of yi by (4) may
-  suffer of innaccuracies and of spurious overflow.  Noting
-  that  s = 2 b / (1 + |z|**2), we have:
-
-  yi = 0.25 log ( (1 + s)/(1 - s) )                 (5)
-
-  3        5
-  yi = 0.25*( 2 * ( s + 1/3 s  + 1/5 s  + ... ))
-
-  yi = 0.25 * lnp1m1(s)    if  |s| < SLIM
-
-  So if |s| is less than SLIM we switch to a special
-  evaluation done by the function lnp1m1. The
-  threshold value SLIM is choosen by experiment
-  (with the Pari-gp software). For |s|
-  "very small" we used a truncated taylor dvp,
-  else a minimax polynome (see lnp1m1).
-
-  To avoid spurious overflows (which result in spurious
-  underflows for s) in computing s with s= 2 b / (1 + |z|**2)
-  when |z|^2 > RMAX (max positive float) we use :
-
-  s = 2d0 / ( (a/b)*a + b )
-
-  but if |b| = Inf  this formula leads to NaN when
-  |a| is also Inf. As we have :
-
-  |s| <= 2 / |b|
-
-  we impose simply : s = 0  when |b| = Inf
-
-  Evaluation for z very near to i or -i:
-  --------------------------------------
-  Floating point numbers of the form a+i or a-i with 0 <
-  a**2 < tiny (approximately 1d-308) may lead to underflow
-  (i.e., a**2 = 0) and the logarithm will break formula (4).
-  So we switch to the following formulas:
-
-  If b = +-1 and |a| < sqrt(tiny) approximately 1d-150 (say)
-  then (by using that a**2 + 4 = 4 in machine for such a):
-
-  yi = 0.5 * log( 2/|a| )   for b=1
-
-  yi = 0.5 * log( |a|/2 )   for b=-1
-
-  finally: yi = 0.5 * sign(b) * log( 2/|a| )
-  yi = 0.5 * sign(b) * (log(2) - log(|a|)) (6)
-
-  The last trick is to avoid overflow for |a|=tiny!  In fact
-  this formula may be used until a**2 + 4 = 4 so that the
-  threshold value may be larger.
-*/
-
-#include <math.h>
-#include "lapack.h"
-#include "atan.h"
-#include "abs.h"
-#include "lnp1m1.h"
-
-#define _sign(a, b) b >=0 ? a : -a
-
-doubleComplex		zatans(doubleComplex z) {
-  static double sSlim	= 0.2;
-  static double sAlim	= 1E-150;
-  static double sTol	= 0.3;
-  static double sLn2	= 0.6931471805599453094172321;
-
-  double RMax				= getOverflowThreshold();
-  double Pi_2				= 2.0 * datans(1);
-
-  double _inReal = zreals(z);
-  double _inImg = zimags(z);
-  double _outReal = 0;
-  double _outImg = 0;
-
-  /* Temporary variables */
-  double R2 = 0;
-  double S = 0;
-
-
-  if(_inImg == 0)
-    {
-      _outReal	= datans(_inReal);
-      _outImg	= 0;
-    }
-  else
-    {
-      R2 = _inReal * _inReal + _inImg * _inImg; /* Oo */
-      if(R2 > RMax)
-	{
-	  if( dabss(_inImg) > RMax)
-	    S = 0;
-	  else
-	    S = 1 / (((0.5 * _inReal) / _inImg) * _inReal + 0.5 * _inImg );
-	}
-      else
-	S = (2 * _inImg) / (1+R2);
-
-      if(dabss(S) < sSlim)
-	{
-	  /*
-	    s is small: |s| < SLIM  <=>  |z| outside the following disks:
-	    D+ = D(center = [0;  1/slim], radius = sqrt(1/slim**2 - 1)) if b > 0
-	    D- = D(center = [0; -1/slim], radius = sqrt(1/slim**2 - 1)) if b < 0
-	    use the special evaluation of log((1+s)/(1-s)) (5)
-	  */
-	  _outImg = dlnp1m1s(S) * 0.25;
-	}
-      else
-	{
-	  if(dabss(S) == 1 && dabss(_inReal) <= sAlim)
-	    {
-	      /* |s| >= SLIM  => |z| is inside D+ or D- */
-	      _outImg = _sign(0.5,_inImg) * ( sLn2 - log(dabss(_inReal)));
-	    }
-	  else
-	    {
-	      _outImg = 0.25 * log((pow(_inReal,2) + pow((_inImg + 1),2)) / (pow(_inReal,2) + pow((_inImg - 1),2)));
-	    }
-	}
-      if(_inReal == 0)
-	{/* z is purely imaginary */
-	  if( dabss(_inImg) > 1)
-	    {/* got sign(b) * pi/2 */
-	      _outReal = _sign(1, _inImg) * Pi_2;
-	    }
-	  else if( dabss(_inImg) == 1)
-	    {/* got a Nan with 0/0 */
-	      _outReal = (_inReal - _inReal) / (_inReal - _inReal); /* Oo */
-	    }
-	  else
-	    _outReal = 0;
-	}
-      else if(R2 > RMax)
-	{/* _outImg is necessarily very near sign(a)* pi/2 */
-	  _outReal = _sign(1, _inReal) * Pi_2;
-	}
-      else if(dabss(1 - R2) + dabss(_inReal) <= sTol)
-	{/* |b| is very near 1 (and a is near 0)  some cancellation occur in the (next) generic formula */
-	  _outReal = 0.5 * atan2(2 * _inReal, (1-_inImg) * (1 + _inImg) - pow(_inReal,2));
-	}
-      else
-	_outReal = 0.5 * atan2(2 * _inReal, 1 - R2);
-    }
-
-  return DoubleComplex(_outReal, _outImg);
-}