build with visual studio (dynamic libraries)

1 files changed, 110 insertions, 105 deletions

diff --git a/src/elementaryFunctions/sqrt/csqrts.c b/src/elementaryFunctions/sqrt/csqrts.c
index f06d2dd1..a2ed819e 100644
--- a/src/elementaryFunctions/sqrt/csqrts.c
+++ b/src/elementaryFunctions/sqrt/csqrts.c
@@ -1,105 +1,110 @@
-/*
- *  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 <math.h>
-#include "sqrt.h"
-#include "lapack.h"
-#include "abs.h"
-#include "sign.h"
-#include "pythag.h"
-
-#define _sign(a, b) b >=0 ? a : -a
-
-floatComplex	csqrts(floatComplex in) {
-  float RMax =  (float) getOverflowThreshold();
-  float BRMin = 2.0f * (float) getUnderflowThreshold();
-
-  float RealIn = creals(in);
-  float ImgIn = cimags(in);
-
-  float RealOut = 0;
-  float ImgOut = 0;
-
-  if(RealIn == 0)
-    {/* pure imaginary case */
-      if(dabss(ImgIn >= BRMin))
-	RealOut = ssqrts(0.5f * sabss(ImgIn));
-      else
-	RealOut = ssqrts(sabss(ImgIn)) * ssqrts(0.5);
-
-      ImgOut = _sign(1, ImgIn) * RealOut;
-    }
-  else if( sabss(RealIn) <= RMax && sabss(ImgIn) <= RMax)
-    {/* standard case : a (not zero) and b are finite */
-      float Temp = ssqrts(2.0f * (sabss(RealIn) + spythags(RealIn, ImgIn)));
-      /* overflow test */
-      if(Temp > RMax)
-	{/*     handle (spurious) overflow by scaling a and b */
-	  float RealTemp = RealIn / 16.0f;
-	  float ImgTemp = ImgIn / 16.0f;
-	  Temp = ssqrts(2.0f * (sabss(RealIn) + spythags(RealIn, ImgTemp)));
-	  if(RealTemp >= 0)
-	    {
-	      RealOut	= 2 * Temp;
-	      ImgOut	= 4 * ImgTemp / Temp;
-	    }
-	  else
-	    {
-	      RealOut	= 4 * sabss(ImgIn) / Temp;
-	      ImgOut	= _sign(2, ImgIn) * Temp;
-	    }
-	}
-      else if(RealIn >= 0) /* classic switch to get the stable formulas */
-	{
-	  RealOut	= 0.5f * Temp;
-	  ImgOut	= ImgIn / Temp;
-	}
-      else
-	{
-	  RealOut	= sabss(ImgIn) / Temp;
-	  ImgOut	= (_sign(0.5f, ImgIn)) * Temp;
-	}
-    }
-  else
-    {
-      /*
-      //Here we treat the special cases where a and b are +- 00 or NaN.
-      //The following is the treatment recommended by the C99 standard
-      //with the simplification of returning NaN + i NaN if the
-      //the real part or the imaginary part is NaN (C99 recommends
-      //something more complicated)
-      */
-
-      if(isnan(RealIn) == 1 || isnan(ImgIn) == 1)
-	{/* got NaN + i NaN */
-	  RealOut	= RealIn + ImgIn;
-	  ImgOut	= RealOut;
-	}
-      else if( dabss(ImgIn) > RMax)
-	{/* case a +- i oo -> result must be +oo +- i oo  for all a (finite or not) */
-	  RealOut	= sabss(ImgIn);
-	  ImgOut	= ImgIn;
-	}
-      else if(RealIn < -RMax)
-	{/* here a is -Inf and b is finite */
-	  RealOut	= 0;
-	  ImgOut	= _sign(1, ImgIn) * sabss(RealIn);
-	}
-      else
-	{/* here a is +Inf and b is finite */
-	  RealOut	= RealIn;
-	  ImgOut	= 0;
-	}
-    }
-
-  return FloatComplex(RealOut, ImgOut);
-}
+/*

+ *  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 <math.h>

+#include "sqrt.h"

+#include "lapack.h"

+#include "abs.h"

+#include "sign.h"

+#include "pythag.h"

+

+#ifdef _MSC_VER

+#include <float.h>

+#define isnan(x) _isnan((double)x)

+#endif

+

+#define _sign(a, b) b >=0 ? a : -a

+

+floatComplex	csqrts(floatComplex in) {

+  float RMax =  (float) getOverflowThreshold();

+  float BRMin = 2.0f * (float) getUnderflowThreshold();

+

+  float RealIn = creals(in);

+  float ImgIn = cimags(in);

+

+  float RealOut = 0;

+  float ImgOut = 0;

+

+  if(RealIn == 0)

+    {/* pure imaginary case */

+      if(dabss(ImgIn >= BRMin))

+	RealOut = ssqrts(0.5f * sabss(ImgIn));

+      else

+	RealOut = ssqrts(sabss(ImgIn)) * ssqrts(0.5);

+

+      ImgOut = _sign(1, ImgIn) * RealOut;

+    }

+  else if( sabss(RealIn) <= RMax && sabss(ImgIn) <= RMax)

+    {/* standard case : a (not zero) and b are finite */

+      float Temp = ssqrts(2.0f * (sabss(RealIn) + spythags(RealIn, ImgIn)));

+      /* overflow test */

+      if(Temp > RMax)

+	{/*     handle (spurious) overflow by scaling a and b */

+	  float RealTemp = RealIn / 16.0f;

+	  float ImgTemp = ImgIn / 16.0f;

+	  Temp = ssqrts(2.0f * (sabss(RealIn) + spythags(RealIn, ImgTemp)));

+	  if(RealTemp >= 0)

+	    {

+	      RealOut	= 2 * Temp;

+	      ImgOut	= 4 * ImgTemp / Temp;

+	    }

+	  else

+	    {

+	      RealOut	= 4 * sabss(ImgIn) / Temp;

+	      ImgOut	= _sign(2, ImgIn) * Temp;

+	    }

+	}

+      else if(RealIn >= 0) /* classic switch to get the stable formulas */

+	{

+	  RealOut	= 0.5f * Temp;

+	  ImgOut	= ImgIn / Temp;

+	}

+      else

+	{

+	  RealOut	= sabss(ImgIn) / Temp;

+	  ImgOut	= (_sign(0.5f, ImgIn)) * Temp;

+	}

+    }

+  else

+    {

+      /*

+      //Here we treat the special cases where a and b are +- 00 or NaN.

+      //The following is the treatment recommended by the C99 standard

+      //with the simplification of returning NaN + i NaN if the

+      //the real part or the imaginary part is NaN (C99 recommends

+      //something more complicated)

+      */

+

+      if(isnan(RealIn) == 1 || isnan(ImgIn) == 1)

+	{/* got NaN + i NaN */

+	  RealOut	= RealIn + ImgIn;

+	  ImgOut	= RealOut;

+	}

+      else if( dabss(ImgIn) > RMax)

+	{/* case a +- i oo -> result must be +oo +- i oo  for all a (finite or not) */

+	  RealOut	= sabss(ImgIn);

+	  ImgOut	= ImgIn;

+	}

+      else if(RealIn < -RMax)

+	{/* here a is -Inf and b is finite */

+	  RealOut	= 0;

+	  ImgOut	= _sign(1, ImgIn) * sabss(RealIn);

+	}

+      else

+	{/* here a is +Inf and b is finite */

+	  RealOut	= RealIn;

+	  ImgOut	= 0;

+	}

+    }

+

+  return FloatComplex(RealOut, ImgOut);

+}