diff --git a/matlab/disp_dr.m b/matlab/disp_dr.m
index 69826e2e1c4c14cd7cc146263ef629ca18852827..26c06ba2bad0cd7f367cf0de09fb9b23f84be7e7 100644
--- a/matlab/disp_dr.m
+++ b/matlab/disp_dr.m
@@ -24,14 +24,18 @@ function disp_dr(dr,order,var_list)
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
 
-global M_
+global M_ options_
 
 nx =size(dr.ghx,2);
 nu =size(dr.ghu,2);
-k = find(dr.kstate(:,2) <= M_.maximum_lag+1);
-klag = dr.kstate(k,[1 2]);
-
-k1 = dr.order_var;
+if options_.block
+    k = dr.npred + dr.nboth;
+    k1 = M_.block_structure.variable_reordered;
+else
+    k = find(dr.kstate(:,2) <= M_.maximum_lag+1);
+    klag = dr.kstate(k,[1 2]);
+    k1 = dr.order_var;
+end;
 
 if size(var_list,1) == 0
     var_list = M_.endo_names(1:M_.orig_endo_nbr, :);
@@ -98,7 +102,11 @@ end
 %
 for k=1:nx
     flag = 0;
-    str1 = subst_auxvar(k1(klag(k,1)),klag(k,2)-M_.maximum_lag-2);
+    if options_.block
+        str1 = subst_auxvar(k1(dr.nstatic+k),-1);
+    else
+        str1 = subst_auxvar(k1(klag(k,1)),klag(k,2)-M_.maximum_lag-2);
+    end;
     str = sprintf('%-20s',str1);
     for i=1:nvar
         x = dr.ghx(ivar(i),k);
diff --git a/matlab/simult_.m b/matlab/simult_.m
index aab1ad97216853fa1cf45aebbe82244ddb9bc059..b48db083e3d8cc51d1673c59f90d81add01dad10 100644
--- a/matlab/simult_.m
+++ b/matlab/simult_.m
@@ -67,32 +67,39 @@ if options_.k_order_solver% Call dynare++ routines.
     mexErrCheck('dynare_simul_', err);
     y_(dr.order_var,:) = y_;
 else
-    k2 = dr.kstate(find(dr.kstate(:,2) <= 1+M_.maximum_lag),[1 2]);
-    k2 = k2(:,1)+(M_.maximum_lag+1-k2(:,2))*M_.endo_nbr;
+    if options_.block
+        k2 = dr.nstatic +1 : dr.nstatic + dr.npred+dr.nboth;
+        order_var = M_.block_structure.variable_reordered;        
+    else
+        k2 = dr.kstate(find(dr.kstate(:,2) <= M_.maximum_lag+1),[1 2]);
+        k2 = k2(:,1)+(M_.maximum_lag+1-k2(:,2))*M_.endo_nbr;
+        order_var = dr.order_var;
+    end;
+    
     switch iorder
       case 1
         if isempty(dr.ghu)% For (linearized) deterministic models.
             for i = 2:iter+M_.maximum_lag
-                yhat = y_(dr.order_var(k2),i-1);
-                y_(dr.order_var,i) = dr.ghx*yhat;
+                yhat = y_(order_var(k2),i-1);
+                y_(order_var,i) = dr.ghx*yhat;
             end
         elseif isempty(dr.ghx)% For (linearized) purely forward variables (no state variables).
             y_(dr.order_var,:) = dr.ghu*transpose(ex_);
         else
             epsilon = dr.ghu*transpose(ex_);
             for i = 2:iter+M_.maximum_lag
-                yhat = y_(dr.order_var(k2),i-1);
-                y_(dr.order_var,i) = dr.ghx*yhat + epsilon(:,i-1);
+                yhat = y_(order_var(k2),i-1);
+                y_(order_var,i) = dr.ghx*yhat + epsilon(:,i-1);
             end
         end
         y_ = bsxfun(@plus,y_,dr.ys);
       case 2
