# HG changeset patch
# User Jaroslav Hajek <highegg@gmail.com>
# Date 1236980006 -3600
# Node ID 0a58c4cd1405860ee2b7615e92e8caa6242a6fff
# Parent  f464119ec1656585320380c73de2656170ebd143
optimize element-wise power

diff -r f464119ec165 -r 0a58c4cd1405 src/ChangeLog
--- a/src/ChangeLog	Fri Mar 13 21:26:13 2009 +0100
+++ b/src/ChangeLog	Fri Mar 13 22:33:26 2009 +0100
@@ -1,3 +1,9 @@
+2009-03-13  Jaroslav Hajek  <highegg@gmail.com>
+
+	* xpow.cc (xpow (const NDArray&, double), xpow (const ComplexNDArray&, double),
+	xpow (const FloatNDArray&, float), xpow (const FloatComplexNDArray&,
+	float)): Use xisint for testing ints. Optimize w.r.t int exponents.
+
 2009-03-13  John W. Eaton  <jwe@octave.org>
 
 	* parse.y (maybe_warn_assign_as_truth_value, make_binary_op,
diff -r f464119ec165 -r 0a58c4cd1405 src/xpow.cc
--- a/src/xpow.cc	Fri Mar 13 21:26:13 2009 +0100
+++ b/src/xpow.cc	Fri Mar 13 22:33:26 2009 +0100
@@ -79,7 +79,7 @@
 {
   double retval;
 
-  if (a < 0.0 && static_cast<int> (b) != b)
+  if (a < 0.0 && ! xisint (b))
     {
       Complex atmp (a);
 
@@ -700,7 +700,7 @@
   octave_idx_type nr = a.rows ();
   octave_idx_type nc = a.cols ();
 
-  if (static_cast<int> (b) != b && a.any_element_is_negative ())
+  if (! xisint (b) && a.any_element_is_negative ())
     {
       ComplexMatrix result (nr, nc);
 
@@ -1084,30 +1084,53 @@
 {
   octave_value retval;
 
-  if (static_cast<int> (b) != b && a.any_element_is_negative ())
+  if (! xisint (b))
     {
-      ComplexNDArray result (a.dims ());
-
-      for (octave_idx_type i = 0; i < a.length (); i++)
-	{
-	  OCTAVE_QUIT;
-
-	  Complex atmp (a (i));
-
-	  result(i) = std::pow (atmp, b);
-	}
-
-      retval = result;
+      if (a.any_element_is_negative ())
+        {
+          ComplexNDArray result (a.dims ());
+
+          for (octave_idx_type i = 0; i < a.length (); i++)
+            {
+              OCTAVE_QUIT;
+
+              Complex atmp (a (i));
+
+              result(i) = std::pow (atmp, b);
+            }
+
+          retval = result;
+        }
+      else
+        {
+          NDArray result (a.dims ());
+          for (octave_idx_type i = 0; i < a.length (); i++)
+            {
+              OCTAVE_QUIT;
+              result(i) = std::pow (a(i), b);
+            }
+
+          retval = result;
+        }
     }
   else
     {
       NDArray result (a.dims ());
 
-      for (octave_idx_type i = 0; i < a.length (); i++)
-	{
-	  OCTAVE_QUIT;
-	  result(i) = std::pow (a(i), b);
-	}
+      int ib = static_cast<int> (b);
+      if (ib == 2)
+        {
+          for (octave_idx_type i = 0; i < a.length (); i++)
+            result.xelem (i) = a(i) * a(i);
+        }
+      else
+        {
+          for (octave_idx_type i = 0; i < a.length (); i++)
+            {
+              OCTAVE_QUIT;
+              result(i) = std::pow (a(i), ib);
+            }
+        }
 
       retval = result;
     }
@@ -1371,7 +1394,7 @@
 {
   float retval;
 
-  if (a < 0.0 && static_cast<int> (b) != b)
+  if (a < 0.0 && ! xisint (b))
     {
       FloatComplex atmp (a);
 
@@ -1979,7 +2002,7 @@
   octave_idx_type nr = a.rows ();
   octave_idx_type nc = a.cols ();
 
-  if (static_cast<int> (b) != b && a.any_element_is_negative ())
+  if (! xisint (b) && a.any_element_is_negative ())
     {
       FloatComplexMatrix result (nr, nc);
 
@@ -2363,30 +2386,53 @@
 {
   octave_value retval;
 
-  if (static_cast<int> (b) != b && a.any_element_is_negative ())
+  if (! xisint (b))
     {
-      FloatComplexNDArray result (a.dims ());
-
-      for (octave_idx_type i = 0; i < a.length (); i++)
-	{
-	  OCTAVE_QUIT;
-
-	  FloatComplex atmp (a (i));
-
-	  result(i) = std::pow (atmp, b);
-	}
-
-      retval = result;
+      if (a.any_element_is_negative ())
+        {
+          FloatComplexNDArray result (a.dims ());
+
+          for (octave_idx_type i = 0; i < a.length (); i++)
+            {
+              OCTAVE_QUIT;
+
+              FloatComplex atmp (a (i));
+
+              result(i) = std::pow (atmp, b);
+            }
+
+          retval = result;
+        }
+      else
+        {
+          FloatNDArray result (a.dims ());
+          for (octave_idx_type i = 0; i < a.length (); i++)
+            {
+              OCTAVE_QUIT;
+              result(i) = std::pow (a(i), b);
+            }
+
+          retval = result;
+        }
     }
   else
     {
       FloatNDArray result (a.dims ());
 
-      for (octave_idx_type i = 0; i < a.length (); i++)
-	{
-	  OCTAVE_QUIT;
-	  result(i) = std::pow (a(i), b);
-	}
+      int ib = static_cast<int> (b);
+      if (ib == 2)
+        {
+          for (octave_idx_type i = 0; i < a.length (); i++)
+            result.xelem (i) = a(i) * a(i);
+        }
+      else
+        {
+          for (octave_idx_type i = 0; i < a.length (); i++)
+            {
+              OCTAVE_QUIT;
+              result(i) = std::pow (a(i), ib);
+            }
+        }
 
       retval = result;
     }