Resolvendo o Set Cover com OpenSolver e Pyomo

A ideia deste exercício é resolver o problema da cobertura de conjuntos, utilizando o Pyomo (http://www.pyomo.org/) e o OpenSolver (www.opensolver.org).

Link para download: https://github.com/asgunzi/SetCover.

O Set cover é o problema de cobertura de conjuntos.

Tenho uma quantidade de elementos, digamos 50.

E tenho uma quantidade de conjuntos, digamos 15.

Cada conjunto tem um preço (abaixo está como peso), e cada conjunto contém os elementos.

Algumas aplicações: imagine que há várias doenças, e uma série de vacinas. Cada vacina cobre um conjunto diferente de doenças, e cada vacina tem um custo. Minimizar o custo total, cobrindo todas as doenças.

Open Solver:

Variável de decisão binária: utilizo ou não o conjunto.

Função objetivo: se utilizo o conjunto, é o preço deste.

Restrição: devo cobrir todos os elementos.

Resultado: FO de 223 e as escolhas pintadas em rosa.

Formulação Pyomo:

Leitura de dados utilizando o módulo pandas.

    df1 = pd.read_excel(io=’SetCover_Instancia01.xlsm’, sheet_name=’Base’, usecols = “B:P”, skiprows = 1, nrows = 1,header = None)

    df2 = pd.read_excel(io=’SetCover_Instancia01.xlsm’, sheet_name=’Base’, usecols = “B:P”, skiprows = 7, header = None)

    Nsets =  len(df1.columns)

    Nelementos =  len(df2)

Crio o model com o primeiro comando abaixo. Defino os índices (conjuntos e elementos).

   #—————- Modelo Pyomo ————#     

    model = ConcreteModel()

    model.S = list(range(0,Nsets)) #Num sets

    model.E = list(range(0,Nelementos)) #Num elementos

Defino a variável de decisão com o comando Var. O nome da variável é AtivaSet, ela é binária e inicializa com o valor zero.

    model.AtivaSet = Var(model.S, within=Binary, initialize=0)

Função objetivo.

Soma do produto entre a variável de decisão e o peso de cada conjunto.

def obj_rule(model):

        return sum(model.AtivaSet[s]*pesos[s] for s in model.S)

model.obj = Objective(rule = obj_rule, sense=minimize)

Restrição.

Para cada elemento, ele deve pertencer a um conjunto, ou seja, a soma do produto da matriz de alocação e a variável de decisão devem ser maior do que 1.

    def constrAlocacao(model, e):

        return sum(model.AtivaSet[s]*Malocacao[e][s] for s in model.S) >= 1

    model.restr1 = Constraint(model.E, rule = constrAlocacao)      

Algumas configurações do solver, e depois mando resolver

    SOLVER_NAME = ‘cbc’

    TIME_LIMIT = 60*1   #Em segundos

    solver = SolverFactory(SOLVER_NAME)

    solver.options[‘seconds’] = TIME_LIMIT

    results  = solver.solve(model)

    results.write() #para mostrar o status da solucao

O resultado obtido foi o mesmo do OpenSolver.

Além do CBC, outro solver free é o GLPK.

Todos os outros são pagos. O Pyomo pode chamar CPLEX, Gurobi, se tiver a licença destes.

Toda a manipulação de entrada de dados, matrizes, e saída de dados, é Python puro.

Vantagens: é free, e muito melhor que o Open Solver para manipular grande quantidade de dados.

Desvantagens: requer um grau maior de conhecimento que o OpenSolver. Não é tão fácil de manipular ou tão poderoso quanto um AIMMS + CPLEX.

Port in a storm

Here’s follow the solution (aimms) to the challenge Port in a Storm (puzzleor).

We need to create the sets docks and groups.

The distances, given in the spreadsheet:

Tip: in the definition, we can fill the matrix with the values.

Another tip: it is easy to exchange indexes. Make sure it does not happen.

The variable has range integer in dimensions (g,d)

Contraints: max docks and the minimum we have to distribute from the groups

Next, the mathematical programming.

We’re minimizing.

FO is the sum of variables times distances.

In the main execution, we write solve port (name of the mathematical programming);

Type F6 to run:

Result is FO = 31, the same as in the excel.

Note that the assignment is different from excel – in this case there’s more than one solution.

The zip file here contains the formulation. And here, the Excel solver version.

Below transcribed the full formulation.

MAIN MODEL Main_Port_in_a_storm

DECLARATION SECTION

SET:

   identifier   :  Docks

   subset of    :  Integers

   index        :  d

   definition   :  {1..3} ;



SET:

   identifier   :  Groups

   subset of    :  Integers

   index        :  g

   definition   :  {1..3} ;



PARAMETER:

   identifier   :  distances

   index domain :  (g,d)

   definition   :  data { ( 1, 1 ) : 1,  ( 1, 2 ) : 3,  ( 1, 3 ) : 2,  ( 2, 1 ) : 2,  ( 2, 2 ) : 2,  ( 2, 3 ) : 1,  ( 3, 1 ) : 3,  ( 3, 2 ) : 2,  ( 3, 3 ) : 1 } ;



PARAMETER:

   identifier   :  max_docks

   index domain :  d

   definition   :  data { 1 : 4,  2 : 7,  3 : 9 } ;



PARAMETER:

   identifier   :  min_groups

   index domain :  (g)

   definition   :  data { 1 : 8,  2 : 5,  3 : 7 } ;



VARIABLE:

   identifier   :  aloc

   index domain :  (g,d)

   range        :  integer ;



VARIABLE:

   identifier   :  fo

   range        :  free

   definition   :  sum((g,d),aloc(g, d)*distances(g, d)) ;



CONSTRAINT:

   identifier   :  c_docks

   index domain :  d

   definition   :  sum(g, aloc(g,d))<=max_docks(d) ;



CONSTRAINT:

   identifier   :  c_groups

   index domain :  (g)

   definition   :  sum(d, aloc(g,d))>=min_groups(g) ;



MATHEMATICAL PROGRAM:

   identifier   :  Port

   objective    :  fo

   direction    :  minimize

   constraints  :  AllConstraints

   variables    :  AllVariables

   type         :  Automatic ;

ENDSECTION ;

PROCEDURE

identifier :  MainInitialization

ENDPROCEDURE ;

PROCEDURE

identifier :  MainExecution

body       : 

  solve port;

ENDPROCEDURE ;

PROCEDURE

identifier :  MainTermination

body       : 

  return DataManagementExit();

ENDPROCEDURE ;

ENDMODEL Main_Port_in_a_storm ;

Ideias técnicas com uma pitada de filosofia:

https://ideiasesquecidas.com/

A simple exercise in AIMMS

I want to model the following Linear program (in the spreadsheet, using solver).

First of all, you’ll need to create a new Project, and one Set for each dimension.

Two dimensions = two sets. Each of this varying from 1 to 2.

We can write the data directly using {1,2} in definition area.

Now, you’ll have to create the arrays of data.

MatrizDados(i,j)

Rhs(i)  !Right hand side

c(j) !Coefficients in the FO

Next, create the variable. In this case, the name given was “var” and the dimension is j.

(You can name the dimensions as you wish, but do not make confusion with these).

Now, the constraint equations.

Note that the index domain is i. There is one equation like this for each line.

Another variable is the objective function:  sum(j,c(j)*var(j));

The mathematical programming is where we indicate the direction, the objective function and the type of problem (LP, MIP).

Last step. In the Main Execution, write down “solve mp;” – it is the command to solve the mathematical program we’ve called mp.

Click F6 to run.

The variable var(j) will contain the results: {12.5, 0} – the same as in the excel solver.

In this link, the aimms model. Perhaps it will not run in an older version of the software.

Regards.

Model Main_Teste01 {

    Set Set1 {

        SubsetOf: Integers;

        Index: i;

        Definition: {

            {1,2}

        }

    }

    Set Set2 {

        SubsetOf: Integers;

        Index: j;

        Definition: {

            {1,2}

        }

    }

   Parameter MatrizDados {

        IndexDomain: (i,j);

        Definition: data { ( 1, 1 ) : 10,  ( 1, 2 ) : 30,  ( 2, 1 ) : 20,  ( 2, 2 ) : 40 };

    }

    Parameter c {

        IndexDomain: j;

        Definition: data { 1 : 12,  2 : 23 };

    }

    Parameter rhs {

        IndexDomain: i;

        Definition: data { 1 : 150,  2 : 250 };

    }

    Variable var {

        IndexDomain: j;

        Range: nonnegative;

    }

    Constraint restr {

        IndexDomain: i;

        Definition: sum(j, var(j)*MatrizDados(i,j))<= rhs(i);

    }

    Variable fo {

        Range: free;

        Definition: sum(j,c(j)*var(j));

    }

    MathematicalProgram mp {

        Objective: fo;

        Direction: minimize;

        Constraints: AllConstraints;

        Variables: AllVariables;

        Type: LP;

    }

    Procedure MainInitialization;

    Procedure MainExecution {

        Body: {

            solve mp;

        }

    }

    Procedure MainTermination {

        Body: {

            return DataManagementExit();

        }

    }

}