userfcn_direction_mll_formulation

 USERFCN_DIRECTION_MLL_FORMULATION adds one variable and as many
 constraints as dispatchable loads to enforce the load increase direction

0001 function om = userfcn_direction_mll_formulation(om,mpopt,args)
0002 % USERFCN_DIRECTION_MLL_FORMULATION adds one variable and as many
0003 % constraints as dispatchable loads to enforce the load increase direction
0004 
0005 %   MATPOWER
0006 %   Copyright (c) 2015-2016, Power Systems Engineering Research Center (PSERC)
0007 %   by Camille Hamon
0008 %
0009 %   This file is part of MATPOWER/mx-maxloadlim.
0010 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0011 %   See https://github.com/MATPOWER/mx-maxloadlim/ for more info.
0012 
0013 define_constants;
0014 mpc = om.get_mpc();
0015 dir_mll = mpc.dir_mll;
0016 
0017 % identify dispatchable loads
0018 idx_vl = isload(mpc.gen);
0019 n_vl = sum(idx_vl);
0020 n_g = size(mpc.gen,1)-n_vl;
0021 if length(dir_mll) ~= n_vl
0022     error_msg = ['The number of dispatchable loads is not equal to the '...
0023         'length of the direction vector'];
0024     error(error_msg);
0025 end
0026 % Add the amount of load increase alpha with 0 <= alpha <= inf
0027 om.add_var('alpha',1,0,0,inf);
0028 
0029 %% Load increase
0030 % Add the constraint for enforcing the direction of load increase
0031 Pl0 = -mpc.gen(idx_vl,PG)/mpc.baseMVA;
0032 % finding the internal indices of the dispatchable loads
0033 int_idx_disp_loads = find(idx_vl); 
0034 idx_A_dirmll_i = [1:n_vl 1:n_vl]';
0035 idx_A_dirmll_j = [int_idx_disp_loads' (n_g+n_vl+1)*ones(1,n_vl)]';
0036 vals_A_dirmll = [ones(n_vl,1);dir_mll];
0037 A_dirmll = sparse(idx_A_dirmll_i,idx_A_dirmll_j,vals_A_dirmll,n_vl,n_g+n_vl+1);
0038 om.add_lin_constraint('dir_mll',A_dirmll,-Pl0,-Pl0,{'Pg','alpha'});
0039 % Add cost of alpha to -1 to maximize loads in the given direction
0040 om.add_quad_cost('alpha_cost',[],-1,0,{'alpha'});
0041 
0042 %% Generator changes
0043 % Add the constraint for enforcing the direction of generation change
0044 idx_var_gen = find(mpc.dir_var_gen_all);
0045 if ~isempty(idx_var_gen)
0046     Pg0 = mpc.gen(idx_var_gen,PG)/mpc.baseMVA;
0047     nb_var_gen = length(idx_var_gen);
0048     idx_A_var_gen_i = [1:nb_var_gen 1:nb_var_gen]'; % Constraint number
0049     idx_A_var_gen_j = [idx_var_gen;(n_g+n_vl+1)*ones(nb_var_gen,1)]; % Generator number and alpha column (column in constraint matrix)
0050     vals_A_var_gen = [ones(nb_var_gen,1);-nonzeros(mpc.dir_var_gen_all)];
0051     A_var_gen = sparse(idx_A_var_gen_i,idx_A_var_gen_j,vals_A_var_gen,nb_var_gen,n_g+n_vl+1);
0052     om.add_constraints('dir_var_gen',A_var_gen,...
0053         Pg0,Pg0,{'Pg','alpha'});
0054 end