Fixed missing copyright notice.

author: baudin 2010-06-04 12:56:48 +0000
committer: baudin 2010-06-04 12:56:48 +0000
commit: d91785f60376859d9dc0be0739b02ce6b8166bd9 (patch)
tree: c44e2f98fccd1f05bedb9ed9ea2f2dc3bb05417f
parent: a2b9321ab9a83532e3e0250fd370b37b6a70520b (diff)
download: scilab2c-d91785f60376859d9dc0be0739b02ce6b8166bd9.tar.gz
scilab2c-d91785f60376859d9dc0be0739b02ce6b8166bd9.tar.bz2
scilab2c-d91785f60376859d9dc0be0739b02ce6b8166bd9.zip
15 files changed, 1536 insertions, 1518 deletions
diff --git a/scilab2c/src/c/elementaryFunctions/acos/cacoss.c b/scilab2c/src/c/elementaryFunctions/acos/cacoss.c
index 1059c8c3..97420313 100644
--- a/scilab2c/src/c/elementaryFunctions/acos/cacoss.c
+++ b/scilab2c/src/c/elementaryFunctions/acos/cacoss.c
@@ -1,146 +1,147 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2007-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
- *
- */
-
-/*
- * This fonction is a translation of fortran wacos write by Bruno Pincon <Bruno.Pincon@iecn.u-nancy.fr>
- *     REFERENCE
- *        This is a Fortran-77 translation of an algorithm by
- *        T.E. Hull, T. F. Fairgrieve and P.T.P. Tang which
- *        appears in their article :
- *          "Implementing the Complex Arcsine and Arccosine
- *           Functions Using Exception Handling", ACM, TOMS,
- *           Vol 23, No. 3, Sept 1997, p. 299-335
- */
-
-#include "acos.h"
-#include "atan.h"
-#include "log.h"
-#include "log1p.h"
-#include "sqrt.h"
-#include "abs.h"
-#include "lapack.h"
-#include "min.h"
-#include "max.h"
-
-#define localSign(x) (x>0 ? 1.0f : -1.0f)
-
-floatComplex		cacoss(floatComplex z) {
-	static float sfltPi		= 3.1415926535897932384626433f;
-	static float sfltPi_2		= 1.5707963267948966192313216f;
-	static float sfltLn2		= 0.6931471805599453094172321f;
-	static float sfltAcross	= 1.5f;
-	static float sfltBcross	= 0.6417f;
-
-	float fltLsup = ssqrts((float) getOverflowThreshold())/8.0f;
-	float fltLinf = 4.0f * ssqrts((float) getUnderflowThreshold());
-	float fltEpsm = ssqrts((float) getRelativeMachinePrecision());
-
-	float fltAbsReal	= sabss(creals(z));
-	float fltAbsImg		= sabss(cimags(z));
-	float fltSignReal	= localSign(creals(z));
-	float fltSignImg	= localSign(cimags(z));
-
-	float fltR = 0, fltS = 0, fltA = 0, fltB = 0;
-
-	float fltTemp = 0;
-
-	float _pfltReal = 0;
-	float _pfltImg = 0;
-
-	if( min(fltAbsReal, fltAbsImg) > fltLinf && max(fltAbsReal, fltAbsImg) <= fltLsup)
-	  {/* we are in the safe region */
-		fltR = ssqrts( (fltAbsReal + 1 )*(fltAbsReal + 1 ) + fltAbsImg*fltAbsImg);
-		fltS = ssqrts( (fltAbsReal - 1 )*(fltAbsReal - 1 ) + fltAbsImg*fltAbsImg);
-		fltA = 0.5f * ( fltR + fltS );
-		fltB = fltAbsReal / fltA;
-
-
-		/* compute the real part */
-		if(fltB <= sfltBcross)
-			_pfltReal = sacoss(fltB);
-		else if( fltAbsReal <= 1)
-			_pfltReal = satans(ssqrts(0.5f * (fltA + fltAbsReal) * (fltAbsImg*fltAbsImg / (fltR + (fltAbsReal + 1)) + (fltS + (1 - fltAbsReal)))) / fltAbsReal);
-		else
-			_pfltReal = satans((fltAbsImg * ssqrts(0.5f * ((fltA + fltAbsReal) / (fltR + (fltAbsReal + 1)) + (fltA + fltAbsReal) / (fltS + (fltAbsReal - 1))))) / fltAbsReal);
-
-		/* compute the imaginary part */
-		if(fltA <= sfltAcross)
-		{
-			float fltImg1 = 0;
-
-			if(fltAbsReal < 1)
-			  /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(y**2)/(S+(1.d0-x))) */
-				fltImg1 = 0.5f * (fltAbsImg*fltAbsImg / (fltR + (fltAbsReal + 1)) + fltAbsImg*fltAbsImg / (fltS + (1 - fltAbsReal)));
-			else
-			  /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(S+(x-1.d0))) */
-				fltImg1 = 0.5f * (fltAbsImg*fltAbsImg / (fltR + (fltAbsReal + 1)) + (fltS + (fltAbsReal - 1)));
-			/* ai = logp1(Am1 + sqrt(Am1*(A+1.d0))) */
-			fltTemp = fltImg1 + ssqrts(fltImg1 *( fltA + 1));
-			_pfltImg = slog1ps(fltTemp);
-		}
-		else
-		  /* ai = log(A + sqrt(A**2 - 1.d0)) */
-			_pfltImg = slogs(fltA + ssqrts(fltA*fltA - 1));
-	}
-	else
-	  {/* evaluation in the special regions ... */
-		if(fltAbsImg <= fltEpsm * sabss(fltAbsReal - 1))
-		{
-			if(fltAbsReal < 1)
-			{
-				_pfltReal	= sacoss(fltAbsReal);
-				_pfltImg	= fltAbsImg / ssqrts((1 + fltAbsReal) * (1 - fltAbsReal));
-			}
-			else
-			{
-				_pfltReal = 0;
-				if(fltAbsReal <= fltLsup)
-				{
-					fltTemp		= (fltAbsReal - 1) + ssqrts((fltAbsReal - 1) * (fltAbsReal + 1));
-					_pfltImg	= slog1ps(fltTemp);
-				}
-				else
-					_pfltImg	= sfltLn2 + slogs(fltAbsReal);
-			}
-		}
-		else if(fltAbsImg < fltLinf)
-		{
-			_pfltReal	= ssqrts(fltAbsImg);
-			_pfltImg	= _pfltReal;
-		}
-		else if((fltEpsm * fltAbsImg - 1 >= fltAbsReal))
-		{
-			_pfltReal	= sfltPi_2;
-			_pfltImg	= sfltLn2 + slogs(fltAbsImg);
-		}
-		else if(fltAbsReal > 1)
-		{
-			_pfltReal	= satans(fltAbsImg / fltAbsReal);
-			fltTemp		= (fltAbsReal / fltAbsImg)*(fltAbsReal / fltAbsImg);
-			_pfltImg	= sfltLn2 + slogs(fltAbsImg) + 0.5f * slog1ps(fltTemp);
-		}
-		else
-		{
-			float fltTemp2 = ssqrts(1 + fltAbsImg*fltAbsImg);
-			_pfltReal	= sfltPi_2;
-			fltTemp		= 2 * fltAbsImg * (fltAbsImg + fltTemp2);
-			_pfltImg	= 0.5f * slog1ps(fltTemp);
-		}
-	}
-	if(fltSignReal < 0)
-		_pfltReal = sfltPi - _pfltReal;
-
-	if(fltAbsImg != 0 || fltSignReal < 0)
-		_pfltImg = - fltSignImg * _pfltImg;
-
-	return FloatComplex(_pfltReal, _pfltImg);
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2007-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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
+ *
+ */
+
+/*
+ * This fonction is a translation of fortran wacos write by Bruno Pincon <Bruno.Pincon@iecn.u-nancy.fr>
+ *     REFERENCE
+ *        This is a Fortran-77 translation of an algorithm by
+ *        T.E. Hull, T. F. Fairgrieve and P.T.P. Tang which
+ *        appears in their article :
+ *          "Implementing the Complex Arcsine and Arccosine
+ *           Functions Using Exception Handling", ACM, TOMS,
+ *           Vol 23, No. 3, Sept 1997, p. 299-335
+ */
+
+#include "acos.h"
+#include "atan.h"
+#include "log.h"
+#include "log1p.h"
+#include "sqrt.h"
+#include "abs.h"
+#include "lapack.h"
+#include "min.h"
+#include "max.h"
+
+#define localSign(x) (x>0 ? 1.0f : -1.0f)
+
+floatComplex		cacoss(floatComplex z) {
+	static float sfltPi		= 3.1415926535897932384626433f;
+	static float sfltPi_2		= 1.5707963267948966192313216f;
+	static float sfltLn2		= 0.6931471805599453094172321f;
+	static float sfltAcross	= 1.5f;
+	static float sfltBcross	= 0.6417f;
+
+	float fltLsup = ssqrts((float) getOverflowThreshold())/8.0f;
+	float fltLinf = 4.0f * ssqrts((float) getUnderflowThreshold());
+	float fltEpsm = ssqrts((float) getRelativeMachinePrecision());
+
+	float fltAbsReal	= sabss(creals(z));
+	float fltAbsImg		= sabss(cimags(z));
+	float fltSignReal	= localSign(creals(z));
+	float fltSignImg	= localSign(cimags(z));
+
+	float fltR = 0, fltS = 0, fltA = 0, fltB = 0;
+
+	float fltTemp = 0;
+
+	float _pfltReal = 0;
+	float _pfltImg = 0;
+
+	if( min(fltAbsReal, fltAbsImg) > fltLinf && max(fltAbsReal, fltAbsImg) <= fltLsup)
+	  {/* we are in the safe region */
+		fltR = ssqrts( (fltAbsReal + 1 )*(fltAbsReal + 1 ) + fltAbsImg*fltAbsImg);
+		fltS = ssqrts( (fltAbsReal - 1 )*(fltAbsReal - 1 ) + fltAbsImg*fltAbsImg);
+		fltA = 0.5f * ( fltR + fltS );
+		fltB = fltAbsReal / fltA;
+
+
+		/* compute the real part */
+		if(fltB <= sfltBcross)
+			_pfltReal = sacoss(fltB);
+		else if( fltAbsReal <= 1)
+			_pfltReal = satans(ssqrts(0.5f * (fltA + fltAbsReal) * (fltAbsImg*fltAbsImg / (fltR + (fltAbsReal + 1)) + (fltS + (1 - fltAbsReal)))) / fltAbsReal);
+		else
+			_pfltReal = satans((fltAbsImg * ssqrts(0.5f * ((fltA + fltAbsReal) / (fltR + (fltAbsReal + 1)) + (fltA + fltAbsReal) / (fltS + (fltAbsReal - 1))))) / fltAbsReal);
+
+		/* compute the imaginary part */
+		if(fltA <= sfltAcross)
+		{
+			float fltImg1 = 0;
+
+			if(fltAbsReal < 1)
+			  /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(y**2)/(S+(1.d0-x))) */
+				fltImg1 = 0.5f * (fltAbsImg*fltAbsImg / (fltR + (fltAbsReal + 1)) + fltAbsImg*fltAbsImg / (fltS + (1 - fltAbsReal)));
+			else
+			  /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(S+(x-1.d0))) */
+				fltImg1 = 0.5f * (fltAbsImg*fltAbsImg / (fltR + (fltAbsReal + 1)) + (fltS + (fltAbsReal - 1)));
+			/* ai = logp1(Am1 + sqrt(Am1*(A+1.d0))) */
+			fltTemp = fltImg1 + ssqrts(fltImg1 *( fltA + 1));
+			_pfltImg = slog1ps(fltTemp);
+		}
+		else
+		  /* ai = log(A + sqrt(A**2 - 1.d0)) */
+			_pfltImg = slogs(fltA + ssqrts(fltA*fltA - 1));
+	}
+	else
+	  {/* evaluation in the special regions ... */
+		if(fltAbsImg <= fltEpsm * sabss(fltAbsReal - 1))
+		{
+			if(fltAbsReal < 1)
+			{
+				_pfltReal	= sacoss(fltAbsReal);
+				_pfltImg	= fltAbsImg / ssqrts((1 + fltAbsReal) * (1 - fltAbsReal));
+			}
+			else
+			{
+				_pfltReal = 0;
+				if(fltAbsReal <= fltLsup)
+				{
+					fltTemp		= (fltAbsReal - 1) + ssqrts((fltAbsReal - 1) * (fltAbsReal + 1));
+					_pfltImg	= slog1ps(fltTemp);
+				}
+				else
+					_pfltImg	= sfltLn2 + slogs(fltAbsReal);
+			}
+		}
+		else if(fltAbsImg < fltLinf)
+		{
+			_pfltReal	= ssqrts(fltAbsImg);
+			_pfltImg	= _pfltReal;
+		}
+		else if((fltEpsm * fltAbsImg - 1 >= fltAbsReal))
+		{
+			_pfltReal	= sfltPi_2;
+			_pfltImg	= sfltLn2 + slogs(fltAbsImg);
+		}
+		else if(fltAbsReal > 1)
+		{
+			_pfltReal	= satans(fltAbsImg / fltAbsReal);
+			fltTemp		= (fltAbsReal / fltAbsImg)*(fltAbsReal / fltAbsImg);
+			_pfltImg	= sfltLn2 + slogs(fltAbsImg) + 0.5f * slog1ps(fltTemp);
+		}
+		else
+		{
+			float fltTemp2 = ssqrts(1 + fltAbsImg*fltAbsImg);
+			_pfltReal	= sfltPi_2;
+			fltTemp		= 2 * fltAbsImg * (fltAbsImg + fltTemp2);
+			_pfltImg	= 0.5f * slog1ps(fltTemp);
+		}
+	}
+	if(fltSignReal < 0)
+		_pfltReal = sfltPi - _pfltReal;
+
+	if(fltAbsImg != 0 || fltSignReal < 0)
+		_pfltImg = - fltSignImg * _pfltImg;
+
+	return FloatComplex(_pfltReal, _pfltImg);
+}
diff --git a/scilab2c/src/c/elementaryFunctions/acos/zacoss.c b/scilab2c/src/c/elementaryFunctions/acos/zacoss.c
index 7758a932..de6f3fe9 100644
--- a/scilab2c/src/c/elementaryFunctions/acos/zacoss.c
+++ b/scilab2c/src/c/elementaryFunctions/acos/zacoss.c
@@ -1,146 +1,147 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2007-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
- *
- */
-
-/*
- * This fonction is a translation of fortran wacos write by Bruno Pincon <Bruno.Pincon@iecn.u-nancy.fr>
- *     REFERENCE
- *        This is a Fortran-77 translation of an algorithm by
- *        T.E. Hull, T. F. Fairgrieve and P.T.P. Tang which
- *        appears in their article :
- *          "Implementing the Complex Arcsine and Arccosine
- *           Functions Using Exception Handling", ACM, TOMS,
- *           Vol 23, No. 3, Sept 1997, p. 299-335
- */
-
-#include "acos.h"
-#include "atan.h"
-#include "log.h"
-#include "log1p.h"
-#include "sqrt.h"
-#include "abs.h"
-#include "lapack.h"
-#include "min.h"
-#include "max.h"
-
-#define localSign(x) (x>0 ? 1 : -1)
-
-doubleComplex		zacoss(doubleComplex z) {
-	static double sdblPi		= 3.1415926535897932384626433;
-	static double sdblPi_2		= 1.5707963267948966192313216;
-	static double sdblLn2		= 0.6931471805599453094172321;
-	static double sdblAcross	= 1.5;
-	static double sdblBcross	= 0.6417;
-
-	double dblLsup = dsqrts(getOverflowThreshold())/8.0;
-	double dblLinf = 4.0 * dsqrts(getUnderflowThreshold());
-	double dblEpsm = dsqrts(getRelativeMachinePrecision());
-
-	double dblAbsReal	= dabss(zreals(z));
-	double dblAbsImg	= dabss(zimags(z));
-	double dblSignReal	= localSign(zreals(z));
-	double dblSignImg	= localSign(zimags(z));
-
-	double dblR = 0, dblS = 0, dblA = 0, dblB = 0;
-
-	double dblTemp = 0;
-
-	double _pdblReal = 0;
-	double _pdblImg = 0;
-
-	if( min(dblAbsReal, dblAbsImg) > dblLinf && max(dblAbsReal, dblAbsImg) <= dblLsup)
-	  {/* we are in the safe region */
-		dblR = dsqrts( (dblAbsReal + 1 )*(dblAbsReal + 1 ) + dblAbsImg*dblAbsImg);
-		dblS = dsqrts( (dblAbsReal - 1 )*(dblAbsReal - 1 ) + dblAbsImg*dblAbsImg);
-		dblA = 0.5 * ( dblR + dblS );
-		dblB = dblAbsReal / dblA;
-
-
-		/* compute the real part */
-		if(dblB <= sdblBcross)
-			_pdblReal = dacoss(dblB);
-		else if( dblAbsReal <= 1)
-			_pdblReal = datans(dsqrts(0.5 * (dblA + dblAbsReal) * (dblAbsImg*dblAbsImg / (dblR + (dblAbsReal + 1)) + (dblS + (1 - dblAbsReal)))) / dblAbsReal);
-		else
-			_pdblReal = datans((dblAbsImg * dsqrts(0.5 * ((dblA + dblAbsReal) / (dblR + (dblAbsReal + 1)) + (dblA + dblAbsReal) / (dblS + (dblAbsReal - 1))))) / dblAbsReal);
-
-		/* compute the imaginary part */
-		if(dblA <= sdblAcross)
-		{
-			double dblImg1 = 0;
-
-			if(dblAbsReal < 1)
-			  /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(y**2)/(S+(1.d0-x))) */
-				dblImg1 = 0.5 * (dblAbsImg*dblAbsImg / (dblR + (dblAbsReal + 1)) + dblAbsImg*dblAbsImg / (dblS + (1 - dblAbsReal)));
-			else
-			  /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(S+(x-1.d0))) */
-				dblImg1 = 0.5 * (dblAbsImg*dblAbsImg / (dblR + (dblAbsReal + 1)) + (dblS + (dblAbsReal - 1)));
-			/* ai = logp1(Am1 + sqrt(Am1*(A+1.d0))) */
-			dblTemp = dblImg1 + dsqrts(dblImg1 *( dblA + 1));
-			_pdblImg = dlog1ps(dblTemp);
-		}
-		else
-		  /* ai = log(A + sqrt(A**2 - 1.d0)) */
-			_pdblImg = dlogs(dblA + dsqrts(dblA*dblA - 1));
-	}
-	else
-	  {/* evaluation in the special regions ... */
-		if(dblAbsImg <= dblEpsm * dabss(dblAbsReal - 1))
-		{
-			if(dblAbsReal < 1)
-			{
-				_pdblReal	= dacoss(dblAbsReal);
-				_pdblImg	= dblAbsImg / dsqrts((1 + dblAbsReal) * (1 - dblAbsReal));
-			}
-			else
-			{
-				_pdblReal = 0;
-				if(dblAbsReal <= dblLsup)
-				{
-					dblTemp		= (dblAbsReal - 1) + dsqrts((dblAbsReal - 1) * (dblAbsReal + 1));
-					_pdblImg	= dlog1ps(dblTemp);
-				}
-				else
-					_pdblImg	= sdblLn2 + dlogs(dblAbsReal);
-			}
-		}
-		else if(dblAbsImg < dblLinf)
-		{
-			_pdblReal	= dsqrts(dblAbsImg);
-			_pdblImg	= _pdblReal;
-		}
-		else if((dblEpsm * dblAbsImg - 1 >= dblAbsReal))
-		{
-			_pdblReal	= sdblPi_2;
-			_pdblImg	= sdblLn2 + dlogs(dblAbsImg);
-		}
-		else if(dblAbsReal > 1)
-		{
-			_pdblReal	= datans(dblAbsImg / dblAbsReal);
-			dblTemp		= (dblAbsReal / dblAbsImg)*(dblAbsReal / dblAbsImg);
-			_pdblImg	= sdblLn2 + dlogs(dblAbsImg) + 0.5 * dlog1ps(dblTemp);
-		}
-		else
-		{
-			double dblTemp2 = dsqrts(1 + dblAbsImg*dblAbsImg);
-			_pdblReal	= sdblPi_2;
-			dblTemp		= 2 * dblAbsImg * (dblAbsImg + dblTemp2);
-			_pdblImg	= 0.5 * dlog1ps(dblTemp);
-		}
-	}
-	if(dblSignReal < 0)
-		_pdblReal = sdblPi - _pdblReal;
-
-	if(dblAbsImg != 0 || dblSignReal < 0)
-		_pdblImg = - dblSignImg * _pdblImg;
-
-	return DoubleComplex(_pdblReal, _pdblImg);
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2007-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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
+ *
+ */
+
+/*
+ * This fonction is a translation of fortran wacos write by Bruno Pincon <Bruno.Pincon@iecn.u-nancy.fr>
+ *     REFERENCE
+ *        This is a Fortran-77 translation of an algorithm by
+ *        T.E. Hull, T. F. Fairgrieve and P.T.P. Tang which
+ *        appears in their article :
+ *          "Implementing the Complex Arcsine and Arccosine
+ *           Functions Using Exception Handling", ACM, TOMS,
+ *           Vol 23, No. 3, Sept 1997, p. 299-335
+ */
+
+#include "acos.h"
+#include "atan.h"
+#include "log.h"
+#include "log1p.h"
+#include "sqrt.h"
+#include "abs.h"
+#include "lapack.h"
+#include "min.h"
+#include "max.h"
+
+#define localSign(x) (x>0 ? 1 : -1)
+
+doubleComplex		zacoss(doubleComplex z) {
+	static double sdblPi		= 3.1415926535897932384626433;
+	static double sdblPi_2		= 1.5707963267948966192313216;
+	static double sdblLn2		= 0.6931471805599453094172321;
+	static double sdblAcross	= 1.5;
+	static double sdblBcross	= 0.6417;
+
+	double dblLsup = dsqrts(getOverflowThreshold())/8.0;
+	double dblLinf = 4.0 * dsqrts(getUnderflowThreshold());
+	double dblEpsm = dsqrts(getRelativeMachinePrecision());
+
+	double dblAbsReal	= dabss(zreals(z));
+	double dblAbsImg	= dabss(zimags(z));
+	double dblSignReal	= localSign(zreals(z));
+	double dblSignImg	= localSign(zimags(z));
+
+	double dblR = 0, dblS = 0, dblA = 0, dblB = 0;
+
+	double dblTemp = 0;
+
+	double _pdblReal = 0;
+	double _pdblImg = 0;
+
+	if( min(dblAbsReal, dblAbsImg) > dblLinf && max(dblAbsReal, dblAbsImg) <= dblLsup)
+	  {/* we are in the safe region */
+		dblR = dsqrts( (dblAbsReal + 1 )*(dblAbsReal + 1 ) + dblAbsImg*dblAbsImg);
+		dblS = dsqrts( (dblAbsReal - 1 )*(dblAbsReal - 1 ) + dblAbsImg*dblAbsImg);
+		dblA = 0.5 * ( dblR + dblS );
+		dblB = dblAbsReal / dblA;
+
+
+		/* compute the real part */
+		if(dblB <= sdblBcross)
+			_pdblReal = dacoss(dblB);
+		else if( dblAbsReal <= 1)
+			_pdblReal = datans(dsqrts(0.5 * (dblA + dblAbsReal) * (dblAbsImg*dblAbsImg / (dblR + (dblAbsReal + 1)) + (dblS + (1 - dblAbsReal)))) / dblAbsReal);
+		else
+			_pdblReal = datans((dblAbsImg * dsqrts(0.5 * ((dblA + dblAbsReal) / (dblR + (dblAbsReal + 1)) + (dblA + dblAbsReal) / (dblS + (dblAbsReal - 1))))) / dblAbsReal);
+
+		/* compute the imaginary part */
+		if(dblA <= sdblAcross)
+		{
+			double dblImg1 = 0;
+
+			if(dblAbsReal < 1)
+			  /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(y**2)/(S+(1.d0-x))) */
+				dblImg1 = 0.5 * (dblAbsImg*dblAbsImg / (dblR + (dblAbsReal + 1)) + dblAbsImg*dblAbsImg / (dblS + (1 - dblAbsReal)));
+			else
+			  /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(S+(x-1.d0))) */
+				dblImg1 = 0.5 * (dblAbsImg*dblAbsImg / (dblR + (dblAbsReal + 1)) + (dblS + (dblAbsReal - 1)));
+			/* ai = logp1(Am1 + sqrt(Am1*(A+1.d0))) */
+			dblTemp = dblImg1 + dsqrts(dblImg1 *( dblA + 1));
+			_pdblImg = dlog1ps(dblTemp);
+		}
+		else
+		  /* ai = log(A + sqrt(A**2 - 1.d0)) */
+			_pdblImg = dlogs(dblA + dsqrts(dblA*dblA - 1));
+	}
+	else
+	  {/* evaluation in the special regions ... */
+		if(dblAbsImg <= dblEpsm * dabss(dblAbsReal - 1))
+		{
+			if(dblAbsReal < 1)
+			{
+				_pdblReal	= dacoss(dblAbsReal);
+				_pdblImg	= dblAbsImg / dsqrts((1 + dblAbsReal) * (1 - dblAbsReal));
+			}
+			else
+			{
+				_pdblReal = 0;
+				if(dblAbsReal <= dblLsup)
+				{
+					dblTemp		= (dblAbsReal - 1) + dsqrts((dblAbsReal - 1) * (dblAbsReal + 1));
+					_pdblImg	= dlog1ps(dblTemp);
+				}
+				else
+					_pdblImg	= sdblLn2 + dlogs(dblAbsReal);
+			}
+		}
+		else if(dblAbsImg < dblLinf)
+		{
+			_pdblReal	= dsqrts(dblAbsImg);
+			_pdblImg	= _pdblReal;
+		}
+		else if((dblEpsm * dblAbsImg - 1 >= dblAbsReal))
+		{
+			_pdblReal	= sdblPi_2;
+			_pdblImg	= sdblLn2 + dlogs(dblAbsImg);
+		}
+		else if(dblAbsReal > 1)
+		{
+			_pdblReal	= datans(dblAbsImg / dblAbsReal);
+			dblTemp		= (dblAbsReal / dblAbsImg)*(dblAbsReal / dblAbsImg);
+			_pdblImg	= sdblLn2 + dlogs(dblAbsImg) + 0.5 * dlog1ps(dblTemp);
+		}
+		else
+		{
+			double dblTemp2 = dsqrts(1 + dblAbsImg*dblAbsImg);
+			_pdblReal	= sdblPi_2;
+			dblTemp		= 2 * dblAbsImg * (dblAbsImg + dblTemp2);
+			_pdblImg	= 0.5 * dlog1ps(dblTemp);
+		}
+	}
+	if(dblSignReal < 0)
+		_pdblReal = sdblPi - _pdblReal;
+
+	if(dblAbsImg != 0 || dblSignReal < 0)
+		_pdblImg = - dblSignImg * _pdblImg;
+
+	return DoubleComplex(_pdblReal, _pdblImg);
+}
diff --git a/scilab2c/src/c/elementaryFunctions/asin/casins.c b/scilab2c/src/c/elementaryFunctions/asin/casins.c
index 9c15ce53..35a4a8d8 100644
--- a/scilab2c/src/c/elementaryFunctions/asin/casins.c
+++ b/scilab2c/src/c/elementaryFunctions/asin/casins.c
@@ -1,144 +1,146 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2008-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
- *
- */
-
-/*
- *     REFERENCE
- *        This is a Fortran-77 translation of an algorithm by
- *        T.E. Hull, T. F. Fairgrieve and P.T.P. Tang which
- *        appears in their article :
- *          "Implementing the Complex Arcsine and Arccosine
- *           Functions Using Exception Handling", ACM, TOMS,
- *           Vol 23, No. 3, Sept 1997, p. 299-335
- */
-
-#include "lapack.h"
-#include "asin.h"
-#include "atan.h"
-#include "sqrt.h"
-#include "abs.h"
-#include "log.h"
-#include "log1p.h"
-#include "min.h"
-#include "max.h"
-
-floatComplex		casins(floatComplex z) {
-  static float sdblPi_2		= 1.5707963267948966192313216f;
-  static float sdblLn2		= 0.6931471805599453094172321f;
-  static float sdblAcross	= 1.5f;
-  static float sdblBcross	= 0.6417f;
-
-  float dblLsup = ssqrts((float) getOverflowThreshold())/ 8.0f;
-  float dblLinf = 4.0f * ssqrts((float) getUnderflowThreshold());
-  float dblEpsm = ssqrts((float) getRelativeMachinePrecision());
-
-  float _dblReal	= creals(z);
-  float _dblImg		= cimags(z);
-
-  float dblAbsReal	= sabss(_dblReal);
-  float dblAbsImg	= sabss(_dblImg);
-  float iSignReal	= _dblReal < 0 ? -1.0f : 1.0f;
-  float iSignImg	= _dblImg < 0 ? -1.0f : 1.0f;
-
-  float dblR = 0, dblS = 0, dblA = 0, dblB = 0;
-
-  float dblTemp = 0;
-
-  float _pdblReal = 0;
-  float _pdblImg = 0;
-
-  if( min(dblAbsReal, dblAbsImg) > dblLinf && max(dblAbsReal, dblAbsImg) <= dblLsup)
-    {
-      /* we are in the safe region */
-      dblR = ssqrts( (dblAbsReal + 1) * (dblAbsReal + 1) + dblAbsImg * dblAbsImg);
-      dblS = ssqrts( (dblAbsReal - 1) * (dblAbsReal - 1) + dblAbsImg * dblAbsImg);
-      dblA = (float) 0.5 * ( dblR + dblS );
-      dblB = dblAbsReal / dblA;
-
-
-      /* compute the real part */
-      if(dblB <= sdblBcross)
-	_pdblReal = sasins(dblB);
-      else if( dblAbsReal <= 1)
-	_pdblReal = satans(dblAbsReal / ssqrts( 0.5f * (dblA + dblAbsReal) * ( (dblAbsImg * dblAbsImg) / (dblR + (dblAbsReal + 1)) + (dblS + (1 - dblAbsReal)))));
-      else
-	_pdblReal = satans(dblAbsReal / (dblAbsImg * ssqrts( 0.5f * ((dblA + dblAbsReal) / (dblR + (dblAbsReal + 1)) + (dblA + dblAbsReal) / (dblS + (dblAbsReal-1))))));
-
-      /* compute the imaginary part */
-      if(dblA <= sdblAcross)
-	{
-	  float dblImg1 = 0;
-
-	  if(dblAbsReal < 1)
-	    /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(y**2)/(S+(1.d0-x))) */
-	    dblImg1 = 0.5f * (dblAbsImg * dblAbsImg / (dblR + (dblAbsReal + 1)) + dblAbsImg * dblAbsImg / (dblS + (1 - dblAbsReal)));
-	  else
-	    /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(S+(x-1.d0))) */
-	    dblImg1 = 0.5f * (dblAbsImg * dblAbsImg / (dblR + (dblAbsReal + 1)) + (dblS + (dblAbsReal - 1)));
-	  /* ai = logp1(Am1 + sqrt(Am1*(A+1.d0))) */
-	  dblTemp = dblImg1 + ssqrts(dblImg1 * (dblA + 1));
-	  _pdblImg = slog1ps(dblTemp);
-	}
-      else
-	/* ai = log(A + sqrt(A**2 - 1.d0)) */
-	_pdblImg = slogs(dblA + ssqrts(dblA * dblA - (float) 1.0));
-    }
-  else
-    {
-      /* evaluation in the special regions ... */
-      if(dblAbsImg <= dblEpsm * dabss(dblAbsReal - 1))
-	{
-	  if(dblAbsReal < 1)
-	    {
-	      _pdblReal	= sasins(dblAbsReal);
-	      _pdblImg	= dblAbsImg / ssqrts((1 + dblAbsReal) * (1 - dblAbsReal));
-	    }
-	  else
-	    {
-	      _pdblReal = sdblPi_2;
-	      if(dblAbsReal <= dblLsup)
-		{
-		  dblTemp		= (dblAbsReal - 1) + ssqrts((dblAbsReal - 1) * (dblAbsReal + 1));
-		  _pdblImg	= slog1ps(dblTemp);
-		}
-	      else
-		_pdblImg	= sdblLn2 + slogs(dblAbsReal);
-	    }
-	}
-      else if(dblAbsImg < dblLinf)
-	{
-	  _pdblReal	= sdblPi_2 - ssqrts(dblAbsImg);
-	  _pdblImg	= ssqrts(dblAbsImg);
-	}
-      else if((dblEpsm * dblAbsImg - 1 >= dblAbsReal))
-	{
-	  _pdblReal	= dblAbsReal * dblAbsImg;
-	  _pdblImg	= sdblLn2 + slogs(dblAbsReal);
-	}
-      else if(dblAbsReal > 1)
-	{
-	  _pdblReal	= satans(dblAbsReal / dblAbsImg);
-	  dblTemp	= (dblAbsReal / dblAbsImg) * (dblAbsReal / dblAbsImg);
-	  _pdblImg	= sdblLn2 + slogs(dblAbsReal) + 0.5f * slog1ps(dblTemp);
-	}
-      else
-	{
-	  float dblTemp2 = ssqrts(1 + dblAbsImg * dblAbsImg);
-	  _pdblReal	= dblAbsReal / dblTemp2;
-	  dblTemp	= 2.0f * dblAbsImg * (dblAbsImg + dblTemp2);
-	  _pdblImg	= 0.5f * slog1ps(dblTemp);
-	}
-    }
-  _pdblReal *= iSignReal;
-  _pdblImg *= iSignImg;
-
-  return (FloatComplex(_pdblReal, _pdblImg));
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2008-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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
+ *
+ */
+
+/*
+ *     REFERENCE
+ *        This is a Fortran-77 translation of an algorithm by
+ *        T.E. Hull, T. F. Fairgrieve and P.T.P. Tang which
+ *        appears in their article :
+ *          "Implementing the Complex Arcsine and Arccosine
+ *           Functions Using Exception Handling", ACM, TOMS,
+ *           Vol 23, No. 3, Sept 1997, p. 299-335
+ *     Thanks to Tom Fairgrieve
+ */
+
+#include "lapack.h"
+#include "asin.h"
+#include "atan.h"
+#include "sqrt.h"
+#include "abs.h"
+#include "log.h"
+#include "log1p.h"
+#include "min.h"
+#include "max.h"
+
+floatComplex		casins(floatComplex z) {
+  static float sdblPi_2		= 1.5707963267948966192313216f;
+  static float sdblLn2		= 0.6931471805599453094172321f;
+  static float sdblAcross	= 1.5f;
+  static float sdblBcross	= 0.6417f;
+
+  float dblLsup = ssqrts((float) getOverflowThreshold())/ 8.0f;
+  float dblLinf = 4.0f * ssqrts((float) getUnderflowThreshold());
+  float dblEpsm = ssqrts((float) getRelativeMachinePrecision());
+
+  float _dblReal	= creals(z);
+  float _dblImg		= cimags(z);
+
+  float dblAbsReal	= sabss(_dblReal);
+  float dblAbsImg	= sabss(_dblImg);
+  float iSignReal	= _dblReal < 0 ? -1.0f : 1.0f;
+  float iSignImg	= _dblImg < 0 ? -1.0f : 1.0f;
+
+  float dblR = 0, dblS = 0, dblA = 0, dblB = 0;
+
+  float dblTemp = 0;
+
+  float _pdblReal = 0;
+  float _pdblImg = 0;
+
+  if( min(dblAbsReal, dblAbsImg) > dblLinf && max(dblAbsReal, dblAbsImg) <= dblLsup)
+    {
+      /* we are in the safe region */
+      dblR = ssqrts( (dblAbsReal + 1) * (dblAbsReal + 1) + dblAbsImg * dblAbsImg);
+      dblS = ssqrts( (dblAbsReal - 1) * (dblAbsReal - 1) + dblAbsImg * dblAbsImg);
+      dblA = (float) 0.5 * ( dblR + dblS );
+      dblB = dblAbsReal / dblA;
+
+
+      /* compute the real part */
+      if(dblB <= sdblBcross)
+	_pdblReal = sasins(dblB);
+      else if( dblAbsReal <= 1)
+	_pdblReal = satans(dblAbsReal / ssqrts( 0.5f * (dblA + dblAbsReal) * ( (dblAbsImg * dblAbsImg) / (dblR + (dblAbsReal + 1)) + (dblS + (1 - dblAbsReal)))));
+      else
+	_pdblReal = satans(dblAbsReal / (dblAbsImg * ssqrts( 0.5f * ((dblA + dblAbsReal) / (dblR + (dblAbsReal + 1)) + (dblA + dblAbsReal) / (dblS + (dblAbsReal-1))))));
+
+      /* compute the imaginary part */
+      if(dblA <= sdblAcross)
+	{
+	  float dblImg1 = 0;
+
+	  if(dblAbsReal < 1)
+	    /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(y**2)/(S+(1.d0-x))) */
+	    dblImg1 = 0.5f * (dblAbsImg * dblAbsImg / (dblR + (dblAbsReal + 1)) + dblAbsImg * dblAbsImg / (dblS + (1 - dblAbsReal)));
+	  else
+	    /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(S+(x-1.d0))) */
+	    dblImg1 = 0.5f * (dblAbsImg * dblAbsImg / (dblR + (dblAbsReal + 1)) + (dblS + (dblAbsReal - 1)));
+	  /* ai = logp1(Am1 + sqrt(Am1*(A+1.d0))) */
+	  dblTemp = dblImg1 + ssqrts(dblImg1 * (dblA + 1));
+	  _pdblImg = slog1ps(dblTemp);
+	}
+      else
+	/* ai = log(A + sqrt(A**2 - 1.d0)) */
+	_pdblImg = slogs(dblA + ssqrts(dblA * dblA - (float) 1.0));
+    }
+  else
+    {
+      /* evaluation in the special regions ... */
+      if(dblAbsImg <= dblEpsm * dabss(dblAbsReal - 1))
+	{
+	  if(dblAbsReal < 1)
+	    {
+	      _pdblReal	= sasins(dblAbsReal);
+	      _pdblImg	= dblAbsImg / ssqrts((1 + dblAbsReal) * (1 - dblAbsReal));
+	    }
+	  else
+	    {
+	      _pdblReal = sdblPi_2;
+	      if(dblAbsReal <= dblLsup)
+		{
+		  dblTemp		= (dblAbsReal - 1) + ssqrts((dblAbsReal - 1) * (dblAbsReal + 1));
+		  _pdblImg	= slog1ps(dblTemp);
+		}
+	      else
+		_pdblImg	= sdblLn2 + slogs(dblAbsReal);
+	    }
+	}
+      else if(dblAbsImg < dblLinf)
+	{
+	  _pdblReal	= sdblPi_2 - ssqrts(dblAbsImg);
+	  _pdblImg	= ssqrts(dblAbsImg);
+	}
+      else if((dblEpsm * dblAbsImg - 1 >= dblAbsReal))
+	{
+	  _pdblReal	= dblAbsReal * dblAbsImg;
+	  _pdblImg	= sdblLn2 + slogs(dblAbsReal);
+	}
+      else if(dblAbsReal > 1)
+	{
+	  _pdblReal	= satans(dblAbsReal / dblAbsImg);
+	  dblTemp	= (dblAbsReal / dblAbsImg) * (dblAbsReal / dblAbsImg);
+	  _pdblImg	= sdblLn2 + slogs(dblAbsReal) + 0.5f * slog1ps(dblTemp);
+	}
+      else
+	{
+	  float dblTemp2 = ssqrts(1 + dblAbsImg * dblAbsImg);
+	  _pdblReal	= dblAbsReal / dblTemp2;
+	  dblTemp	= 2.0f * dblAbsImg * (dblAbsImg + dblTemp2);
+	  _pdblImg	= 0.5f * slog1ps(dblTemp);
+	}
+    }
+  _pdblReal *= iSignReal;
+  _pdblImg *= iSignImg;
+
+  return (FloatComplex(_pdblReal, _pdblImg));
+}
diff --git a/scilab2c/src/c/elementaryFunctions/asin/zasins.c b/scilab2c/src/c/elementaryFunctions/asin/zasins.c
index 5d2110ff..5bd586a8 100644
--- a/scilab2c/src/c/elementaryFunctions/asin/zasins.c
+++ b/scilab2c/src/c/elementaryFunctions/asin/zasins.c
@@ -1,144 +1,146 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2007-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
- *
- */
-
-/*
- *     REFERENCE
- *        This is a Fortran-77 translation of an algorithm by
- *        T.E. Hull, T. F. Fairgrieve and P.T.P. Tang which
- *        appears in their article :
- *          "Implementing the Complex Arcsine and Arccosine
- *           Functions Using Exception Handling", ACM, TOMS,
- *           Vol 23, No. 3, Sept 1997, p. 299-335
- */
-
-#include "lapack.h"
-#include "asin.h"
-#include "atan.h"
-#include "sqrt.h"
-#include "abs.h"
-#include "log.h"
-#include "log1p.h"
-#include "min.h"
-#include "max.h"
-
-doubleComplex		zasins(doubleComplex z) {
-  static double sdblPi_2	= 1.5707963267948966192313216;
-  static double sdblLn2		= 0.6931471805599453094172321;
-  static double sdblAcross	= 1.5;
-  static double sdblBcross	= 0.6417;
-
-  double dblLsup = dsqrts(getOverflowThreshold())/8.0;
-  double dblLinf = 4 * dsqrts(getUnderflowThreshold());
-  double dblEpsm = dsqrts(getRelativeMachinePrecision());
-
-  double _dblReal	= zreals(z);
-  double _dblImg	= zimags(z);
-
-  double dblAbsReal	= dabss(_dblReal);
-  double dblAbsImg	= dabss(_dblImg);
-  int iSignReal		= _dblReal < 0 ? -1 : 1;
-  int iSignImg		= _dblImg < 0 ? -1 : 1;
-
-  double dblR = 0, dblS = 0, dblA = 0, dblB = 0;
-
-  double dblTemp = 0;
-
-  double _pdblReal = 0;
-  double _pdblImg = 0;
-
-  if( min(dblAbsReal, dblAbsImg) > dblLinf && max(dblAbsReal, dblAbsImg) <= dblLsup)
-    {
-      /* we are in the safe region */
-      dblR = dsqrts( (dblAbsReal + 1) * (dblAbsReal + 1) + dblAbsImg * dblAbsImg);
-      dblS = dsqrts( (dblAbsReal - 1) * (dblAbsReal - 1) + dblAbsImg * dblAbsImg);
-      dblA = 0.5 * ( dblR + dblS );
-      dblB = dblAbsReal / dblA;
-
-
-      /* compute the real part */
-      if(dblB <= sdblBcross)
-	_pdblReal = dasins(dblB);
-      else if( dblAbsReal <= 1)
-	_pdblReal = datans(dblAbsReal / dsqrts( 0.5 * (dblA + dblAbsReal) * ( (dblAbsImg * dblAbsImg) / (dblR + (dblAbsReal + 1)) + (dblS + (1 - dblAbsReal)))));
-      else
-	_pdblReal = datans(dblAbsReal / (dblAbsImg * dsqrts(0.5 * ((dblA + dblAbsReal) / (dblR + (dblAbsReal + 1)) + (dblA + dblAbsReal) / (dblS + (dblAbsReal-1))))));
-
-      /* compute the imaginary part */
-      if(dblA <= sdblAcross)
-	{
-	  double dblImg1 = 0;
-
-	  if(dblAbsReal < 1)
-	    /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(y**2)/(S+(1.d0-x))) */
-	    dblImg1 = 0.5 * (dblAbsImg * dblAbsImg / (dblR + (dblAbsReal + 1)) + dblAbsImg * dblAbsImg / (dblS + (1 - dblAbsReal)));
-	  else
-	    /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(S+(x-1.d0))) */
-	    dblImg1 = 0.5 * (dblAbsImg * dblAbsImg / (dblR + (dblAbsReal + 1)) + (dblS + (dblAbsReal - 1)));
-	  /* ai = logp1(Am1 + sqrt(Am1*(A+1.d0))) */
-	  dblTemp = dblImg1 + dsqrts(dblImg1 * (dblA + 1));
-	  _pdblImg = dlog1ps(dblTemp);
-	}
-      else
-	/* ai = log(A + sqrt(A**2 - 1.d0)) */
-	_pdblImg = dlogs(dblA + dsqrts(dblA * dblA - 1));
-    }
-  else
-    {
-      /* evaluation in the special regions ... */
-      if(dblAbsImg <= dblEpsm * dabss(dblAbsReal - 1))
-	{
-	  if(dblAbsReal < 1)
-	    {
-	      _pdblReal	= dasins(dblAbsReal);
-	      _pdblImg	= dblAbsImg / dsqrts((1 + dblAbsReal) * (1 - dblAbsReal));
-	    }
-	  else
-	    {
-	      _pdblReal = sdblPi_2;
-	      if(dblAbsReal <= dblLsup)
-		{
-		  dblTemp		= (dblAbsReal - 1) + dsqrts((dblAbsReal - 1) * (dblAbsReal + 1));
-		  _pdblImg	= dlog1ps(dblTemp);
-		}
-	      else
-		_pdblImg	= sdblLn2 + dlogs(dblAbsReal);
-	    }
-	}
-      else if(dblAbsImg < dblLinf)
-	{
-	  _pdblReal	= sdblPi_2 - dsqrts(dblAbsImg);
-	  _pdblImg	= dsqrts(dblAbsImg);
-	}
-      else if((dblEpsm * dblAbsImg - 1 >= dblAbsReal))
-	{
-	  _pdblReal	= dblAbsReal * dblAbsImg;
-	  _pdblImg	= sdblLn2 + dlogs(dblAbsReal);
-	}
-      else if(dblAbsReal > 1)
-	{
-	  _pdblReal	= datans(dblAbsReal / dblAbsImg);
-	  dblTemp		= (dblAbsReal / dblAbsImg) * (dblAbsReal / dblAbsImg);
-	  _pdblImg	= sdblLn2 + dlogs(dblAbsReal) + 0.5 * dlog1ps(dblTemp);
-	}
-      else
-	{
-	  double dblTemp2 = dsqrts(1 + dblAbsImg * dblAbsImg);
-	  _pdblReal	= dblAbsReal / dblTemp2;
-	  dblTemp		= 2 * dblAbsImg * (dblAbsImg + dblTemp2);
-	  _pdblImg	= 0.5 * dlog1ps(dblTemp);
-	}
-    }
-  _pdblReal *= iSignReal;
-  _pdblImg *= iSignImg;
-
-  return (DoubleComplex(_pdblReal, _pdblImg));
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2007-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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
+ *
+ */
+
+/*
+ *     REFERENCE
+ *        This is a Fortran-77 translation of an algorithm by
+ *        T.E. Hull, T. F. Fairgrieve and P.T.P. Tang which
+ *        appears in their article :
+ *          "Implementing the Complex Arcsine and Arccosine
+ *           Functions Using Exception Handling", ACM, TOMS,
+ *           Vol 23, No. 3, Sept 1997, p. 299-335
+ *     Thanks to Tom Fairgrieve
+ */
+
+#include "lapack.h"
+#include "asin.h"
+#include "atan.h"
+#include "sqrt.h"
+#include "abs.h"
+#include "log.h"
+#include "log1p.h"
+#include "min.h"
+#include "max.h"
+
+doubleComplex		zasins(doubleComplex z) {
+  static double sdblPi_2	= 1.5707963267948966192313216;
+  static double sdblLn2		= 0.6931471805599453094172321;
+  static double sdblAcross	= 1.5;
+  static double sdblBcross	= 0.6417;
+
+  double dblLsup = dsqrts(getOverflowThreshold())/8.0;
+  double dblLinf = 4 * dsqrts(getUnderflowThreshold());
+  double dblEpsm = dsqrts(getRelativeMachinePrecision());
+
+  double _dblReal	= zreals(z);
+  double _dblImg	= zimags(z);
+
+  double dblAbsReal	= dabss(_dblReal);
+  double dblAbsImg	= dabss(_dblImg);
+  int iSignReal		= _dblReal < 0 ? -1 : 1;
+  int iSignImg		= _dblImg < 0 ? -1 : 1;
+
+  double dblR = 0, dblS = 0, dblA = 0, dblB = 0;
+
+  double dblTemp = 0;
+
+  double _pdblReal = 0;
+  double _pdblImg = 0;
+
+  if( min(dblAbsReal, dblAbsImg) > dblLinf && max(dblAbsReal, dblAbsImg) <= dblLsup)
+    {
+      /* we are in the safe region */
+      dblR = dsqrts( (dblAbsReal + 1) * (dblAbsReal + 1) + dblAbsImg * dblAbsImg);
+      dblS = dsqrts( (dblAbsReal - 1) * (dblAbsReal - 1) + dblAbsImg * dblAbsImg);
+      dblA = 0.5 * ( dblR + dblS );
+      dblB = dblAbsReal / dblA;
+
+
+      /* compute the real part */
+      if(dblB <= sdblBcross)
+	_pdblReal = dasins(dblB);
+      else if( dblAbsReal <= 1)
+	_pdblReal = datans(dblAbsReal / dsqrts( 0.5 * (dblA + dblAbsReal) * ( (dblAbsImg * dblAbsImg) / (dblR + (dblAbsReal + 1)) + (dblS + (1 - dblAbsReal)))));
+      else
+	_pdblReal = datans(dblAbsReal / (dblAbsImg * dsqrts(0.5 * ((dblA + dblAbsReal) / (dblR + (dblAbsReal + 1)) + (dblA + dblAbsReal) / (dblS + (dblAbsReal-1))))));
+
+      /* compute the imaginary part */
+      if(dblA <= sdblAcross)
+	{
+	  double dblImg1 = 0;
+
+	  if(dblAbsReal < 1)
+	    /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(y**2)/(S+(1.d0-x))) */
+	    dblImg1 = 0.5 * (dblAbsImg * dblAbsImg / (dblR + (dblAbsReal + 1)) + dblAbsImg * dblAbsImg / (dblS + (1 - dblAbsReal)));
+	  else
+	    /* Am1 = 0.5d0*((y**2)/(R+(x+1.d0))+(S+(x-1.d0))) */
+	    dblImg1 = 0.5 * (dblAbsImg * dblAbsImg / (dblR + (dblAbsReal + 1)) + (dblS + (dblAbsReal - 1)));
+	  /* ai = logp1(Am1 + sqrt(Am1*(A+1.d0))) */
+	  dblTemp = dblImg1 + dsqrts(dblImg1 * (dblA + 1));
+	  _pdblImg = dlog1ps(dblTemp);
+	}
+      else
+	/* ai = log(A + sqrt(A**2 - 1.d0)) */
+	_pdblImg = dlogs(dblA + dsqrts(dblA * dblA - 1));
+    }
+  else
+    {
+      /* evaluation in the special regions ... */
+      if(dblAbsImg <= dblEpsm * dabss(dblAbsReal - 1))
+	{
+	  if(dblAbsReal < 1)
+	    {
+	      _pdblReal	= dasins(dblAbsReal);
+	      _pdblImg	= dblAbsImg / dsqrts((1 + dblAbsReal) * (1 - dblAbsReal));
+	    }
+	  else
+	    {
+	      _pdblReal = sdblPi_2;
+	      if(dblAbsReal <= dblLsup)
+		{
+		  dblTemp		= (dblAbsReal - 1) + dsqrts((dblAbsReal - 1) * (dblAbsReal + 1));
+		  _pdblImg	= dlog1ps(dblTemp);
+		}
+	      else
+		_pdblImg	= sdblLn2 + dlogs(dblAbsReal);
+	    }
+	}
+      else if(dblAbsImg < dblLinf)
+	{
+	  _pdblReal	= sdblPi_2 - dsqrts(dblAbsImg);
+	  _pdblImg	= dsqrts(dblAbsImg);
+	}
+      else if((dblEpsm * dblAbsImg - 1 >= dblAbsReal))
+	{
+	  _pdblReal	= dblAbsReal * dblAbsImg;
+	  _pdblImg	= sdblLn2 + dlogs(dblAbsReal);
+	}
+      else if(dblAbsReal > 1)
+	{
+	  _pdblReal	= datans(dblAbsReal / dblAbsImg);
+	  dblTemp		= (dblAbsReal / dblAbsImg) * (dblAbsReal / dblAbsImg);
+	  _pdblImg	= sdblLn2 + dlogs(dblAbsReal) + 0.5 * dlog1ps(dblTemp);
+	}
+      else
+	{
+	  double dblTemp2 = dsqrts(1 + dblAbsImg * dblAbsImg);
+	  _pdblReal	= dblAbsReal / dblTemp2;
+	  dblTemp		= 2 * dblAbsImg * (dblAbsImg + dblTemp2);
+	  _pdblImg	= 0.5 * dlog1ps(dblTemp);
+	}
+    }
+  _pdblReal *= iSignReal;
+  _pdblImg *= iSignImg;
+
+  return (DoubleComplex(_pdblReal, _pdblImg));
+}
diff --git a/scilab2c/src/c/elementaryFunctions/atan/catans.c b/scilab2c/src/c/elementaryFunctions/atan/catans.c
index d3a232fb..d2081181 100644
--- a/scilab2c/src/c/elementaryFunctions/atan/catans.c
+++ b/scilab2c/src/c/elementaryFunctions/atan/catans.c
@@ -1,248 +1,249 @@
-/*
- *  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 "atan.h"
-#include "abs.h"
-#include "lnp1m1.h"
-#include "lapack.h"
-
-#define _sign(a, b) b >=0 ? a : -a
-
-floatComplex		catans(floatComplex z) {
-  static float sSlim	= 0.2f;
-  /*    .
-  **   / \	WARNING : this algorithm was based on double precision
-  **  / ! \	using float truncate the value to 0.
-  ** `----'
-  **
-  ** static float sAlim	= 1E-150f;
-  */
-  static float sAlim	= 0.0f;
-  static float sTol	= 0.3f;
-  static float sLn2	= 0.6931471805599453094172321f;
-
-  float RMax				= (float) getOverflowThreshold();
-  float Pi_2				= 2.0f * satans(1);
-
-  float _inReal = creals(z);
-  float _inImg = cimags(z);
-  float _outReal = 0;
-  float _outImg = 0;
-
-  /* Temporary variables */
-  float R2 = 0;
-  float S = 0;
-
-
-  if(_inImg == 0)
-    {
-      _outReal	= satans(_inReal);
-      _outImg	= 0;
-    }
-  else
-    {
-      R2 = _inReal * _inReal + _inImg * _inImg; /* Oo */
-      if(R2 > RMax)
-	{
-	  if( dabss(_inImg) > RMax)
-	    S = 0;
-	  else
-	    S = 1.0f / (((0.5f * _inReal) / _inImg) * _inReal + 0.5f * _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 = slnp1m1s(S) * 0.25f;
-	}
-      else
-	{
-	  if(sabss(S) == 1 && sabss(_inReal) <= sAlim)
-	    {
-	      /* |s| >= SLIM  => |z| is inside D+ or D- */
-	      _outImg = _sign(0.5f,_inImg) * ( sLn2 - logf(sabss(_inReal)));
-	    }
-	  else
-	    {
-	      _outImg = 0.25f * logf((powf(_inReal,2) + powf((_inImg + 1.0f),2)) / (powf(_inReal,2) + powf((_inImg - 1.0f),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(sabss(1 - R2) + sabss(_inReal) <= sTol)
-	{/* |b| is very near 1 (and a is near 0)  some cancellation occur in the (next) generic formula */
-	  _outReal = 0.5f * atan2f(2.0f * _inReal, (1.0f - _inImg) * (1.0f + _inImg) - powf(_inReal,2.0f));
-	}
-      else
-	_outReal = 0.5f * atan2f(2.0f * _inReal, 1.0f - R2);
-    }
-
-  return FloatComplex(_outReal, _outImg);
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2006-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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 "atan.h"
+#include "abs.h"
+#include "lnp1m1.h"
+#include "lapack.h"
+
+#define _sign(a, b) b >=0 ? a : -a
+
+floatComplex		catans(floatComplex z) {
+  static float sSlim	= 0.2f;
+  /*    .
+  **   / \	WARNING : this algorithm was based on double precision
+  **  / ! \	using float truncate the value to 0.
+  ** `----'
+  **
+  ** static float sAlim	= 1E-150f;
+  */
+  static float sAlim	= 0.0f;
+  static float sTol	= 0.3f;
+  static float sLn2	= 0.6931471805599453094172321f;
+
+  float RMax				= (float) getOverflowThreshold();
+  float Pi_2				= 2.0f * satans(1);
+
+  float _inReal = creals(z);
+  float _inImg = cimags(z);
+  float _outReal = 0;
+  float _outImg = 0;
+
+  /* Temporary variables */
+  float R2 = 0;
+  float S = 0;
+
+
+  if(_inImg == 0)
+    {
+      _outReal	= satans(_inReal);
+      _outImg	= 0;
+    }
+  else
+    {
+      R2 = _inReal * _inReal + _inImg * _inImg; /* Oo */
+      if(R2 > RMax)
+	{
+	  if( dabss(_inImg) > RMax)
+	    S = 0;
+	  else
+	    S = 1.0f / (((0.5f * _inReal) / _inImg) * _inReal + 0.5f * _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 = slnp1m1s(S) * 0.25f;
+	}
+      else
+	{
+	  if(sabss(S) == 1 && sabss(_inReal) <= sAlim)
+	    {
+	      /* |s| >= SLIM  => |z| is inside D+ or D- */
+	      _outImg = _sign(0.5f,_inImg) * ( sLn2 - logf(sabss(_inReal)));
+	    }
+	  else
+	    {
+	      _outImg = 0.25f * logf((powf(_inReal,2) + powf((_inImg + 1.0f),2)) / (powf(_inReal,2) + powf((_inImg - 1.0f),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(sabss(1 - R2) + sabss(_inReal) <= sTol)
+	{/* |b| is very near 1 (and a is near 0)  some cancellation occur in the (next) generic formula */
+	  _outReal = 0.5f * atan2f(2.0f * _inReal, (1.0f - _inImg) * (1.0f + _inImg) - powf(_inReal,2.0f));
+	}
+      else
+	_outReal = 0.5f * atan2f(2.0f * _inReal, 1.0f - R2);
+    }
+
+  return FloatComplex(_outReal, _outImg);
+}
diff --git a/scilab2c/src/c/elementaryFunctions/atan/zatans.c b/scilab2c/src/c/elementaryFunctions/atan/zatans.c
index c511d790..4b8e9640 100644
--- a/scilab2c/src/c/elementaryFunctions/atan/zatans.c
+++ b/scilab2c/src/c/elementaryFunctions/atan/zatans.c
@@ -1,241 +1,242 @@
-/*
- *  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);
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2006-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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);
+}
diff --git a/scilab2c/src/c/elementaryFunctions/lnp1m1/dlnp1m1s.c b/scilab2c/src/c/elementaryFunctions/lnp1m1/dlnp1m1s.c
index 28cdbe0c..7e1759be 100644
--- a/scilab2c/src/c/elementaryFunctions/lnp1m1/dlnp1m1s.c
+++ b/scilab2c/src/c/elementaryFunctions/lnp1m1/dlnp1m1s.c
@@ -1,76 +1,77 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2008-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
- *
- */
-
-#include "lnp1m1.h"
-#include "abs.h"
-
-/*
-	PURPOSE :  Compute   v = log ( (1 + s)/(1 - s) )
-	for small s, this is for |s| < SLIM = 0.20
-
-	ALGORITHM :
-	1/ if |s| is "very small" we use a truncated
-	taylor dvp (by keeping 3 terms) from :
-	2        4          6
-	t = 2 * s * ( 1 + 1/3 s  + 1/5 s  + [ 1/7 s  + ....] )
-	2        4
-	t = 2 * s * ( 1 + 1/3 s  + 1/5 s  + er)
-
-	The limit E until we use this formula may be simply
-	gotten so that the negliged part er is such that :
-	2        4
-	(#) er <= epsm * ( 1 + 1/3 s  + 1/5 s )   for all |s|<= E
-
-	As  er  = 1/7 s^6 + 1/9 s^8 + ...
-	er <= 1/7 * s^6 ( 1 + s^2 + s^4 + ...) = 1/7  s^6/(1-s^2)
-
-	the inequality (#) is forced if :
-
-	1/7  s^6 / (1-s^2)  <= epsm * ( 1 + 1/3 s^2  + 1/5 s^4 )
-
-	s^6 <= 7 epsm * (1 - 2/3 s^2 - 3/15 s^4 - 1/5 s^6)
-
-	So that E is very near (7 epsm)^(1/6) (approximately 3.032d-3):
-
-	2/ For larger |s| we used a minimax polynome :
-
-	yi = s * (2  + d3 s^3 + d5 s^5 .... + d13 s^13 + d15 s^15)
-
-	This polynome was computed (by some remes algorithm) following
-	(*) the sin(x) example (p 39) of the book :
-
-	"ELEMENTARY FUNCTIONS"
-	"Algorithms and implementation"
-	J.M. Muller (Birkhauser)
-
-	(*) without the additionnal raffinement to get the first coefs
-	very near floating point numbers)
-*/
-double dlnp1m1s(double Var)
-{
-	static double D3	= 0.66666666666672679472;
-	static double D5	= 0.39999999996176889299;
-	static double D7	= 0.28571429392829380980;
-	static double D9	= 0.22222138684562683797;
-	static double D11	= 0.18186349187499222459;
-	static double D13	= 0.15250315884469364710;
-	static double D15	= 0.15367270224757008114;
-	static double E		= 3.032E-3;
-	static double C3	= 2.0/3.0;
-	static double C5	= 2.0/5.0;
-
-	double S2 = Var * Var;
-	if( dabss(Var) <= E)
-		return Var * (2 + S2 * (C3 + C5 * S2));
-	else
-		return Var * (2 + S2 * (D3 + S2 * (D5 + S2 * (D7 + S2 * (D9 + S2 * (D11 + S2 * (D13 + S2 * D15)))))));
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2008-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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
+ *
+ */
+
+#include "lnp1m1.h"
+#include "abs.h"
+
+/*
+	PURPOSE :  Compute   v = log ( (1 + s)/(1 - s) )
+	for small s, this is for |s| < SLIM = 0.20
+
+	ALGORITHM :
+	1/ if |s| is "very small" we use a truncated
+	taylor dvp (by keeping 3 terms) from :
+	2        4          6
+	t = 2 * s * ( 1 + 1/3 s  + 1/5 s  + [ 1/7 s  + ....] )
+	2        4
+	t = 2 * s * ( 1 + 1/3 s  + 1/5 s  + er)
+
+	The limit E until we use this formula may be simply
+	gotten so that the negliged part er is such that :
+	2        4
+	(#) er <= epsm * ( 1 + 1/3 s  + 1/5 s )   for all |s|<= E
+
+	As  er  = 1/7 s^6 + 1/9 s^8 + ...
+	er <= 1/7 * s^6 ( 1 + s^2 + s^4 + ...) = 1/7  s^6/(1-s^2)
+
+	the inequality (#) is forced if :
+
+	1/7  s^6 / (1-s^2)  <= epsm * ( 1 + 1/3 s^2  + 1/5 s^4 )
+
+	s^6 <= 7 epsm * (1 - 2/3 s^2 - 3/15 s^4 - 1/5 s^6)
+
+	So that E is very near (7 epsm)^(1/6) (approximately 3.032d-3):
+
+	2/ For larger |s| we used a minimax polynome :
+
+	yi = s * (2  + d3 s^3 + d5 s^5 .... + d13 s^13 + d15 s^15)
+
+	This polynome was computed (by some remes algorithm) following
+	(*) the sin(x) example (p 39) of the book :
+
+	"ELEMENTARY FUNCTIONS"
+	"Algorithms and implementation"
+	J.M. Muller (Birkhauser)
+
+	(*) without the additionnal raffinement to get the first coefs
+	very near floating point numbers)
+*/
+double dlnp1m1s(double Var)
+{
+	static double D3	= 0.66666666666672679472;
+	static double D5	= 0.39999999996176889299;
+	static double D7	= 0.28571429392829380980;
+	static double D9	= 0.22222138684562683797;
+	static double D11	= 0.18186349187499222459;
+	static double D13	= 0.15250315884469364710;
+	static double D15	= 0.15367270224757008114;
+	static double E		= 3.032E-3;
+	static double C3	= 2.0/3.0;
+	static double C5	= 2.0/5.0;
+
+	double S2 = Var * Var;
+	if( dabss(Var) <= E)
+		return Var * (2 + S2 * (C3 + C5 * S2));
+	else
+		return Var * (2 + S2 * (D3 + S2 * (D5 + S2 * (D7 + S2 * (D9 + S2 * (D11 + S2 * (D13 + S2 * D15)))))));
+}
diff --git a/scilab2c/src/c/elementaryFunctions/lnp1m1/slnp1m1s.c b/scilab2c/src/c/elementaryFunctions/lnp1m1/slnp1m1s.c
index 6c991cc0..9940810c 100644
--- a/scilab2c/src/c/elementaryFunctions/lnp1m1/slnp1m1s.c
+++ b/scilab2c/src/c/elementaryFunctions/lnp1m1/slnp1m1s.c
@@ -1,76 +1,77 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2007-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
- *
- */
-
-#include "lnp1m1.h"
-#include "abs.h"
-
-/*
-  PURPOSE :  Compute   v = log ( (1 + s)/(1 - s) )
-  for small s, this is for |s| < SLIM = 0.20
-
-  ALGORITHM :
-  1/ if |s| is "very small" we use a truncated
-  taylor dvp (by keeping 3 terms) from :
-  2        4          6
-  t = 2 * s * ( 1 + 1/3 s  + 1/5 s  + [ 1/7 s  + ....] )
-  2        4
-  t = 2 * s * ( 1 + 1/3 s  + 1/5 s  + er)
-
-  The limit E until we use this formula may be simply
-  gotten so that the negliged part er is such that :
-  2        4
-  (#) er <= epsm * ( 1 + 1/3 s  + 1/5 s )   for all |s|<= E
-
-  As  er  = 1/7 s^6 + 1/9 s^8 + ...
-  er <= 1/7 * s^6 ( 1 + s^2 + s^4 + ...) = 1/7  s^6/(1-s^2)
-
-  the inequality (#) is forced if :
-
-  1/7  s^6 / (1-s^2)  <= epsm * ( 1 + 1/3 s^2  + 1/5 s^4 )
-
-  s^6 <= 7 epsm * (1 - 2/3 s^2 - 3/15 s^4 - 1/5 s^6)
-
-  So that E is very near (7 epsm)^(1/6) (approximately 3.032d-3):
-
-  2/ For larger |s| we used a minimax polynome :
-
-  yi = s * (2  + d3 s^3 + d5 s^5 .... + d13 s^13 + d15 s^15)
-
-  This polynome was computed (by some remes algorithm) following
-  (*) the sin(x) example (p 39) of the book :
-
-  "ELEMENTARY FUNCTIONS"
-  "Algorithms and implementation"
-  J.M. Muller (Birkhauser)
-
-  (*) without the additionnal raffinement to get the first coefs
-  very near floating point numbers)
-*/
-float slnp1m1s(float Var)
-{
-  static float D3	= 0.66666666666672679472f;
-  static float D5	= 0.39999999996176889299f;
-  static float D7	= 0.28571429392829380980f;
-  static float D9	= 0.22222138684562683797f;
-  static float D11	= 0.18186349187499222459f;
-  static float D13	= 0.15250315884469364710f;
-  static float D15	= 0.15367270224757008114f;
-  static float E	= 3.032E-3f;
-  static float C3	= 2.0f/3.0f;
-  static float C5	= 2.0f/5.0f;
-
-  float S2 = Var * Var;
-  if( sabss(Var) <= E)
-    return Var * (2 + S2 * (C3 + C5 * S2));
-  else
-    return Var * (2 + S2 * (D3 + S2 * (D5 + S2 * (D7 + S2 * (D9 + S2 * (D11 + S2 * (D13 + S2 * D15)))))));
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2007-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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
+ *
+ */
+
+#include "lnp1m1.h"
+#include "abs.h"
+
+/*
+  PURPOSE :  Compute   v = log ( (1 + s)/(1 - s) )
+  for small s, this is for |s| < SLIM = 0.20
+
+  ALGORITHM :
+  1/ if |s| is "very small" we use a truncated
+  taylor dvp (by keeping 3 terms) from :
+  2        4          6
+  t = 2 * s * ( 1 + 1/3 s  + 1/5 s  + [ 1/7 s  + ....] )
+  2        4
+  t = 2 * s * ( 1 + 1/3 s  + 1/5 s  + er)
+
+  The limit E until we use this formula may be simply
+  gotten so that the negliged part er is such that :
+  2        4
+  (#) er <= epsm * ( 1 + 1/3 s  + 1/5 s )   for all |s|<= E
+
+  As  er  = 1/7 s^6 + 1/9 s^8 + ...
+  er <= 1/7 * s^6 ( 1 + s^2 + s^4 + ...) = 1/7  s^6/(1-s^2)
+
+  the inequality (#) is forced if :
+
+  1/7  s^6 / (1-s^2)  <= epsm * ( 1 + 1/3 s^2  + 1/5 s^4 )
+
+  s^6 <= 7 epsm * (1 - 2/3 s^2 - 3/15 s^4 - 1/5 s^6)
+
+  So that E is very near (7 epsm)^(1/6) (approximately 3.032d-3):
+
+  2/ For larger |s| we used a minimax polynome :
+
+  yi = s * (2  + d3 s^3 + d5 s^5 .... + d13 s^13 + d15 s^15)
+
+  This polynome was computed (by some remes algorithm) following
+  (*) the sin(x) example (p 39) of the book :
+
+  "ELEMENTARY FUNCTIONS"
+  "Algorithms and implementation"
+  J.M. Muller (Birkhauser)
+
+  (*) without the additionnal raffinement to get the first coefs
+  very near floating point numbers)
+*/
+float slnp1m1s(float Var)
+{
+  static float D3	= 0.66666666666672679472f;
+  static float D5	= 0.39999999996176889299f;
+  static float D7	= 0.28571429392829380980f;
+  static float D9	= 0.22222138684562683797f;
+  static float D11	= 0.18186349187499222459f;
+  static float D13	= 0.15250315884469364710f;
+  static float D15	= 0.15367270224757008114f;
+  static float E	= 3.032E-3f;
+  static float C3	= 2.0f/3.0f;
+  static float C5	= 2.0f/5.0f;
+
+  float S2 = Var * Var;
+  if( sabss(Var) <= E)
+    return Var * (2 + S2 * (C3 + C5 * S2));
+  else
+    return Var * (2 + S2 * (D3 + S2 * (D5 + S2 * (D7 + S2 * (D9 + S2 * (D11 + S2 * (D13 + S2 * D15)))))));
+}
diff --git a/scilab2c/src/c/elementaryFunctions/log/clogs.c b/scilab2c/src/c/elementaryFunctions/log/clogs.c
index d402d1d1..3b9f6911 100644
--- a/scilab2c/src/c/elementaryFunctions/log/clogs.c
+++ b/scilab2c/src/c/elementaryFunctions/log/clogs.c
@@ -1,64 +1,65 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2007-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
- *
- */
-
-#include <math.h>
-#include "log.h"
-#include "lapack.h"
-#include "log1p.h"
-#include "pythag.h"
-
-floatComplex	clogs(floatComplex in) {
-  static float sR2	= 1.41421356237309504f;
-
-  float _RealIn = creals(in);
-  float _ImgIn = cimags(in);
-
-  float _RealOut = 0;
-  float _ImgOut = 0;
-
-  float RMax = (float) getOverflowThreshold();
-  float LInf = sqrtf((float) getUnderflowThreshold());
-  float LSup = sqrtf(0.5f * RMax);
-
-  float AbsReal = fabsf(_RealIn);
-  float AbsImg	= fabsf(_ImgIn);
-
-  _ImgOut = atan2f(_ImgIn, _RealIn);
-
-  if(_ImgIn > _RealIn)
-    {/* switch Real part and Imaginary part */
-      float Temp	= AbsReal;
-      AbsReal		= AbsImg;
-      AbsImg		= Temp;
-    }
-
-  if((0.5 <= AbsReal) && (AbsReal <= sR2))
-    _RealOut = 0.5f * slog1ps((AbsReal - 1.0f) * (AbsReal + 1.0f) + AbsImg * AbsImg);
-  else if(LInf < AbsImg && AbsReal < LSup)
-    _RealOut = 0.5f * slogs(AbsReal * AbsReal + AbsImg * AbsImg);
-  else if(AbsReal > RMax)
-    _RealOut = AbsReal;
-  else
-    {
-      float Temp = spythags(AbsReal, AbsImg);
-      if(Temp <= RMax)
-	{
-	  _RealOut = slogs(Temp);
-	}
-      else /* handle rare spurious overflow with : */
-	{
-	  float Temp2 = AbsImg/AbsReal;
-	  _RealOut = slogs(AbsReal) + 0.5f * slog1ps(Temp2 * Temp2);
-	}
-    }
-  return FloatComplex(_RealOut, _ImgOut);
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2007-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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
+ *
+ */
+
+#include <math.h>
+#include "log.h"
+#include "lapack.h"
+#include "log1p.h"
+#include "pythag.h"
+
+floatComplex	clogs(floatComplex in) {
+  static float sR2	= 1.41421356237309504f;
+
+  float _RealIn = creals(in);
+  float _ImgIn = cimags(in);
+
+  float _RealOut = 0;
+  float _ImgOut = 0;
+
+  float RMax = (float) getOverflowThreshold();
+  float LInf = sqrtf((float) getUnderflowThreshold());
+  float LSup = sqrtf(0.5f * RMax);
+
+  float AbsReal = fabsf(_RealIn);
+  float AbsImg	= fabsf(_ImgIn);
+
+  _ImgOut = atan2f(_ImgIn, _RealIn);
+
+  if(_ImgIn > _RealIn)
+    {/* switch Real part and Imaginary part */
+      float Temp	= AbsReal;
+      AbsReal		= AbsImg;
+      AbsImg		= Temp;
+    }
+
+  if((0.5 <= AbsReal) && (AbsReal <= sR2))
+    _RealOut = 0.5f * slog1ps((AbsReal - 1.0f) * (AbsReal + 1.0f) + AbsImg * AbsImg);
+  else if(LInf < AbsImg && AbsReal < LSup)
+    _RealOut = 0.5f * slogs(AbsReal * AbsReal + AbsImg * AbsImg);
+  else if(AbsReal > RMax)
+    _RealOut = AbsReal;
+  else
+    {
+      float Temp = spythags(AbsReal, AbsImg);
+      if(Temp <= RMax)
+	{
+	  _RealOut = slogs(Temp);
+	}
+      else /* handle rare spurious overflow with : */
+	{
+	  float Temp2 = AbsImg/AbsReal;
+	  _RealOut = slogs(AbsReal) + 0.5f * slog1ps(Temp2 * Temp2);
+	}
+    }
+  return FloatComplex(_RealOut, _ImgOut);
+}
diff --git a/scilab2c/src/c/elementaryFunctions/log/zlogs.c b/scilab2c/src/c/elementaryFunctions/log/zlogs.c
index dd0792a9..e5e9ded1 100644
--- a/scilab2c/src/c/elementaryFunctions/log/zlogs.c
+++ b/scilab2c/src/c/elementaryFunctions/log/zlogs.c
@@ -1,64 +1,65 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2007-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
- *
- */
-
-#include <math.h>
-#include "log.h"
-#include "lapack.h"
-#include "log1p.h"
-#include "pythag.h"
-
-doubleComplex	zlogs(doubleComplex in) {
-  static double sR2	= 1.41421356237309504;
-
-  double _RealIn = zreals(in);
-  double _ImgIn = zimags(in);
-
-  double _RealOut = 0;
-  double _ImgOut = 0;
-
-  double RMax = getOverflowThreshold();
-  double LInf = sqrt(getUnderflowThreshold());
-  double LSup = sqrt(0.5 * RMax);
-
-  double AbsReal = fabs(_RealIn);
-  double AbsImg	= fabs(_ImgIn);
-
-  _ImgOut = atan2(_ImgIn, _RealIn);
-
-  if(_ImgIn > _RealIn)
-    {/* switch Real part and Imaginary part */
-      double Temp	= AbsReal;
-      AbsReal		= AbsImg;
-      AbsImg		= Temp;
-    }
-
-  if((0.5 <= AbsReal) && (AbsReal <= sR2))
-    _RealOut = 0.5 * dlog1ps((AbsReal - 1) * (AbsReal + 1) + AbsImg * AbsImg);
-  else if(LInf < AbsImg && AbsReal < LSup)
-    _RealOut = 0.5 * dlogs(AbsReal * AbsReal + AbsImg * AbsImg);
-  else if(AbsReal > RMax)
-    _RealOut = AbsReal;
-  else
-    {
-      double Temp = dpythags(AbsReal, AbsImg);
-      if(Temp <= RMax)
-	{
-	  _RealOut = dlogs(Temp);
-	}
-      else /* handle rare spurious overflow with : */
-	{
-	  double Temp2 = AbsImg/AbsReal;
-	  _RealOut = dlogs(AbsReal) + 0.5 * dlog1ps(Temp2 * Temp2);
-	}
-    }
-  return DoubleComplex(_RealOut, _ImgOut);
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2007-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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
+ *
+ */
+
+#include <math.h>
+#include "log.h"
+#include "lapack.h"
+#include "log1p.h"
+#include "pythag.h"
+
+doubleComplex	zlogs(doubleComplex in) {
+  static double sR2	= 1.41421356237309504;
+
+  double _RealIn = zreals(in);
+  double _ImgIn = zimags(in);
+
+  double _RealOut = 0;
+  double _ImgOut = 0;
+
+  double RMax = getOverflowThreshold();
+  double LInf = sqrt(getUnderflowThreshold());
+  double LSup = sqrt(0.5 * RMax);
+
+  double AbsReal = fabs(_RealIn);
+  double AbsImg	= fabs(_ImgIn);
+
+  _ImgOut = atan2(_ImgIn, _RealIn);
+
+  if(_ImgIn > _RealIn)
+    {/* switch Real part and Imaginary part */
+      double Temp	= AbsReal;
+      AbsReal		= AbsImg;
+      AbsImg		= Temp;
+    }
+
+  if((0.5 <= AbsReal) && (AbsReal <= sR2))
+    _RealOut = 0.5 * dlog1ps((AbsReal - 1) * (AbsReal + 1) + AbsImg * AbsImg);
+  else if(LInf < AbsImg && AbsReal < LSup)
+    _RealOut = 0.5 * dlogs(AbsReal * AbsReal + AbsImg * AbsImg);
+  else if(AbsReal > RMax)
+    _RealOut = AbsReal;
+  else
+    {
+      double Temp = dpythags(AbsReal, AbsImg);
+      if(Temp <= RMax)
+	{
+	  _RealOut = dlogs(Temp);
+	}
+      else /* handle rare spurious overflow with : */
+	{
+	  double Temp2 = AbsImg/AbsReal;
+	  _RealOut = dlogs(AbsReal) + 0.5 * dlog1ps(Temp2 * Temp2);
+	}
+    }
+  return DoubleComplex(_RealOut, _ImgOut);
+}
diff --git a/scilab2c/src/c/elementaryFunctions/log1p/dlog1ps.c b/scilab2c/src/c/elementaryFunctions/log1p/dlog1ps.c
index 4c055c86..e75a05fd 100644
--- a/scilab2c/src/c/elementaryFunctions/log1p/dlog1ps.c
+++ b/scilab2c/src/c/elementaryFunctions/log1p/dlog1ps.c
@@ -1,33 +1,34 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2008-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
- *
- */
-
-#include "log1p.h"
-#include "log.h"
-#include "lnp1m1.h"
-
-double		dlog1ps(double in) {
-  static double A = -1.0/3.0;
-  static double B = 0.5;
-
-  if(in < -1)
-    {/* got NaN */
-      return (in - in) / (in - in); /* NaN */
-    }
-  else if(A <= in && in <= B)
-    {/* use the function log((1+g)/(1-g)) with g = x/(x + 2) */
-      return dlnp1m1s(in / (in + 2));
-    }
-  else
-    {/* use the standard formula */
-      return dlogs(in + 1);
-    }
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2008-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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
+ *
+ */
+
+#include "log1p.h"
+#include "log.h"
+#include "lnp1m1.h"
+
+double		dlog1ps(double in) {
+  static double A = -1.0/3.0;
+  static double B = 0.5;
+
+  if(in < -1)
+    {/* got NaN */
+      return (in - in) / (in - in); /* NaN */
+    }
+  else if(A <= in && in <= B)
+    {/* use the function log((1+g)/(1-g)) with g = x/(x + 2) */
+      return dlnp1m1s(in / (in + 2));
+    }
+  else
+    {/* use the standard formula */
+      return dlogs(in + 1);
+    }
+}
diff --git a/scilab2c/src/c/elementaryFunctions/log1p/slog1ps.c b/scilab2c/src/c/elementaryFunctions/log1p/slog1ps.c
index 501de6b0..04786524 100644
--- a/scilab2c/src/c/elementaryFunctions/log1p/slog1ps.c
+++ b/scilab2c/src/c/elementaryFunctions/log1p/slog1ps.c
@@ -1,33 +1,35 @@
-/*
- *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
- *  Copyright (C) 2008-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
- *
- */
-
-#include "log1p.h"
-#include "log.h"
-#include "lnp1m1.h"
-
-float		slog1ps(float in) {
-  static double A = -1.0/3.0;
-  static double B = 0.5;
-
-  if(in < -1)
-    {/* got NaN */
-      return (in - in) / (in - in); /* NaN */
-    }
-  else if(A <= in && in <= B)
-    {/* use the function log((1+g)/(1-g)) with g = x/(x + 2) */
-      return slnp1m1s(in / (in + 2));
-    }
-  else
-    {/* use the standard formula */
-      return slogs(in + 1);
-    }
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2008-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+
+ *
+ *  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
+ *
+ */
+
+#include "log1p.h"
+#include "log.h"
+#include "lnp1m1.h"
+
+float		slog1ps(float in) {
+  static double A = -1.0/3.0;
+  static double B = 0.5;
+
+  if(in < -1)
+    {/* got NaN */
+      return (in - in) / (in - in); /* NaN */
+    }
+  else if(A <= in && in <= B)
+    {/* use the function log((1+g)/(1-g)) with g = x/(x + 2) */
+      return slnp1m1s(in / (in + 2));
+    }
+  else
+    {/* use the standard formula */
+      return slogs(in + 1);
+    }
+}
diff --git a/scilab2c/src/c/elementaryFunctions/sqrt/csqrts.c b/scilab2c/src/c/elementaryFunctions/sqrt/csqrts.c
index a2ed819e..a24f9558 100644
--- a/scilab2c/src/c/elementaryFunctions/sqrt/csqrts.c
+++ b/scilab2c/src/c/elementaryFunctions/sqrt/csqrts.c
@@ -1,6 +1,7 @@
 /*
  *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
  *  Copyright (C) 2008-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
  *
  *  This file must be used under the terms of the CeCILL.
  *  This source file is licensed as described in the file COPYING, which
diff --git a/scilab2c/src/c/elementaryFunctions/sqrt/zsqrts.c b/scilab2c/src/c/elementaryFunctions/sqrt/zsqrts.c
index f1a2b05c..3637ddd6 100644
--- a/scilab2c/src/c/elementaryFunctions/sqrt/zsqrts.c
+++ b/scilab2c/src/c/elementaryFunctions/sqrt/zsqrts.c
@@ -1,6 +1,7 @@
 /*
  *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
  *  Copyright (C) 2008-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
  *
  *  This file must be used under the terms of the CeCILL.
  *  This source file is licensed as described in the file COPYING, which
diff --git a/scilab2c/src/c/elementaryFunctions/tan/ztans.c b/scilab2c/src/c/elementaryFunctions/tan/ztans.c
index ff130dcc..761da36b 100644
--- a/scilab2c/src/c/elementaryFunctions/tan/ztans.c
+++ b/scilab2c/src/c/elementaryFunctions/tan/ztans.c
@@ -1,103 +1,104 @@
-/*
- *  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
- *
- */
-
-/*
-  ALGORITHM
-   based on the formula :
-
-                    0.5 sin(2 xr) +  i 0.5 sinh(2 xi)
-   tan(xr + i xi) = ---------------------------------
-                        cos(xr)^2  + sinh(xi)^2
-
-   noting  d = cos(xr)^2 + sinh(xi)^2, we have :
-
-           yr = 0.5 * sin(2 * xr) / d       (1)
-
-           yi = 0.5 * sinh(2 * xi) / d      (2)
-
-   to avoid spurious overflows in computing yi with
-   formula (2) (which results in NaN for yi)
-   we use also the following formula :
-
-           yi = sign(xi)   when |xi| > LIM  (3)
-
-   Explanations for (3) :
-
-   we have d = sinh(xi)^2 ( 1 + (cos(xr)/sinh(xi))^2 ),
-   so when :
-
-      (cos(xr)/sinh(xi))^2 < epsm   ( epsm = max relative error
-                                     for coding a real in a f.p.
-                                     number set F(b,p,emin,emax)
-                                        epsm = 0.5 b^(1-p) )
-   which is forced  when :
-
-       1/sinh(xi)^2 < epsm   (4)
-   <=> |xi| > asinh(1/sqrt(epsm)) (= 19.06... in ieee 754 double)
-
-   sinh(xi)^2 is a good approximation for d (relative to the f.p.
-   arithmetic used) and then yr may be approximate with :
-
-    yr = cosh(xi)/sinh(xi)
-       = sign(xi) (1 + exp(-2 |xi|))/(1 - exp(-2|xi|))
-       = sign(xi) (1 + 2 u + 2 u^2 + 2 u^3 + ...)
-
-   with u = exp(-2 |xi|)). Now when :
-
-       2 exp(-2|xi|) < epsm  (2)
-   <=> |xi| > 0.5 * log(2/epsm) (= 18.71... in ieee 754 double)
-
-   sign(xi)  is a good approximation for yr.
-
-   Constraint (1) is stronger than (2) and we take finaly
-
-   LIM = 1 + log(2/sqrt(epsm))
-
-   (log(2/sqrt(epsm)) being very near asinh(1/sqrt(epsm))
-
-AUTHOR
-   Bruno Pincon <Bruno.Pincon@iecn.u-nancy.fr>
-*/
-
-#include <math.h>
-#include "lapack.h"
-#include "tan.h"
-#include "sqrt.h"
-#include "cos.h"
-#include "sinh.h"
-#include "sin.h"
-#include "log.h"
-#include "abs.h"
-
-#define localSign(x) x >= 0 ? 1.0 : -1.0
-
-doubleComplex		ztans(doubleComplex z) {
-  double Temp = 0;
-  double Lim = 1 + dlogs(2.0 / dsqrts( getRelativeMachinePrecision()));
-  double RealIn = zreals(z);
-  double ImagIn = zimags(z);
-  double RealOut = 0;
-  double ImagOut = 0;
-
-  Temp = pow(dcoss(RealIn), 2) + pow(dsinhs(ImagIn), 2);
-  RealOut = 0.5 * dsins(2 * RealIn) / Temp;
-  if(dabss(ImagIn) < Lim)
-    {
-      ImagOut = 0.5 * dsinhs(2 * ImagIn) / Temp;
-    }
-  else
-    {
-      ImagOut = localSign(ImagIn);
-    }
-
-  return DoubleComplex(RealOut, ImagOut);
-}
+/*
+ *  Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+ *  Copyright (C) 2006-2008 - INRIA - Bruno JOFRET
+ *  Copyright (C) Bruno Pincon
+ *
+ *  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
+ *
+ */
+
+/*
+  ALGORITHM
+   based on the formula :
+
+                    0.5 sin(2 xr) +  i 0.5 sinh(2 xi)
+   tan(xr + i xi) = ---------------------------------
+                        cos(xr)^2  + sinh(xi)^2
+
+   noting  d = cos(xr)^2 + sinh(xi)^2, we have :
+
+           yr = 0.5 * sin(2 * xr) / d       (1)
+
+           yi = 0.5 * sinh(2 * xi) / d      (2)
+
+   to avoid spurious overflows in computing yi with
+   formula (2) (which results in NaN for yi)
+   we use also the following formula :
+
+           yi = sign(xi)   when |xi| > LIM  (3)
+
+   Explanations for (3) :
+
+   we have d = sinh(xi)^2 ( 1 + (cos(xr)/sinh(xi))^2 ),
+   so when :
+
+      (cos(xr)/sinh(xi))^2 < epsm   ( epsm = max relative error
+                                     for coding a real in a f.p.
+                                     number set F(b,p,emin,emax)
+                                        epsm = 0.5 b^(1-p) )
+   which is forced  when :
+
+       1/sinh(xi)^2 < epsm   (4)
+   <=> |xi| > asinh(1/sqrt(epsm)) (= 19.06... in ieee 754 double)
+
+   sinh(xi)^2 is a good approximation for d (relative to the f.p.
+   arithmetic used) and then yr may be approximate with :
+
+    yr = cosh(xi)/sinh(xi)
+       = sign(xi) (1 + exp(-2 |xi|))/(1 - exp(-2|xi|))
+       = sign(xi) (1 + 2 u + 2 u^2 + 2 u^3 + ...)
+
+   with u = exp(-2 |xi|)). Now when :
+
+       2 exp(-2|xi|) < epsm  (2)
+   <=> |xi| > 0.5 * log(2/epsm) (= 18.71... in ieee 754 double)
+
+   sign(xi)  is a good approximation for yr.
+
+   Constraint (1) is stronger than (2) and we take finaly
+
+   LIM = 1 + log(2/sqrt(epsm))
+
+   (log(2/sqrt(epsm)) being very near asinh(1/sqrt(epsm))
+
+AUTHOR
+   Bruno Pincon <Bruno.Pincon@iecn.u-nancy.fr>
+*/
+
+#include <math.h>
+#include "lapack.h"
+#include "tan.h"
+#include "sqrt.h"
+#include "cos.h"
+#include "sinh.h"
+#include "sin.h"
+#include "log.h"
+#include "abs.h"
+
+#define localSign(x) x >= 0 ? 1.0 : -1.0
+
+doubleComplex		ztans(doubleComplex z) {
+  double Temp = 0;
+  double Lim = 1 + dlogs(2.0 / dsqrts( getRelativeMachinePrecision()));
+  double RealIn = zreals(z);
+  double ImagIn = zimags(z);
+  double RealOut = 0;
+  double ImagOut = 0;
+
+  Temp = pow(dcoss(RealIn), 2) + pow(dsinhs(ImagIn), 2);
+  RealOut = 0.5 * dsins(2 * RealIn) / Temp;
+  if(dabss(ImagIn) < Lim)
+    {
+      ImagOut = 0.5 * dsinhs(2 * ImagIn) / Temp;
+    }
+  else
+    {
+      ImagOut = localSign(ImagIn);
+    }
+
+  return DoubleComplex(RealOut, ImagOut);
+}
author	baudin	2010-06-04 12:56:48 +0000
committer	baudin	2010-06-04 12:56:48 +0000
commit	d91785f60376859d9dc0be0739b02ce6b8166bd9 (patch)
tree	c44e2f98fccd1f05bedb9ed9ea2f2dc3bb05417f
parent	a2b9321ab9a83532e3e0250fd370b37b6a70520b (diff)
download	scilab2c-d91785f60376859d9dc0be0739b02ce6b8166bd9.tar.gz scilab2c-d91785f60376859d9dc0be0739b02ce6b8166bd9.tar.bz2 scilab2c-d91785f60376859d9dc0be0739b02ce6b8166bd9.zip