simrsmpc.m 4.53 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124
function [experiment] = simcdc(name, mpcmode, zetaVariant)
% Runs named experiments. First run init.m.
%
% The named experiments are:
%   'singlestep'    Performs a single step 0->0.50 and runs 5 seconds. Drops are random.
%   'repstep'       Repeated step changes.
%   'largerepstep'  Performs large steps close to the constraint. Note: Keeps u, does not use OL data.
%   'largerepstepol'  Same as above but uses open loop when switching.
%   'disturb'    A noisy disturbance is applied through the input u.
%
% mpc_mode is one of:
%   0: Local (LF) control only
%   1: Mode A
%   2: Mode B
%   3: Mode C
%   4: Disregard drop and run HF controller only (local HF control)
%
% zetaVariant is one of:
%   0: LQ
%   1: u(0)
%   2: Open loop, u(k-i)
%   3: u(1)
%
% After finished experiment use savesim(experiment) to store the result on disk.
% The result is stored in sim/{name}/{mode}/.

h_hf = evalin('base', 'h_hf');
refscale = 10/0.55;
rng(1);

assignin('base', 'zetaVariant', zetaVariant);

if strcmp(name, 'singlestep')
  sp=0.35;
  xinit = [-0.401 0 0]; % Small difference from reference allows Matlab to solve algebraic loop for LQc
  T = h_hf * [0:1/h_hf * 2] - 1/60;
  % Delay the reference so LQ doesn't kick in for all modes at the first step!
  simin_ref = [T; ones(1, 5)*-0.4*refscale ones(1, length(T)-5) * sp*refscale]'; 
  simin_drop = [T; randi(2, 1, length(T))-1]';
  simin_au = [T; zeros(3,length(T))]';
  stopTime=T(end);
elseif strcmp(name, 'repstep')
  xinit = [-0.401 0 0]; % Small difference from reference allows Matlab to solve algebraic loop for LQc
  T = h_hf * [0:1/h_hf * 2] - 1/60;
  % Delay the reference so LQ doesn't kick in for all modes at the first step!
  simin_ref = [ones(1, 5)*-0.4*refscale ones(1, length(T)-5) * 0.4*refscale]; 
  simin_ref = kron(repmat([1 -1], 1, 90), simin_ref);
  simin_drop = randi(2, 1, length(simin_ref))-1;
  simin_au = zeros(3,length(simin_ref));
  T = h_hf * [0:length(simin_ref)-1];
  stopTime=T(end);
  simin_ref = [T; simin_ref]';
  simin_drop = [T; simin_drop]';
  simin_au = [T; simin_au]';
elseif strcmp(name, 'largerepstep') || strcmp(name, 'largerepstepol')
  xinit = [-0.501 0 0]; % Small difference from reference allows Matlab to solve algebraic loop for LQc
  T = h_hf * [0:1/h_hf * 2] - 1/60;
  % Delay the reference so LQ doesn't kick in for all modes at the first step!
  simin_ref = [ones(1, 5)*-0.5*refscale ones(1, length(T)-5) * 0.5*refscale]; 
  simin_ref = kron(repmat([1 -1], 1, 90), simin_ref);
  simin_drop = randi(2, 1, length(simin_ref))-1;
  simin_au = zeros(3,length(simin_ref));
  T = h_hf * [0:length(simin_ref)-1];
  stopTime=T(end);
  simin_ref = [T; simin_ref]';
  simin_drop = [T; simin_drop]';
  simin_au = [T; simin_au]';
  if strcmp(name, 'largerepstepol')
    assignin('base', 'zetaVariant', 2);
  end
%  else
%    assignin('base', 'zetaVariant', 1);
%  end
elseif strcmp(name, 'disturb')
%  assignin('base', 'zetaVariant', 2);
  t = 360;
  xinit = [0.541 0 0]; % Small difference from reference allows Matlab to solve algebraic loop for LQc
  T = h_hf * [0:1/h_hf * t] - 1/60;
  % Delay the reference so LQ doesn't kick in for all modes at the first step!
  simin_ref = [ones(1, length(T))*0.54*refscale]; 
  simin_drop = randi(2, 1, length(simin_ref))-1;
  simin_au = 0 + 1.*randn(1, length(T));
  T = h_hf * [0:length(simin_ref)-1];
  stopTime=T(end);
  simin_ref = [T; simin_ref]';
  simin_drop = [T; simin_drop]';
  simin_au = [T; simin_au; zeros(2, length(simin_au))]';

elseif strcmp(name, 'basic')
  xinit = [-0.401 0 0]; % Small difference from reference allows Matlab to solve algebraic loop for LQc
  stopTime = 5;
  T = h_hf * [0:1/h_hf * stopTime] - 1/60;
  ref=(10/0.55)*0.5;
  simin_ref = [];
  steplen=3;
  for sec = 1:steplen:stopTime
    l=min(steplen, stopTime-sec);
    simin_ref = [simin_ref ref*ones(1, 1/h_hf*l)];
    ref = -ref;
  end
  simin_ref = [simin_ref ones(1, length(T)-length(simin_ref))*simin_ref(end)];
  simin_ref = [T; simin_ref]';
  simin_drop = [T; randi(2, 1, length(T))-1]';
  simin_au = [T; zeros(1,length(T))]';
else
  ERROR
end

SimStopTime=sprintf('%d', stopTime);
assignin('base', 'xinit', xinit);
assignin('base', 'simin_drop', simin_drop);
assignin('base', 'simin_au', simin_au);
assignin('base', 'simin_ref', simin_ref);
assignin('base', 'mpc_mode', mpcmode);
simout = sim('rsmpc', 'ReturnWorkspaceOutputs', 'on', 'StopTime', SimStopTime);

experiment = {};
experiment.sim = simout;
modes = {'LF', 'A', 'B', 'C', 'HF'};
experiment.mpc_mode = modes{evalin('base', 'mpc_mode')+1};
experiment.zeta_variant = sprintf('%d', evalin('base', 'zetaVariant'));
experiment.name = name;

end