cdf2mpc

CDF2MPC  Converts an IEEE CDF data file into a MATPOWER case struct.
   MPC = CDF2MPC(CDF_FILE_NAME)
   MPC = CDF2MPC(CDF_FILE_NAME, VERBOSE)
   MPC = CDF2MPC(CDF_FILE_NAME, MPC_NAME)
   MPC = CDF2MPC(CDF_FILE_NAME, MPC_NAME, VERBOSE)
   [MPC, WARNINGS] = CDF2MPC(CDF_FILE_NAME, ...)

   Converts an IEEE Common Data Format (CDF) data file into a MATPOWER case
   struct.

   Input:
       CDF_FILE_NAME :  name of the IEEE CDF file to be converted
       MPC_NAME      :  (optional) file name to use to save the resulting
                         MATPOWER case
       VERBOSE       :  1 (default) to display progress info, 0 otherwise

   Output(s):
       MPC      : resulting MATPOWER case struct
       WARNINGS : (optional) cell array of strings containing warning
                  messages (included by default in comments of MPC_NAME).

   The IEEE CDF does not include some data need to run an optimal power
   flow. This script creates default values for some of this data as
   follows:

       Bus data:
           Vmin = 0.94 p.u.
           Vmax = 1.06 p.u.
       Gen data:
           Pmin = 0 MW
           Pmax = Pg + baseMVA
       Gen cost data:
           Quadratic costs with:
               c2 = 10 / Pg, c1 = 20, c0 = 0, if Pg is non-zero, and
               c2 = 0.01,    c1 = 40, c0 = 0, if Pg is zero
           This should yield an OPF solution "close" to the
           existing solution (assuming it is a solved case)
           with lambdas near $40/MWh. See 'help caseformat'
           for details on the cost curve format.

   CDF2MPC may modify some of the data which are "infeasible" for
   running optimal power flow. If so, warning information will be
   printed out on screen.

   Note: Since our code can not handle transformers with variable tap,
   you may not expect to get exactly the same power flow solution
   using converted data. This is the case when we converted ieee300.cdf.

