MATLAB-YALMIP optimal trendy program cannot converge

I use MATLAB-YALMIP to write the optimal trend program based on the IEEEE14 system, but I can't converge. Please help to see if the program is wrong?

code:

% Optimal Power Flow Example
% Based on the IEEE 14-bus standard case (i.e., case14 in Matpower)
% Implemented in MATLAB with YALMIP
clear; clc;

% Load power system data
mpc = loadcase('case14'); % You can also load other data files, such as case9, case30, etc.

% Define variables and parameters
nb = size(mpc.bus, 1); % Number of buses (nodes)
nl = size(mpc.branch, 1); % Number of branches (lines)
ng = size(mpc.gen, 1); % Number of generators
nr = nb - ng; % Number of load buses

% Extract node data
bus_id = mpc.bus(:, 1); % Bus IDs
bus_type = mpc.bus(:, 2); % Bus types (1-PQ, 2-PV, 3-Slack)
Pd = mpc.bus(:, 3); % Active power demand (MW)
Qd = mpc.bus(:, 4); % Reactive power demand (MVAr)
Vmin = mpc.bus(:, 13); % Minimum bus voltage (p.u.)
Vmax = mpc.bus(:, 12); % Maximum bus voltage (p.u.)

% Extract generator data
gen_id = mpc.gen(:, 1); % Bus IDs where generators are connected
Pgmin = mpc.gen(:, 10); % Minimum generator active power output (MW)
Pgmax = mpc.gen(:, 9); % Maximum generator active power output (MW)
Qgmin = mpc.gen(:, 5); % Minimum generator reactive power output (MVAr)
Qgmax = mpc.gen(:, 4); % Maximum generator reactive power output (MVAr)
a2 = mpc.gencost(:, 5); % Generator cost curve coefficient for quadratic term
a1 = mpc.gencost(:, 6); % Generator cost curve coefficient for linear term
a0 = mpc.gencost(:, 7); % Generator cost curve constant term

% Extract branch data
fbus = mpc.branch(:, 1); % "From" bus (starting node) of each branch
tbus = mpc.branch(:, 2); % "To" bus (ending node) of each branch
Plmax = ones(nl, 1) * 4000; % Maximum branch active power flow limit (MW)
% Plmax = mpc.branch(:, 6) .* mpc.branch(:, 3) * 1e-3; % Maximum branch active power flow limit (MW)

% Calculate the bus admittance matrix
[Ybus, ~, ~] = makeYbus(mpc); % Use a Matpower function to calculate the bus admittance matrix
Gbus = real(Ybus); % Real part of bus admittance matrix (conductance)
Bbus = imag(Ybus); % Imaginary part of bus admittance matrix (susceptance)

% Define YALMIP variables
Pg = sdpvar(ng, 1, 'full'); % Generator active power output variables
Qg = sdpvar(ng, 1, 'full'); % Generator reactive power output variables
V = sdpvar(nb, 1, 'full'); % Bus voltage magnitude variables
theta = sdpvar(nb, 1, 'full'); % Bus voltage angle variables
Pl = sdpvar(nl, 1, 'full'); % Branch active power flow variables

% Define the objective function
objective = sum(a2 .* Pg.^2 + a1 .* Pg + a0); % Minimize the total operating cost

% Define constraints
constraints = [];
% Node power balance constraints
for i = 1:nb
    % Active power balance equation for node i
    if ismember(i, gen_id)
        record = find(gen_id == i);
        constraints = [constraints, Pg(record) - Pd(i) - V(i) * sum(V * (Gbus(i, :) * cos(theta(i) - theta) + Bbus(i, :) * sin(theta(i) - theta))) == 0];
        constraints = [constraints, Qg(record) - Qd(i) + V(i) * sum(V * (Gbus(i, :) * sin(theta(i) - theta) - Bbus(i, :) * cos(theta(i) - theta))) == 0];
    else
        constraints = [constraints, -Pd(i) - V(i) * sum(V * (Gbus(i, :) * cos(theta(i) - theta) + Bbus(i, :) * sin(theta(i) - theta))) == 0];
        constraints = [constraints, -Qd(i) + V(i) * sum(V * (Gbus(i, :) * sin(theta(i) - theta) - Bbus(i, :) * cos(theta(i) - theta))) == 0];
    end
end

% Generator output limit constraints
for i = 1:ng
    % Minimum active power output constraint for generator i
    constraints = [constraints, Pgmin(i) <= Pg(i)];
    % Maximum active power output constraint for generator i
    constraints = [constraints, Pg(i) <= Pgmax(i)];
    % Minimum reactive power output constraint for generator i
    constraints = [constraints, Qgmin(i) <= Qg(i)];
    % Maximum reactive power output constraint for generator i
    constraints = [constraints, Qg(i) <= Qgmax(i)];
end
% Bus voltage limit constraints
for i = 1:nb
    % Minimum bus voltage magnitude constraint for bus i
    constraints = [constraints, Vmin(i) <= V(i)];
    % Maximum bus voltage magnitude constraint for bus i
    constraints = [constraints, V(i) <= Vmax(i)];
end
% Branch power flow limit constraints
for i = 1:nl
    % Calculate active power flow on branch i (assuming the "from" bus is the positive direction)
    Pl(i) = V(fbus(i))^2 * Gbus(fbus(i), tbus(i)) - V(fbus(i)) * V(tbus(i)) * (Gbus(fbus(i), tbus(i)) * cos(theta(fbus(i) - theta(tbus(i))) + Bbus(fbus(i), tbus(i)) * sin(theta(fbus(i) - theta(tbus(i)))));
    % Absolute value constraint for active power flow on branch i
    constraints = [constraints, abs(Pl(i)) <= Plmax(i)];
end

% Solve the optimal power flow problem
result = solvesdp(constraints, objective); % Use YALMIP to solve the optimal power flow problem

% Output and analyze the results
if result.problem == 0 % If the solution is successful
    disp('Optimal Power Flow Solved Successfully!'); % Display success message
    disp('Total System Operating Cost:'); % Display total system operating cost
    value(objective) % Display the objective function value
    disp('Generator Active Power Output:'); % Display generator active power output
    value(Pg) % Display the generator active power output variable values
    disp('Generator Reactive Power Output:'); % Display generator reactive power output
    value(Qg) % Display the generator reactive power output variable values
    disp('Bus Voltage Magnitude:'); % Display bus voltage magnitude
    value(V) % Display the bus voltage magnitude variable values
    disp('Bus Voltage Angle (in degrees):'); % Display bus voltage angle (in degrees)
    rad2deg(value(theta)) % Display bus voltage angle variable values and convert to degrees
    disp('Branch Active Power Flow:'); % Display branch active power flow
    value(Pl)
end

[Torsten]

You should ask your question in the YALMIP forum:

https://groups.google.com/g/yalmip

https://github.com/yalmip/YALMIP/discussions

[Younes]

Hi,

Could you please share the version without the error, if possible. ?

[Kirk]

this link is error https://github.com/yalmip/YALMIP/discussionsSnow Rider 3D??