Run E4ST with full model

[Ethan-Russell]

Gurobi settings trial and error

All of this is running with Gurobi 9.51

Suboptimal Solution

The first time I was able to successfully run, Gurobi returned a suboptimal solution from the barrier method. I'm not sure why this would be, b/c the objective value was pretty low, at least compared to HiGHS solution. Adjusting the BarConvTol gave more iterations, but ultimately had the same result. I wonder if there are more tolerances that I should be looking into adjusting? With the following settings, it ran for a while, and eventually returned that the problem was infeasible or unbounded.

LogToConsole    = false,
NumericFocus    = 3,
BarHomogeneous  = 1,
method          = 2,
BarIterLimit    = 1000,
Crossover       = 1,
FeasibilityTol  = 1e-2,
OptimalityTol   = 1e-6,
BarConvTol      = 1e-8,
LogFile         = log_file_full,
Threads         = 1,

To determine whether infeasible or unbounded, I set DualReductions to zero and re-solved, as per this.

Out of Memory

When Threads was not set, I sometimes received an out of memory exception. I think that is because it was allocating a bunch to each thread, as per this site?

I tried setting Threads to 1, and that seemed to be ok

Numerical Trouble Encountered

However, when I got to the barrier problem, I received a note that "Numerical trouble encountered" after only 6 iterations of the barrier algorithm, so (since I had crossover turned on) it went to the simplex method and started ticking away. The log file gave a note that it may be a good idea to run with the BarHomogeneous algorithm (which is set to -1, or off for LP's by default)

[Ethan-Russell]

Gurobi confirmed that the problem is INFEASIBLE. Now to figure out why...

I checked the point where all variables are zero, which in my mind should be feasible. The only thing that would make that infeasible are pcap_min constraining the lower bound of capacity, and cf_min constraining the lower bound of generation (if pcap_min is provided too). However, this problem has neither of those.

I tested out this theory, as follows:

function compute_feasibility_of_zeros(model)
    d = (;
        θ_bus = 0.0,
        pgen_gen = 0.0,
        pserv_bus=0.0,
        pcap_gen = 0.0,
    )

    point =  Dict(((model[k][i] => v) for (k,v) in pairs(d) for i in eachindex(model[k]))...)

    primal_feasibility_report(model, point)
end

This returned an empty Dict, which signifies that it is a feasible point!

Currently running compute_conflict! on the model to see what Gurobi says is the infeasibility.

[Ethan-Russell]

Variable	Description	Lowest Lower Bound	Highest Upper Bound
theta_bus	The voltage angle of the bus	-1.00E+04	1.00E+04
pcap_gen	The power capacity of each generator	0	9.00E+03
pgen_gen	The power generation of each generator	0	9.00E+03
pserv_bus	The power served at each bus	0	8.00E+03

Constraint	Description	Lowest LHS Coefficient Magnitude	Highest LHS Coefficient Magnitude	Comments
cons_pflow	Constrain the total power flowing in and out to equal the generated power less the power served	1.00E-01	3.00E+05	This is adding power flows, which are a function of theta_bus values, often with large coefficients (up to 3e5), plus pgen_gen variable(s) minus the pserv_bus variable
cons_pgen_max	Constrain the generation to be less than capacity	1	1
cons_ref_bus	Constrain reference buses to have 0 theta_bus	1	1
cons_branch_pflow_pos	Constrain the power flow to be less than limit	1.00E-01	3.00E+05	RHS between 0 and 1e4

[Ethan-Russell]

Here are some combinations that worked:

• With pcurt_bus_kw as variable, set all the reactances to 0.1 and disable line limits. This is in order to reduce the range of coefficients.
• With pcurt_bus_kw as variable, simply disable line limits. This still worked
• I tried combining parallel lines and re-running with the line limits, did NOT work.
• Tried re-running with pcurt_bus as variable (instead of kw scaling) to reduce the coefficient matrix, and it converged! Only modification is the objective function is in terms of millions of dollars.

[Ethan-Russell]

I noticed that in the reduced model I was running, the demand was far greater than the demand to be matched in the load_match sheet of the esc file. I added demand_match to the config, and now I am seeing:

• 10E-5 (~0) curtailment
• 7 minutes solve time (with Gurobi 10.0.0)

I'm very excited to see what solve time might look like once we add more generators.

[Ethan-Russell]

Closed in #76

[Ethan-Russell]

Reopening - unable to converge with full set of buildable generators, unfortunately.

[Ethan-Russell]

Running without expansion.

After adding in some annual adjustments as well as the capex/fom fix, now seeing longer runtime despite better problem statistics being slightly better. I wonder if it is because previously, pcap was pretty much able to be neglected since its coefficients would have been so small compared to generation.

Optimize a model with 2002728 rows, 1347635 columns and 6343740 nonzeros
Model fingerprint: 0xa4226789
Coefficient statistics:
  Matrix range     [2e-01, 1e+05]
  Objective range  [3e-06, 2e+00]
  Bounds range     [8e-04, 1e+06]
  RHS range        [9e-04, 2e+04]
...
Barrier solved model in 58 iterations and 6016.58 seconds (3844.77 work units)
Optimal objective 2.95350800e+04

Model Summary:
9×4 DataFrame
 Row │ variable               type           dimensions     length
     │ Symbol                 Symbol         Any            Any
─────┼─────────────────────────────────────────────────────────────
   1 │ cons_branch_pflow_neg  ConstraintRef  irregular      559884
   2 │ cons_branch_pflow_pos  ConstraintRef  irregular      559884
   3 │ cons_pbal              ConstraintRef  (8761, 1, 52)  455572
   4 │ cons_pgen_max          ConstraintRef  (8219, 1, 52)  427388
   5 │ pcap_gen               VariableRef    (8219, 1)      8219
   6 │ pflow_dc               VariableRef    (17, 1, 52)    884
   7 │ pgen_gen               VariableRef    (8219, 1, 52)  427388
   8 │ plcurt_bus             VariableRef    (8761, 1, 52)  455572
   9 │ θ_bus                  VariableRef    (8761, 1, 52)  455572

[Ethan-Russell]

Ok, after adding the availability factors and making sure the regions are right, the full model is now solving with and without buildable generators.

• no-expansion.log
• expansion.log
• old.log

Specific questions:

• Do you recommend changing any of our default Gurobi settings?
• Do you know why the memory is higher in the rewritten model?
• Do you know why each barrier iteration would be taking longer?
• Do you think it is worthwhile to reformulate to decrease coefficient ranges?

[Ethan-Russell]

• Ray: often simplex is faster (on smaller problems)
• Ray: when running many smaller cases, could get better results with fewer threads
• Ray: number of threads could reduce the time per iteration.
• Ray: number of iterations can be very problem-specific, even separate from problem size.
• Ray has used CPLEX and Mosek, sometimes one of the other solvers will be better for different problems
• Could look into concurrent method (does JuMP have this to try multiple solvers???)
• Worth looking into how dual simplex works on this problem
• With barrier method, Gurobi's default is to have crossover turned on, but we turned it off because crossover was taking a very long time.
• Gurobi tuning mode for finding best parameters?
• May be worth contacting Gurobi to set up a call for tuning
• Ray has not done much with reformulation to decrease coefficient ranges.

[Ethan-Russell]

Tried Gurobi tuning with the default set of parameters. It was not able to find a method that would improve on the default set of parameters. One observation. When leaving the Method parameter unspecified, the solver chooses to use a concurrent LP optimizer, with primal simplex, dual simplex, and barrier. It displays the barrier log, and barrier converges really fast, 35 minutes! But then it crosses over, and takes hours (haven't had one converge yet after 24 hours). If I leave the method as barrier with the rest of the default settings, it takes much longer and does not converge. Why would that be? Is there something with the concurrent method that helps the barrier method? It is not possible to run the concurrent method without crossover. Only other thought is that the concurrent method removes a couple of threads from barrier, so maybe 6-threaded barrier has a better starting seed than 8-threaded.

[Ethan-Russell]

Closing, though this is a helpful issue for parameter tuning, etc.