diff --git a/prism/src/mtbdd/PM_JOR.cc b/prism/src/mtbdd/PM_JOR.cc
index 8ef964ab..44a33fef 100644
--- a/prism/src/mtbdd/PM_JOR.cc
+++ b/prism/src/mtbdd/PM_JOR.cc
@@ -34,6 +34,8 @@
 #include "PrismMTBDDGlob.h"
 #include "jnipointer.h"
 #include "prism.h"
+#include "ExportIterations.h"
+#include <memory>
 
 //------------------------------------------------------------------------------
 
@@ -122,6 +124,16 @@ jdouble omega		// omega (jor parameter)
 		sol = DD_PermuteVariables(ddman, sol, rvars, cvars, num_rvars);
 	}
 	
+	std::unique_ptr<ExportIterations> iterationExport;
+	if (PM_GetFlagExportIterations()) {
+		std::string title("PM_JOR (");
+		title += (omega == 1.0)?"Jacobi": ("JOR omega=" + std::to_string(omega));
+		title += ")";
+
+		iterationExport.reset(new ExportIterations(title.c_str()));
+		iterationExport->exportVector(sol, rvars, num_rvars, odd, 0);
+	}
+
 	// get setup time
 	stop = util_cpu_time();
 	time_for_setup = (double)(stop - start2)/1000;
@@ -149,7 +161,10 @@ jdouble omega		// omega (jor parameter)
 			Cudd_Ref(sol);
 			tmp = DD_Apply(ddman, APPLY_PLUS, tmp, DD_Apply(ddman, APPLY_TIMES, sol, DD_Constant(ddman, 1-omega)));
 		}
-		
+
+		if (iterationExport)
+			iterationExport->exportVector(tmp, rvars, num_rvars, odd, 0);
+
 		// check convergence
 		switch (term_crit) {
 		case TERM_CRIT_ABSOLUTE:
diff --git a/prism/src/mtbdd/PM_NondetReachReward.cc b/prism/src/mtbdd/PM_NondetReachReward.cc
index f2f72bcb..05c10e94 100644
--- a/prism/src/mtbdd/PM_NondetReachReward.cc
+++ b/prism/src/mtbdd/PM_NondetReachReward.cc
@@ -34,6 +34,8 @@
 #include "PrismMTBDDGlob.h"
 #include "jnipointer.h"
 #include "prism.h"
+#include "ExportIterations.h"
+#include <memory>
 
 //------------------------------------------------------------------------------
 
@@ -121,7 +123,13 @@ jboolean min		// min or max probabilities (true = min, false = max)
 	// print memory usage
 	i = DD_GetNumNodes(ddman, a);
 	PM_PrintToMainLog(env, "[nodes=%d] [%.1f Kb]\n", i, i*20.0/1024.0);
-	
+
+	std::unique_ptr<ExportIterations> iterationExport;
+	if (PM_GetFlagExportIterations()) {
+		iterationExport.reset(new ExportIterations("PM_NondetReachReward"));
+		iterationExport->exportVector(sol, rvars, num_rvars, odd, 0);
+	}
+
 	// get setup time
 	stop = util_cpu_time();
 	time_for_setup = (double)(stop - start2)/1000;
@@ -163,7 +171,10 @@ jboolean min		// min or max probabilities (true = min, false = max)
 		// put infinities (for 'inf' states) back into into solution vector		
 		Cudd_Ref(inf);
 		tmp = DD_ITE(ddman, inf, DD_PlusInfinity(ddman), tmp);
-		
+
+		if (iterationExport)
+			iterationExport->exportVector(tmp, rvars, num_rvars, odd, 0);
+
 		// check convergence
 		switch (term_crit) {
 		case TERM_CRIT_ABSOLUTE:
diff --git a/prism/src/mtbdd/PM_NondetUntil.cc b/prism/src/mtbdd/PM_NondetUntil.cc
index 099d3a57..e540a544 100644
--- a/prism/src/mtbdd/PM_NondetUntil.cc
+++ b/prism/src/mtbdd/PM_NondetUntil.cc
@@ -34,6 +34,8 @@
 #include "PrismMTBDDGlob.h"
 #include "jnipointer.h"
 #include "prism.h"
+#include "ExportIterations.h"
+#include <memory>
 
 //------------------------------------------------------------------------------
 
@@ -90,6 +92,12 @@ jboolean min		// min or max probabilities (true = min, false = max)
 	Cudd_Ref(yes);
 	sol = yes;
 
+	std::unique_ptr<ExportIterations> iterationExport;
+	if (PM_GetFlagExportIterations()) {
+		iterationExport.reset(new ExportIterations("PM_NondetUntil"));
+		iterationExport->exportVector(sol, rvars, num_rvars, odd, 0);
+	}
+
 	// get setup time
 	stop = util_cpu_time();
 	time_for_setup = (double)(stop - start2)/1000;
@@ -136,7 +144,10 @@ jboolean min		// min or max probabilities (true = min, false = max)
 		// put 1s (for 'yes' states) back into into solution vector		
 		Cudd_Ref(yes);
 		tmp = DD_Apply(ddman, APPLY_MAX, tmp, yes);
-		
+
+		if (iterationExport)
+			iterationExport->exportVector(tmp, rvars, num_rvars, odd, 0);
+
 		// check convergence
 		switch (term_crit) {
 		case TERM_CRIT_ABSOLUTE:
diff --git a/prism/src/mtbdd/PM_Power.cc b/prism/src/mtbdd/PM_Power.cc
index 1190c059..222e64d2 100644
--- a/prism/src/mtbdd/PM_Power.cc
+++ b/prism/src/mtbdd/PM_Power.cc
@@ -34,6 +34,8 @@
 #include "PrismMTBDDGlob.h"
 #include "jnipointer.h"
 #include "prism.h"
+#include "ExportIterations.h"
+#include <memory>
 
 //------------------------------------------------------------------------------
 
@@ -94,7 +96,13 @@ jboolean transpose	// transpose A? (i.e. solve xA=b not Ax=b?)
 	if (transpose) {
 		sol = DD_PermuteVariables(ddman, sol, rvars, cvars, num_rvars);
 	}
-	
+
+	std::unique_ptr<ExportIterations> iterationExport;
+	if (PM_GetFlagExportIterations()) {
+		iterationExport.reset(new ExportIterations("PM_Power"));
+		iterationExport->exportVector(sol, rvars, num_rvars, odd, 0);
+	}
+
 	// get setup time
 	stop = util_cpu_time();
 	time_for_setup = (double)(stop - start2)/1000;
@@ -118,6 +126,9 @@ jboolean transpose	// transpose A? (i.e. solve xA=b not Ax=b?)
 		Cudd_Ref(b);
 		tmp = DD_Apply(ddman, APPLY_PLUS, tmp, b);
 		
+		if (iterationExport)
+			iterationExport->exportVector(tmp, rvars, num_rvars, odd, 0);
+
 		// check convergence
 		switch (term_crit) {
 		case TERM_CRIT_ABSOLUTE: