IEA-ETSAP Forum

Hi TIMES modelling community, After the conversation at the Bonn meeting around open-source tools + solvers I was curious what the time and memory penalty of moving away from GAMS / CPLEX would be for the very small TIMES-NZ model. I am with those who argued there would be big benefits if we could migrate to more open and modern tools - if we can get them to work for us. I had a bit of a play - I'm coming at this without much experience of optimization problems at this scale, and I know the performance of different algorithms can depend on details of the problem, so I'm not sure how this would generalize to other / larger models. In any case the results, although a little confusing, seem to corroborate what was said in the conversation. I did a simple experiment for TIMES-NZ (kea scenario) and got the following results (not attempting to be scientific, just ballpark): ===================================  Running time for TIMES-NZ scenario kea-v2_1_2  VEDA (GAMS/CPLEX solver) running kea-v2_1_2    :  76 seconds * 7 threads Julia (HiGHS solver) running kea-v2_1_2        :  242 seconds * 6 threads Python (HiGHS solver) running kea-v2_1_2      :  405 seconds Open-source solver relative slowdown :    3 to 5 times longer (approx) ===================================  Memory usage for TIMES-NZ scenario kea-v2_1_2 On-disk kea-v2_1_2, saved as an mps file        :  85 MB VEDA (GAMS/CPLEX solver) running kea-v2_1_2    :    8 MB  # used by GAMS process. How can this be? Am I missing some of the memory usage? Julia (HiGHS solver) running kea-v2_1_2        : 2426 MB  # total memory used by julia process. Not clear how much of this is due to the model. The size of the model object in memory is 159 MB according to Julia. Python (HiGHS solver) running kea-v2_1_2      :  400 MB  # total memory used by python process. Open-source memory usage increase          :  50 to 300 times more (approx) =================================== I'm a bit confused by how much less memory is being used by GAMS / CPLEX while solving the problem. First thought was perhaps the CPLEX solver does a more effective presolve? But from the logging it looks like Highs is just as successful. Presolved Problem Specifications: GAMS/CPLEX: Rows: 166,606 Columns: 139,519 Nonzeros: 1,192,044 HiGHS: Rows: 166,922 Columns: 140,319 Nonzeros: 1,196,177 A few things I tried made little difference to the results: - Pulled out the Matrix representation of the LP and create a new JuMP model, attempting to suppress variable naming (result: same size and solve time) - Implemented the model as a Convex.jl optimization problem and tried solving it with SCS (Bad idea. Way too slow!) Will PS - in case anyone else wants to try profiling their TIMES model in other languages, here's a crude way: the following steps are one way to ensure a "writelp" line gets added at the bottom of the cplex.opt file when you run TIMES using VEDA. The "writelp" line tells cplex to dump the LP to disk. In my case the line that gets added is "writelp TIMES_kea-v2_1_2.lp", because my scenario is called kea-v2_1_2.     1) Go to the latest TIMES model file directory (GAMS_SrcTIMES.vA.B.C) present in your VEDA directory     2) Open the file solve.mod using a text editor like VSCode.     3) Find a line (e.g. for TIMES v4.8.1, it's line 33) that looks like:     $  IF    %METHOD% == P    $SETLOCAL METHOD LP     4) Immediately below it, insert the following six lines:     DISPLAY "Configure cpex to export an MPS file for use with other solvers.";     $onecho >> cplex.opt     * Keep the line empty above this one - it puts a newline where it's needed!     writelp %MODEL_NAME%_%RUN_NAME%.lp     $offecho     5) Run your TIMES scenarios using VEDA as usual.     6) After running TIMES, in the relevant scenario directory e.g. GAMS_WrkTIMES\kea-v2_1_2, you should find a new ".lp" file that contains the full model problem. In my case it's called TIMES_kea-v2_1_2.lp

Thanks for reporting your findings from the experiments. However, am I correct that the experiments were basically only about open source solver performance? > I am with those who argued there would be big benefits if we could migrate to more open and modern tools - if we can get them to work for us. I have understood those expected benefits to be referring to migrating the whole TIMES code to Julia (or some other open source modeling environment).  Should that not be accomplished, GAMS would be needed anyway, and many open source solvers like HiGHS are already available under GAMS as a free solver options (no solver cost). > VEDA (GAMS/CPLEX solver) running kea-v2_1_2    :  76 seconds * 7 threads > Julia (HiGHS solver) running kea-v2_1_2        :  242 seconds * 6 threads > Python (HiGHS solver) running kea-v2_1_2      :  405 seconds I am a bit confused about this, because that would seem to suggest that you have already converted the full TIMES code into Julia and Python, but that would sound too good to be true!  So, do you mean that the role of Julia and Python is here limited to reading the model from the file TIMES_kea-v2_1_2.lp (written out by GAMS/Cplex) into a Julia/Python model object and then simply solving that with HiGHS?  If so, do you in fact see some notable benefit of using Julia/Python for just that purpose (only the solver interface)?  > VEDA (GAMS/CPLEX solver) running kea-v2_1_2    :    8 MB  # used by GAMS process. How can this be? Am I missing some of the memory usage? I suspect that you may be missing the memory use of the Cplex solver process (which is not "gams.exe" but I guess would be gmsgennx.exe). >  $onecho >> cplex.opt [...] I think editing the solve.mod file is in general not a very good idea (and that location is not even relevant, because $onecho is a compile time directive), but of course could be done for a quick experiment. A better option would be to modify the options file either by the RUN file or by vtrun.cmd (with user commands that can be added into them automatically), if the writelp file should have the casename, or by just creating a new version of the cplex.opt having the writelp option (if a generic writelp file name is sufficient when creating that file). Note also that GAMS convert can also convert the model into a Cplex lp format file. Anyway, in my view the main issue would be about migrating the TIMES code, while the open source solvers are in fact available to users regardless of that, unless I am missing something?  Could you perhaps consider submitting a proposal to the related ETSAP call about a migration evaluation?

Hi Antti-L, thanks for your response! > Thanks for reporting your findings form the experiments. However, am I correct that the experiments were basically only about open source solver performance? Yes! I should have paused to explain my thinking. As you say, all I have done here is to intercept the LP that is produced by TIMES / GAMS, loaded it into other languages, and checked how an open-source solver performs when called from these open languages. The results are indicative of overall performance only if the model generation part, currently done by GAMS, is a relatively small part of the total computation. (Dr Parzen pointed out at the ETSAP meeting that most of the computational work typically happens in the solver.) In the case of TIMES-NZ (GAMS), model-generation is about the first 6s and roughly the next 70s is spent with the solver. So I should have reported the CPLEX time as 70s instead of 76. And hypothetically, if we did re-implement TIMES in Julia, the slowdown will be greater than reported here if model generation in Julia turns out to be slower than in GAMS. The assumption here is that since it's a small fraction (under 10 percent) of the time, the impact is less important - the goal here being a quick / rough estimate. (And translating TIMES looks harder!) For what it's worth, re-generating the model from its sparse matrix representation in Julia (attached script) took 2 seconds. Probably not that informative, but maybe slightly reassuring. > I have understood those expected benefits to be referring to migrating the whole TIMES code to Julia (or some other open source modeling environment).  Should that not be accomplished, GAMS would be needed anyway, and many open source solvers like HiGHS are already available under GAMS as a free solver options (no solver cost). Yes, the cheap experiment I performed here doesn't address the much more challenging task of porting TIMES to Julia, and doesn't yield any of the benefits. But I was curious what the order of magnitude of performance change would be. Based on the evidence to hand, I suspect TIMES-NZ would take not much more than about 3x longer if we used Julia to do the model generation and HiGHs to do the solve. I can't speak for other modelling groups, especially those running much bigger models, but for TIMES-NZ this performance penalty would be outweighed by the benefits. > ...I am a bit confused about this, because that would seem to suggest that you have already converted the full TIMES code into Julia and Python, but that would sound too good to be true!  So, do you mean that the role of Julia and Python is here limited to reading the model from the file IMES_kea-v2_1_2.lp (written out by GAMS/Cplex) into a Julia/Python model object and then simply solving that with HiGHS?  If so, do you in fact see some notable benefit of using Julia/Python for just that purpose (only the solver interface)? Sorry, no, I haven't done the conversion, and like you, I can't see any benefit of what I've done here, except as a way of getting a feel for what the performance might be (subject to the assumptions discussed) if we did move the whole thing across to Julia/Python. > I suspect that you may be missing the memory use of the Cplex solver process (which is not "gams.exe" but I guess would be gmsgennx.exe). Ah - thanks! With the gmsgennx.exe process it looks like memory usage by GAMS / CPLEX may be in the same ballpark (gmsgennx.exe was using 645 MB for most of the solve) as what is seen with HiGHS called from the open-source languages. > I think editing the solve.mod file is in general not a very good idea, but of course could be done for a quick experiment. A better option would be to modify the options file by vtrun.cmd (with user commands that can be added there automatically), if the writelp file should have the casename, or create a new version of the cplex.opt having the writelp option (if a generic writelp file name is sufficient when creating that file). Note also that GAMS convert can also convert the model into a Cplex lp format file. Thanks! (Still learning) > Anyway, in my view the main issue would be about migrating the TIMES code, while the open source solvers are in fact available to users regardless of that, unless I am missing something?  Could you perhaps consider submitting a proposal to the related ETSAP call about a migration evaluation? I'd be happy to submit this information to the migration evaluation (sorry, I don't see a thread for this), but as you point out all that I have done here is to get a rough estimate (in the case of TIMES-NZ) of the performance penalty of going open source. If others did similar exercises, we could build up a picture of how the community would be affected by a migration to Julia, possibly this has already been done. If so I'd be interested in the results! Sharing my (hacky) approach here was aimed to make that easier if others want to do the same experiment for their models. Hope that makes my original post a bit clearer!

