如何限制 GAMS 混合整数非线性规划中非零变量的数量？
How to limit the count of nonzero variables in GAMS mixed-integer nonlinear programming?
问题背景：
Neo 区块链上有 25 位候选人获得投票。除了我以外的每个选民都投票了。
排名第 1 至第 7 的候选人每人将给予 200 美元，按比例分配给投票给他们的人。 8日至21日每人派发100元。比如我投了100票给目前第8位的候选人，把him/her推到了第7位，其他人投了900票给他，所以我得到候选人200*100/(900+100)=20美金.
我知道所有其他选民的投票结果，我最多可以额外投N票给n个候选人，因为我只有n个代理人代表们。也就是说，我的投票向量不应超过 n 个非零值，并且我的投票向量的总和应等于 N。
我应该如何投票获得最高利润？
("Dollar"在下面的代码中写成"GAS"。我可以投的票数在代码中写成"NEO")
我的 GAMS 代码：
这些代码出错：Endogenous $ operation not allowed
（因为我的投票，我还没有想好如何处理立场的变化。这可能很难。）
（LINDO 求解器给出了目前最好的结果）
Set i "candidate index" /
0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24
/;
Parameters
V(i) "other people's votes" /
0 11409980
1 3024833
2 2200734
3 2040928
4 2032100
5 2017067
6 2010927
7 2007339
8 2001625
9 1008817
10 1000025
11 1000000
12 641913
13 396926
14 392171
15 390364
16 383919
17 383722
18 382447
19 382310
20 381850
21 923
22 272
23 180
24 0
/
NEO "our NEOs" /1162/
num_agents "our num of agents" /2/
reward_factor(i) "how much gas is given by the candidate at each rank" /
0 200
1 200
2 200
3 200
4 200
5 200
6 200
7 100
8 100
9 100
10 100
11 100
12 100
13 100
14 100
15 100
16 100
17 100
18 100
19 100
20 100
21 0
22 0
23 0
24 0
/;
Integer Variable vote(i) "voting vector that needs to be optimized";
vote.up(i) = NEO;
Free Variable GAS_reward "total GAS reward";
Binary Variable is_nonzero(i) "whether an element in the voting vector is zero";
Integer Variable count_nonzero "count of nonzero values in vote(i)";
Equations
GAS_reward_eqn "the reward according to our voting vector v"
sum_v_eqn "sum of voting vector equals num of our NEOs"
is_nonzero_eqn(i) "whether the element vote(i) is not zero"
num_agents_eqn "sum of nonzero element in voting vector does not exceed num of our agents";
GAS_reward_eqn.. GAS_reward =e= sum(i, reward_factor(i)*vote(i)/(V(i)+vote(i)+0.000001));
sum_v_eqn.. NEO =e= sum(i, vote(i));
is_nonzero_eqn(i).. is_nonzero(i) =e= 1 $ (vote(i));
num_agents_eqn.. sum(i, is_nonzero(i)) =l= num_agents;
Model NEOBurger /all/;
Option MINLP=lindo;
Solve NEOBurger using MINLP maximizing GAS_reward;
（其中代码 is_nonzero_eqn(i).. is_nonzero(i) =e= 1 $ (vote(i)); 是非线性的！）
如果您需要额外的测试用例和参考，这里是问题的 Python（3.7 或 3.8）代码
strategy_multi_agents.py
# Good profit
# Can somehow deal with sorting problem
# Cannot restrict num of agents
# Cannot restrict to integer
from typing import Union, List
import numpy as np
np.set_printoptions(suppress=True)
from scipy.optimize import minimize
class BurgerNEOMultiAgentStrategy:
"""
To find the best distribution of our GAS assuming that we do not affect the ranking
"""
consensus_rank, consensus_gain = 7, 200
council_rank, council_gain = 21, 100
def __init__(self, candidate_votes: List[int], our_NEOs: int, num_of_agents: int = 2):
"""
:param candidate_votes: should be sorted in descending order. All the results are based on the descending order.
"""
self.len_candidate_votes = len(candidate_votes)
if self.len_candidate_votes >= self.council_rank:
self.candidate_votes = sorted(candidate_votes, reverse=True)
if self.candidate_votes != candidate_votes:
raise ValueError('You did not input a descending list of candidate_votes.')
self.candidate_votes = np.array(self.candidate_votes, dtype=np.int_) # descending order
else:
raise ValueError(f'len(candidate_votes) should be >= {self.council_rank}. Got {self.len_candidate_votes}.')
self.our_NEOs = np.int_(our_NEOs)
reward_factor = np.zeros(self.len_candidate_votes, dtype=np.int_)
reward_factor[0:self.consensus_rank] = self.consensus_gain
reward_factor[self.consensus_rank:self.council_rank] = self.council_gain
self.reward_factor = reward_factor
self.num_of_agents = num_of_agents
def calc_GAS_gain(self, voting_vector: np.ndarray) -> Union[np.float_, np.int_]:
new_voting_result = voting_vector + self.candidate_votes
argsort = np.flip(new_voting_result.argsort())
voting_vector = voting_vector[argsort]
new_voting_result = new_voting_result[argsort]
return np.sum(self.reward_factor * voting_vector / (new_voting_result + 1e-5))
def calc_voting_vector(self):
def sum_is_our_NEO(voting_vector: np.ndarray):
return np.sum(voting_vector) - self.our_NEOs
starting_point = self.our_NEOs / self.len_candidate_votes * np.ones(self.len_candidate_votes)
constraints = [
{'type': 'ineq', 'fun': lambda voting_vector: voting_vector}, # all greater than or equals 0
{'type': 'eq','fun': sum_is_our_NEO}, # sum is our_NEO
# {'type': 'ineq', 'fun': lambda x: -np.abs(x - np.round(x))}, # all integers
# {'type': 'ineq', 'fun':
# lambda x: - np.count_nonzero(np.invert(np.isclose(x, np.zeros(self.len_candidate_votes), atol=1e-3/self.len_candidate_votes))) + self.num_of_agents}, # num of agents
]
result = minimize(lambda voting_vector: (-self.calc_GAS_gain(voting_vector)),
starting_point, constraints=constraints, tol=1e-15, options={'maxiter': 1e9})
# return np.round(result.x).astype(np.int_), -result.fun, np.sum(np.round(result.x).astype(np.int_))
return result.x, -result.fun, np.sum(result.x)
if __name__ == '__main__':
print(
BurgerNEOMultiAgentStrategy(
[11409980, 3024833, 2200734, 2040928, 2032100, 2017067, 2010927, 2007339, 2001625, 1008817, 1000025, 1000000, 641913, 396926, 392171, 390364, 383919, 383722, 382447, 382310, 381850, 923, 272, 180, 0],
1162, 2).calc_voting_vector()
)
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 0,0,0,0, 0,0], 50, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 0,0], 25, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,1], 25, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,1], 25, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,0], 25, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,0], 1000, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 0,0], 26, 2).calc_voting_vector())
print()
# Restricted to integer
# Restricted num of agents
# Not robust (sometimes constraints are broken)
# Not good profit
# Cannot deal with sorting
from typing import List
import random
from gekko import GEKKO
from strategy_multi_agents import BurgerNEOMultiAgentStrategy as ScipyStrategy
class BurgerNEOMultiAgentStrategy:
"""
To find the best distribution of our GAS assuming that we do not affect the ranking
"""
consensus_rank, consensus_gain = 7, 200
council_rank, council_gain = 21, 100
def __init__(self, candidate_votes: List[int], our_NEOs: int, num_of_agents: int = 2):
"""
:param candidate_votes: should be sorted in descending order. All the results are based on the descending order.
"""
self.model = GEKKO(remote=False, name='neoburger-strategy')
self.model.options.IMODE = 3
self.model.options.MAX_ITER = 10000
self.model.options.COLDSTART = 1
self.model.options.AUTO_COLD = 1
self.len_candidate_votes = len(candidate_votes)
if self.len_candidate_votes >= self.council_rank:
self.candidate_votes = sorted(candidate_votes, reverse=True)
if self.candidate_votes != candidate_votes:
raise ValueError('You did not input a descending list of candidate_votes.')
self.candidate_votes = [self.model.Param(can) for can in self.candidate_votes] # descending order
else:
raise ValueError(f'len(candidate_votes) should be >= {self.council_rank}. Got {self.len_candidate_votes}.')
self.our_NEOs = self.model.Param(our_NEOs, name='our_NEOs')
reward_factor = [self.consensus_gain] * self.consensus_rank \
+ [self.council_gain] * (self.council_rank - self.consensus_rank) \
+ [0] * (self.len_candidate_votes - self.council_rank)
self.reward_factor = [self.model.Param(r) for r in reward_factor]
self.num_of_agents = self.model.Param(num_of_agents, name='num_of_agents')
# voting_vector = [self.our_NEOs / self.len_candidate_votes] * self.len_candidate_votes
voting_vector = ScipyStrategy(candidate_votes, our_NEOs, num_of_agents).calc_voting_vector()[0]
voting_vector = [max(0, round(v)) for v in voting_vector]
self.voting_vector = [self.model.Var(v, lb=0, ub=our_NEOs, integer=True) for v in voting_vector]
def calc_GAS_gain(self):
# new_vote_result = [v+V for v, V in zip(self.voting_vector, self.candidate_votes)]
# sorted_vote_result = [(v, V) for v, V in zip(self.voting_vector, new_vote_result)] # sorted(, key=lambda x: x[1].value)
#
# sorted_voting_vector = [v for v, V in sorted_vote_result]
# sorted_vote_result = [V for v, V in sorted_vote_result]
# return sum([v*r/(V+v+0.0001) for v, V, r in zip(sorted_voting_vector, sorted_vote_result, self.reward_factor)])
return sum([v*r/(V+v+1) for v, V, r in zip(self.voting_vector, self.candidate_votes, self.reward_factor)])
def calc_voting_vector(self):
# constraints
# sum(voting_vector) == our_NEOs
self.model.Equation(sum(self.voting_vector) == self.our_NEOs)
# agent number
nonzero_vote = sum([self.model.if3(v - 0.1, 0, 1) for v in self.voting_vector])
# self.model.Equation(nonzero_vote >= 0)
self.model.Equation(nonzero_vote <= self.num_of_agents)
# objective
self.model.Maximize(self.calc_GAS_gain())
# self.model.options.RTOL = 1e-14
self.model.options.OTOL = 1e-14
self.model.options.REQCTRLMODE = 1
self.model.solve(disp=False, debug=False, GUI=False)
return self.voting_vector, -self.model.options.OBJFCNVAL
if __name__ == '__main__':
print(
BurgerNEOMultiAgentStrategy(
[11409980, 3024833, 2200734, 2040928, 2032100, 2017067, 2010927, 2007339, 2001625, 1008817, 1000025, 1000000, 641913, 396926, 392171, 390364, 383919, 383722, 382447, 382310, 381850, 923, 272, 180, 0],
1162, 6).calc_voting_vector()
)
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 0,0,0,0, 0,0], 50, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 0,0], 25, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,1], 25, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,1], 25, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,0], 25, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,0], 1000, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 0,0], 26, 6).calc_voting_vector())
print()
而且我还有一些代码（此处未显示），它可以让 1 或 2 个代理的一切正确。这不是通过优化实现的，而是通过手动数学推导实现的。
我会这样写：
* This is not using proper GAMS syntax:
*is_nonzero_eqn(i).. is_nonzero(i) =e= 1 $ (vote(i));
* a linear version is:
is_nonzero_eqn(i).. vote(i) =l= is_nonzero(i)*vote.up(i);
这是一个标准的 MIP 公式。使用：
Model NEOBurger /all/;
Option MINLP=baron, optcr = 0;
Solve NEOBurger using MINLP maximizing GAS_reward;
我们看到：
GAMS/BARON 36.1.0 r2c0a44a Released Aug 02, 2021 WEI x86 64bit/MS Window
===========================================================================
BARON version 21.1.13. Built: WIN-64 Wed Jan 13 16:04:12 EST 2021
BARON is a product of The Optimization Firm.
For information on BARON, see https://minlp.com/about-baron
If you use this software, please cite publications from
https://minlp.com/baron-publications, such as:
Kilinc, M. and N. V. Sahinidis, Exploiting integrality in the global
optimization of mixed-integer nonlinear programming problems in BARON,
Optimization Methods and Software, 33, 540-562, 2018.
===========================================================================
This BARON run may utilize the following subsolver(s)
For LP/MIP/QP: ILOG CPLEX
For NLP: MINOS, SNOPT, GAMS external NLP, IPOPT, FILTERSD, FILTERSQP
===========================================================================
Doing local search
Solving bounding LP
Preprocessing found feasible solution with value 0.203660597686E-001
Starting multi-start local search
Done with local search
===========================================================================
Iteration Open nodes Time (s) Lower bound Upper bound
* 1 1 0.50 0.303390 3.71066
1 1 0.56 0.303390 0.303737
* 2 2 0.62 0.303405 0.303737
* 2 2 0.69 0.303428 0.303737
* 4 2 0.83 0.303473 0.303737
* 4 1 0.88 0.303639 0.303737
* 5 1 0.92 0.303681 0.303737
5 0 0.92 0.303681 0.303681
Calculating duals
*** Normal completion ***
Wall clock time: 4.05
Total CPU time used: 0.94
Total no. of BaR iterations: 5
Best solution found at node: 5
Max. no. of nodes in memory: 2
All done
===========================================================================
Solution = 0.30368112407161 found at node 5
Best possible = 0.303681125072
Absolute gap = 1.00039027062238E-9 optca = 1E-9
Relative gap = 3.29421287011233E-9 optcr = 1E-9
(Note that BARON uses a different formula to compute the relative gap as
was used for the above reported value.)
如果我们显示我们看到的相关变量：
---- 90 VARIABLE vote.L voting vector that needs to be optimized
19 466.000, 20 696.000
---- 90 VARIABLE is_nonzero.L whether an element in the voting vector is zero
19 1.000, 20 1.000
这是正确的行为。
请注意 MINLP = NLP + MIP，因此我们可以（并且应该）在编写 MINLP 模型时使用我们了解的有关 MIP 建模的所有工具和技巧。
问题背景：
Neo 区块链上有 25 位候选人获得投票。除了我以外的每个选民都投票了。 排名第 1 至第 7 的候选人每人将给予 200 美元，按比例分配给投票给他们的人。 8日至21日每人派发100元。比如我投了100票给目前第8位的候选人，把him/her推到了第7位，其他人投了900票给他，所以我得到候选人200*100/(900+100)=20美金.
我知道所有其他选民的投票结果，我最多可以额外投N票给n个候选人，因为我只有n个代理人代表们。也就是说，我的投票向量不应超过 n 个非零值，并且我的投票向量的总和应等于 N。
我应该如何投票获得最高利润？
("Dollar"在下面的代码中写成"GAS"。我可以投的票数在代码中写成"NEO")
我的 GAMS 代码：
这些代码出错：Endogenous $ operation not allowed
（因为我的投票，我还没有想好如何处理立场的变化。这可能很难。）
（LINDO 求解器给出了目前最好的结果）
Set i "candidate index" /
0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24
/;
Parameters
V(i) "other people's votes" /
0 11409980
1 3024833
2 2200734
3 2040928
4 2032100
5 2017067
6 2010927
7 2007339
8 2001625
9 1008817
10 1000025
11 1000000
12 641913
13 396926
14 392171
15 390364
16 383919
17 383722
18 382447
19 382310
20 381850
21 923
22 272
23 180
24 0
/
NEO "our NEOs" /1162/
num_agents "our num of agents" /2/
reward_factor(i) "how much gas is given by the candidate at each rank" /
0 200
1 200
2 200
3 200
4 200
5 200
6 200
7 100
8 100
9 100
10 100
11 100
12 100
13 100
14 100
15 100
16 100
17 100
18 100
19 100
20 100
21 0
22 0
23 0
24 0
/;
Integer Variable vote(i) "voting vector that needs to be optimized";
vote.up(i) = NEO;
Free Variable GAS_reward "total GAS reward";
Binary Variable is_nonzero(i) "whether an element in the voting vector is zero";
Integer Variable count_nonzero "count of nonzero values in vote(i)";
Equations
GAS_reward_eqn "the reward according to our voting vector v"
sum_v_eqn "sum of voting vector equals num of our NEOs"
is_nonzero_eqn(i) "whether the element vote(i) is not zero"
num_agents_eqn "sum of nonzero element in voting vector does not exceed num of our agents";
GAS_reward_eqn.. GAS_reward =e= sum(i, reward_factor(i)*vote(i)/(V(i)+vote(i)+0.000001));
sum_v_eqn.. NEO =e= sum(i, vote(i));
is_nonzero_eqn(i).. is_nonzero(i) =e= 1 $ (vote(i));
num_agents_eqn.. sum(i, is_nonzero(i)) =l= num_agents;
Model NEOBurger /all/;
Option MINLP=lindo;
Solve NEOBurger using MINLP maximizing GAS_reward;
（其中代码 is_nonzero_eqn(i).. is_nonzero(i) =e= 1 $ (vote(i)); 是非线性的！）
如果您需要额外的测试用例和参考，这里是问题的 Python（3.7 或 3.8）代码
strategy_multi_agents.py
# Good profit
# Can somehow deal with sorting problem
# Cannot restrict num of agents
# Cannot restrict to integer
from typing import Union, List
import numpy as np
np.set_printoptions(suppress=True)
from scipy.optimize import minimize
class BurgerNEOMultiAgentStrategy:
"""
To find the best distribution of our GAS assuming that we do not affect the ranking
"""
consensus_rank, consensus_gain = 7, 200
council_rank, council_gain = 21, 100
def __init__(self, candidate_votes: List[int], our_NEOs: int, num_of_agents: int = 2):
"""
:param candidate_votes: should be sorted in descending order. All the results are based on the descending order.
"""
self.len_candidate_votes = len(candidate_votes)
if self.len_candidate_votes >= self.council_rank:
self.candidate_votes = sorted(candidate_votes, reverse=True)
if self.candidate_votes != candidate_votes:
raise ValueError('You did not input a descending list of candidate_votes.')
self.candidate_votes = np.array(self.candidate_votes, dtype=np.int_) # descending order
else:
raise ValueError(f'len(candidate_votes) should be >= {self.council_rank}. Got {self.len_candidate_votes}.')
self.our_NEOs = np.int_(our_NEOs)
reward_factor = np.zeros(self.len_candidate_votes, dtype=np.int_)
reward_factor[0:self.consensus_rank] = self.consensus_gain
reward_factor[self.consensus_rank:self.council_rank] = self.council_gain
self.reward_factor = reward_factor
self.num_of_agents = num_of_agents
def calc_GAS_gain(self, voting_vector: np.ndarray) -> Union[np.float_, np.int_]:
new_voting_result = voting_vector + self.candidate_votes
argsort = np.flip(new_voting_result.argsort())
voting_vector = voting_vector[argsort]
new_voting_result = new_voting_result[argsort]
return np.sum(self.reward_factor * voting_vector / (new_voting_result + 1e-5))
def calc_voting_vector(self):
def sum_is_our_NEO(voting_vector: np.ndarray):
return np.sum(voting_vector) - self.our_NEOs
starting_point = self.our_NEOs / self.len_candidate_votes * np.ones(self.len_candidate_votes)
constraints = [
{'type': 'ineq', 'fun': lambda voting_vector: voting_vector}, # all greater than or equals 0
{'type': 'eq','fun': sum_is_our_NEO}, # sum is our_NEO
# {'type': 'ineq', 'fun': lambda x: -np.abs(x - np.round(x))}, # all integers
# {'type': 'ineq', 'fun':
# lambda x: - np.count_nonzero(np.invert(np.isclose(x, np.zeros(self.len_candidate_votes), atol=1e-3/self.len_candidate_votes))) + self.num_of_agents}, # num of agents
]
result = minimize(lambda voting_vector: (-self.calc_GAS_gain(voting_vector)),
starting_point, constraints=constraints, tol=1e-15, options={'maxiter': 1e9})
# return np.round(result.x).astype(np.int_), -result.fun, np.sum(np.round(result.x).astype(np.int_))
return result.x, -result.fun, np.sum(result.x)
if __name__ == '__main__':
print(
BurgerNEOMultiAgentStrategy(
[11409980, 3024833, 2200734, 2040928, 2032100, 2017067, 2010927, 2007339, 2001625, 1008817, 1000025, 1000000, 641913, 396926, 392171, 390364, 383919, 383722, 382447, 382310, 381850, 923, 272, 180, 0],
1162, 2).calc_voting_vector()
)
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 0,0,0,0, 0,0], 50, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 0,0], 25, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,1], 25, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,1], 25, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,0], 25, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,0], 1000, 2).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 0,0], 26, 2).calc_voting_vector())
print()
# Restricted to integer
# Restricted num of agents
# Not robust (sometimes constraints are broken)
# Not good profit
# Cannot deal with sorting
from typing import List
import random
from gekko import GEKKO
from strategy_multi_agents import BurgerNEOMultiAgentStrategy as ScipyStrategy
class BurgerNEOMultiAgentStrategy:
"""
To find the best distribution of our GAS assuming that we do not affect the ranking
"""
consensus_rank, consensus_gain = 7, 200
council_rank, council_gain = 21, 100
def __init__(self, candidate_votes: List[int], our_NEOs: int, num_of_agents: int = 2):
"""
:param candidate_votes: should be sorted in descending order. All the results are based on the descending order.
"""
self.model = GEKKO(remote=False, name='neoburger-strategy')
self.model.options.IMODE = 3
self.model.options.MAX_ITER = 10000
self.model.options.COLDSTART = 1
self.model.options.AUTO_COLD = 1
self.len_candidate_votes = len(candidate_votes)
if self.len_candidate_votes >= self.council_rank:
self.candidate_votes = sorted(candidate_votes, reverse=True)
if self.candidate_votes != candidate_votes:
raise ValueError('You did not input a descending list of candidate_votes.')
self.candidate_votes = [self.model.Param(can) for can in self.candidate_votes] # descending order
else:
raise ValueError(f'len(candidate_votes) should be >= {self.council_rank}. Got {self.len_candidate_votes}.')
self.our_NEOs = self.model.Param(our_NEOs, name='our_NEOs')
reward_factor = [self.consensus_gain] * self.consensus_rank \
+ [self.council_gain] * (self.council_rank - self.consensus_rank) \
+ [0] * (self.len_candidate_votes - self.council_rank)
self.reward_factor = [self.model.Param(r) for r in reward_factor]
self.num_of_agents = self.model.Param(num_of_agents, name='num_of_agents')
# voting_vector = [self.our_NEOs / self.len_candidate_votes] * self.len_candidate_votes
voting_vector = ScipyStrategy(candidate_votes, our_NEOs, num_of_agents).calc_voting_vector()[0]
voting_vector = [max(0, round(v)) for v in voting_vector]
self.voting_vector = [self.model.Var(v, lb=0, ub=our_NEOs, integer=True) for v in voting_vector]
def calc_GAS_gain(self):
# new_vote_result = [v+V for v, V in zip(self.voting_vector, self.candidate_votes)]
# sorted_vote_result = [(v, V) for v, V in zip(self.voting_vector, new_vote_result)] # sorted(, key=lambda x: x[1].value)
#
# sorted_voting_vector = [v for v, V in sorted_vote_result]
# sorted_vote_result = [V for v, V in sorted_vote_result]
# return sum([v*r/(V+v+0.0001) for v, V, r in zip(sorted_voting_vector, sorted_vote_result, self.reward_factor)])
return sum([v*r/(V+v+1) for v, V, r in zip(self.voting_vector, self.candidate_votes, self.reward_factor)])
def calc_voting_vector(self):
# constraints
# sum(voting_vector) == our_NEOs
self.model.Equation(sum(self.voting_vector) == self.our_NEOs)
# agent number
nonzero_vote = sum([self.model.if3(v - 0.1, 0, 1) for v in self.voting_vector])
# self.model.Equation(nonzero_vote >= 0)
self.model.Equation(nonzero_vote <= self.num_of_agents)
# objective
self.model.Maximize(self.calc_GAS_gain())
# self.model.options.RTOL = 1e-14
self.model.options.OTOL = 1e-14
self.model.options.REQCTRLMODE = 1
self.model.solve(disp=False, debug=False, GUI=False)
return self.voting_vector, -self.model.options.OBJFCNVAL
if __name__ == '__main__':
print(
BurgerNEOMultiAgentStrategy(
[11409980, 3024833, 2200734, 2040928, 2032100, 2017067, 2010927, 2007339, 2001625, 1008817, 1000025, 1000000, 641913, 396926, 392171, 390364, 383919, 383722, 382447, 382310, 381850, 923, 272, 180, 0],
1162, 6).calc_voting_vector()
)
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 0,0,0,0, 0,0], 50, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 0,0], 25, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,99, 75,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,1], 25, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,1], 25, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,0], 25, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 1,0], 1000, 6).calc_voting_vector())
print(BurgerNEOMultiAgentStrategy([100,100,100,100,100, 100,5, 5,5,5, 5,5,5,5,5, 5,5,5,5,5, 5, 1,1,1,1, 0,0], 26, 6).calc_voting_vector())
print()
而且我还有一些代码（此处未显示），它可以让 1 或 2 个代理的一切正确。这不是通过优化实现的，而是通过手动数学推导实现的。
我会这样写：
* This is not using proper GAMS syntax:
*is_nonzero_eqn(i).. is_nonzero(i) =e= 1 $ (vote(i));
* a linear version is:
is_nonzero_eqn(i).. vote(i) =l= is_nonzero(i)*vote.up(i);
这是一个标准的 MIP 公式。使用：
Model NEOBurger /all/;
Option MINLP=baron, optcr = 0;
Solve NEOBurger using MINLP maximizing GAS_reward;
我们看到：
GAMS/BARON 36.1.0 r2c0a44a Released Aug 02, 2021 WEI x86 64bit/MS Window
===========================================================================
BARON version 21.1.13. Built: WIN-64 Wed Jan 13 16:04:12 EST 2021
BARON is a product of The Optimization Firm.
For information on BARON, see https://minlp.com/about-baron
If you use this software, please cite publications from
https://minlp.com/baron-publications, such as:
Kilinc, M. and N. V. Sahinidis, Exploiting integrality in the global
optimization of mixed-integer nonlinear programming problems in BARON,
Optimization Methods and Software, 33, 540-562, 2018.
===========================================================================
This BARON run may utilize the following subsolver(s)
For LP/MIP/QP: ILOG CPLEX
For NLP: MINOS, SNOPT, GAMS external NLP, IPOPT, FILTERSD, FILTERSQP
===========================================================================
Doing local search
Solving bounding LP
Preprocessing found feasible solution with value 0.203660597686E-001
Starting multi-start local search
Done with local search
===========================================================================
Iteration Open nodes Time (s) Lower bound Upper bound
* 1 1 0.50 0.303390 3.71066
1 1 0.56 0.303390 0.303737
* 2 2 0.62 0.303405 0.303737
* 2 2 0.69 0.303428 0.303737
* 4 2 0.83 0.303473 0.303737
* 4 1 0.88 0.303639 0.303737
* 5 1 0.92 0.303681 0.303737
5 0 0.92 0.303681 0.303681
Calculating duals
*** Normal completion ***
Wall clock time: 4.05
Total CPU time used: 0.94
Total no. of BaR iterations: 5
Best solution found at node: 5
Max. no. of nodes in memory: 2
All done
===========================================================================
Solution = 0.30368112407161 found at node 5
Best possible = 0.303681125072
Absolute gap = 1.00039027062238E-9 optca = 1E-9
Relative gap = 3.29421287011233E-9 optcr = 1E-9
(Note that BARON uses a different formula to compute the relative gap as
was used for the above reported value.)
如果我们显示我们看到的相关变量：
---- 90 VARIABLE vote.L voting vector that needs to be optimized
19 466.000, 20 696.000
---- 90 VARIABLE is_nonzero.L whether an element in the voting vector is zero
19 1.000, 20 1.000
这是正确的行为。
请注意 MINLP = NLP + MIP，因此我们可以（并且应该）在编写 MINLP 模型时使用我们了解的有关 MIP 建模的所有工具和技巧。