output profile constraint changed to equality

Author: SinaHKazemi

docs/src/reference/model.md

@@ -546,11 +546,9 @@ In all constraints, if a parameter is not defined, the default value is used. If
 
 - The output profile is used to model demands
 
-In theory, these constraints must be equalities, but in practice, when combined with the total\_energy\_out constraint, one of these constraints becomes redundant, which causes numerical issues for the solver.
- One solution is to remove one of these constraints; another is to keep them all and change them to '≥' or '≤' (this works the same).
 ```math
 \text{var\_energy\_out\_time}[p,y,t] \geq \text{var\_total\_energy\_out}[p,y] 
-* \text{output\_profile}[p,t] \quad \forall p \in P, \forall y \in Y, \forall t \in T
+* \text{output\_profile}[p,t] = \forall p \in P, \forall y \in Y, \forall t \in T
 ```
 
 ### Storage

src/core/model.jl

@@ -21,12 +21,11 @@ function build_model(input::Input, model::Union{JuMP.Model,Nothing}=nothing)
     if model === nothing
         # Gurobi
         model = JuMP.Model(Gurobi.Optimizer)
-        # set_attribute(model, "Crossover", 0)
-        # set_attribute(model, "Method", 2)
+        set_attribute(model, "Crossover", 0)
+        set_attribute(model, "Method", 2)
 
         # HiGHS
         # model = JuMP.Model(HiGHS.Optimizer)
-        model = JuMP.Model(Gurobi.Optimizer)
     end
     vars = add_vars!(model, input)
     constraints = add_constraints!(model, vars, input)
@@ -515,24 +514,21 @@ function  add_constraints!(model, vars, input::Input)::Dict
         end
     end
 
-    # In theory, these constraints should be equalities, but in practice, when combined with the energy_out constraints, one of these constraints becomes redundant, which causes numerical issues for the solver.
-    # One solution is to remove one of these constraints; another is to keep them all and change them to '≥' or '≤' (this works the same).
     constrs["load_shape"] = Dict()
     for p in processes
         !has_param("output_profile",p) && continue
         for y in years
             for t in timesteps
                 constrs["load_shape"][p,y,t] = @constraint(
                     model,
-                    vars["energy_out_time"][p,y,t] >= get_param("output_profile",(p,t)) * vars["total_energy_out"][p,y],
+                    vars["energy_out_time"][p,y,t] == get_param("output_profile",(p,t)) * vars["total_energy_out"][p,y],
                     base_name = "load_shape_$(p)_$(y)_$(t)"
                 )
             end
         end
     end
 
     # Storage
-
     constrs["c_rate_relation"] = Dict()
     for p in processes
         !(get_param("is_storage", p)) && continue