PS - my attempt to include an attachment didn't work. In case useful, here are the contents of the Julia script I used, saved as profile_optimize_lp.jl: using JuMP using HiGHS using SparseArrays using LinearAlgebra using DelimitedFiles import MathOptInterface as MOI using Base.Filesystem: filesize #### Constants const MB = 1024 * 1024 # 1 MB in bytes const lp_path = "TIMES_kea-v2_1_2.lp" #### Functions """ recreate_matrix_model(m0::Model) -> Model Recreate a linear programming model from its matrix representation. # Arguments - `m0::Model`: Original JuMP model to extract matrix data from. # Returns - `Model`: New JuMP model built from the extracted matrix data. This function extracts matrix data such as constraints and objective coefficients from `m0` and constructs a new model `m1` using the HiGHS optimizer. """ function recreate_matrix_model(m0::Model) lp_data = try # Extract matrix data lp_matrix_data(m0); catch e error("Failed to extract matrix data: ", e) end # Access the matrix data A = lp_data.A b_lower = lp_data.b_lower b_upper = lp_data.b_upper x_lower = lp_data.x_lower x_upper = lp_data.x_upper c = lp_data.c c_offset = lp_data.c_offset # Recreate model from the matrix data m1 = Model(HiGHS.Optimizer) # Add variables @variable(m1, x_lower[i] <= x[i=1ize(A, 2)] <= x_upper[i]) # Add all constraints at once using matrix-vector multiplication @constraint(m1, b_lower .<= A * x .<= b_upper) # Set the objective @objective(m1, Min, dot(c, x)) return m1 end #### Main function main() if !isfile(lp_path) error("LP file does not exist: ", lp_path) end # Load and prepare the model println("Load LP from disk") @time m0 = read_from_file(lp_path) set_optimizer(m0, HiGHS.Optimizer) println("Regenerate model from its matrix representation") @time m1 = recreate_matrix_model(m0) println("Size of LP file on disk: ", filesize(lp_path) / MB, " MB") println("Memory used by original LP model before solving: ", Base.summarysize(m0) / MB, " MB") println("Memory used by matrix-form version before solving: ", Base.summarysize(m1) / MB, " MB") # Solve the model @time optimize!(m0) @time optimize!(m1) println("Memory used by original LP model after solving: ", Base.summarysize(m0) / MB, " MB") println("Memory used by matrix-form version after solving: ", Base.summarysize(m1) / MB, " MB") @assert termination_status(m0) == MOI.OPTIMAL @assert termination_status(m1) == MOI.OPTIMAL end main()

Hello, Thanks for bringing this topic. Here are the results of a quick test I have made on an instance of my times model. Only on solver though. My model: 2,706,769 rows  2,793,159 columns  19,279,758 non-zeroes Computer and software: 2 socket(s), 80 core(s), and 80 thread(s) available, Linux virtual machine with around 500GB of RAM with highs 1.7.1, gurobi 11, cplex 22, gams43, times 479 With the following solver parameters: barrier algorithm with crossover on 16 threads The results are: Gurobi: 23min Cplex: 48min Highs: hit the 50000s time limit, no output I do not easily have the stats on RAM but I would say roughly: -cplex < 30GB -gurobi ~60GB max during the presolve phase I switched to gurobi a couple months ago for multiple reasons: -faster than cplex on few tests -my colleagues ready to help on fixing numerical issues where more familiar with gurobi -this is more personal but I found the blog and documentation clearer and more complete for gurobi than cplex some materials related to the topic: https://www.iea-etsap.org/workshop/gothe...on_v03.pdf https://www.linkedin.com/posts/bakytzhan...29249--baA Thanks Victor