-        constant = dr.ys(dr.order_var)+.5*dr.ghs2;
+        constant = dr.ys(order_var)+.5*dr.ghs2;
         if options_.pruning
             y__ = y0;
             for i = 2:iter+M_.maximum_lag
-                yhat1 = y__(dr.order_var(k2))-dr.ys(dr.order_var(k2));
-                yhat2 = y_(dr.order_var(k2),i-1)-dr.ys(dr.order_var(k2));
+                yhat1 = y__(order_var(k2))-dr.ys(order_var(k2));
+                yhat2 = y_(order_var(k2),i-1)-dr.ys(order_var(k2));
                 epsilon = ex_(i-1,:)';
                 [err, abcOut1] = A_times_B_kronecker_C(.5*dr.ghxx,yhat1,options_.threads.kronecker.A_times_B_kronecker_C);
                 mexErrCheck('A_times_B_kronecker_C', err);
@@ -100,13 +107,13 @@ else
                 mexErrCheck('A_times_B_kronecker_C', err);
                 [err, abcOut3] = A_times_B_kronecker_C(dr.ghxu,yhat1,epsilon,options_.threads.kronecker.A_times_B_kronecker_C);
                 mexErrCheck('A_times_B_kronecker_C', err);
-                y_(dr.order_var,i) = constant + dr.ghx*yhat2 + dr.ghu*epsilon ...
+                y_(order_var,i) = constant + dr.ghx*yhat2 + dr.ghu*epsilon ...
                     + abcOut1 + abcOut2 + abcOut3;
-                y__(dr.order_var) = dr.ys(dr.order_var) + dr.ghx*yhat1 + dr.ghu*epsilon;
+                y__(order_var) = dr.ys(order_var) + dr.ghx*yhat1 + dr.ghu*epsilon;
             end
         else
             for i = 2:iter+M_.maximum_lag
-                yhat = y_(dr.order_var(k2),i-1)-dr.ys(dr.order_var(k2));
+                yhat = y_(order_var(k2),i-1)-dr.ys(order_var(k2));
                 epsilon = ex_(i-1,:)';
                 [err, abcOut1] = A_times_B_kronecker_C(.5*dr.ghxx,yhat,options_.threads.kronecker.A_times_B_kronecker_C);
                 mexErrCheck('A_times_B_kronecker_C', err);
diff --git a/matlab/stoch_simul.m b/matlab/stoch_simul.m
index b13ee6a74b051c380a54e081d0ed5c0a514096e2..615b7cd256274670887bae3ca0bcf2f4d4fe849a 100644
--- a/matlab/stoch_simul.m
+++ b/matlab/stoch_simul.m
@@ -74,10 +74,15 @@ if ~options_.noprint
     disp(' ')
     disp(['  Number of variables:         ' int2str(M_.endo_nbr)])
     disp(['  Number of stochastic shocks: ' int2str(M_.exo_nbr)])
-    disp(['  Number of state variables:   ' ...
-          int2str(length(find(oo_.dr.kstate(:,2) <= M_.maximum_lag+1)))])
-    disp(['  Number of jumpers:           ' ...
-          int2str(length(find(oo_.dr.kstate(:,2) == M_.maximum_lag+2)))])
+    if (options_.block)
+        disp(['  Number of state variables:   ' int2str(oo_.dr.npred+oo_.dr.nboth)])
+        disp(['  Number of jumpers:           ' int2str(oo_.dr.nfwrd+oo_.dr.nboth)])
+    else
+        disp(['  Number of state variables:   ' ...
+              int2str(length(find(oo_.dr.kstate(:,2) <= M_.maximum_lag+1)))])
+        disp(['  Number of jumpers:           ' ...
+              int2str(length(find(oo_.dr.kstate(:,2) == M_.maximum_lag+2)))])
+    end;
     disp(['  Number of static variables:  ' int2str(oo_.dr.nstatic)])
     my_title='MATRIX OF COVARIANCE OF EXOGENOUS SHOCKS';
     labels = deblank(M_.exo_names);