0001 function [mpc, warnings] = cdf2mpc(cdf_file_name, mpc_name, verbose)
0002 %CDF2MPC  Converts an IEEE CDF data file into a MATPOWER case struct.
0003 %   MPC = CDF2MPC(CDF_FILE_NAME)
0004 %   MPC = CDF2MPC(CDF_FILE_NAME, VERBOSE)
0005 %   MPC = CDF2MPC(CDF_FILE_NAME, MPC_NAME)
0006 %   MPC = CDF2MPC(CDF_FILE_NAME, MPC_NAME, VERBOSE)
0007 %   [MPC, WARNINGS] = CDF2MPC(CDF_FILE_NAME, ...)
0008 %
0009 %   Converts an IEEE Common Data Format (CDF) data file into a MATPOWER case
0010 %   struct.
0011 %
0012 %   Input:
0013 %       CDF_FILE_NAME :  name of the IEEE CDF file to be converted
0014 %       MPC_NAME      :  (optional) file name to use to save the resulting
0015 %                         MATPOWER case
0016 %       VERBOSE       :  1 (default) to display progress info, 0 otherwise
0017 %
0018 %   Output(s):
0019 %       MPC      : resulting MATPOWER case struct
0020 %       WARNINGS : (optional) cell array of strings containing warning
0021 %                  messages (included by default in comments of MPC_NAME).
0022 %
0023 %   The IEEE CDF does not include some data need to run an optimal power
0024 %   flow. This script creates default values for some of this data as
0025 %   follows:
0026 %
0027 %       Bus data:
0028 %           Vmin = 0.94 p.u.
0029 %           Vmax = 1.06 p.u.
0030 %       Gen data:
0031 %           Pmin = 0 MW
0032 %           Pmax = Pg + baseMVA
0033 %       Gen cost data:
0034 %           Quadratic costs with:
0035 %               c2 = 10 / Pg, c1 = 20, c0 = 0, if Pg is non-zero, and
0036 %               c2 = 0.01,    c1 = 40, c0 = 0, if Pg is zero
0037 %           This should yield an OPF solution "close" to the
0038 %           existing solution (assuming it is a solved case)
0039 %           with lambdas near $40/MWh. See 'help caseformat'
0040 %           for details on the cost curve format.
0041 %
0042 %   CDF2MPC may modify some of the data which are "infeasible" for
0043 %   running optimal power flow. If so, warning information will be
0044 %   printed out on screen.
0045 %
0046 %   Note: Since our code can not handle transformers with variable tap,
0047 %   you may not expect to get exactly the same power flow solution
0048 %   using converted data. This is the case when we converted ieee300.cdf.
0049 
0050 %   MATPOWER
0051 %   Copyright (c) 1996-2016, Power Systems Engineering Research Center (PSERC)
0052 %   by Deqiang (David) Gan, PSERC Cornell & Zhejiang University
0053 %   and Ray Zimmerman, PSERC Cornell
0054 %
0055 %   This file is part of MATPOWER.
0056 %   Covered by the 3-clause BSD License (see LICENSE file for details).
0057 %   See https://matpower.org for more info.
0058 
0059 %% define named indices into bus, gen, branch matrices
0060 [PQ, PV, REF, NONE, BUS_I, BUS_TYPE, PD, QD, GS, BS, BUS_AREA, VM, ...
0061     VA, BASE_KV, ZONE, VMAX, VMIN, LAM_P, LAM_Q, MU_VMAX, MU_VMIN] = idx_bus;
0062 [GEN_BUS, PG, QG, QMAX, QMIN, VG, MBASE, GEN_STATUS, PMAX, PMIN, ...
0063     MU_PMAX, MU_PMIN, MU_QMAX, MU_QMIN, PC1, PC2, QC1MIN, QC1MAX, ...
0064     QC2MIN, QC2MAX, RAMP_AGC, RAMP_10, RAMP_30, RAMP_Q, APF] = idx_gen;
0065 [F_BUS, T_BUS, BR_R, BR_X, BR_B, RATE_A, RATE_B, RATE_C, ...
0066     TAP, SHIFT, BR_STATUS, PF, QF, PT, QT, MU_SF, MU_ST, ...
0067     ANGMIN, ANGMAX, MU_ANGMIN, MU_ANGMAX] = idx_brch;
0068 [PW_LINEAR, POLYNOMIAL, MODEL, STARTUP, SHUTDOWN, NCOST, COST] = idx_cost;
0069 
0070 %% handle input args
0071 if nargin < 2
0072     verbose = 1;
0073     mpc_name = '';
0074 elseif ischar(mpc_name)     %% save the file
0075     if nargin < 3
0076         verbose = 1;
0077     end
0078 else                        %% don't save the file
0079     verbose = mpc_name;
0080     mpc_name = '';
0081 end
0082 
0083 %%-----  read data from CDF file into mpc fields
0084 %% verify valid input filename
0085 [cdf_path cdf_name cdf_ext] = fileparts(cdf_file_name);
0086 if isempty(cdf_ext)
0087     cdf_ext = '.cdf';
0088     cdf_file_name = strcat(cdf_name , cdf_ext);
0089 end
0090 
0091 %% open input file
0092 [fid, msg] = fopen(cdf_file_name, 'r');
0093 if fid < 0
0094     disp(msg);
0095     error('cdf2mpc: Can not read the input file:  %s', cdf_file_name )
0096 end
0097 if verbose
0098     fprintf('Converting file ''%s''\n', cdf_file_name);
0099     fprintf('  WARNINGS:\n');
0100 end
0101 
0102 %% initialize list of warnings
0103 warnings = {};
0104 
0105 % get baseMVA
0106 title_cdf = fgetl(fid);
0107 if isnumeric(str2num(title_cdf(32:37))) && not(isempty(str2num(title_cdf(32:37))))
0108     baseMVA = str2num(title_cdf(32:37));
0109     if length(strfind(title_cdf(2:9), '/')) == 2    %% date in the file
0110         warnings{end+1} = sprintf('check the title format in the first line of the cdf file.');
0111     end
0112 else
0113     error('cdf2mpc: Error getting the Base MVA, check the title format in the first line of the file.')
0114 end
0115     
0116 % find string 'BUS DATA FOLLOWS'
0117 while 1
0118     line = fgetl(fid);
0119     if line(1:16) == 'BUS DATA FOLLOWS', break, end
0120 end
0121 
0122 % ----- get bus data, feed them into matrix bus, gen, gencost
0123 ibus = 0;
0124 igen = 0;
0125 iarea = 0;
0126 
0127 while 1
0128     line = fgetl(fid);
0129     if line(1:4) == '-999', break, end
0130 
0131     % feed bus data
0132     ibus = ibus + 1;
0133     bus(ibus, BUS_I) = str2num(line(1:4));          % bus number
0134     bus_name{ibus} = strtrim(line(6:17));           % bus names
0135     bus(ibus, BUS_TYPE) = str2num(line(25:26));
0136     if bus(ibus, BUS_TYPE) == 0                     % bus type
0137         bus(ibus, BUS_TYPE) = 1;
0138     end
0139     if (bus(ibus, BUS_TYPE) < 2)                    % Pd
0140         bus(ibus, PD) = str2num(line(41:49)) - str2num(line(60:67));
0141     elseif (bus(ibus, BUS_TYPE) >= 2)
0142         bus(ibus, PD) = str2num(line(41:49));
0143     end
0144     bus(ibus, QD) = str2num(line(50:59));           % Qd
0145     bus(ibus, GS) = baseMVA*str2num(line(107:114)); % Gs
0146     bus(ibus, BS) = baseMVA*str2num(line(115:122)); % Bs
0147     bus(ibus, BUS_AREA) = str2num(line(19:20));     % area
0148     bus(ibus, VM) = str2num(line(28:33));           % Vm
0149     bus(ibus, VA) = str2num(line(34:40));           % Va
0150     bus(ibus, BASE_KV) = str2num(line(77:83));      % baseKV
0151     bus(ibus, ZONE) = str2num(line(21:23));         % zone
0152     bus(ibus, VMAX) = 1.06;                         % default voltage upper limit
0153     bus(ibus, VMIN) = 0.94;                         % default voltage lower limit
0154 
0155     % feed gen and gencost
0156     Pg = str2num(line(60:67));
0157     Qg = str2num(line(68:75));
0158     Qmax = str2num(line(91:98));
0159     Qmin = str2num(line(99:106));
0160     if bus(ibus, BUS_TYPE) >= 2
0161         igen = igen + 1;
0162         if bus(ibus, BUS_TYPE) == 3, refgen = igen; end
0163         gen(igen, GEN_BUS) = bus(ibus, BUS_I);      % bus number
0164         gen(igen, PG) = Pg;                         % Pg
0165         if gen(igen, PG) < 0                        % negative Pg is transformed as load
0166             bus(ibus, PD) = bus(ibus, PD) - gen(igen, PG);
0167             warnings{end+1} = sprintf('negative Pg at bus %g treated as Pd', bus(ibus, BUS_I));
0168             if verbose
0169                 fprintf('    %s\n', warnings{end});
0170             end
0171             gen(igen, PG) = 0;
0172         end
0173         gen(igen, QG)   = Qg;                       % Qg
0174         gen(igen, QMAX) = Qmax;                     % Qmax
0175         gen(igen, QMIN) = Qmin;                     % Qmin
0176         if Qmax - Qmin < 0.01                       % Qmax is modified
0177             gen(igen, QMAX) = Qmin + 0.1 * baseMVA;
0178             warnings{end+1} = sprintf('Qmax = Qmin at generator at bus %4i (Qmax set to Qmin + %g)', bus(ibus, BUS_I), baseMVA/10);
0179             if verbose
0180                 fprintf('    %s\n', warnings{end});
0181             end
0182         end
0183         gen(igen, VG)    = str2num(line(85:90));    % specified voltage
0184         gen(igen, MBASE) = baseMVA;                 % baseMVA
0185         gen(igen, GEN_STATUS) = 1;                  % default status is 'on'
0186         gen(igen, PMAX)  = gen(igen, 2) + baseMVA;  % Pmax
0187         gen(igen, PMIN)  = 0;                       % Pmin = 0 by default
0188 
0189         gencost(igen, MODEL)    = POLYNOMIAL;       % by default, sets the model as polynomial
0190         gencost(igen, STARTUP)  = 0;                % start up cost is zero by default
0191         gencost(igen, SHUTDOWN) = 0;                % shut down cost is zero by default
0192         gencost(igen, NCOST)    = 3;                % number of coefficients in polynomial cost
0193 %       gencost(igen, COST)     = 0.01;             % default c2
0194 %       gencost(igen, COST+1)   = 0.3;              % default c1
0195 %       gencost(igen, COST+2)   = 0.2;              % default c0
0196     end
0197 end
0198 
0199 totload = sum(bus(:, PD));
0200 totgen = sum(gen(:, PG));
0201 if totgen < 1.04 * totload
0202     gen(refgen, PMAX) = gen(refgen, PG) + 1.1 * totload - totgen;   % Pg at slack bus is modified
0203     warnings{end+1} = sprintf('Insufficient generation, setting Pmax at slack bus (bus %d) to %g', gen(refgen, [GEN_BUS, PMAX]));
0204     if verbose
0205         fprintf('    %s\n', warnings{end});
0206     end
0207 end
0208 
0209 % ----- set up the cost coefficients of generators
0210 ng = size(gen, 1);
0211 % gencost(:, COST)    = zeros(ng, 1);
0212 % gencost(:, COST+1)  = 100*ones(ng, 1) ./ (gen(:, PG) + 10*ones(ng, 1));
0213 % gencost(:, COST+2)  = 100*ones(ng, 1) ./ (gen(:, PG) + 10*ones(ng, 1));
0214 zg  = find(gen(:, PG) == 0);                %% for Pg = 0
0215 gencost(zg, COST)  = 0.01 * ones(size(zg));
0216 gencost(zg, COST+1) = 40 * ones(size(zg));
0217 nzg = find(gen(:, PG) ~= 0);                %% Pg non-zero
0218 gencost(nzg, COST) = 10 * ones(size(nzg)) ./ gen(nzg, PG);
0219 gencost(nzg, COST+1) = 20 * ones(size(nzg));
0220 gencost(:, COST+2) = zeros(ng, 1);
0221 
0222 % find string 'BRANCH DATA FOLLOWS'
0223 while 1
0224     line = fgetl(fid);
0225     if line(1:19) == 'BRANCH DATA FOLLOWS', break, end
0226 end
0227 
0228 % ----- get branch data, feed them into matrix branch
0229 k = 0;
0230 while 1
0231     line = fgetl(fid);
0232     if line(1:4) == '-999', break, end
0233 
0234     k = k + 1;
0235     branch(k, F_BUS)  = str2num(line(1:4));     % fbus (also the tap bus)
0236     branch(k, T_BUS)  = str2num(line(6:9));     % tbus
0237     branch(k, BR_R)   = str2num(line(20:29));   % R
0238     branch(k, BR_X)   = str2num(line(30:40));   % X
0239     branch(k, BR_B)   = str2num(line(41:50));   % B
0240     branch(k, RATE_A) = str2num(line(51:55));   % RATE A
0241     if branch(k, RATE_A) < 0.000001
0242         branch(k, RATE_A) = 0;                  % RATE A is modified
0243         warnings{end+1} = sprintf('MVA limit of branch %d - %d not given, set to %g', branch(k, [F_BUS, T_BUS, RATE_A]));
0244         if verbose
0245             fprintf('    %s\n', warnings{end});
0246         end
0247     end
0248     branch(k, RATE_B) = str2num(line(57:61));   % RATE B
0249     branch(k, RATE_C) = str2num(line(63:67));   % RATE C
0250     branch(k, TAP)    = str2num(line(77:82));   % transformer turns ratio
0251     branch(k, SHIFT)  = 0;                      % phase shifter can not be modelled
0252     branch(k, BR_STATUS) = 1;                   % by default, branch is on
0253 end
0254 if verbose
0255     fprintf('Done.\n');
0256 end
0257 fclose(fid);
0258 
0259 % put in struct
0260 mpc.baseMVA = baseMVA;
0261 mpc.bus     = bus;
0262 mpc.branch  = branch;
0263 mpc.gen     = gen;
0264 mpc.gencost = gencost;
0265 mpc.bus_name = bus_name;
0266 mpc = loadcase(mpc);    %% convert to internal (e.g. v. '2') case format
0267 
0268 %% (optionally) save MATPOWER case file
0269 if ~isempty(mpc_name)
0270     comments = {''};
0271     if ~isempty(title_cdf)
0272         comments{end+1} = sprintf('   %s', title_cdf);
0273     end
0274     comments{end+1} = '';
0275     comments{end+1} = sprintf('   Converted by MATPOWER %s using CDF2MPC on %s', mpver, date);
0276     comments{end+1} = sprintf('   from ''%s''.', cdf_file_name);
0277 
0278     %% warnings
0279     comments{end+1} = '';
0280     comments{end+1} = '   WARNINGS:';
0281     for k = 1:length(warnings)
0282         comments{end+1} = sprintf('       %s', warnings{k});
0283     end
0284     comments{end+1} = '';
0285     comments{end+1} = sprintf('   See CASEFORMAT for details on the MATPOWER case file format.');
0286 
0287     if verbose
0288         spacers = repmat('.', 1, 45-length(mpc_name));
0289         fprintf('Saving to MATPOWER case ''%s'' %s', mpc_name, spacers);
0290     end
0291     savecase(mpc_name, comments, mpc);
0292     if verbose
0293         fprintf(' done.\n');
0294     end
0295 end