@willcattoneeca : Ok, many thanks for the clarifications. Indeed, if users would like to give up using high-performance solvers like Cplex, there could be a notable performance penalty. And your results of about 3.5 times longer solution time for your model instance still do look reasonably good in that respect.  Interestingly, Amit tried HiGHS last February with a comparable-sized model, and got less promising results: AK> "encouraged by your enthusiasm for HiGHS, I also tried it on a small model that CPLEX solves in under a minute. It took one hour to solve it" I am not sure which solver options he used, but it looks like the performance gap can be also much higher, depending on the model, like also VictorG confirmed. Nevertheless, I think many users thus far happy with Cplex might want to work with it also under Julia, and in that case one might have a cost penalty instead (but not for academic users, as I think IBM offers free Cplex for academic use). Concerning the ETSAP call for proposals on the feasibility evaluation, such was decided upon in the Bonn ExCo, and I thought you might be interested to submit such a proposal, if you have expertise on Julia.

@VictorG: Thanks a lot for the comparison with your bigger model instance. Indeed, many TIMES models have already been in the 2–10 million range in terms of number of constraints, and the TIMES GAMS code usually still performs quite well with such larger models.  For example, I would expect your test model to be generated by GAMS in less than 2 minutes, correct?  With such larger TIMES models, I suspect that under Julia/Python the performance of the model generator might not scale up so well and the overhead could become notable, but that is something the feasibility evaluation should be able to estimate.

Indeed, less than two minutes, before gurobi starts, the log file says: --- Executing GUROBI (Solvelink=0): elapsed 0:01:32.021 I feel I could wait a bit longer if it were in python/julia though

@Antti-L and @VictorG: thank you both for sharing the further context and information. Much appreciated! Based on my somewhat limited experience, the conversation linked by Bakha appears very relevant. It's often the case that non-performant Julia code can be significantly improved with subtle tweaks. While a clever implementation of TIMES in Julia/JuMP could probably avoid these issues, ensuring the success of a migration or feasibility project might require input from someone with fairly deep expertise in Julia profiling and optimization. No doubt there are many capable people who could do justice to this effort!

Very interesting discussion! Thanks for initiating it @willcattoneeca. Out of curiosity I've also done some benchmarking on the version of TIMES-DK available here: https://github.com/olejandro/times-dk_gams Total execution time was 5927 and 288 seconds for running with HiGHS and CPLEX respectively. Regarding the prototype, I could probably try converting this pyomo example (https://github.com/olejandro/times-gams-pyomo-demo) from the previous feasibility study to Julia to kick start the process. What do you think @Antti-L?

I am all in support of such a voluntary endeavour.  It would certainly represent a valuable contribution in view of the forthcoming renewed ETSAP feasibility evaluation project, for which you might then subsequently have a good position of submitting a proposal?

Cool, thanks! Well, let's see how this conversion goes. :-)

I've converted the pyomo example to Julia/JuMP. It is available here: https://github.com/olejandro/TIMES.jl. It can be run by executing the following in e.g. cmd: It takes quite a while to run, so please be patient. :-) I'll add some timing outputs during the execution, as well as include some info in the README soon. Please feel free to open any issues in the repository on GitHub.

@Olex : That's a great start!  Which version of Julia should work with it?  When I did my first experiments in 2018 (essentially generating just 300,000 EQ_ACTFLO equations), it worked quite well at that time, but when I tried it again in 2023 with Julia 1.6.7, hardly anything in my dataframe manipulations worked any longer, and I did not have the time and patience to investigate how I would be able to make my code work again. So, it seems different versions of Julia may have notable compatibility issues.

Thanks @Antii-L! Interesting! As far as I remember we ended the project in summer 2018. Julia v1.0 was released in August that year, so before that the API would have been unstable. This is my first experience with Julia, so I installed the latest version available 2 weeks ago - v1.10.4. If we decide to proceed further with this we could e.g. decide on the versions of Julia to support and just run automated testing on GitHub to discover any incompatibilities.