diff --git a/python/modules/py_vector.cpp b/python/modules/py_vector.cpp
index 017591715..78b481b5e 100644
--- a/python/modules/py_vector.cpp
+++ b/python/modules/py_vector.cpp
@@ -149,6 +149,7 @@ py::class_<Type> bind_eigen_2(py::module &m, const char *name,
         .def("mean", [](const Type &m) {return m.mean();})
 
         .def("sum", [](const Type &m) {return m.sum();})
+        .def("square", [](const Type &m) {return m.array().square();})
         .def("prod", [](const Type &m) {return m.prod();})
         .def("trace", [](const Type &m) {return m.trace();})
         .def("norm", [](const Type &m) {return m.norm();})
diff --git a/tutorial/701_Statistics/main.cpp b/tutorial/701_Statistics/main.cpp
index e7f87b84d..ebbf03251 100755
--- a/tutorial/701_Statistics/main.cpp
+++ b/tutorial/701_Statistics/main.cpp
@@ -36,7 +36,12 @@ int main(int argc, char *argv[])
   double area_avg   = area.mean();
   double area_min   = area.minCoeff() / area_avg;
   double area_max   = area.maxCoeff() / area_avg;
-  double area_sigma = sqrt(((area.array()-area_avg)/area_avg).square().mean());
+  auto area_ns = ((area.array()-area_avg)/area_avg);
+
+  std::cout << area_ns(0) << std::endl;
+//  area_ns.square();
+  std::cout << area_ns(0) << std::endl;
+  double area_sigma = sqrt(area_ns.square().mean());
 
   printf("Areas (Min/Max)/Avg_Area Sigma: \n%.2f/%.2f (%.2f)\n",
     area_min,area_max,area_sigma);