OPF_HESSFCN Evaluates Hessian of Lagrangian for AC OPF.
LXX = OPF_HESSFCN(X, LAMBDA, COST_MULT, OM, YBUS, YF, YT, MPOPT, IL)
Hessian evaluation function for AC optimal power flow, suitable
for use with MIPS or FMINCON's interior-point algorithm.
Inputs:
X : optimization vector
LAMBDA (struct)
.eqnonlin : Lagrange multipliers on power balance equations
.ineqnonlin : Kuhn-Tucker multipliers on constrained branch flows
COST_MULT : (optional) Scale factor to be applied to the cost
(default = 1).
OM : OPF model object
YBUS : bus admittance matrix
YF : admittance matrix for "from" end of constrained branches
YT : admittance matrix for "to" end of constrained branches
MPOPT : MATPOWER options struct
IL : (optional) vector of branch indices corresponding to
branches with flow limits (all others are assumed to be
unconstrained). The default is [1:nl] (all branches).
YF and YT contain only the rows corresponding to IL.
Outputs:
LXX : Hessian of the Lagrangian.
Examples:
Lxx = opf_hessfcn(x, lambda, cost_mult, om, Ybus, Yf, Yt, mpopt);
Lxx = opf_hessfcn(x, lambda, cost_mult, om, Ybus, Yf, Yt, mpopt, il);
See also OPF_COSTFCN, OPF_CONSFCN.0001 function Lxx = opf_hessfcn(x, lambda, cost_mult, om, Ybus, Yf, Yt, mpopt, il) 0002 %OPF_HESSFCN Evaluates Hessian of Lagrangian for AC OPF. 0003 % LXX = OPF_HESSFCN(X, LAMBDA, COST_MULT, OM, YBUS, YF, YT, MPOPT, IL) 0004 % 0005 % Hessian evaluation function for AC optimal power flow, suitable 0006 % for use with MIPS or FMINCON's interior-point algorithm. 0007 % 0008 % Inputs: 0009 % X : optimization vector 0010 % LAMBDA (struct) 0011 % .eqnonlin : Lagrange multipliers on power balance equations 0012 % .ineqnonlin : Kuhn-Tucker multipliers on constrained branch flows 0013 % COST_MULT : (optional) Scale factor to be applied to the cost 0014 % (default = 1). 0015 % OM : OPF model object 0016 % YBUS : bus admittance matrix 0017 % YF : admittance matrix for "from" end of constrained branches 0018 % YT : admittance matrix for "to" end of constrained branches 0019 % MPOPT : MATPOWER options struct 0020 % IL : (optional) vector of branch indices corresponding to 0021 % branches with flow limits (all others are assumed to be 0022 % unconstrained). The default is [1:nl] (all branches). 0023 % YF and YT contain only the rows corresponding to IL. 0024 % 0025 % Outputs: 0026 % LXX : Hessian of the Lagrangian. 0027 % 0028 % Examples: 0029 % Lxx = opf_hessfcn(x, lambda, cost_mult, om, Ybus, Yf, Yt, mpopt); 0030 % Lxx = opf_hessfcn(x, lambda, cost_mult, om, Ybus, Yf, Yt, mpopt, il); 0031 % 0032 % See also OPF_COSTFCN, OPF_CONSFCN. 0033 0034 % MATPOWER 0035 % Copyright (c) 1996-2017, Power Systems Engineering Research Center (PSERC) 0036 % by Ray Zimmerman, PSERC Cornell 0037 % and Carlos E. Murillo-Sanchez, PSERC Cornell & Universidad Nacional de Colombia 0038 % 0039 % This file is part of MATPOWER. 0040 % Covered by the 3-clause BSD License (see LICENSE file for details). 0041 % See https://matpower.org for more info. 0042 0043 %% ----- evaluate d2f ----- 0044 [f, df, d2f] = opf_costfcn(x, om); 0045 d2f = d2f * cost_mult; 0046 0047 %%----- evaluate Hessian of power balance constraints ----- 0048 d2G = om.eval_nln_constraint_hess(x, lambda.eqnonlin, 1); 0049 0050 %%----- evaluate Hessian of flow constraints ----- 0051 d2H = om.eval_nln_constraint_hess(x, lambda.ineqnonlin, 0); 0052 0053 %%----- do numerical check using (central) finite differences ----- 0054 if 0 0055 mpc = om.get_mpc(); 0056 nl = size(mpc.branch, 1); %% number of branches 0057 if nargin < 9 0058 il = (1:nl); %% all lines have limits by default 0059 end 0060 0061 nx = length(x); 0062 step = 1e-5; 0063 num_d2f = sparse(nx, nx); 0064 num_d2G = sparse(nx, nx); 0065 num_d2H = sparse(nx, nx); 0066 for i = 1:nx 0067 xp = x; 0068 xm = x; 0069 xp(i) = x(i) + step/2; 0070 xm(i) = x(i) - step/2; 0071 % evaluate cost & gradients 0072 [fp, dfp] = opf_costfcn(xp, om); 0073 [fm, dfm] = opf_costfcn(xm, om); 0074 % evaluate constraints & gradients 0075 [Hp, Gp, dHp, dGp] = opf_consfcn(xp, om, Ybus, Yf, Yt, mpopt, il); 0076 [Hm, Gm, dHm, dGm] = opf_consfcn(xm, om, Ybus, Yf, Yt, mpopt, il); 0077 num_d2f(:, i) = cost_mult * (dfp - dfm) / step; 0078 num_d2G(:, i) = (dGp - dGm) * lambda.eqnonlin / step; 0079 num_d2H(:, i) = (dHp - dHm) * lambda.ineqnonlin / step; 0080 end 0081 d2f_err = full(max(max(abs(d2f - num_d2f)))); 0082 d2G_err = full(max(max(abs(d2G - num_d2G)))); 0083 d2H_err = full(max(max(abs(d2H - num_d2H)))); 0084 if d2f_err > 1e-6 0085 fprintf('Max difference in d2f: %g\n', d2f_err); 0086 end 0087 if d2G_err > 1e-5 0088 fprintf('Max difference in d2G: %g\n', d2G_err); 0089 end 0090 if d2H_err > 1e-6 0091 fprintf('Max difference in d2H: %g\n', d2H_err); 0092 end 0093 end 0094 0095 Lxx = d2f + d2G + d2H;