在 pyomo 中重新定义变量域

Question

我对如何在 pyomo 中使用实例更改模型有点困惑。例如，我有一个模型，然后我创建了一个变量

    #non negative variable
    SS.lambda_DA_E=pe.Var(SS.T, domain=pe.PositiveReals)

然后我创建一个实例并求解模型

    instance=SS.create_instance()
    solver = po.SolverFactory('mpec_nlp')
    results = solver.solve(instance, load_solutions=True, tee=True)

现在我想把这个变量重新定义为一个自由变量，给它赋上我从前面的解中得到的初始值，然后重新求解，我尝试使用

    def _init_rule_L_E2(m, t):
        return (pe.value(instance.lambda_DA_E[t]))
    instance.lambda_DA_E=pe.Var(SS.T, domain=pe.Reals,  initialize=_init_rule_L_E2)

然后我得到一个错误 ValueError：未初始化的 NumericValue 对象没有值 lambda_DA_E[t1]

我是否正确地重新定义了变量？因为我也收到警告

WARNING: Implicitly replacing the Component attribute lambda_DA_E (type=<class
    'pyomo.core.base.var.IndexedVar'>) on block unknown with a new Component
    (type=<class 'pyomo.core.base.var.IndexedVar'>). This is usually
    indicative of a modelling error. To avoid this warning, use
    block.del_component() and block.add_component().

这是完整的代码 ```

    import pandas as pd
    import pyomo.environ as pe
    import pyomo.opt as po
    import pyomo.mpec as mpc
    
    #DATA
    T=2;
    G=2;
    
    time = ['t{0}'.format(t+1) for t in range(T)]
    gen=['G{0}'.format(g+1) for g in range(G)]
    
    
    #Technical characteristic
    
    elec_maxprod= {'G1': 200, 'G2': 200,} #Only for Generators
    
    D_E={'t1':350, 't2': 150} #El demand
    
    C_fuel = {'G1': 15.0, 'G2': 50.0}
    #MODEL
    SS=pe.ConcreteModel()
    #Duals are desired
    SS.dual = pe.Suffix(direction=pe.Suffix.IMPORT) 
    
    
    ### SETS
    SS.T = pe.Set(initialize = time)
    SS.G = pe.Set(initialize = gen ) 
    
    ### PARAMETERS
    
    SS.P_max = pe.Param(SS.G, initialize = elec_maxprod) #Only for Generators
    SS.D_E = pe.Param(SS.T, initialize = D_E) #El demand
    SS.C_fuel = pe.Param(SS.G, initialize = C_fuel)
    
    ### VARIABLES
    
    SS.p_DA=pe.Var(SS.G, SS.T, domain=pe.PositiveReals)
    
    ##Dual variables
    #positive variables
    
    SS.mu_min_G=pe.Var(SS.G, SS.T, domain=pe.PositiveReals)
    SS.mu_max_G=pe.Var(SS.G, SS.T, domain=pe.PositiveReals)
    
    #free variables
    SS.lambda_DA_E=pe.Var(SS.T, domain=pe.PositiveReals) # day-ahead electricity price in period t [$ per MWh]
    
    ##CONSTRAINTS
    
    ##Maximum and Minimum Generation 
    SS.comp1 = mpc.ComplementarityList()
    for t in SS.T:
        for g in SS.G:
           SS.comp1.add(expr=mpc.complements( SS.p_DA[g,t] >= 0, SS.mu_min_G[g,t] >= 0))
            
    SS.comp2 = mpc.ComplementarityList()
    for t in SS.T:
        for g in SS.G:
            SS.comp2.add(expr=mpc.complements( SS.P_max[g] - SS.p_DA[g,t] >= 0, SS.mu_max_G[g,t] >= 0))
            
    ##Stationarity         
    SS.L1_p_DA = pe.ConstraintList()
    for t in SS.T:
        for g in SS.G:
            SS.L1_p_DA.add( SS.C_fuel[g] - SS.mu_min_G[g,t] + SS.mu_max_G[g,t] - SS.lambda_DA_E[t]== 0)       
            
    #El Balance
    SS.El_DA_bal=pe.ConstraintList()
    for t in SS.T:
        SS.El_DA_bal.add( sum(SS.p_DA[g,t] for g in SS.G) - SS.D_E[t]  == 0 )
            
    
    instance=SS.create_instance() 
    solver = po.SolverFactory('mpec_nlp')
    results = solver.solve(instance, load_solutions=True, tee=True)
    
    print("\nDisplaying Solution\n" + '-'*60)     
    ##Printing the variables results
    for t in SS.T:
        print('price_DA[',t,']', pe.value(instance.lambda_DA_E[t]))
    
    ##Redefine variable with a new feasible region and new initial value    
    def _init_rule_L_E2(m, t):
        return (pe.value(instance.lambda_DA_E[t]))
    instance.lambda_DA_E=pe.Var(SS.T, domain=pe.Reals,  initialize=_init_rule_L_E2) # day-ahead electricity price in period t [$ per MWh]
    
    ##Solve it again
    results = solver.solve(instance, load_solutions=True, tee=True)
    ##Printing the variables results
    for t in SS.T:
        print('price_DA[',t,']', pe.value(instance.lambda_DA_E[t]))

Answer 1

错误是因为您没有修改当前组件 lambda_DA_E 而是用新组件替换它。

如果你想改变组件的任何属性，你可以使用这样的赋值：

#After solving the whole model.
>>>instance.lambda_DA_E[t].display()
lambda_DA_E : Size=2, Index=T
Key : Lower : Value : Upper : Fixed : Stale : Domain
 t1 :     0 :  None :  None : False :  True : PositiveReals
 t2 :     0 :  None :  None : False :  True : PositiveReals
#After solving the first time
>>>print('changing Domain:')
>>>instance.lambda_DA_E[t].domain=pe.Reals
>>>results = solver.solve(instance, load_solutions=True, tee=True)
>>>instance.lambda_DA_E[t].display()
lambda_DA_E : Size=2, Index=T
Key : Lower : Value              : Upper : Fixed : Stale : Domain
 t1 :  None :  50.00000000097778 :  None : False : False :  Reals
 t2 :  None : 15.000000000977778 :  None : False : False :  Reals

然后你可以用solver.solve(instance, ...)

再次求解模型

只是一些注释：

因为 SS 是一个 ConcreteModel() 你不能调用 instance=SS.create_instance()。根据 pyomo 版本，这可能会引发 WARNING:
```
 WARNING: DEPRECATED: Cannot call Model.create_instance() on a constructed model; returning a clone of the current model instance.
```
如果你想创建 SS ConcreteModel 的副本，你可以使用 clone():
```
 instance = SS.clone()
```

我从未使用过 pyomo.environ.mpc 的 ComplementarityList，但尝试解决这个问题会 return 警告我：

 WARNING: Implicitly replacing the Component attribute v (type=<class 'pyomo.core.base.var.SimpleVar'>) on block comp2[2] with a new Component (type=<class 'pyomo.core.base.var.SimpleVar'>). This is usually indicative of a modelling error. To avoid this warning, use block.del_component() and block.add_component().

我认为 mpec_nlp 求解器使用了某种在解析时再次调用的重新表述。这可能不是问题，但仅供参考。

在 pyomo 中重新定义变量域

Redefine variable domain in pyomo

python

mathematical-optimization

pyomo