Added a routine for setting the random number generator and the seed,...

Added a routine for setting the random number generator and the seed, independantly of the version of matlab. Removed the calls to rand('state',x) in simult.m and added a test file. Still need to adapt the routines related to mcmc and swz.

matlab/global_initialization.m

@@ -307,5 +307,9 @@ options_.ms.draws_nbr_mean_var_estimate = 200000;
 options_.ms.draws_nbr_modified_harmonic_mean = 1000000;
 options_.ms.thinning_factor = 1;
 options_.ms.dirichlet_scale = [1.0 1.5 2.0];
+
 % initialize persistent variables in priordens()
 priordens([],[],[],[],[],[],1);
+
+% Set dynare random generator and seed.
+set_dynare_seed('default');
\ No newline at end of file

matlab/set_dynare_seed.m

+function set_dynare_seed(a,b,c)
+% Set seeds depending on matlab (octave) version. This routine is called in dynare_config and can be called by the 
+% user in the mod file.
+%    
+% Copyright (C) 2010 Dynare Team
+%
+% This file is part of Dynare.
+%
+% Dynare is free software: you can redistribute it and/or modify
+% it under the terms of the GNU General Public License as published by
+% the Free Software Foundation, either version 3 of the License, or
+% (at your option) any later version.
+%
+% Dynare is distributed in the hope that it will be useful,
+% but WITHOUT ANY WARRANTY; without even the implied warranty of
+% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+% GNU General Public License for more details.
+%
+% You should have received a copy of the GNU General Public License
+% along with Dynare.  If not, see <http://www.gnu.org/licenses/>.
+global options_
+
+if ~nargin
+    error('set_dynare_seed:: I need at least one input argument!')
+end
+
+matlab_random_streams = 1;
+if  exist('OCTAVE_VERSION') || matlab_ver_less_than('7.7')
+    matlab_random_streams = 1;
+end
+
+if matlab_random_streams% Use new matlab interface.
+    if nargin==1
+        if ischar(a) && strcmpi(a,'default')
+            options_.DynareRandomStreams.algo = 'mt19937ar';
+            options_.DynareRandomStreams.seed = 0;
+            s = RandStream(options_.DynareRandomStreams.algo,'Seed',options_.DynareRandomStreams.seed);
+            reset(RandStream.setDefaultStream(s));
+            return
+        end
+        if ischar(a) && strcmpi(a,'reset')
+            s = RandStream(options_.DynareRandomStreams.algo,'Seed',options_.DynareRandomStreams.seed);
+            reset(RandStream.setDefaultStream(s));
+            return
+        end
+        if ischar(a)
+            error('set_dynare_seed:: something is wrong in the calling sequence!')
+        end
+        if ~ischar(a)
+            options_.DynareRandomStreams.algo = 'mt19937ar';
+            options_.DynareRandomStreams.seed = a;
+            s = RandStream(options_.DynareRandomStreams.algo,'Seed',options_.DynareRandomStreams.seed);
+            reset(RandStream.setDefaultStream(s));
+            return
+        end
+    elseif nargin==2
+        if ~ischar(a) || ~( strcmpi(a,'mcg16807') || ...
+              strcmpi(a,'mlfg6331_64') || ...
+              strcmpi(a,'mrg32k3a') || ...
+              strcmpi(a,'mt19937ar') || ...
+              strcmpi(a,'shr3cong') || ...
+              strcmpi(a,'swb2712') )
+            disp('set_dynare_seed:: First argument must be string designing the uniform random number algorithm!')
+            RandStream.list
+            disp(' ')
+            disp('set_dynare_seed:: Change the first input accordingly...')
+            disp(' ')
+            error(' ')
+        end
+        if ~isint(b)
+            error('set_dynare_seed:: The second input argument must be an integer!')
+        end
+        options_.DynareRandomStreams.algo = a;
+        options_.DynareRandomStreams.seed = b;
+        s = RandStream(options_.DynareRandomStreams.algo,'Seed',options_.DynareRandomStreams.seed);
+        reset(RandStream.setDefaultStream(s));
+    end
+else% Use old matlab interface.
+    if nargin==1
+        if ischar(a) && strcmpi(a,'default')
+            options_.DynareRandomStreams.algo = 'twister';
+            options_.DynareRandomStreams.seed = 0;
+            rand(options_.DynareRandomStreams.algo,options_.DynareRandomStreams.seed);
+            randn('state',options_.DynareRandomStreams.seed);
+            return
+        end
+        if ischar(a) && strcmpi(a,'reset')
+            rand(options_.DynareRandomStreams.algo,options_.DynareRandomStreams.seed);
+            randn('state',options_.DynareRandomStreams.seed);
+            return
+        end
+        if ~ischar(a) && isint(a)
+            options_.DynareRandomStreams.seed = a;
+            rand(options_.DynareRandomStreams.algo,options_.DynareRandomStreams.seed);
+            randn('state',options_.DynareRandomStreams.seed);
+        else
+            error('set_dynare_seed:: Something is wrong in the calling sequence!')
+        end
+    else
+        error('set_dynare_seed:: Cannot use more than one input argument with your version of Matlab/Octave!')
+    end
+end
\ No newline at end of file

matlab/simult.m

@@ -34,10 +34,10 @@ global M_ options_ oo_
 
 order = options_.order;
 replic = options_.replic;
+
 if replic == 0
     replic = 1;
 end
