diff --git a/Code/Mantid/Framework/CurveFitting/test/EvaluateFunctionTest.h b/Code/Mantid/Framework/CurveFitting/test/EvaluateFunctionTest.h
index 6815af212f12f4b6c658681899dcc617d6078ec9..6e83db7e8d7a48353b02368779116c2378f751a9 100644
--- a/Code/Mantid/Framework/CurveFitting/test/EvaluateFunctionTest.h
+++ b/Code/Mantid/Framework/CurveFitting/test/EvaluateFunctionTest.h
@@ -4,11 +4,14 @@
 #include <cxxtest/TestSuite.h>
 
 #include "MantidCurveFitting/EvaluateFunction.h"
+#include "MantidAPI/AlgorithmManager.h"
 #include "MantidAPI/AnalysisDataService.h"
 #include "MantidAPI/FunctionFactory.h"
 #include "MantidAPI/FunctionDomain1D.h"
 #include "MantidAPI/FunctionValues.h"
 #include "MantidAPI/IFunction1D.h"
+#include "MantidAPI/IMDHistoWorkspace.h"
+#include "MantidAPI/IMDIterator.h"
 #include "MantidAPI/MatrixWorkspace.h"
 #include "MantidKernel/EmptyValues.h"
 
@@ -78,6 +81,61 @@ public:
     tester.checkResult();
   }
 
+  void test_MD_histo() {
+    int nx = 5;
+    int ny = 6;
+    std::vector<double> signal(nx*ny);
+    std::vector<double> extents(4);
+    extents[0] = -3; extents[1] = 3;
+    extents[2] = -3; extents[3] = 3;
+    std::vector<int> nBins(2);
+    nBins[0] = nx;
+    nBins[1] = ny;
+
+    auto alg = AlgorithmManager::Instance().create("CreateMDHistoWorkspace");
+    alg->initialize();
+    alg->setProperty("Dimensionality",2);
+    alg->setProperty("SignalInput",signal);
+    alg->setProperty("ErrorInput",signal);
+    alg->setProperty("Extents",extents);
+    alg->setProperty("NumberOfBins",nBins);
+    alg->setProperty("Names","x,y");
+    alg->setProperty("Units","U,V");
+    alg->setProperty("OutputWorkspace","EvaluateFunction_inWS");
+    alg->execute();
+    auto inWS = AnalysisDataService::Instance().retrieve("EvaluateFunction_inWS");
+    TS_ASSERT(inWS);
+
+    alg = AlgorithmManager::Instance().create("EvaluateFunction");
+    alg->initialize();
+    alg->setPropertyValue("Function", "name=UserFunctionMD,Formula=sin(x)*sin(y)");
+    alg->setProperty("InputWorkspace", "EvaluateFunction_inWS");
+    alg->setProperty("OutputWorkspace", "EvaluateFunction_outWS");
+    alg->execute();
+    auto outWS = AnalysisDataService::Instance().retrieve("EvaluateFunction_outWS");
+    TS_ASSERT(outWS);
+
+    auto mdws =  AnalysisDataService::Instance().retrieveWS<IMDHistoWorkspace>("EvaluateFunction_outWS");
+    TS_ASSERT(mdws);
+
+    auto iter = mdws->createIterator();
+    TS_ASSERT(iter);
+    if (!iter)
+      return;
+    do {
+      auto xy = iter->getCenter();
+      auto signal = iter->getSignal();
+      double value = sin(xy[0])*sin(xy[1]);
+      if (value == 0.0) {
+        TS_ASSERT_DELTA(signal, 0.0, 1e-14);
+      } else {
+        // Precision is lost due to the use of floats in MD workspaces.
+        TS_ASSERT_DELTA((signal-value)/value, 0.0, 1e-6);
+      }
+    } while (iter->next());
+    delete iter;
+  }
+
 private:
 
   class Tester1D {
@@ -200,7 +258,6 @@ private:
       setDefaultXRange();
       TS_ASSERT_DIFFERS( nData, 0 );
       auto &Y = outputWorkspace->readY(1);
-      //TS_ASSERT_EQUALS(y.size(), Y.size());
       size_t j = 0;
       for(size_t i = 0; i < nData; ++i) {
         if (xValues[i] <= StartX || xValues[i] >= EndX ) continue;