Trying to build a pyomo model to solve "SEND+MORE=MONEY" task.
C4 C3 C2 C1
S E N D
+ M O R E
______________
M O N E Y
Model building:
from pyomo.environ import *
m = ConcreteModel()
m.letters = Set(initialize = ['S', 'E', 'N', 'D', 'M', 'O', 'R', 'Y'])
m.tens = Set(initialize = ['C1', 'C2', 'C3', 'C4'])
m.let_val = Var(m.letters, bounds=(0,9), initialize=9, within=Integers)
m.ten_val = Var(m.tens, within=Binary)
m.C1 = Constraint(expr = m.ten_val['C3'] + m.let_val['S'] + m.let_val['M'] == m.let_val['O'] + m.ten_val['C4']*10)
m.C2 = Constraint(expr = m.ten_val['C2'] + m.let_val['E'] + m.let_val['O'] == m.let_val['N'] + m.ten_val['C3']*10)
m.C3 = Constraint(expr = m.ten_val['C1'] + m.let_val['N'] + m.let_val['R'] == m.let_val['E'] + m.ten_val['C2']*10)
m.C4 = Constraint(expr = m.let_val['D'] + m.let_val['E'] == m.let_val['Y'] + m.ten_val['C1']*10)
m.C5 = Constraint(expr = m.let_val['M'] == m.ten_val['C4'])
m.C6 = Constraint(expr = m.let_val['M'] == 1)
m.f1 = Objective(expr = m.let_val['S']*1000 + m.let_val['E']*100 + m.let_val['N']*10 + m.let_val['D'] +
m.let_val['M']*1000 + m.let_val['O']*100 + m.let_val['R']*10 + m.let_val['E'], sense=minimize)
solver = SolverFactory('glpk')
results = solver.solve(m)
print(' ', int((value(m.let_val['S'])*1000 + value(m.let_val['E'])*100 + value(m.let_val['N'])*10 + value(m.let_val['D']))))
print('+ ', int((value(m.let_val['M'])*1000 + value(m.let_val['O'])*100 + value(m.let_val['R'])*10 + value(m.let_val['E']))))
print('=', int((value(m.let_val['M'])*10000 + value(m.let_val['O'])*1000 + value(m.let_val['N'])*100 + value(m.let_val['E'])*10 + value(m.let_val['Y']))))
'glpk' solution:
9000
+ 1000
= 10000
'ipopt' solution:
8942
+ 1057
= 10000
Ipopt is closer but not good enough. I had to help solver and add to constraints (C5 and C6). Without it the answer is 0.
I couldn't find any information about how to declare instructions that all letter values should be different from each other. How it can be done?
How to solve such a task with pyomo correctly?
Here is one solution:
from pyomo.environ import *
m = ConcreteModel()
m.letters = Set(initialize = ['S', 'E', 'N', 'D', 'M', 'O', 'R', 'Y'])
m.vals = RangeSet(0, 9)
m.tens = Set(initialize = ['C1', 'C2', 'C3', 'C4'])
m.ten_val = Var(m.tens, within=Binary, initialize=0)
m.x = Var(m.letters, m.vals, within=Binary, initialize=1)
def c_rule(m, i):
return sum(m.x[i, j] for j in m.vals) <= 1
m.col = Constraint(m.letters, rule=c_rule) # not more than one value for the letter
def r_rule(m, j):
return sum(m.x[i, j] for i in m.letters) <= 1
m.row = Constraint(m.vals, rule=r_rule) # not more than one value for the digit
m.C1 = Constraint(expr = m.ten_val['C4'] == m.x['M', 1])
m.C2 = Constraint(expr = m.x['M', 1] == 1)
m.C3 = Constraint(expr = m.ten_val['C3'] + sum(m.x['S', j]*j for j in m.vals) + m.x['M', 1] == sum(m.x['O', j]*j for j in m.vals) + m.ten_val['C4']*10)
m.C4 = Constraint(expr = m.ten_val['C2'] + sum(m.x['E', j]*j for j in m.vals) + sum(m.x['O', j]*j for j in m.vals) == sum(m.x['N', j]*j for j in m.vals) + m.ten_val['C3']*10)
m.C5 = Constraint(expr = m.ten_val['C1'] + sum(m.x['N', j]*j for j in m.vals) + sum(m.x['R', j]*j for j in m.vals) == sum(m.x['E', j]*j for j in m.vals) + m.ten_val['C2']*10)
m.C6 = Constraint(expr = sum(m.x['D', j]*j for j in m.vals) + sum(m.x['E', j]*j for j in m.vals) == sum(m.x['Y', j]*j for j in m.vals) + m.ten_val['C1']*10)
m.f1 = Objective(expr = sum(m.x['S', j]*j for j in m.vals), sense=minimize) # we are searching for working solution, no matter what to minimize
solver = SolverFactory('glpk')
results = solver.solve(m);
x_leters = [i for i, v in m.x.get_values().items() if v == 1]
digits = {}
for i in m.letters:
val = 0
for j in x_leters:
if i in j:
val = j[1]
digits[i] = val
print(' ',digits['S'],digits['E'],digits['N'],digits['D'])
print('+',digits['M'],digits['O'],digits['R'],digits['E'])
print(digits['M'],digits['O'],digits['N'],digits['E'],digits['Y'])
9 5 6 7
+ 1 0 8 5
1 0 6 5 2
Great job!