-seed = options_.simul_seed;
 
 if replic > 1
     fname = [M_.fname,'_simul'];
@@ -51,11 +51,6 @@ oo_.exo_simul = zeros(M_.maximum_lag+M_.maximum_lead+options_.periods,M_.exo_nbr
 chol_S = chol(M_.Sigma_e(i_exo_var,i_exo_var));
 
 for i=1:replic
-    if isempty(seed)
-        randn('state',sum(100*clock));
-    else
-        randn('state',seed+i-1);
-    end
     if ~isempty(M_.Sigma_e)
         oo_.exo_simul(:,i_exo_var) = randn(M_.maximum_lag+M_.maximum_lead+options_.periods,nxs)*chol_S;
     end
@@ -67,15 +62,4 @@ end
 
 if replic > 1
     fclose(fh);
-end
-
-
-% 02/20/01 MJ replaced ys by dr.ys
-% 02/22/01 MJ removed commented out lines
-%             removed useless temps
-%             stderr_ replaced by M_.Sigma_e
-% 02/28/01 MJ changed expression for M_.Sigma_e
-% 02/18/03 MJ added ys in the calling sequence for arbitrary initial values
-%             suppressed useless calling parameter istoch
-% 05/10/03 MJ removed repmat() in call to simult_() for lag > 1
-% 05/29/03 MJ test for 0 variances
+end
\ No newline at end of file

tests/seeds.mod

+// Example 1 from Collard's guide to Dynare
+var y, c, k, a, h, b;
+varexo e, u;
+
+parameters beta, rho, alpha, delta, theta, psi, tau;
+
+alpha = 0.36;
+rho   = 0.95;
+tau   = 0.025;
+beta  = 0.99;
+delta = 0.025;
+psi   = 0;
+theta = 2.95;
+
+phi   = 0.1;
+
+model;
+c*theta*h^(1+psi)=(1-alpha)*y;
+k = beta*(((exp(b)*c)/(exp(b(+1))*c(+1)))
+    *(exp(b(+1))*alpha*y(+1)+(1-delta)*k));
+y = exp(a)*(k(-1)^alpha)*(h^(1-alpha));
+k = exp(b)*(y-c)+(1-delta)*k(-1);
+a = rho*a(-1)+tau*b(-1) + e;
+b = tau*a(-1)+rho*b(-1) + u;
+end;
+
+initval;
+y = 1.08068253095672;
+c = 0.80359242014163;
+h = 0.29175631001732;
+k = 11.08360443260358;
+a = 0;
+b = 0;
+e = 0;
+u = 0;
+end;
+
+shocks;
+var e; stderr 0.009;
+var u; stderr 0.009;
+var e, u = phi*0.009*0.009;
+end;
+
+stoch_simul(periods=1000,irf=0,nomoments);
+endo_simul_1 = oo_.endo_simul;
+
+set_dynare_seed('reset')
+stoch_simul(periods=1000,irf=0,nomoments);
+endo_simul_2 = oo_.endo_simul;
+
+set_dynare_seed('default')
+stoch_simul(periods=1000,irf=0,nomoments);
+endo_simul_3 = oo_.endo_simul;
+
+t1 = endo_simul_1-endo_simul_2;
+t2 = endo_simul_1-endo_simul_3;
+
+if any(abs(t1(:))>1e-12) || any(abs(t2(:))>1e-12)
+   disp('Test 1.')
+   error('Test failure:: Problem with the seed of the random number algorithm')
+end
+
+set_dynare_seed(57)
+stoch_simul(periods=1000,irf=0,nomoments);
+endo_simul_4 = oo_.endo_simul;
+
+set_dynare_seed('reset')
+stoch_simul(periods=1000,irf=0,nomoments);
+endo_simul_5 = oo_.endo_simul;
+
+t3 = endo_simul_4-endo_simul_5;
+
+if any(abs(t3(:))>1e-12)
+   disp('Test 2.')
+   error('Test failure:: Problem with the seed of the random number algorithm')
+end
+
+set_dynare_seed('default')
+stoch_simul(periods=1000,irf=0,nomoments);
+endo_simul_6 = oo_.endo_simul;
+
+t4 = endo_simul_6-endo_simul_1;
+
+if any(abs(t4(:))>1e-12)
+   disp('Test 3.')
+   error('Test failure:: Problem with the seed of the random number algorithm')
+end
+
+if ~exist('OCTAVE_VERSION') || ~matlab_ver_less_than(7.7)
+
+    set_dynare_seed('mlfg6331_64',0)
+    stoch_simul(periods=1000,irf=0,nomoments);
+    endo_simul_7 = oo_.endo_simul;
+
+    set_dynare_seed('reset')
+    stoch_simul(periods=1000,irf=0,nomoments);
+    endo_simul_8 = oo_.endo_simul;
+
+    t5 = endo_simul_7-endo_simul_8;
+
+    if any(abs(t5(:))>1e-12)
+        disp('Test 4.')
+        error('Test failure:: Problem with the seed of the random number algorithm')
+    end
+
+    set_dynare_seed('default')
+    stoch_simul(periods=1000,irf=0,nomoments);
+    endo_simul_9 = oo_.endo_simul;
+    
+    t6 = endo_simul_9-endo_simul_1;
+
+    if any(abs(t6(:))>1e-12)
+        disp('Test 5.')
+        error('Test failure:: Problem with the seed of the random number algorithm')
+    end
+
+end