adding set_dynare_random_generator_state() and

get_dynare_random_generator_state(). Use now different seeds for
different Metropolis chains. Fixed handling of random generator state
thourghout the code.

matlab/cmaes.m

@@ -653,7 +653,7 @@ else % flgresume
 %  load(opts.SaveFilename, 'startseed');
 %  randn('state', startseed);
 %  disp(['SEED RELOADED FROM ' opts.SaveFilename]);
-  startseed = randn('state');         % for retrieving in saved variables
+%  startseed = randn('state');         % for retrieving in saved variables
 
   % Initialize further constants
   chiN=N^0.5*(1-1/(4*N)+1/(21*N^2));  % expectation of 

matlab/get_dynare_random_generator_state.m

+function [state_u,state_n] = get_dynare_random_generator_state()
+% Get state of Matlab/Octave random generator depending on matlab
+% (octave) version. 
+% In older versions, Matlab kept one generator for uniformly distributed numbers and
+% one for normally distributed numbers. 
+% For backward compatibility, we return two vectors, but, in recent
+% versions of Matlab and in Octave, we return two identical vectors.
+%    
+% Copyright (C) 2010-2012 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/>.
+
+    matlab_random_streams = ~(exist('OCTAVE_VERSION') || matlab_ver_less_than('7.7'));
+
+    if matlab_random_streams% Use new matlab interface.
+        if matlab_ver_less_than('7.12')
+            s = RandStream.getDefaultStream();
+        else
+            s = RandStream.getGlobalStream();
+        end
+        if isequal(s.Type,'legacy')
+            state_u = rand('state')
+            state_n = randn('state')
+        else            
+            state_u = s.State;
+            state_n = state_u;
+        end
+    else% Use old matlab interface.
+        state_u = rand('state')
+        state_n = randn('state')
+    end
\ No newline at end of file

matlab/independent_metropolis_hastings_core.m

@@ -101,8 +101,7 @@ jloop=0;
 
 for b = fblck:nblck,
     jloop=jloop+1;
-    randn('state',record.Seeds(b).Normal);
-    rand('state',record.Seeds(b).Unifor);
+    set_dynare_random_generator_state(record.Seeds(b).Unifor,record.Seeds(b).Normal);
     if (options_.load_mh_file~=0)  && (fline(b)>1) && OpenOldFile(b)
         load(['./' MhDirectoryName '/' ModelName '_mh' int2str(NewFile(b)) ...
               '_blck' int2str(b) '.mat'])
@@ -236,8 +235,7 @@ for b = fblck:nblck,
     elseif ~whoiam
         close(hh);
     end
-    record.Seeds(b).Normal = randn('state');
-    record.Seeds(b).Unifor = rand('state');
+    [record.Seeds(b).Unifor, record.Seeds(b).Normal] = get_dynare_random_generator_state();
     OutputFileName(jloop,:) = {[MhDirectoryName,filesep], [ModelName '_mh*_blck' int2str(b) '.mat']};
 end% End of the loop over the mh-blocks.
 

matlab/metropolis_hastings_initialization.m

@@ -184,8 +184,12 @@ if ~options_.load_mh_file && ~options_.mh_recover
     % separate initializaton for each chain
     JSUM = 0;
     for j=1:nblck,
-        record.Seeds(j).Normal = randn('state');
-        record.Seeds(j).Unifor = rand('state');
+        % we set a different seed for the random generator for each block
+        % then we record the corresponding random generator state (vector)
+        set_dynare_seed(options_.DynareRandomStreams.seed+j);
+        % record.Seeds keeps a vector of the random generator state and
+        % not the scalar seed despite its name
+        [record.Seeds(j).Unifor,record.Seeds(j).Normal] = get_dynare_random_generator_state();
     end
     record.InitialParameters = ix2;
     record.InitialLogLiK = ilogpo2;

matlab/random_walk_metropolis_hastings_core.m

@@ -135,14 +135,7 @@ jloop=0;
 JSUM = 0;
 for b = fblck:nblck,
     jloop=jloop+1;
-    try  % Trap in the case matlab slave is called from an octave master.
-        randn('state',record.Seeds(b).Normal);
-        rand('state',record.Seeds(b).Unifor);
-    catch
-        JSUM  = JSUM + sum(100*clock);
-        randn('state',JSUM);
-        rand('state',JSUM);
-    end
+    set_dynare_seed(options_.DynareRandomStreams.Seed+b);
     if (options_.load_mh_file~=0) && (fline(b)>1) && OpenOldFile(b)
         load(['./' MhDirectoryName '/' ModelName '_mh' int2str(NewFile(b)) ...
             '_blck' int2str(b) '.mat'])
@@ -276,8 +269,7 @@ for b = fblck:nblck,
 %         close(hh);
 %     end
     dyn_waitbar_close(hh);
-    record.Seeds(b).Normal = randn('state');
-    record.Seeds(b).Unifor = rand('state');
+    [record.Seeds(b).Unifor, record.Seeds(b).Normal] = get_dynare_random_generator_state();
     OutputFileName(jloop,:) = {[MhDirectoryName,filesep], [ModelName '_mh*_blck' int2str(b) '.mat']};
 end% End of the loop over the mh-blocks.
 

matlab/set_dynare_random_generator_state.m

+function [state_u,state_n] = get_dynare_random_generator_state(state_u,state_n)
+% Wet state of Matlab/Octave random generator depending on matlab
+% (octave) version. 
+% In older versions, Matlab kept one generator for uniformly distributed numbers and
+% one for normally distributed numbers. 
+% For backward compatibility, we return two vectors, but, in recent
+% versions of Matlab and in Octave, we return two identical vectors.
+%    
+% Copyright (C) 2010-2012 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/>.
+
+    matlab_random_streams = ~(exist('OCTAVE_VERSION') || matlab_ver_less_than('7.7'));
+
+    if matlab_random_streams% Use new matlab interface.
+        if matlab_ver_less_than('7.12')
+            s = RandStream.getDefaultStream();
+        else
+            s = RandStream.getGlobalStream();
+        end
+        if isequal(s.Type,'legacy')
+            rand('state',state_u);
+            randn('state',state_n);
+        else            
+            if ~isequal(state_u,state_n)
+                error(['You are using the new Matlab RandStream mechanism ' ...
+                       'with a single random generator, but the values ' ...
+                       'of the state of the uniformly ' ...
+                       'distributed numbers and of the state of the ' ...
+                       'normally distributed numbers are different. Something must be ' ...
+                       'wrong, such as reloading old Metropolis runs, ' ...
+                       'computed on a different version of Matlab. If you ' ...
+                       'don''t understand the origin of the problem, ' ...
+                       'please, contact Dynare''s development team.'])
+            end
+            s_new = RandStream(s.type,'State',state_u);
+            if matlab_ver_less_than('7.12')
+                RandStream.setDefaultStream(s);
+            else
+                RandStream.setGlobalStream(s);
+            end
+        end
+    else% Use old matlab interface.
+        rand('state',state_u);
+        randn('state',state_n);
+    end
\ No newline at end of file