diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..fc991d4 --- /dev/null +++ b/.gitignore @@ -0,0 +1,5 @@ +solveurGLPK/ +*.stdout_glpk +*.stdout_lp +*.lp +*.sol diff --git a/distrib.py b/distrib.py new file mode 100644 index 0000000..353421f --- /dev/null +++ b/distrib.py @@ -0,0 +1,29 @@ +import os + +result = [os.path.join(dp, f) for dp, dn, filenames in os.walk("data") for f in filenames if os.path.splitext(f)[1] == '.opb'] + +with open("machines.txt") as file: + machines = [machine.rstrip() for machine in file.readlines()] + +# i = 0 +# while len(result) > 0: +# knap = result.pop()[:-4] +# victime = machines[i % len(machines)] +# i+=1 +# cmd = f"""ssh lfainsin@{victime} -o "ConnectTimeout 3" -o "StrictHostKeyChecking no" -o "UserKnownHostsFile /dev/null" "cd 2A/RO/tp2/ && tmux new-session -d 'julia test.jl {knap}.opb > {knap}.stdout_jl'" """ +# print(cmd) +# os.system(cmd) + +# i = 0 +# while len(result) > 0: +# knap = result.pop()[:-4] +# victime = machines[i % len(machines)] +# i+=1 +# cmd = f"""ssh lfainsin@{victime} -o "ConnectTimeout 3" -o "StrictHostKeyChecking no" -o "UserKnownHostsFile /dev/null" "cd 2A/RO/tp2/ && tmux new-session -d 'glpsol --lp {knap}.lp -o {knap}.sol > {knap}.stdout_lp'" """ +# print(cmd) +# os.system(cmd) + +# for machine in machines: +# cmd = f"""ssh lfainsin@{machine} -o "ConnectTimeout 3" -o "StrictHostKeyChecking no" -o "UserKnownHostsFile /dev/null" "tmux kill-server" """ +# print(cmd) +# os.system(cmd) diff --git a/install_glpk.sh b/install_glpk.sh new file mode 100755 index 0000000..beb7526 --- /dev/null +++ b/install_glpk.sh @@ -0,0 +1,22 @@ +mkdir solveurGLPK +cd solveurGLPK +mkdir sources +cd sources + +wget ftp://ftp.gnu.org/gnu/glpk/glpk-5.0.tar.gz +tar -xzvf glpk-5.0.tar.gz +cd glpk-5.0 + +./configure --disable-shared +make +#make --jobs=4 +make check + +cd ../.. +mkdir executables +cp sources/glpk-5.0/examples/glpsol executables/. +cp sources/glpk-5.0/doc/glpk.pdf executables/. +cp sources/glpk-5.0/doc/gmpl.pdf executables/. + +cd executables +./glpsol --help diff --git a/machines.txt b/machines.txt new file mode 100644 index 0000000..4d9414e --- /dev/null +++ b/machines.txt @@ -0,0 +1,252 @@ +ablette.enseeiht.fr +acdc.enseeiht.fr +ackbar.enseeiht.fr +actinium.enseeiht.fr +ader.enseeiht.fr +aerosmith.enseeiht.fr +alose.enseeiht.fr +anguille.enseeiht.fr +antimoine.enseeiht.fr +apollinaire.enseeiht.fr +aragorn.enseeiht.fr +archimede.enseeiht.fr +arryn.enseeiht.fr +arwen.enseeiht.fr +aspicot.enseeiht.fr +aston.enseeiht.fr +atanasoff.enseeiht.fr +atkinson.enseeiht.fr +azote.enseeiht.fr +babagge.enseeiht.fr +backus.enseeiht.fr +banshee.enseeiht.fr +baratheon.enseeiht.fr +basilic.enseeiht.fr +bastie.enseeiht.fr +baudelaire.enseeiht.fr +behemot.enseeiht.fr +bilbo.enseeiht.fr +bishop.enseeiht.fr +bleriot.enseeiht.fr +bobafett.enseeiht.fr +boeing.enseeiht.fr +bonite.enseeiht.fr +boole.enseeiht.fr +boromir.enseeiht.fr +bouba.enseeiht.fr +boucher.enseeiht.fr +brassens.enseeiht.fr +bravoos.enseeiht.fr +brochet.enseeiht.fr +bulbizarre.enseeiht.fr +cable.enseeiht.fr +calimero.enseeiht.fr +candy.enseeiht.fr +carapuce.enseeiht.fr +carbone.enseeiht.fr +carmack.enseeiht.fr +carpe.enseeiht.fr +casimir.enseeiht.fr +chenipan.enseeiht.fr +chevesne.enseeiht.fr +chewie.enseeiht.fr +clapton.enseeiht.fr +clash.enseeiht.fr +clementine.enseeiht.fr +cobalt.enseeiht.fr +colossus.enseeiht.fr +cooper.enseeiht.fr +copernic.enseeiht.fr +corb.enseeiht.fr +cyclope.enseeiht.fr +dagobah.enseeiht.fr +darwin.enseeiht.fr +daurade.enseeiht.fr +dazzler.enseeiht.fr +deeppurple.enseeiht.fr +demusset.enseeiht.fr +descartes.enseeiht.fr +diabolo.enseeiht.fr +dijkstra.enseeiht.fr +doors.enseeiht.fr +dorne.enseeiht.fr +dragon.enseeiht.fr +dylan.enseeiht.fr +eagles.enseeiht.fr +edison.enseeiht.fr +einsten.enseeiht.fr +elrond.enseeiht.fr +eomer.enseeiht.fr +eowyn.enseeiht.fr +eperlan.enseeiht.fr +epica.enseeiht.fr +esteban.enseeiht.fr +ewok.enseeiht.fr +farman.enseeiht.fr +fermat.enseeiht.fr +ferre.enseeiht.fr +fletan.enseeiht.fr +fluor.enseeiht.fr +forge.enseeiht.fr +frodon.enseeiht.fr +galilee.enseeiht.fr +galium.enseeiht.fr +gambit.enseeiht.fr +gandalf.enseeiht.fr +gardon.enseeiht.fr +garros.enseeiht.fr +gautier.enseeiht.fr +gimli.enseeiht.fr +gobelin.enseeiht.fr +gobie.enseeiht.fr +goldorak.enseeiht.fr +gorgone.enseeiht.fr +gosling.enseeiht.fr +goujon.enseeiht.fr +greyjoy.enseeiht.fr +griffon.enseeiht.fr +grove.enseeiht.fr +guynemer.enseeiht.fr +havok.enseeiht.fr +hawking.enseeiht.fr +heidi.enseeiht.fr +hendrix.enseeiht.fr +hippogriffe.enseeiht.fr +hopper.enseeiht.fr +hugo.enseeiht.fr +hydre.enseeiht.fr +hypotrempe.enseeiht.fr +iceberg.enseeiht.fr +iode.enseeiht.fr +jabba.enseeiht.fr +jeangrey.enseeiht.fr +jobs.enseeiht.fr +jubele.enseeiht.fr +kenobi.enseeiht.fr +kepler.enseeiht.fr +kernighan.enseeiht.fr +kiss.enseeiht.fr +knuth.enseeiht.fr +ladyoscar.enseeiht.fr +lafontaine.enseeiht.fr +lamartine.enseeiht.fr +lando.enseeiht.fr +lannister.enseeiht.fr +latecoere.enseeiht.fr +legolas.enseeiht.fr +leia.enseeiht.fr +lindbergh.enseeiht.fr +liskov.enseeiht.fr +loureed.enseeiht.fr +lovelace.enseeiht.fr +luke.enseeiht.fr +magicarpe.enseeiht.fr +magma.enseeiht.fr +magneto.enseeiht.fr +malicia.enseeiht.fr +mallarme.enseeiht.fr +manticore.enseeiht.fr +marcheursblancs.enseeiht.fr +martell.enseeiht.fr +maupassant.enseeiht.fr +maya.enseeiht.fr +meeren.enseeiht.fr +melofee.enseeiht.fr +mermoz.enseeiht.fr +merry.enseeiht.fr +messerschmitt.enseeiht.fr +metallica.enseeiht.fr +minotaure.enseeiht.fr +mordocet.enseeiht.fr +motorhead.enseeiht.fr +murdock.enseeiht.fr +muse.enseeiht.fr +mystique.enseeiht.fr +neon.enseeiht.fr +nerophis.enseeiht.fr +neumann.enseeiht.fr +newton.enseeiht.fr +nickel.enseeiht.fr +nymphe.enseeiht.fr +ohm.enseeiht.fr +omble.enseeiht.fr +orphie.enseeiht.fr +oxygene.enseeiht.fr +palomete.enseeiht.fr +palpatine.enseeiht.fr +pascal.enseeiht.fr +pegase.enseeiht.fr +phenix.enseeiht.fr +phosphore.enseeiht.fr +piafabec.enseeiht.fr +pikachu.enseeiht.fr +pinkfloyd.enseeiht.fr +pippin.enseeiht.fr +poe.enseeiht.fr +polaris.enseeiht.fr +polonium.enseeiht.fr +portreal.enseeiht.fr +prevert.enseeiht.fr +psykokwak.enseeiht.fr +psylocke.enseeiht.fr +queen.enseeiht.fr +r2d2.enseeiht.fr +rattata.enseeiht.fr +rimbaud.enseeiht.fr +rocket.enseeiht.fr +rondoudou.enseeiht.fr +roucool.enseeiht.fr +rouget.enseeiht.fr +saintexupery.enseeiht.fr +salameche.enseeiht.fr +sam.enseeiht.fr +sand.enseeiht.fr +sandre.enseeiht.fr +sar.enseeiht.fr +sauvageons.enseeiht.fr +scorpion.enseeiht.fr +scoubidou.enseeiht.fr +shadowcat.enseeiht.fr +shannon.enseeiht.fr +silure.enseeiht.fr +snorki.enseeiht.fr +snow.enseeiht.fr +socrate.enseeiht.fr +sodium.enseeiht.fr +solo.enseeiht.fr +stallman.enseeiht.fr +stark.enseeiht.fr +stones.enseeiht.fr +succube.enseeiht.fr +sunfire.enseeiht.fr +superbus.enseeiht.fr +swartz.enseeiht.fr +tacaud.enseeiht.fr +tanche.enseeiht.fr +tao.enseeiht.fr +targaryen.enseeiht.fr +taupiqueur.enseeiht.fr +tesla.enseeiht.fr +titane.enseeiht.fr +tornade.enseeiht.fr +torvalds.enseeiht.fr +truite.enseeiht.fr +turing.enseeiht.fr +tyrell.enseeiht.fr +uranium.enseeiht.fr +vador.enseeiht.fr +vairon.enseeiht.fr +verlaine.enseeiht.fr +voisin.enseeiht.fr +volantis.enseeiht.fr +volta.enseeiht.fr +winterfell.enseeiht.fr +wright.enseeiht.fr +wyvern.enseeiht.fr +xorn.enseeiht.fr +yoda.enseeiht.fr +z6po.enseeiht.fr +zemanek.enseeiht.fr +zia.enseeiht.fr +zirconium.enseeiht.fr +zztop.enseeiht.fr diff --git a/notebook-exemple.ipynb b/notebook-exemple.ipynb new file mode 100644 index 0000000..0ca813d --- /dev/null +++ b/notebook-exemple.ipynb @@ -0,0 +1,366 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# TP 2-3 : Branch-and-bound applied to a knapsack problem" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Initialisation (à faire une seule fois)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import Pkg; \n", + "Pkg.add(\"GraphRecipes\");\n", + "Pkg.add(\"Plots\"); \n", + "using GraphRecipes, Plots #only used to visualize the search tree at the end of the branch-and-bound" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Récupération des données" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "function readKnaptxtInstance(filename)\n", + " price=[]\n", + " weight=[]\n", + " KnapCap=[]\n", + " open(filename) do f\n", + " for i in 1:3\n", + " tok = split(readline(f))\n", + " if (tok[1] == \"ListPrices=\")\n", + " for i in 2:(length(tok)-1)\n", + " push!(price,parse(Int64, tok[i]))\n", + " end\n", + " elseif(tok[1] == \"ListWeights=\")\n", + " for i in 2:(length(tok)-1)\n", + " push!(weight,parse(Int64, tok[i]))\n", + " end\n", + " elseif(tok[1] == \"Capacity=\")\n", + " push!(KnapCap, parse(Int64, tok[2]))\n", + " else\n", + " println(\"Unknown read :\", tok)\n", + " end \n", + " end\n", + " end\n", + " capacity=KnapCap[1]\n", + " return price, weight, capacity\n", + "end" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Procédure d'application des tests de sondabilités TA, TO et TR pour le cas de la relaxation linéaire" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "function TestsSondabilite_relaxlin(model2, x, varsbin, BestProfit, Bestsol)\n", + " TA, TO, TR = false, false, false\n", + " if (termination_status(model2) == MOI.INFEASIBLE)#Test de faisabilite\n", + " TA=true\n", + " println(\"TA\")\n", + " elseif (objective_value(model2) <= BestProfit) #Test d'optimalite\n", + " TO=true\n", + " println(\"TO\")\n", + " elseif ( prod(abs.([round.(v, digits=0) for v in value.(varsbin)]-value.(varsbin)) .<= fill(10^-5, size(varsbin))) \n", + " ) #Test de resolution\n", + " TR=true\n", + " println(\"TR\")\n", + " #if (value(benef) >= BestProfit)\n", + " if (objective_value(model2) >= BestProfit)\n", + " Bestsol = value.(x)\n", + " #BestProfit=value(benef)\n", + " BestProfit=objective_value(model2)\n", + " end\n", + " else\n", + " println(\"non sondable\")\n", + " end\n", + " TA, TO, TR, Bestsol, BestProfit\n", + "end" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Procédure de séparation et stratégie d'exploration permettant de se placer au prochain noeud à traiter" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "\n", + "function SeparerNoeud_relaxlin(varsshouldbebinary, listvars, listvals)\n", + " # le noeud est non-sondable. Appliquer le critère de séparation pour le séparer en sous-noeuds \n", + " # et choisir un noeud-fils le plus à gauche \n", + " \n", + " #find a fractionnal variable\n", + " i, var = 1, 0\n", + " while ((i <= length(varsshouldbebinary)) && (var==0))\n", + " #if (varsshouldbebinary[i] ∉ listvars)\n", + " if (abs(round(value(varsshouldbebinary[i]), digits=0) - value(varsshouldbebinary[i]) ) >= 10^-5)\n", + " var=varsshouldbebinary[i]\n", + " end\n", + " i+=1\n", + " end\n", + " \n", + " #=\n", + " #find most fractionnal variable ?\n", + " i, var, maxfrac = -1, 0, 0.0\n", + " for i in 1:length(varsshouldbebinary)\n", + " if (abs(round(value(varsshouldbebinary[i]), digits=0) - value(varsshouldbebinary[i]) ) >= maxfrac) \n", + " #if a variable is more fractinonal\n", + " var=varsshouldbebinary[i]\n", + " maxfrac=abs(round(value(varsshouldbebinary[i]), digits=0) - value(varsshouldbebinary[i]) )\n", + " #println(i, \" \", var, \" \", maxfrac)\n", + " end\n", + " end\n", + " =#\n", + " \n", + "\n", + " set_lower_bound(var,1.0)\n", + " set_upper_bound(var,1.0)\n", + "\n", + " push!(listvars,var) #stocker l'identite de la variable choisie pour la séparation\n", + " push!(listvals,1.0) #stocker la branche choisie, identifiee par la valeur de la variable choisie\n", + " listvars, listvals\n", + "end\n", + "\n", + "\n", + "function ExplorerAutreNoeud_relaxlin(listvars, listvals, listnodes)\n", + " #this node is sondable, go back to parent node then right child if possible\n", + " \n", + " stop=false\n", + " #check if we are not at the root node\n", + " if (length(listvars)>= 1)\n", + " #go back to parent node\n", + " var=pop!(listvars)\n", + " theval=pop!(listvals)\n", + " tmp=pop!(listnodes)\n", + " set_lower_bound(var,0.0)\n", + " set_upper_bound(var,1.0)\n", + "\n", + " #go to right child if possible, otherwise go back to parent\n", + " while ( (theval==0.0) && (length(listvars)>= 1))\n", + " var=pop!(listvars)\n", + " theval=pop!(listvals)\n", + " tmp=pop!(listnodes)\n", + " set_lower_bound(var,0.0) \n", + " set_upper_bound(var,1.0)\n", + " end\n", + " if theval==1.0\n", + " set_lower_bound(var,0.0)\n", + " set_upper_bound(var,0.0)\n", + " push!(listvars,var)\n", + " push!(listvals,0.0)\n", + " else\n", + " println(\"\\nFINISHED\")\n", + " stop=true\n", + " end\n", + " else\n", + " #the root node was sondable\n", + " println(\"\\nFINISHED\")\n", + " stop=true\n", + " end\n", + " listvars, listvals, listnodes, stop \n", + "end" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Création de la relaxation linéaire (= modèle associé au noeud 0): SECTION A SUPPRIMER !!!! \n", + "\n", + " Cette section est à commenter/supprimer et remplacer par vos propres calculs de bornes supérieures et autres, par exemple basées sur les bornes 1 et 2 vues en cours, ou d'autres calculs de bornes de votre choix/conception validés au préalable par votre encadrant/e de TP " + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "Pkg.add(\"Clp\");\n", + "using JuMP, Clp" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "function CreationModeleLP(price, weight, capacity)\n", + "# ROOT NODE\n", + " \n", + " model2 = Model(Clp.Optimizer) # set optimizer\n", + " set_optimizer_attribute(model2, \"LogLevel\", 0) #don't display anything during solve\n", + " set_optimizer_attribute(model2, \"Algorithm\", 4) #LP solver chosen is simplex\n", + "\n", + " # define x variables as CONTINUOUS (recall that it is not possible to define binary variables in Clp)\n", + " @variable(model2, 0 <= x[i in 1:4] <= 1)\n", + " varsshouldbebinary=[x[1] x[2] x[3] x[4]]\n", + "\n", + " # define objective function\n", + " @objective(model2, Max, sum(price[i]*x[i] for i in 1:4))\n", + "\n", + " # define the capacity constraint \n", + " @constraint(model2, sum(weight[i]*x[i] for i in 1:4) <= capacity)\n", + "\n", + " println(model2)\n", + "\n", + " return model2, x, varsshouldbebinary\n", + "end\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Boucle principale : résoudre la relaxation linéaire, appliquer les tests de sondabilité, identifier le prochain noeud, répéter." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "\n", + "function SolveKnapInstance(filename)\n", + "\n", + " if (split(filename,\"/\")[end] != \"test.opb\")\n", + " println(\"This version of the code works only for the test instance !!!!\")\n", + " else\n", + " price, weight, capacity = readKnaptxtInstance(filename)\n", + " model2, x, varsshouldbebinary = CreationModeleLP(price, weight, capacity)\n", + " \n", + " #create the structure to memorize the search tree for visualization at the end\n", + " trParentnodes=Int64[] #will store orig node of arc in search tree\n", + " trChildnodes=Int64[] #will store destination node of arc in search tree\n", + " trNamenodes=[] #will store names of nodes in search tree\n", + " \n", + " #intermediate structure to navigate in the search tree\n", + " listvars=[]\n", + " listvals=[]\n", + " listnodes=[]\n", + "\n", + " BestProfit=-1\n", + " Bestsol=[]\n", + "\n", + " current_node_number=0\n", + " stop = false\n", + "\n", + " while (!stop)\n", + "\n", + " println(\"\\nNode number \", current_node_number, \": \\n-----\\n\", model2)\n", + "\n", + " #Update the search tree\n", + " push!(trNamenodes,current_node_number+1) \n", + " if (length(trNamenodes)>=2)\n", + " push!(trParentnodes,listnodes[end]+1) # +1 because the 1st node is \"node 0\"\n", + " push!(trChildnodes, current_node_number+1) # +1 because the 1st node is \"node 0\"\n", + " end\n", + " push!(listnodes, current_node_number)\n", + "\n", + "\n", + " print(\"Solve model2 to compute the bounds of the current node: start ... \")\n", + " status = optimize!(model2)\n", + " println(\"... end\")\n", + "\n", + " print(\"\\nSolution relax lin\"); \n", + " if (termination_status(model2) == MOI.INFEASIBLE)#(has_values(model2))\n", + " print(\" : NOT AVAILABLE (probably infeasible or ressources limit reached)\")\n", + " else\n", + " [print(\"\\t\", name(v),\"=\",value(v)) for v in all_variables(model2)] \n", + " end\n", + " println(\" \"); println(\"\\nPrevious Solution memorized \", Bestsol, \" with bestprofit \", BestProfit, \"\\n\")\n", + "\n", + " TA, TO, TR, Bestsol, BestProfit = TestsSondabilite_relaxlin(model2, x, varsshouldbebinary, BestProfit, Bestsol)\n", + "\n", + " is_node_sondable = TA || TO || TR\n", + "\n", + " if (!is_node_sondable)\n", + " listvars, listvals = SeparerNoeud_relaxlin(varsshouldbebinary, listvars, listvals)\n", + " else\n", + " listvars, listvals, listnodes, stop = ExplorerAutreNoeud_relaxlin(listvars, listvals, listnodes)\n", + " end\n", + "\n", + " current_node_number = current_node_number + 1\n", + " end\n", + "\n", + " println(\"\\n******\\n\\nOptimal value = \", BestProfit, \"\\n\\nOptimal x=\", Bestsol)\n", + "\n", + " return BestProfit, Bestsol, trParentnodes, trChildnodes, trNamenodes\n", + " end\n", + "\n", + "end\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Affichage du résultat final" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "BestProfit, Bestsol, trParentnodes, trChildnodes, trNamenodes = SolveKnapInstance(\"data/test.opb\")\n", + "println(\"\\n******\\n\\nOptimal value = \", BestProfit, \"\\n\\nOptimal x=\", Bestsol)\n", + "graphplot(trParentnodes, trChildnodes, names=trNamenodes, method=:tree)\n", + "println(trParentnodes)\n", + "println(trChildnodes)\n", + "println(trNamenodes)" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Julia 1.6.3", + "language": "julia", + "name": "julia-1.6" + }, + "language_info": { + "file_extension": ".jl", + "mimetype": "application/julia", + "name": "julia", + "version": "1.6.3" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/notebook.ipynb b/notebook.ipynb new file mode 100644 index 0000000..2ddbd64 --- /dev/null +++ b/notebook.ipynb @@ -0,0 +1,967 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# TP 2-3 : Branch-and-bound applied to a knapsack problem" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Initialisation (à faire une seule fois)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import Pkg;\n", + "Pkg.add(\"GraphRecipes\");\n", + "Pkg.add(\"Plots\");\n", + "using GraphRecipes, Plots #only used to visualize the search tree at the end of the branch-and-bound" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Récupération des données" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "readKnapInstance" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "\"\"\"Open and read a KnapFile.\n", + "\n", + "Args: \\\\\n", + " - filename (String): the name of the file to read.\n", + "\n", + "Returns: \\\\\n", + " - price (Vector{Integer}): prices of items to put in the KnapSack. \\\\\n", + " - weight (Vector{Integer}): weights of items to put in the KnapSack. \\\\\n", + " - capacity (Integer): the maximum capacity of the KnapSack.\n", + "\"\"\"\n", + "function readKnapInstance(filename)\n", + " price = []\n", + " weight = []\n", + " capacity = -1\n", + " open(filename) do f\n", + " for i = 1:3\n", + " tok = split(readline(f))\n", + " if (tok[1] == \"ListPrices=\")\n", + " for i = 2:(length(tok)-1)\n", + " push!(price, parse(Int64, tok[i]))\n", + " end\n", + " elseif (tok[1] == \"ListWeights=\")\n", + " for i = 2:(length(tok)-1)\n", + " push!(weight, parse(Int64, tok[i]))\n", + " end\n", + " elseif (tok[1] == \"Capacity=\")\n", + " capacity = parse(Int64, tok[2])\n", + " else\n", + " println(\"Unknown read :\", tok)\n", + " end\n", + " end\n", + " end\n", + " return price, weight, capacity\n", + "end" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "# # on convert tous les .opb en .lp pour les utiliser dans glpk\n", + "# for (root, dirs, files) in walkdir(\"data\")\n", + "# for file in files\n", + "# if endswith(file, \".opb\")\n", + "# price, weight, capacity = readKnapInstance(root * \"/\" * file)\n", + "# filename = splitext(file)[1]\n", + "# f = open(root * \"/\" * filename * \".lp\", \"w\");\n", + "# write(f, \"Maximize\\n\")\n", + "# write(f, \" Knap: \")\n", + "# for (i, p) in enumerate(price)\n", + "# write(f, \"+ \" * string(p) * \" obj\" * string(i) * \" \")\n", + "# end\n", + "# write(f, \"\\n\")\n", + "\n", + "# write(f, \"\\nSubject To\\n\")\n", + "# write(f, \" MaxKnap: \")\n", + "# for (i, w) in enumerate(weight)\n", + "# write(f, \"+ \" * string(w) * \" obj\" * string(i) * \" \")\n", + "# end\n", + "# write(f, \"<= \" * string(capacity) * \"\\n\")\n", + "\n", + "# write(f, \"\\nBinary\\n\")\n", + "# for (i, p) in enumerate(price)\n", + "# write(f, \" obj\" * string(i) * \"\\n\")\n", + "# end\n", + "\n", + "# write(f, \"\\nEnd\\n\")\n", + "# close(f)\n", + "# end\n", + "# end\n", + "# end" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Procédure d'application des tests de sondabilités TA, TO et TR pour le cas de la relaxation linéaire" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "TestsSondabilite" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "\"\"\"Test if a node should be pruned.\n", + "\n", + "Args: \\\\\n", + " - x (Vector{Integer}): the node to be tested. \\\\\n", + " - price (Vector{Integer}): prices of items to put in the KnapSack. \\\\\n", + " - weight (Vector{Integer}): weights of items to put in the KnapSack. \\\\\n", + " - capacity (Integer): the maximum capacity of the KnapSack. \\\\\n", + " - BestProfit (Integer): the current BestProfit value. \\\\\n", + " - Bestsol (Integer): the current BestSol values. \\\\\n", + " - affich (Bool): determine if the function should print to stdout.\n", + "\n", + "Returns: \\\\\n", + " - TA (Bool): true if the node is feasible. \\\\\n", + " - TO (Bool): true if the node is optimal. \\\\\n", + " - TR (Bool): true if the node is resolvable. \\\\\n", + " - BestProfit (Integer): the updated value of BestProfit. \\\\\n", + " - Bestsol (Vector{Integer}): the updated values of BestSol.\n", + "\"\"\"\n", + "function TestsSondabilite(x, price, weight, capacity, BestProfit, Bestsol, affich)\n", + " TA, TO, TR = false, false, false\n", + "\n", + " if (!Constraints(x, weight, capacity)) # Test de faisabilite\n", + " TA = true\n", + " if affich\n", + " println(\"TA\\n\")\n", + " end\n", + "\n", + " elseif (Objective(x, price) <= BestProfit) # Test d'optimalite\n", + " TO = true\n", + " if affich\n", + " println(\"TO\\n\")\n", + " end\n", + "\n", + " elseif (AllDef(x)) # Test de resolution\n", + " TR = true\n", + " if affich\n", + " println(\"TR : solution \", \" de profit \", Objective(x, price), \"\\n\")\n", + " end\n", + " if (Objective(x, price) >= BestProfit) # Le profit de la solution trouvée est meilleur que les autres\n", + " if affich\n", + " println(\"\\t-> Cette solution a un meilleur profit.\\n\")\n", + " end\n", + " # On remplace la solution et le profit par les nouvelles valeurs\n", + " Bestsol = x\n", + " BestProfit = Objective(x, price)\n", + " else\n", + " if affich\n", + " println(\"\\t-> Cette solution est moins bonne.\\n\")\n", + " end\n", + " end\n", + "\n", + " elseif affich\n", + " println(\"non sondable\\n\")\n", + " end\n", + "\n", + " return TA, TO, TR, Bestsol, BestProfit\n", + "end" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Procédure de séparation et stratégie d'exploration permettant de se placer au prochain noeud à traiter" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "SeparerNoeud" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "\"\"\"Split a node in two.\n", + "\n", + "Args: \\\\\n", + " - price (Vector{Integer}): prices of items to put in the KnapSack. \\\\\n", + " - listvars (Vector{Vector{Integer}}): the current values of listvars. \\\\\n", + " - listvals (Vector{Integer}): the current values of listvals.\n", + "\n", + "Returns: \\\\\n", + " - listvars (Vector{Vector{Integer}}): the updated values of listvars. \\\\\n", + " - listvals (Vector{Integer}): the updated values of listvals.\n", + "\"\"\"\n", + "function SeparerNoeud(price, listvars, listvals)\n", + " # Le noeud est non-sondable. Appliquer le critère de séparation pour le séparer en sous-noeuds \n", + "\n", + " # Cas du noeud le plus à gauche\n", + "\n", + " # On sépare le noeud actuel en 2 sous-noeuds\n", + " predX = pop!(listvars)\n", + " nextX0 = copy(predX)\n", + " nextX1 = copy(predX)\n", + "\n", + " # On initialise leurs valeurs à zéro\n", + " val0 = 0\n", + " val1 = 0\n", + "\n", + " # On fixe la nouvelle variable des deux sous-noeuds\n", + " n = length(predX)\n", + " for i = 1:n\n", + " if predX[i] == -1\n", + " # L'un a zéro\n", + " nextX0[i] = 0\n", + " # L'autre a un\n", + " nextX1[i] = 1\n", + "\n", + " # On calcule leurs valeurs\n", + " val0 = Objective(nextX0, price)\n", + " val1 = Objective(nextX1, price)\n", + " break\n", + " end\n", + " end\n", + " \n", + " # On ajoute les sous-noeuds a la pile des noeuds a explorer\n", + " push!(listvars, nextX0)\n", + " push!(listvars, nextX1)\n", + "\n", + " # On ajoute aussi leurs valeurs\n", + " push!(listvals, val0)\n", + " push!(listvals, val1)\n", + "\n", + " return listvars, listvals\n", + "end" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "ExplorerAutreNoeud" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "\"\"\"Pop node fom the list to explore another node.\n", + "\n", + "Args: \\\\\n", + " - price (Vector{Integer}): prices of items to put in the KnapSack. \\\\\n", + " - listvars (Vector{Vector{Integer}}): the current values of listvars. \\\\\n", + " - listvals (Vector{Integer}): the current values of listvals.\n", + "\n", + "Returns: \\\\\n", + " - listvars (Vector{Vector{Integer}}): the updated values of listvars. \\\\\n", + " - listvals (Vector{Integer}): the updated values of listvals. \\\\\n", + " - stop (Bool): true if the tree search is finished.\n", + "\"\"\"\n", + "function ExplorerAutreNoeud(listvars, listvals)\n", + " # Le noeud est sondable, on l'enlève de la pile des noeuds à sonder\n", + "\n", + " stop = false\n", + " if (length(listvars) > 1)\n", + " # On passe au noeud suivant\n", + " var = pop!(listvars)\n", + " val = pop!(listvals)\n", + " else\n", + " # Il n'y a pas d'autre noeud\n", + " stop = true\n", + " end\n", + "\n", + " return listvars, listvals, stop\n", + "end" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Fonctions décrivant l'objectif et les contraintes" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "AllDef (generic function with 1 method)" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "# Fonction objectif que l'on souhaite maximiser/minimiser (évalué dans le meilleur des cas)\n", + "Objective(x, price) =\n", + " sum(\n", + " if x[i] < 0\n", + " price[i]\n", + " else\n", + " price[i] * x[i]\n", + " end\n", + " for i = 1:length(x)\n", + " )\n", + "\n", + "# Fonction permettant de vérfier toutes les contraintes du modèle (dans le meilleur des cas)\n", + "Constraints(x, weight, capacity) =\n", + " sum(\n", + " if x[i] < 0\n", + " 0\n", + " else\n", + " weight[i] * x[i]\n", + " end\n", + " for i = 1:length(x)\n", + " ) <= capacity\n", + "\n", + "# Fonction qui nous dis si toutes les variables de x sont fixées\n", + "function AllDef(x)\n", + " for i = 1:length(x)\n", + " if x[i] < 0\n", + " return false\n", + " end\n", + " end\n", + " return true\n", + "end" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Résolution du problème KnapSack" + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "SolveKnapInstance" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "\"\"\"Solve the KnapSack problem for the data contained in `filename`.\n", + "\n", + "Args: \\\\\n", + " - filename (String): the name of the file to read.\n", + "\n", + "Returns: \\\\\n", + " - trParentnodes (Vector{Integer}): the parents nodes, to plot the tree.\n", + " - trChildnodes (Vector{Integer}): the child nodes, to plot the tree.\n", + " - trNamenodes (Vector{Integer}): the name of the nodes, to plot the tree.\n", + "\"\"\"\n", + "function SolveKnapInstance(filename)\n", + "\n", + " stop = false\n", + " affich = true\n", + "\n", + " # Extraction des données\n", + " price, weight, capacity = readKnapInstance(filename)\n", + "\n", + " tri = true\n", + " if tri\n", + " couples = zip(price, weight)\n", + " couples = sort(collect(couples), by = x -> x[1] / x[2])\n", + " unzip(a) = map(x->getfield.(a, x), fieldnames(eltype(a)))\n", + " price, weight = unzip(couples)\n", + " end\n", + "\n", + " expected = split(split(filename, \"-\")[2], \".\")[1]\n", + " nodes_nb = length(price)\n", + " nodes_max = 2^nodes_nb\n", + "\n", + " if affich\n", + " println(\"---\", filename, \"---\")\n", + " println(\"Capacity = \", capacity)\n", + " println(\"Number of objects = \", nodes_nb)\n", + " println(\"Expected optimal value = \", expected)\n", + " println(\"Maximum number of nodes = \", nodes_max)\n", + " println() end\n", + "\n", + " # Pour dessiner le graph\n", + " trParentnodes = Int64[]\n", + " trChildnodes = Int64[]\n", + " trNamenodes = []\n", + "\n", + " # Liste des variable pour naviguer de noeuds en noeuds\n", + " listvars = []\n", + " listvals = []\n", + " listnodes = []\n", + "\n", + " # La meilleur solution et sa valeur\n", + " BestProfit = -1\n", + " Bestsol = []\n", + "\n", + " # Compter le nombre de noeud explorés\n", + " current_node_number = 0\n", + "\n", + " # On ajoute le premier noeud à explorer (la racine)\n", + " push!(listvars, [-1 for p in price])\n", + " push!(listvals, Objective(last(listvars), price))\n", + " push!(listnodes, 1)\n", + " push!(trNamenodes, 0)\n", + " newnodeid = 2\n", + "\n", + " while (!stop)\n", + "\n", + " # Le noeud actuel\n", + " x = last(listvars)\n", + "\n", + " if affich && current_node_number % 10000 == 0\n", + " println(\"Node n°\", current_node_number, \": \\tBestProfit = \", BestProfit)\n", + " LastBestProfit = BestProfit\n", + " end\n", + "\n", + " # Test de sondabilité du noeud actuel\n", + " # -> On mets a jour la solution et sa valeur si besoin\n", + " TA, TO, TR, Bestsol, BestProfit = TestsSondabilite(x, price, weight, capacity, BestProfit, Bestsol, false)\n", + "\n", + " is_node_sondable = TA || TO || TR\n", + " if (!is_node_sondable)\n", + " # Le noeud n'est pas sondable, on le sépare en 2 sous-noeuds\n", + " listvars, listvals = SeparerNoeud(price, listvars, listvals)\n", + "\n", + " curnode = pop!(listnodes)\n", + "\n", + " push!(trParentnodes, curnode)\n", + " push!(trParentnodes, curnode)\n", + "\n", + " push!(listnodes, newnodeid + 1)\n", + " push!(listnodes, newnodeid)\n", + "\n", + " push!(trChildnodes, newnodeid)\n", + " push!(trChildnodes, newnodeid + 1)\n", + "\n", + " push!(trNamenodes, newnodeid - 1)\n", + " push!(trNamenodes, newnodeid)\n", + "\n", + " newnodeid += 2\n", + "\n", + " else\n", + " # Le noeud est sondable, on passe au noeud suivant\n", + " listvars, listvals, stop = ExplorerAutreNoeud(listvars, listvals)\n", + "\n", + " pop!(listnodes)\n", + " end\n", + "\n", + " current_node_number += 1\n", + " end\n", + "\n", + " println(\"\\n--------------------------------------------------\")\n", + " println(\"Expected optimal value = \", expected)\n", + " println(\"Optimal value = \", BestProfit)\n", + " println(\"Optimal x = \", Bestsol)\n", + " println(\"Maximum number of nodes = \", nodes_max)\n", + " println(\"Nodes explored = \", current_node_number) \n", + " \n", + " return trParentnodes, trChildnodes, trNamenodes\n", + "end" + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "---data/test-65.opb---\n", + "Capacity = 10\n", + "Number of objects = 4\n", + "Expected optimal value = 65\n", + "Maximum number of nodes = 16\n", + "\n", + "Node n°0: \tBestProfit = -1\n", + "\n", + "--------------------------------------------------\n", + "Expected optimal value = 65\n", + "Optimal value = 65\n", + "Optimal x = [0, 1, 0, 1]\n", + "Maximum number of nodes = 16\n", + "Nodes explored = 23\n" + ] + }, + { + "data": { + "text/plain": [ + "([1, 1, 2, 2, 4, 4, 7, 7, 5, 5 … 3, 3, 14, 14, 17, 17, 15, 15, 20, 20], [2, 3, 4, 5, 6, 7, 8, 9, 10, 11 … 14, 15, 16, 17, 18, 19, 20, 21, 22, 23], Any[0, 1, 2, 3, 4, 5, 6, 7, 8, 9 … 13, 14, 15, 16, 17, 18, 19, 20, 21, 22])" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "trParentnodes, trChildnodes, trNamenodes = SolveKnapInstance(\"data/test-65.opb\")" + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdeTxU+/8H8PeMfcm+70uIlDUUrZRUVKpbifY97ft6q9v+bd+1qlSiVSrtCSnZKSFLoWTNbhgzvz/cnyskzJk5w7yfj/t4XM6c8/m8H92ul885n/P5UJhMJiCEEEK8ikp2AQghhBCZMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIEQIIcTTMAgRQgjxNAxChBBCPA2DECGEEE/DIESIK+Tk5AQFBeXm5pJdCEI8B4MQIfJNnzm7u77BzLnztXV1l69aTXY5CPEWCpPJJLsGhHjahw8fepuYMsbvAdWekPqG79GetM+pmpqaZNeFEK/AIESIZLaDh74JfQ0ny4BCgdpqWCwzasy4+zd9ya4LIV6Bt0YRIlN8fPz79+9BTAYoFAAAAWGg8j9/9ryiooLs0hDiFRiECJHJ2NhYW1sHqkr+/b6uFupqexj2EBUVJbUuhHgIBiFCZKJSqWdPHKbU1cLHZ5CfAVF3qFTqhVPHKPUDRIQQ++EzQoTI17f/wMiId3Sg8jPrhjo6Pbx/l+yKEOIhOCJEiHzhIcGrV64Q5aec8TqNKYgQh2EQIsQVqFSqpqamhIQE2YUgxHMwCBHiChEREQYGBunp6WQXghDPwSBEiHx0Oj0iIsLBwSEuLo7sWhDiORiECJHv3bt3Ojo6dnZ2MTExZNeCEM/BIESIfLdu3Ro9erSRkdGXL1/Ky8vJLgch3oJBiBDJmEzmnTt3XF1dBQQELC0tQ0NDya4IId6CQYgQyUJCQkRFRXv16gUA9vb2z58/J7sihHgLBiFCJPPy8po7d2791w4ODs+ePSO3HoR4Da4sgxCZ8vPz9fX109LSZGRkAKCurk5FRSUiIgK3YUKIY3BEiBCZjhw5MmnSpPoUBAA+Pj4XF5fbt2+TWxVCPAWDECHSlJWVeXl5rVy5svHBsWPHYhAixEkYhAiR5tixY46Ojt27d2980MHBISUlJS0tjayqEOI1+IwQIXIUFxf36NEjLCysSRACwNKlS6Wlpbdu3UpGXQjxHAxChMixatWqqqqqEydONP8oJiZm3Lhxnz9/plLxng1CbMdPdgEI8aKUlJTLly8nJCS0+KmZmZmMjMzjx4+dnJw4XBhCPAh/30SIBMuXL1+/fr2iouLvTvD09Dx+/DgnS0KIZ+GtUYQ47c6dO5s2bYqNjRUQEPjdOdXV1Zqami0+QUQIEQtHhAhxVGlp6dKlS48dO9ZKCgKAsLDwnDlzDh8+zLHCEOJZOCJEiKMWLFjAZDJPnz79xzPz8vJ69Ojx8eNHJSUlDhSGEM/CESFCnPP8+fMHDx7s3bu3LScrKChMnDjx5MmT7K4KIR6HI0KEOKS0tLR3796nTp1q+1zQ9PR0Gxub1NRUSUlJttaGEC/DIESIQ2bOnCkkJHTq1Kn2XqWlpbVlyxY2VYUQwiBEiBNu3bq1YcOGqKgocXHxdl2YlpZWPyiUkpJiU20I8TgMQoTYLicnx8LC4t69e9bW1h24fNasWaqqqtu3bye8MIQQYBAixG51dXX29vaOjo7r16/vWAs5OTmmpqYJCQk4fRQhdsBZowix1/bt26lU6po1azrcgqqqqoeHx65duwisCiHUAEeECLHRy5cv3d3do6KiWBzMFRQUGBoavn37VldXl6jaEEL1cESIELvk5uZ6eHhcvnyZ9VuacnJyK1as6PDNVYRQK3BEiBBbMBiM4cOH29ra/v3334Q0WF1d3aNHDx8fHzs7O0IaRAjVwxEhQmxR/+bf5s2biWpQWFj4n3/+Wb16Nf7yihCxMAgRIt7Dhw+vXLly9epVYnfWdXd3ZzAYV69eJbBNhBDeGkWIYJmZmTY2Nnfu3Onbty/hjb99+3bChAmfPn0SExMjvHGEeBOOCBEiUnV19YQJE9avX8+OFAQAGxubgQMH7tmzhx2NI8SbcESIEJHmzp1bWlrq6+vLvi5yc3N79+4dFhamp6fHvl4Q4h04IkSIMN7e3qGhoefOnWNrL0pKSitWrGDlDX2EUGM4IkSIGPHx8UOHDn316pWhoSG7+6LRaL169Tp27JijoyO7+0Koy8MRIUIEKCsrmzhx4oEDBziQggAgJCR0+PDhJUuW0Gg0DnSHUNeGI0KEWMVkMsePH6+ionLs2DFO9jt69Oi+ffuuW7eOk50i1PVgECLEqkOHDl2/fj0kJERISIiT/WZkZFhZWUVHR6urq3OyX4S6GAxChFjy9u3bMWPGvH37VktLi/O9b9u2LTEx0d/fn/NdI9Rl4DNChDquuLh48uTJZ86cISUFAWDdunVxcXFBQUGk9I5Q14AjQoQ6ztXVVVNT89ChQyTWEBgYuGrVqvj4eEFBQRLLQKjzwhEhQh108uTJrKysvXv3klvGqFGj9PT0yA1jhDo1HBEi1BEfP34cNGhQaGiovr4+2bVAenp6nz59oqOjNTU1ya4Foc4HR4QItRuNRpsyZcq+ffu4IQUBQEdHx9PTE9eaQahjcESIULutXr06MzOTq+ZqVlVV9ezZ8/z584MHDya7FoQ6GQxChNrnzZs348aNi4uLU1BQILuWX9y5c2fTpk2xsbECAgJk14JQZ4K3RhFqh8rKymnTpnl5eXFbCgLA2LFjNTQ0Tp8+TXYhCHUyOCJEqB1Wr179/ft3Hx8fsgtpWVJS0sCBAz9+/CgnJ0d2LQh1GhiECLVVXFyco6NjfHw8Fw4HGyxbtoxGo506dYrsQhDqNDAIEWoTBoNhZWW1dOlSDw8PsmtpTXFxsaGh4ZMnT3r37k12LQh1DviMEKE2OXv2rIiIiLu7O9mF/IG0tPSWLVtWrFhBdiEIdRo4IkToz4qKioyMjDrLMKuurs7ExGTv3r0jR44kuxaEOgEMQoT+bOXKldXV1SdOnCC7kLZ68ODBqlWrEhIS+Pn5ya4FIW6HQYjQH+Tk5JiZmcXHxyspKZFdSzs4ODhMmDBh3rx5ZBeCELfDIEToD2bNmqWqqrp9+3ayC2mfuLi4kSNHpqSkiIqKkl0LQlwNgxCh1mRlZZmZmaWmpkpLS5NdS7tNnDjRzMxs3bp1ZBeCEFfDIESoNUuXLhUWFiZ9r6WOSUlJsbOzS0lJkZKSIrsWhLgXBiFCv1VSUqKtrf3hwwdlZWWya+mg2bNnq6mpbd26lexCEOJeGIQI/daJEydev35948YNsgvpuC9fvvTp0yc5Obkz3tpFiDPwhXqEfuvMmTOdfdalpqami4vLkSNHyC4EIe6FI0KEWvbhw4cRI0ZkZmZSKBSya2FJWlpa375909PTxcXFya4FIW6EI0KEWubv7z9+/PjOnoIAoKurO2TIkLNnz5JdCEJcCkeECLXMwsLiyJEjdnZ2ZBdCgOjo6DFjxqSlpeGevQg1hyNChFpQXFz8+fNnGxsbsgshhrm5uZ6e3s2bN8kuBCFuhEGIUAvCwsJsbGy60kKdnp6enWitVIQ4CYMQoRa8f//e2tqa7CqI5OLikpOTExkZSXYhCHEdDEKEmsrMzHz9+rWuri7ZhRCJj49v9uzZ58+fJ7sQhLgOTpZB6Bc9jXtlfvlSWV7evbfF1Ml/bV63muyKCJOTk9O7d+/s7GwRERGya0GIi+CIEKH/JCcnZ+T8qJx1FYTFP1sv3bNv/8ePH8kuijCqqqoWFhYBAQFkF4IQd8EgROg/cz2X1ziuAePhwCcAhoOrh6+duXAJ2UURycPD49q1a2RXgRB3wSBE6F+pqanhoa/rBv63phpjwNz4+ISioiISqyKWi4vLq1evysrKyC4EIS6CQYjQv/T09HqbW1Lf+jQcoby9qte9e1darlpSUtLW1vbRo0dkF4IQF8EgROg/504cFg7cDt+SgFEHPz6LBG67cOpoF1hlrbHRo0cHBgaSXQVCXARnjSL0C+ex454+eUKrLBftJmnTt+/zx11t8JSWljZgwICcnByyC0GIW3SdhTMQIsTWTRuqykpev37tff6svLw82eUQT1dXV1BQ8NOnTz169CC7FoS4Ao4IEfoFg8FQVVWl0WhdaY5MEzNnzrS2tu7sWy0iRBR8RojQLx48eKCsrNy1t+6zsLCIjo4muwqEuAUGIUK/2Lt379SpU8XExMguhI3Mzc0xCBFqgEGI0H8CAgJ+/vxpa2srKipKdi1sZGJi8uHDBwaDQXYhCHEFDEKE/lVVVbVixYojR45UVFR069aN7HLYSFRUVFJSMjc3l+xCEOIKGIQI/WvNmjXW1tb29vZlZWVdOwg/ffokLi7u7+9Pp9PJrgUh8mEQIgQAEBgYeP/+/fqta0tKSiQlJcmuiF2OHj1qaWX99Vvuuo2bJKVlMjMzya4IIZJhECIEycnJs2bN8vX1lZKSAoDS0lIJCQmyi2KXoBeva2xn16x4Vr38WRVVZObcBWRXhBDJMAgRr8vPz3dxcdmzZ4+NjU39kZ8/f9YnYtcTGhoaHP6+1nkLaJqDqjGTXyg8MiYiIoLsuhAiEwYh4mnl5eUjRoyYOHHijBkzGg4WFxd3pYW2G9u1/3DlgAUgKAoA4Lcaeg2n2c6+GfCQ7LoQIhMGIeJdFRUVI0eOtLS03L59e+PjXXhEuMJzvuibc0Cnwd2tUJQF43eLvvcZNWwI2XUhRCZcaxTxqPLy8lGjRunr69dPkGmsqKhIRkaGlKrYzcHBwcrY4PUJVwatCjxO8gftG2xnM2DAALLrQohMOCJEvKioqMjBwaFHjx5eXl5UatP/C7rwrVEAGOfsxPzwTCArGrZbioSdmz/Dg+yKECIZBiHiOVlZWQMGDBg4cOCpU6eapyAAFBYWysrKcr4wzhg+fPiNGzcm/zXepKehopyMvr4+2RUhRDK8NYp4S3x8/KhRo5YtW7ZixYrfnVNUVNSFR4Tdu3fv3r17UVGRiIjIjx8/9uzZc/bs2RZ/IUCIR+DffsRDHj16NHTo0P3797eSgtDVR4T1xMXFy8rKrly58vnz5+nTp1dXV5NdEUKkwSBEvOLYsWOzZs26c+fOX3/91cpplZWVANC1d58AAAkJibKyMnFx8YcPH9bU1NjZ2X369InsohAiBwYh6vpqamrmzZt39uzZN2/e9OvXr/WTeWE4CACSkpI/f/4EADExMV9f3zlz5gwYMGDbtm0VFRVkl4YQp2EQoi4uLy/PwcEhLy/vzZs3Wlpafzy/oKBATk6O/XWRTFpauri4uOHbefPmRUVFJScnGxgYnDhxAu+UIp6CQYi6sqioKCsrq8GDB9+6dauNm84XFhZ21ZcIG5ORkSkqKmp8RF1d/dq1a/fu3Xvy5Imuru6uXbsKCgrIKg8hTsIgRF2Wj4+Pk5PToUOHtm3b1vZZkfn5+fLy8mwtjBvIysoWFhY2P25hYXHv3r3Hjx9nZGTo6el5eHiEhoZyvjyEOAmDEHVBdDp9xYoV27dvf/ny5dixY9t1bVFRES88IxQWFhYQECgtLW3xU2Nj47Nnz6alpZmbm8+dO9fAwGDXrl1fv37lcJEIcQYGIepqCgsLhw8fnpSU9O7du549e7b38vz8fF54RggACgoKeXl5rZwgIyOzfPnyjx8/Xr58OScnx9LSctCgQV5eXi0OJRHqvDAIUZeSmJhoZWVlaWkZGBjYsZfiu/BCo038MQgbWFtbnzhxIjs7e8WKFa9evdLV1R01apSPj09ZWRm7i0SIAzAIUddx7949e3v7HTt27Nmzh4+Pr2ONFBYW8siIUElJqY1BWE9QUNDFxeX69es5OTlubm5+fn7q6uoTJ068d+9eTU0N++pEiN0wCFEXsW/fPk9PzwcPHkyePJmVdnjkPUIAUFRU/PHjRwcuFBMTc3NzCwgISE9Pd3BwOHz4sKqq6qJFi3CDX9RJYRCiTo9Op8+ePdvX1zc8PNzS0pLF1ngnCJWUlL5//85KCzIyMnPmzHn58mV0dLSKioq7u3uvXr2OHz9eUlJCVJEIcQAGIercysvLnZ2dc3NzQ0JC1NTUWG+QR16oBwBlZWUWg7CBurr6xo0bk5OTjx07FhYWpqOjs2TJkrS0NEIaR4jdMAhRJ1ZQUDB48GB1dfW7d+8StTpo196MsDECg7AehUIZNGjQ9evXExMTJSQk+vbtO3XqVIxDxP0wCFFnlZOTM3DgwGHDhp05c4afn5gNxeh0ekVFhaSkJCGtcTkVFRVig7CBsrLyjh07Pn/+rK+vb2Njs3LlyvLycnZ0hBAhMAhRp5STkzN48OBp06bt3LmTwGaLi4slJSUpFAqBbXItFRWVnJwc9rUvISGxadOmDx8+FBcX9+rV6+XLl+zrCyFWYBCizic3N3fIkCFz585ds2YNsS0XFxfzyEuEAKCoqFhYWEin09nai4KCwoULF7y8vKZMmXLw4EG29oVQx+AO9aiTKSkpcXJy8vDwWLVqFTsa55H7ogDAx8cnJyf348cPVVVVdvc1bNiwiIgIR0fHvLy8tWvXNhwXEhISFRVld+8ItQ5HhKgzodPp48aNs7Oz27RpEzva//nzp5SUFDta5k4qKirfvn3jTF+ysrKFpZX7Dh9XUNdS1tT99x91rSabYCDEeRiEqDNZuXIlPz//4cOH2dQ+T40IgVNByGAwQkJCXMeN/ymswDxeTD9cQDv4g+YZQLOdXSGqMMLZ5XdrfyPEGRiEqNPw9/cPCgq6ceNGh5dP+6OysjIJCQk2Nc6FOBOE5eXlO3fufPb8OY3O+O9oaR5omtW5/B0VFe3k5MTuGhBqBT4jRJ1Ddnb24sWLHzx4wNYRW2lpabdu3djXPrfhTBBKSEiISslBn7/gR/p/R02d6//NGLkxImAru2tAqBU4IkSdw7JlyxYuXGhhYcHWXsrLy9u4kX3XwJkgfPv27eMXwXTLiU0/yIyCV6cZGRF8EnL37t1jdxkI/Q4GIeoEnj9/HhsbS/jLEs1VVFTwWhCy6Z36xp49f8Ho3g8EhZt+UFsFFcVQXVErpZ6Q/JndZSD0OxiEqBPYunXrzp07hYWb/SQlWkVFBU/N5id8lbUWLV2yWDgtFHJTmn6gZwcj14PHSWZGpJP9IHaXgdDvYBAibhceHp6bmzt+/HgO9FVdXS0iIsKBjrgE6xtQtEW3bt32/rNVJPjEL0fptPp/i9xcKSAgYGBgwO4yEPodnCyDuN3ly5fnzJnDvpmijVVXVyclJX38+NHIyIgD3ZFOXl6+uLi4rq6O3X+8fn436n6kAQNghw0MWw5WE+HYWKBQoKaSlv7u4IH9PHVHGnEbCpPJJLsGhFqWlpaWkZExadKkiIgIHR0dtvb1+PHj6fMWFebn0asqFJRV0lI+EbWdBZdTVlaOjo5WVlZmay8/f/4MCwubMGNB1ci/gV8Q+ASAyYTSHyIhp7z2bPZwd2dr7wi1DoMQcam7d+/OmTufT1AoP++HlKREQlysiooKm/qqqanRNjD+NmAlvD4HqsYCX96vnOa6e8d2NnXHVUxMTC5dumRqasqBvjZu/SciJgEAQkJey8nJGRoa6WmpnzxygANdI9QKDELEpTxmzr2RUllrPAoAKD4LRzgMCrx7m0197d67b9fN0HKKCAxZAAlBQK8RjbjyIea9lpYWm3rkHg4ODmvXrh06dCgnO+3Vq5e4uHh4eDgnO0Xod3CyDOJGMTExt+7dr51wAIyGgKw6U0TqRXBoQkICm7o7ff5SuVIvEJOG7rYAAMLdGPoD3r57x6buuIq8vHx+fj6HO1VTU0tOTuZwpwj9DgYh4kZbd+2jDVwIZflw0AmOjgbrSbQB8y/53mRTdxtXLaM+PgiyWvD6HHyNhdRQSHo+aOBANnXHVeTl5QsKCjjcqZqaGgCkpqZyuF+EWoRBiLiRx8RxotG+oKALm97C/zIh6blg+KVBffuwqbvZs2fJyslSUoIhMxJKvvNlxXq4TVJSUmJTd1xFRkamsLCQw53Ky8traWk9f/6cw/0i1CIMQsSNxo8fb6iuQH11GiqLoaaSwmSoyoiPGjWKTd1RqdSA2zctpOl8cfdFy3J0NVSOHz3Kpr64jYyMTHFxMYc7lZeXV1RUfPz4MYf7RahFGISIS40dMUwo4G/KSjXKSjXKl+izp46ztTsbG5s3zx5Sq36OHTn81bMngoKCbO2Oe0hJSf38+ZPDnSorK4uKir548aKsrIzDXSPUHAYh4lJubm779uwSFxHuZ22Vm5s7ePBgdvf4+fNnRUXFvLw8dr9Ux1UkJSVLSko43KmysnJBQUHfvn2DgoI43DVCzWEQIi6lqanp6elJoVD69u0rLy/PgR6joqLMzMw+fvzIgb64wevXr/39/aOjo9PS0vz9/f39/ZOSkjjTdf2uFxMnTvTx8eFMjwi1ApdYQ9zr58+fNBpNWlqaM91FREQMGDAgJCSkoKBATk6OM52SJTo62mnMeCHDgfQ6ehVVdd5hPwAmI3VJZkqSlJQUu3tXV1fPycmZMGHCypUrc3NzeWReEuJaOCJE3Cs2NlZeXp5K5dDf0pcvXw4cONDS0vJdl36DMCMjY9OmTUMchlYZOxfPuFY224++JLB4xrXiGdcrFXr26m1SWVnJ7hqEhYXFxMQqKyvHjBlz+fJldneHUOswCBH3ioyMVFRUpFAoHOgrNzc3JyfH3Ny8b9++XXvFk5KSku/fv1fUMpnSar9+kFtLq8rOymLfwgWNaWhofP361dPT8+TJk3V1dRzoEaHfwSBE3Cs4OFhdXZ0zQRgYGOjo6MjHxzdw4MCXL19yoEeyGBgY3H/ygm7sBE3+YH2Xw8QDwMe/ZNV6DpShpaWVmZlpbm6uoqISEBDAgR4R+h0MQsSl6urqwsLCVFVVOROE9+/fr39P0dbWNiEhobS0lAOdkmL3vv2V6lYgo/HL0fd+IKMBWhZA5U9Mz3r16hW7y6gPQgBYunTpwYMH2d0dQq3AIERcKiwsTEdHR0REhANBWFpaGhwc7OzsDADCwsK2trZPnz5ld6dkKSgqrhOV+eVQeQE8PQIum+u/o4hK1jEY7C5DW1s7PT0dAMaPH//jx4/Q0FB294jQ72AQIi51//59Z2dnBoPBgckyt27dsre3l5CQqP921KhRgYGB7O6ULDu3bhaIvQNljRbajn8ERVnwPwfYYQP06rqcJH09PXaXoaOjUx+EfHx8q1at2r17N7t7ROh3MAgRl7p7966LiwtngvDy5cseHh4N344ZMyYwMJBOp7O7X1JIS0vv+HuTYPKz/w5ZTYTtcbDiESy6CQC+16+pqqqyuwxdXd20tLT6r6dNm5aYmBgREcHuThFqEQYh4kb1LzCYmZlxIAjT0tI+fvw4YsSIhiOqqqq6urpddcrMhw8fjhw+xCj+Di9PwwYDyE8HfkEQlQZRacGQMxISEvb29hz45UNbWzs7O7u2thYAhISE1q9fv3XrVnZ3ilCLMAgRN7p27Zq7uzsAcCAIz549O3Xq1CaLi/711183btxga79k6dmzZ1pa2qOHD8REhIWtx/O/9oIbK4XurBe+u0Eg5NzHDx/ExcU5UIagoKCysvKXL1/qv505c2ZSUhI+KUSkwCBEXIdGo/n6+tbfq2QymWydLEOj0by9vefOndvk+MSJE+/evUuj0djXNbkcHBxO7N221U7WrDBUNem2q+S3rf3lAm7d4MBN0QZ6enoNWxIKCgpu37591apVTCaTYwUgVA+DEHGd27dvm5qa6ujoAEBdXR0fHx/7+vL39zc1NdVrNjdEVVXVxMTkwYMH7OuadNOmTVu7dm1dXZ2jo2P37t3Xrl07ZMgQThagp6eXkpLS8O2UKVOqq6vv3bvHyRoQAgxCxIXOnDkzZ86c+q/ZfWv06NGjixcvbvGjqVOnXrp0iX1dc4Pq6urk5GQHB4fY2FjO966vr984CKlU6t69e9etW1f/4BAhjsEgRNwlISHh8+fPo0ePrv+WrUEYGhpaWlrq5OTU4qfjx48PCQn5/v07m3rnBhEREb169erXr19kZCTne28ShADg6Oioo6Nz8uRJzheDeBkGIeIuBw4cWLx4sYCAQP23bA3C/fv3L1269Hfti4mJubq6XrlyhU29c4NXr14NGDBAU1OTj4+v/q0+TmoehACwf//+Xbt2FRYWcrgYxMswCBEXycnJuX///uzZsxuOsC8IU1JSwsPDp02b1so5s2fPPn/+fBeevhEUFDRs2DAAsLOzCw4O5nDvmpqa+fn5TTa7MDIymjhx4pYtWzhcDOJlGISIixw8eHDatGkyMv8tAMa+IPzf//63aNEiUVHRVs6xsbERFBTkfEJwRkFBQVJSkp2dHQDY29s/f/6cwwXw8fHp6Og0TBxtsG3bttu3b8fFxXG4HsSzMAgRtygqKrp06dKKFSsaH2TTrNHv37/fvn174cKFfzxz7ty5Z86cIbwAbhAQEDBs2DAhISEAGDp06NOnTxnsX2K0iR49enz69KnJQWlp6a1bty5evLgLj8URV8EgRNzi0KFD48aNU1P7ZZM8No0IDx8+7OHh0ZZt6D08PIKCgvLy8givgXS3bt0aO3Zs/deampqKioqc38OZL3YAACAASURBVJHYwMAgOTm5+fE5c+bQaLSrV69yuB7EmzAIEVf4+fOnl5fXunXrmhxnx4jw58+f58+fbzL0/B0pKakxY8ZcvHiR2BpIV1BQEBYWVr/zVD0XFxfOv8P3uyCkUqlHjx5du3ZtF94PC3EPDELEFfbv3z927Fhtbe0mx9kRhCdPnhw1apSGhsafTwUAgAULFpw5c4bztw3Zys/Pb8SIEY1XU3N1dfX39+dwGb8LQgCwtrYeOXLk33//zeGSEA/CIETky8/PP3PmzKZNm5p/RPit0aqqqmPHjq1Zs6btl/Tp00dGRubJkycElkG6S5cuNZkxa25uDgAxMTGcLMPAwCAlJeV3zwL37Nnj6+tLysv+iKdgECLy7dy5093dXV1dvflHhI8IL1y4YGNjY2Rk1K6rFixY0JXe8k5MTPz27ZuDg0OT45MmTbp+/TonK5GSkhIVFf3dqgUyMjI7duxYtGgRzppBbIVBiEiWkZHh4+PT/OlgPWKDkE6nHzhwYO3ate29cPLkyeHh4Q1bJXR2Z8+enTlzZvM/WDc3t+vXr3P4JnCLr9U3mDFjBgBcuHCBgxUhnoNBiEi2adOmJUuWKCgotPgpsbdG/fz8NDU1bWxs2nuhiIiIm5vb+fPniaqERFVVVdeuXZs5c2bzj3r27KmgoPDixQtO1qOvr/+7x4QAQKVST58+vXHjxoKCAk5WhXgKBiEiU0xMzMuXL1uZwEnsiHD//v2rVq3q2LXz588/d+5cF1gP2tfXt2/fvpqami1+On36dG9vb07W03gzphb16tVr8uTJGzZs4FhJiNdgECIyrVmzZsuWLa3sBEtgED579qyqqup3S2z/kaGhoa6ubhfYmOnEiRMLFiz43adubm4PHjwoLi7mWD2tTBxtsH379kePHoWHh3OmJMRrMAgRaR48eJCdnd14ZdHmCLw1euTIkVWrVrHS2pw5c86dO0dIMWQJDw8vLS11dHT83QmysrLDhw+/du0ax0rq3r17Wlpa6+d069atfkm8uro6zlSFeAoGISIHnU5fu3btvn37+Pn5WzmNqBFhamrq+/fvp0yZwkojEyZMePv2bXZ2Nuv1kOXo0aOLFi1q/beB2bNnc3JVOV1d3czMzD8m3KRJk6SlpU+dOsWZqhBPwSBE5Dh79qyCgoKzs3PrpxEVhMePH589e7awsDArjYiIiEyYMKHz7tabnZ399OnT+nmYrRgyZEhVVRXH7kOKiIjIycllZWX98cxjx45t3bq1a+8QiUiBQYhIUFpaun379gMHDvzxTEJujZaXl1+9enXevHkstgMAM2bMuHTpUid9re3EiRMeHh4SEhKtn0ahUObNm+fl5cWZqgBAR0enLbshGhkZzZgxo8WFFxBiBQYhIsHu3budnJzMzMz+eCadTm/93mlb3Lhxo3///i2+sN9eVlZWgoKCnF+cmnWVlZXnz59fvHhxW06eOXNmQEBAfn4+u6uq15bHhPW2bNny+PHjzvjnj7gZBiHitIyMjHPnzu3cubMtJxMyIjx79uycOXNYbKSBm5ubj48PUa1xzKVLl+zs7HR0dNpysrS0tIuLC8duAmtra2dkZLTlzG7duu3atWvp0qWddFCOuBMGIeK0devWLV26VFlZuS0ns/6MMCkpKTs7u5V5ku01adIkf39/Op1OVIMcwGAwjhw5snz58rZfsmTJklOnTnFmlRltbe223Bqt5+HhQaFQcIcmRCAMQsRRb968CQ8Pb+MWSEBEEF65cmXKlCkEvpWvo6OjpaX16tUrohrkgEePHnXr1q1///5tv8Tc3FxWVjYoKIh9VTXQ0tJq+/J1FArl8OHD69evr6ioYGtViHdgECLOYTKZa9as2bVrl6ioaBsvYf0Zob+/v5ubGystNOfq6nrnzh1i22SrQ4cOLVu2rL1XLVq06Pjx4+yopwlNTc12reNqbW1tZ2d38OBB9pWEeAoGIeKcGzdu0Gi0dsUSi88IY2JiKBSKiYlJh1to0dixY+/du9dZHlMlJiZ++vTpr7/+au+FEydOjIqKauM0FlYoKysXFRXV1NS0/ZLdu3cfOXLkx48f7KsK8Q4MQsQhNTU1mzZt2rt3b7uCjcVbo/fu3RszZkyHL/8dfX19MTGxzrJP3tGjRxcsWCAgINDeC4WFhWfMmHH69Gl2VNUYlUpVUlL69u1b2y/R0tJyd3dv45QrhFqHQYg45NixY8bGxkOGDGnXVSwG4aNHj0aMGNHhy1vh5OT0+PFjdrRMrKKiops3b86dO7djl8+bN+/SpUtVVVXEVtWcmppae5fs2bhx4/Xr1zMzM9lTEeIhGISIE4qLi/fu3btnz572XsjKrdGSkpKkpKS+fft27PLW2dvbc3i7oo45f/68s7OzvLx8xy7X1ta2srK6ceMGsVU1p6Ki0q4RIQDIy8vPnz8fB4WIdRiEiBN27drl6urao0eP9l7IyogwLCzM2tpaSEioY5e3bsCAAW/fvuXylygYDIaXl5enpycrjSxYsIADK3wqKyt3YO20lStX3rt3jwNPMVHXhkGI2C4rK+vixYt///13B65lJQgjIiI6sAdvG0lKSmppacXFxbGpfUI8fvxYWlq6T58+rDTi5OSUl5cXHR1NVFUtkpeXz8vLa+9VUlJSCxYs2Lt3LztKQrwDgxCx3T///DN//vw2vkHfBCu3RuPi4kxNTTt2bVv06dMnMjKSfe2z7uzZs/Pnz2exESqVyoH9KBQUFDoQhACwdOnS27dv5+TkEF4S4h0YhIi9UlNT79y5s3Llyo5d3rH3CJlM5suXL9+/fy8nJ9exftvCxMQkPj6efe2zKDc3Nzg4eOLEiaw3NWPGDH9/f7a+wC4rK1tUVNSBC2VkZDw8PI4cOUJ4SYh3YBAi9tq2bduKFSukpaU7dnnHbo3qdNd3Gev6La9w2HCnnbvbPUPnjwYPHyUqKbNmw+bT571FJWVEJWVUNHU79nOcfXx8fMaOHSsuLs56UyoqKnZ2dn5+fqw39TvS0tId/gNctmzZxYsXcaEZ1GGUzvJSMOqMUlNT+/fvn5qa2q1bt461ICkpmZWV9cedgxo7c+bM/IWLmOteA4UK764LBHt9/5YtKyvbsQKa+Pbt2/Hjxw9736ha+wYolH+P/kjlv7V2iK50UGAApeEg2Xr37n38+PEBAwYQ0tr9+/f37t0bGhpKSGuNvXnzJicnJz093cvLq/5pn6amppWVVbsaGTt2rKOjI+v3gRFvwiBEbDRlypTevXuvXbu2wy2Ii4v/+PFDTEysjefX1NSo6+jl0aiwOwUAIOER5dy0GR5u50+f6HANDUpKSmxtbVPT0mtkdWBrzL9Hv8aA1xQYPJ/vzsZp7lPOnz/Pekesi4+Pd3Z2zszMJCqY6XS6hobGy5cvDQwMCGmwXmpqqqlVP0qPQUwGo7q6WqR+7b1Prz4lxKipqbW9nWfPnq1Zs4bdM3pQV4W3RhG7pKamPn36dNGiRaw00t7JMk+ePKmgigKTAe/9oPALvDjFBObtgEBCfuGTlJQcN2EiVccaRCQbdXkYnNaAwxLmmO2Xr1wpKSlhvSPW+fn5TZo0icDhKT8//6RJk65du0ZUg5WVldu3b7ftP6BaqWfFrGuVc3wZi+9WGAyrKCmuFOhmbtmnXXtd2dvbl5SUdJa1fhC3wSBE7LJv376FCxey+Iyqvc8Ihw0bJsFXB86bIekF3FwPpqMoVP5RTsMJiYScnJz9h49W9/91lZbseNC2BACG/RImlf/s2bOsd8S6W7duTZgwgdg23dzcrl69StQ9JDqdnpub+7OyhkFp9N83Px30bJmL71bU8cnIyLS9NQqF4uHh4e3tTUhtiNdgECK2+P79++3bt9u4H3or2huEgoKCxw/sEX91BAbOheGrIOM9P73ywO5/WCyj3u3bt0HbCiQUfzlaUQzC3QAAqHx1kqrPQ8IJ6YsVnz59qqiosLCwILZZS0tLfn5+ot4YERUVffQypHawJ1B+/SkkIASColUuO5at3dSu3RDd3Nz8/Pw4s4Ei6mIwCBFbnDx50t3dnfUpKkwms72DOVdXV41ufOInnWHvQKEo/x07digoKLBYRr1p06YJ5sRCbsovR7vJQUUxAEBNJbXo63gXtixt2i4PHz4cOXIkO6btTJgwwd/fn5CmTp46XSCkDNq/vuzPxw+xgeA9h3l3c1Yp7ebNW21vUF9fX1ZW9u3bt4SUh3gKBiEiHo1GO3fu3MKFC1lsh8lkMpnMDrxQf/Wyt+soJz4GvbAgb83Ktm4C/EcSEhJ7tv8tEvzrvBsda0gNAQD++9sE+KmTJ08mqrsOe/r0qaOjIztaHj9+/K1b7QinVpSVVzCEm80ldt4M64Jh5ROYdqYmN6Ounbdhx44dGxAQQEh5iKdgECLiXb9+3dLSkvXphR1eX83U1NTCwoLJZLZ9umkbvX8fQSn6AtkJ4DUZ4h8CAAxbAS9Pg89ixtNj69eta/uew2xCp9PDw8OJemuiCRMTEyaT+eHDB9abWrl8qfi3GPj66/SWhlFsdTmjtnqMi3O72hw6dOizZ89Yrw3xGgxCRLwzZ87MmzeP9XZYWV8tJyeHSqVWV1ezXkZjixcvPnzwgIhYN+g9CvgF4Us0VJWA+wnBb3Fz583v2HqqxIqNjdXQ0GjXTJN2GTly5P3791lvR1hYeN8/fwuHegG9BiqLgU4DAMhOACYTamlC1xcZ6OuLiIi0q00bG5vPnz8XFhayXh7iKe1evAqh1iUkJGRnZzs5ObHeFCtB+OXLF1FR0ZKSEmFhYdYraWBiYmJiYvI+7sPz4OONjyupyR7cv4/AjjosMjKSxVW2W+fo6Hjo0KF169ax3tT58+f5qkuhMBMOOsGIdWA+Bp4cgpTXQKvoJsofGBjY3gYFBAQsLS0jIiII+euHeAe+UI8Itnr1akFBQUJ2iauoqFBUVCwvL+/AtWZmZkVFRU+ePCH2BXCudffu3fpFq69du6asrDx48GB9ff1BgwYR3lF5ebmKikpubi4hN4Hj4uL6DXGsM/nlFig15m7suzB9ff0ONLh27VoJCYmNGzeyXhviHRiEiEgMBkNTUzMoKKhnz56st1ZeXq6srFxWVtaBMiQkJHr27Hnw4EFbW1vWK+Fyb968GTp6Ql2vEQBQW1vLx8dHpVL5Y+4mJ8aqqqoS3l3fvn13795NVMo+ePCgyd4RWlpaw4YNa2870+YsDH0bUV5RXlFRoaigCABqyopPA+8ICgoSUifqwvDWKCJSWFiYnJwcISkILEyWSUtLk5eXV1JSKigoIKQSrhUREbFmzZo34W9rxeRAywb6TQUARl0t+K6oYVBNzS3ehIbo6ekR26mtre2bN2+ICsKRI0ey2AKDwfD29va//7Bqtm/9dJt0AGAyvt1cNWfuvH+2b9PQ0CCgUNR1YRAiIt29e9fV1ZWo1jr8jDAqKsrCwkJKSio/P5+oYrgTPz//oEGDonJptdMvwy5b6N4PFLrDixNAK2fuTinfqO/s7Pzp0ydiO7WysvL19SW2zQ7Lzc01NTUtKi6uVeoJWv+/hgCtAg45VQsIX/e9waijX7lyhdQaEbfDWaOISPfv33dxcSGqtQ68TV8vOjrazMxMSUnpx48fRBXDnfT09I6ePlc+7gDIaYGMOpT+AACIvAX9Z4GgKM11b3p6xufPn4nt1NTUlHtW9VRQUFi9Zi2fiiFINFo24flx0LSAlU8oQxbW8QmRVx3qHDAIEWHS09MrKytNTEyIarDDt0bDw8NtbGyUlZW/f/9OVDHcafuuPdUG9qDdBz48gdpq0OoDAFCUBTIaAMC0cWPwC929e5eo7oqLi/39/WNjY3Nycq5everv7+/v79+Bh7gE+vnz57Zde6v7zfrlaGwAaFmC35qaurq79wPfv39PUnWoc8Bbo4gwr169GjhwIIENdmxEWFtbGxcXZ2VlVVJS8vz5cwLr4ULioqLUmkr4GgPXlsKi28AvCABA5QNGHQAAkwGMug6/gtLcgiUrAiKSqDLq9B4Oc47eplKpjILMyc9fEbLLVcecO3++VtcWZH99ClicDcFnYNxOyIqvotVc9r3J1ldKUGeHI0JEmJCQEGKDsGPPCCMjI7t3796tWzc1NbWsrCwC6+FCa1evFMoMh+OuMN8XVAz/PSqjAYWZAEANOcvHrBs/fjzrHd27d8/FxcX/xvUq97MVs64xFt+tmuNboW1bVVvnff7ciROkBeGkiROpKa+h/NdZUaJSYO8JenYwYDa1plJHTYmk6lDngEGICBMZGdnejcVb17EgfPHixZAhQwBAXV09OzubwHq4UEZGBtRUCQvyg5AY5GcArRwAwHoivDwFZXkCd7cYGRkSMmeytrY2MfkzQ1AUqP9/szoxCEK9YeUTcNm8ctXqL1++sN5LB2hoaHgumCcccfWXo+omUFsFAHwvjokICRD7+xnqejAIETFoNFpaWpqhoeGfT22zjt0affny5eDBgwFAUVHx58+fhK+yxlUSExPNzc0Fq0spx13hykLIjAIAGDgXlA0pW80V5aQJWQ4NAKqqq/PrhEFU+r9DqWFgNhpEJBmOqxgCIps3byakow6gMOiM9HeQFQ+3NkJmJADAsOXw5BC8OE4J2N7DwMDU1JSs2lCngEGIiJGamqqtrS0kROQMvQ6MCCsrKyMiIuqXnKZSqWpqal+/fiWwJG4zYcKEp0+fBr96IVxVKCYpLRbqJXbeTeyCu1hhqhg/M+p9hJqaGuu9VFdXL1+zoXzC4V+OSihCfhoAAJWvVlrjzv0HZK3OYWVlNWfWTBFRMVDWBz4BqCwGOS2YflYgwnfwoEGhoaEEPiVFXRJOlkEEYDAYfn5+dXV1ly9fdnNz4+cn5u9VB0aEL1++tLS07Nbt3/19tLS0MjMzO7ZYVydiamr66K5//RJrDXR118vJyRHSfl1dHZ1OB4FfV23tNxXCLsHpSVBTBTWVQN4SVa6urmPHjo1KGBp3aw0AVFdX8/HxCQgIyMrJ+994Quxis6hLwiBErIqOjrbpZ1tHpwOVb6HnEnNzc2NjY0Ja7kAQBgUFNd6KT1tbOz09nZBiuBxbH4OJiYn9b+f2FQeX/bLqq4gEbAqHnI8gLiOw09p1jDM7tgJuIwqFEv7q3w2YLCwsFBUVHz58SFYxqNPBOwaIJUwmc8KUqbXSGoyDOYz9X2lDlqzauJXAxtv7szUwMNDZ+b8VnHV0dHgkCNnNzW2yikA1VJdBVSlUlQIA0GugogjUe1NeeUFV6ZHDh//UBodoaGiQNXMHdVIYhIgl1319s3PzwXAIXFkIJ/+iSyiHRcU9ffqUkMbbG4Tx8fFUKtXIyKjhSPfu3QlfV4U3XbhwgVJTSa0qgSsL4MpCAIDaKjg8CtbrU58e9jp9SkpKiuwa/2VkZNRkFW+EWodBiFhywce3RkoDYu/D2B0w8yI8PVypah76LpKQxts7WSYgIKDxcBAwCImzePHiT58+jRw9VlBCTkhcguI9G3xXCupYCiloTZjsNmPGDLIL/I+5uXlFRUVxcTHZhaBOA4MQsWT7xrUCeclg7AgKuiCjBr1HUOMD3SYQsO729OnTR48enZ2dbWFhERwc3JZLbt++3WTJbz09vbS0NAaDwXo9CADOnTx6bPGEo5P6iMTcNCh8N8tQ4OjSyScO7Se7rl9oaWkJCwvHxMSQXQjqNHA/QsQqa9v+71OymVNPAQCfz8K+PTRDgl+y2GZRUZGCopKgmhGNQRGm1JnqKIW9eNL6JV++fLGysvr27VuT5UnV1NTevHmDG/EQqLCwUFdXd+fOnZGRkRcvXiS7nKaKiopUVFQ2bNiwZcsWsmtBnQOOCBGrbly9IlFXSjk6uttFd9Hqoms+l1lvc92mv5niclVjdjIW+leuD49P/xYYGNj6Jf7+/mPHjm2+SLeBgUFKSgrrJaEGHz586Nmzp6ura0BAQGVlJdnlNCUjIyMgIPDq1SuyC0GdBgYhYpWWlpaT4zAdTfWVSxfn5/1QV1dnscFPnz75+PoxVHvB8xNwxh129C23X7lkzabWr/L3958wYULz4z169CB8Qz4el5CQ0KtXL2Vl5b59+/r7+5NdTgt69OgRERFBo9HILgR1DhiEiFWvXr0KCwtTVVXt378/ISvL1NTUUKh8MNcHltyD9SFgNhoifFu/g//ly5fMzMwW36UzMDBITk5mvSrUIDIy0tLSEgA8PT2PHDlCdjkt6Nmzp5KS0ps3b8guBHUOGISIJbm5uVOnTj137lxJSYm0tPSfL2iD3r17u7qMEgra/e/3un35Poft/6e1pSx9fX3HjRvX4oo2hoaGHz9+JKQwVC8iIqI+CB0dHWk0GlFvyxDI0NBQUVExKCiI7EJQ54CTZVDHlZSU2Nvbu7q6btiwQVlZOSYmRkmJmP1uvn79atjbVKinPZNRV54UJtdN+HtWa69IW1hYHDhwYNCgQc0/ys7Orp9EQ0hhqKioSFtbu6ioqP5x7PXr148fPx4WFkZ2Xb8IDAzcvXt3UVFRUlIS2bWgTgBHhKiDCgsLhw0b1q9fvw0bNtTV1RUVFRG1siUASElJ6Wlp8KcGi3x569DXPPTVi1ZOTk5Ozs3NrV9ouzk1NbXKykp8q4wooaGhNjY2DZOS/vrrr7KysoCAAHKrasLY2Pjr16/l5eU4Twq1BQYh6ojPnz/b2toOHjy4/hFRbm6urKwsUWttA4CEhERsbGxeXt63b98ePXqkq6vbysnXrl2bPHlyK6/e491RAgUHBzd+FsvHx7d3795169bV1taSWFUTmpqaZWVlw4cPv3XrFtm1oE4AgxC124MHD+zs7JYvX75nz576JdC+fv2qqalJVj3Xr1+fOHFiKycYGRlhEBLl2bNn9fseN3ByctLW1j527BhZJTVHoVB69erVq1cv7pzUirgNBiFqBxqNtmrVqoULF96+fXvevHkNx79+/cr6WxMd8+7dOyaT2adPn1bO6dmz54cPHzhWUhf248ePrKys+pkyjR0+fHjPnj3Z2dmkVNUiExOTurq6Hz9+4N1R9EcYhKitYmNj+/Tpk5mZGRMT069fv8YfpaWltX73kn2uXr3q4eHR+jkYhER58uSJvb1983vgenp6np6eixcvJqWqFpmZmcXGxk6aNOnatWtk14K4HQYh+jMajbZ582ZHR8c1a9bcvHlTRkamyQnp6ena2tqcL4xOp/v5+bm5ubV+Ws+ePRMTEzlTUtf26NGjxts9NrZu3bqUlJQbN25wuKTfMTMzi46OnjJlio+PD86NR63DIER/EBwcbGJi8vHjx7i4OHd39xbPSU5OJmUX+GfPnmlra3fv3r3109TU1GpqagoKCjhTVVdFp9OfPHni5OTU4qeCgoIXLlxYtmxZXl4ehwtrkbGxcXp6uoGBgbCwcGhoKNnlIK6GQYh+q7CwcObMmR4eHnv37r1161Yr7wimpKSQEoRXr17943CwnpGREQ4KWVS/drmqqurvTrC2tp4+ffr8+fM5WdXvCAoKGhkZxcfHT5s2zdvbm+xyEFfDIEQtYDKZ3t7exsbGkpKSHz58GD16dCsnFxcX02g0FRUVjpVXr7KyMjAwsPX5og2MjY0xCFkUGBjYZLvH5rZu3fr58+fLlwlYeJ11lpaWUVFRU6dOvXPnTnl5OdnlIO6FQYiaSk5OHjx48MmTJx8+fHjo0KFu3bq1fn5iYmLPnj05U1tj9+7d69u3r4KCQltOxiBkXUBAQOu/EgGAkJDQlStXVq1alZmZyZGiWmNhYREZGamoqDho0CBfX1+yy0HcC4MQ/YdOp+/atat///6urq7h4eFmZmZtuerDhw9GRkbsrq25a9eutfG+KAD06tULg5AVycnJFRUVbfkrYWJisnbtWg8Pj7q6Og4U1gpLS8vIyEgAmDNnztmzZ8ktBnEzDEL0r8TERGtr69DQ0KioqCVLljTf2O93EhISevfuzdbamisqKgoNDR0zZkwbz68fEeLswQ67c+fOmDFj6tdP+KPly5eLiIjs3r37z6eyU/1Ca2VlZY6Ojnl5ebhnPfodDEIETCbz0KFD9vb2np6eDx8+bO+r8QkJCcbGxmyq7Xf8/PyGDx8uLi7exvNlZGTExMS46o3vzuXu3btt/7WDSqV6e3ufOHHi7du3bK2qdfz8/MbGxjExMVQqdc6cOadOnSKxGMTNMAh5XXFxsbOz882bN9+9ezdjxoz2Xs5kMuPj4zk/Irxx48akSZPadUmvXr0SEhLYVE/X9u3bt5SUlN8ta94iFRWV06dPu7u7l5WVsa+wP+rTp0/93dFZs2bdvHmzpKSExGIQ18Ig5GlJSUlWVlYGBgbBwcFaWlodaCEjI0NCQkJWVpbo0lrz/fv3hISE4cOHt+sqfK2+w+7fvz9y5EgBAYF2XTV69GgHBwdyl5tpeEyoqKg4fPhwfI8CtQiDkHeFhYUNHjx48+bNBw4c6PDGEdHR0ebm5sQW9kd+fn7Ozs5CQkLtugonjnbYvXv3/jhftEUHDx589+6dn58f4SW1kaWl5fv37+u/XrRo0fHjxxkMBlnFIK6FQcijwsLCxo4de+nSpalTp7LSTkxMTBsnlxLI39+/vfdFAW+NdlRZWVlYWFh7x9/1REVFr1y5smTJkpycHMILawsDA4O8vLyioiIAsLW1FRcXf/bsGSmVIG6GQciLUlNTx48f7+Pj87t1I9suOjqaw0GYk5Pz6dOnwYMHt/dCIyOjlJQU0uf0dzpBQUH1EdKxyy0tLT09PWfOnEnKlF0qlWpmZtYwX9TT05OrtotCXAKDkOfQaLQJEyZs27Zt2LBhrLcWFRXVfFMetrp9+/bo0aMFBQXbe6GoqKiysvLnz5/ZUVUXFhAQ4OLiwkoL69atKykpOX36NFEltYu5uXn9Y0IAcHNze/fuHW7MhJrAIOQ5+/bt09HRmTt3LutNZWVl8fHxcXhxtZs3r1t3GgAAElpJREFUb7q6unbsWnxM2F51dXVBQUEjR45kpRF+fn5vb+8tW7aQstyMubl5dHR0/dciIiIzZ84kK5IR18Ig5BVBQUH+/v7e3t7/+9//HBwcCLlPFRkZ2fqOuITLy8uLj4+3t7fv2OU9e/bErerb5e3bt2pqaqzvutyjR4/Vq1c33syZY8zNzRu/Sr9w4UJvb+/S0lLOV4K4FgYhrzh4+OiM2XPnLVxcWU1bt35DG5cIaV1UVJSFhQXr7bRdYGDgsGHDhIWFO3Y57tDbXkFBQSNGjCCkqZUrVxYUFFy9epWQ1trOwMDg27dvDa8zamhoDB48mPNlIG6GQcgT0tLSQt+9r1j6pGZzVJ3ZWBrwP336lPVmIyIiODwiZPF5FQZhewUFBXVsvmhzfHx8Xl5eq1ev5vBb7Xx8fMbGxnFxcQ1HlixZcvz4cVxvDzXAIOQJ85eurHVYDhqmIKUEH5/XuO6ZvXBpbW0tK20ymczo6GhOzpSprq5+9eoVKz+XDQwM0tLS6HQ6gVV1YYWFhampqTY2NkQ1aGlp6ezsvG3bNqIabCNTU9PGQThw4EB+fv7nz59zuAzEtTAIu77CwsIXQQ/otjMBACJvg4419JtaUFrB4sKbaWlp4uLibdwFiRAhISHGxsasrGIjLCysqqqalpZGYFVd2MuXL+3s7Nq7oEzrduzYceXKFQ7P3e3du3eTV0g9PT2PHz/OyRoQN8Mg7PpkZWXHT3ITenEEACD0IthNp7y9qq4kr6mpyUqz79+/5/B90cePH7P+4qORkVFSUhIh9XR5r169GjJkCLFtysvLr1ixYuPGjcQ227rms4Xd3d3fvHmTkZHByTIQ18Ig5AmH9u3mCzkHwWchJwEEhETub/H2Ok6lsvRfn/NTRp88ecJ6EPbo0ePTp0+E1NPlBQcHDxw4kPBmlyxZ8vr16/j4eMJb/h1jY+Mmz4ZFREQ8PDy8vLw4VgPiZhiEPEFFRWX/3t1Kwf9TlJYwjz66cM4M1h/8cHhEmJeXl5WVxfokVQMDAwzCtvj58+eXL19MTEwIb1lMTGz58uX/+9//CG/5d2RkZAQFBb9//9744IIFCy5evEij0ThWBuJaHVxqGXU6C+bOHjKwv6KiopSUFOut1dXVxcbGcnK57eDg4P79+7d9u+DfMTQ0xHFAW0RERFhaWnZ4NfbWzZs3T0dH59u3bxxbjcHQ0DApKUlZWbnhSPfu3c3Nzf39/d3d3TlTA+JaOCLkIQYGBoSkIAB8+vRJRUVFUlKSkNbaIiQkpF374f2Onp5ecnIy6+10ee/fv7eysmJT45KSkn/99RcnN0XS19dvvrLaggULcLdeBBiEqGM4/wZhaGiora0t6+3Iycnx8fHl5+ez3lTXFhcXx477og2mT59++fJl9rXfhJ6eXmpqapODI0eOzMrKwj1JEAYh6oioqChO3hetrKxMSUkhqscWBweoicTExF69erGvfWtr65qaGo6teKetrd18jigfH9+MGTPOnTvHmRoQ18IgRB0RGRnJyVfpY2NjDQ0N27sT7+/o6uriq4Stq6ury8jI0NPTY1P7qzds1jexLK6iDxzurG9iqW9iOXDYiMrKSjZ1BwA6OjotLvk9a9asq1evVlVVsa9rxP0wCFG70en0xMREU1NTjvUYExND4AC0e/fuGIQtevjwoaKGrryGroKGbi2fsJqekYKG7vXrvgR2UVNT4+3tfezEqVSHnT9n+xdMu5bqciJ1yLZ3WWWLly5n37JnGhoaX79+bfG4hYXF3bt32dQv6hQwCFG7JSUlqampdevWjWM9JiQkEPi8SldXF3clbKK6ujowMNBt6ow8M/eCBY+KPJ/UbXxXMC8wf8DyuYsWFxQUENWRnp7e8lWraYoGcGEmMBmgaQ55aeC3hiah6u19cdCgQWxaAE9WVraioqK6urr5RzNmzLh48SI7OkWdBQYhajdix2dtkZSUZGRkRFRr2trapGyMx80WL178999bq0Rk4cNTeLwf5LVBXhvubob3fuUlxRv+Jmx10E2bN9Ml1WD1Cxi+Cp4fBwAwdoStMTDnCt+odTGxceHh4UT11RiFQlFSUsrNzW3+0ZgxY6Kjo3NyctjRL+oUMAhRu7F7PmFzycnJBgYGRLWmpaWFQdjEtm3bkjO+1sy/CUsDIOwSMOgAAPOuw+rnICZ9xeda8ymXHVBVVbVh647y8QeBQoWqEhCTAQAQkQAKBQBqHZaVV1Skp6ez3lFjT58+tRk01HrQ0Pwy2qjxk60HDR3kOKrxMnvCwsKurq7Xr18ntl/UiWAQonaLiYkxMzPjWHeVlZWlpaVKSkpENaiiolJUVNTiXTKetXL9Jnq/6aCoBz9SQVoVqI3eo6fw1drNXrRiLeu9XL16tUpaB/TsIDsBQs7D8FW/fHx/F1Pd9GnYe9Y7alBaWjplxux3aqMjLFZWTLvwYeDfEWbLXlcpTnSf3vgerLu7u4+PD4H9os4FgxC1G7FP7H4nIyNj5tz502fPmzJtpoCI+Iw58+cu8CRkN0EqlaqiooK3whpUVVX5X/OhDVgAlcXgPQcmH25yQp2Ne2hoCOurkTk7O0N2Anx8CqcmwnxfkGq0rEzQfvgaI87PHD6wL4u9NJg5c6aqqmr+tyywHAeG9mBoDxpm4L+OWfkzKTlZV1e3qKio/sz+/fsXFxfjauw8C4MQtU/9go2Kiors7mj2wiWXU+iXakzvigwsHbntUo3puTR+txlzCZlYqK6unpWVxXo7XYOIiMiUaTOEnu6HIy7gsAR6OTU5gT/ce8iQIay/vqKoqDh35jTqqYkw9RRoN1q25vkxiHtAsRinJSkwxc2NxV4arF69mskvDMKNZnUlPgYFXZjvS5/v96Og8N69e/WHKRSKq6vrzZs3ieoadS4YhKh9EhMTe/fuzdYufH19dXR0XgQ9qDMYAgNm//tPcjDzvV9CdMTVa9dY7wKDsIl/tmykvzoLPYeC3YymnzHrBMK8jx/Yy3ovtbW1/n43BKhMeHsNgv4H73wBAGLvw80N0L0v/531A20sCByWbd21l243CyiN1qdVMoCCTKj8CRpmtcAfERnV8Mm4ceNu375NVNeoc8EgRO2TmJhobGzM1i4sLCxqqYKg3Qf4Gm0JO2w5rH7BFBJfumpdRUUFi10oKys32YuAx23ZskVcTJQ//DJcWQhXFkJVKQBAxA24shAqinU0VP38/FjvhUKhbN68ecZUD6GPQSAi+e9YTV4bJh/iz0kwMTbu3bu3oKAg6x0BwJcvX27f9KcNW/3LUU1z0LOFjYaw0ZDR0/Fu0IuGT/r165ebm4uzqHgT7j6B2icpKYnd707cD3zwU0wNKCK/HNU0h5LvwCdQrdt/5559u/5haUK/srIyjggbW7Zs2aRJk6bOmpsvLgdKBlCSC+WFIC4HMurdpGT27N5FyM1wfn7+uXPnAkBiSnqoj2fjj4Qlpe99TCRwMwp1dfXuBkbJiY9/uZP+5jLkZ8D+rwBA2WqqqCbd8AmVSnVycgoMDPT09GzaFurqMAhR+3z69MmNuKc4zVVXV2/eur1y5Uu4vanFEypH7zp5YMD2vzezskOQoqLi+/dEzk7s7OrXCbrqfd5j9oK6hP+ygwJwwuvkqFGjiO0u5PljYhtsjkqlnj95ZOiYSZXQKAq/JYGOFfAJQHEO/88cK+dBjS9xcnK6fPkyBiEPwiBE7UPsK33NCQkJ6RkYxme+/+2UmK8xkhKsLmojLy+PG1A0N3To0NwvXWfNHTqd3l1dMT4uF8J9QEEPTEaCsSNcnAUqRkIvjgoLCUyaNKnx+Q4ODrNnz6bRaEStaos6C3xGiNqhpKSksrKSwFf6mqNQKBdOHxMJ2PLvO91NMJnid9edOXaQxQ1j5eTkMAi7vMrKypGOQwUEBaEsH2jlAAA9BsGcKxAbSM1JuHXr1pAhQxqfLy0t3aNHj4iICHLKReTBESFqh/T0dF1dXQqFwtZelJWVrc17B0fGM7ITQEQCDIcAAGRGQUEG0Mr1VPQ1NTVZ7EJWVrbhHTLUVQ0fPnz48OF8gsJ79+3j4xcE3yX1x+m0qguXve3t7ZtfYmdnFxIS0r9/f85WikiGQYjaISMjQ1tbm929fPnyRU1JgVpZxMj9BKU//g3C5Ffwf+3de0yVdRjA8QcOnGNcFEUNQ51JMix3wlKZqFGal2Ta0FJRu2imrmk6TXOmaZKXOQwtU2tZ3hp/5MpSFDXFNa/omCGSdgPkqhyPlwLiiOftj8xLYyEivMDz/YyN7X0P73n+YPvunL3v75ef4eUpoQ+3S0tLCwsLq8lbBAQEXL58+f6Mi/otbuG7M6ZNvf2I1Wr19fWt9MW9evVav359ncyFesSj9vY9QeOTkJBw7ty5hISEOnivJcuWL17xoa3FrZsVXVcuThg7YmX88ppf3DAMb2/v8vJyi8VS9auhRn5+fnh4OF+ba8MnQlRDbm5uu3bt6ua95sx+67mB/d1u980jXl5eNd+Dot/goft3bb95QRHx8W92NjOjbdu2NbwyGoHg4GCLxZKfnx8cHGz2LKg7hBDVUFBQEBERUTfv5enpeX+X9j5//vyaNWsOHT0u60puPKpfeklO7f7r4Prnh484/EMK9wpCROx2e3p6OiFUhbtGUQ1FRUW1estorZo4cWL8Bwnll4puLVizb7VkHXf3GJWent6jR4///WtoERoael/2nEIDwidCVENxcXHLli3NnuIeRfbusy+vQtJ23jo0ZP4/v697+2RsGF9YWNimTRtzhoPZDMOYMn1mdm5ednb2jr0puw8cFJGo3pGzZ0w3ezTUOkKIanA4HA00hBcuXIhburxk+vd3hPBfRvvHxcOSlLRzwoTX6n421AcbN23auH1/yYA50lpEJMdVKoY7Zcn8buH2rl27Nm/evKoLoAHjq1HcFZfLlZiY6HQ6jxw50hDvNE5OTpagUAnqVMm50kvy2ctGr1f2HOJJao1yc3Ptdvv4ceNLrAHSbfiNn5PfSsraMrENjo4OCQlpiP/zuHuEEHelQ8eQcRMmVnh6j4gd3SfqGbPHqbaRI0f6l56XzH3/PVF2VVYNlahJfr+mvD76BTNGg8kCAwOferqvNbSn+N/2bccbW2XeUVl61vAN7Na9e20vIgFzEUJUbceOHYUFBeWzDsjCk9diPz5y+HCDu5vAZrOtXrHMb9vsO466SmX1MOk5xqvU0btbeP/+/U2aDma6ePHiF5u/LH/ixcpOeriuS2rajzXf+Qv1GSFEFVwu1+RpM+WxZ+X4V3KlUDKSxT5o8pszzZ6r2rKysmylThGRT2Jl70oRka/nS0GmnNptbF/s4SphLzqd5r33vitynPgFVnLuTIr4tXB1iEhOrvXtMmAiQogqpKamOv8olYhY+Wm/HFgnjix35KvHUlMrKipbFLsei4mJ+WhVgq1pC+kySFqFSE6aPNpPRsZbjfLYMWPj4uJat25t9owwwYC+UU3O7BXDXcm5Qxuk50se2Sfat6+jdSRgCpZYQxUMw7B375mRcVric6SJvxRkeiztMzp21JbPPzV7tHsxd8Girdu+u/3Ig61a7tq21c/Pz6yRYC7DMJ6MjDppCzNKnDIp8daJsisyJ9QaOXrYQ2WJG1mAtDEjhKhaUlLSkCFDjanfiJdVin/32DKlsLDwvmxZDpjO7XbHx8e/s2BRRVBnGThTOkZIi7YiIilr5fQe35xjv2Sm84Bp48ZzhKhadHR0WOfO+ZvH+TcLuOQoHhATQwXRaBiG4XQ6Oz3S8ecyy/WcNAkKvRFCq88D1/6c+/YsKtjo8YkQd+vEiRN5eXlhYWE13AIJqIeKioo6de7i6dPs5hHDMJraLL+dyWAR2kaPEAKAiIjD4bh69ertR4KCgnx8fMyaB3WGEAIAVOPxCQCAaoQQAKAaIQQAqEYIAQCqEUIAgGqEEACgGiEEAKhGCAEAqhFCAIBqhBAAoBohBACoRggBAKoRQgCAaoQQAKAaIQQAqEYIAQCqEUIAgGqEEACgGiEEAKhGCAEAqhFCAIBqhBAAoBohBACoRggBAKoRQgCAaoQQAKAaIQQAqEYIAQCqEUIAgGqEEACgGiEEAKhGCAEAqhFCAIBqhBAAoBohBACoRggBAKoRQgCAaoQQAKAaIQQAqEYIAQCqEUIAgGqEEACgGiEEAKhGCAEAqhFCAIBqhBAAoBohBACoRggBAKoRQgCAaoQQAKAaIQQAqEYIAQCqEUIAgGqEEACgGiEEAKhGCAEAqhFCAIBqhBAAoBohBACo9jfc823R5tCevgAAAABJRU5ErkJggg==", + "image/svg+xml": "\n\n\n \n \n \n\n\n\n \n \n \n\n\n\n \n \n \n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n", + "text/html": [ + "\n", + "\n", + "\n", + " \n", + " \n", + " \n", + "\n", + "\n", + "\n", + " \n", + " \n", + " \n", + "\n", + "\n", + "\n", + " \n", + " \n", + " \n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n", + "\n" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "# graphplot(trParentnodes, trChildnodes, names = trNamenodes, method = :tree)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Rapport\n", + "\n", + "- Quelques détails sur les points clés et non triviaux de votre implémentation\n", + "- Une courte argumentation de l’adéquation du résultat avec l’instance résolue\n", + "- Quelques éléments d’analyse, par exemple :\n", + "- Analyser l’impact des différents paramètres de la PSE sur la performance en terme de vitesse/délais (en temps cpu et nombre de noeuds explorés) pour l’obtention d’une solution réalisable, pour l’obtention de la meilleure solution, pour l’obtention de bornes supérieures de qualité et pour la complétion de l’algorithme\n", + "\n", + "## Notre implémentation\n", + "\n", + "Pour l'implémentation de notre `branch and bound`, nous avons gardé la même architecture que celle du sujet. Nous avons simplement ré-écrit les fonctions qui calculent la valeur objectif et la valeur de la capacité. \\\n", + "Notre vecteur de décision $x$ est tel que $x\\in\\{-1,0,1\\}^n$ avec $n$ le nombre d'objets disponibles. $1$ et $0$ représentent, respectivement, un objet pris ou non-pris. La valeur $-1$ représente un objet ou la décision n'a pas encore été faite. \\\n", + "Dans la fonction `Objective`, si la valeur est a $-1$, on fait comme-ci elle était à $1$ pour calculer une borne supérieur. Dans la fonction `Constraints` les $-1$ sont concidérés comme des $0$ pour avoir une borne inférieur de notre problème. \\\n", + "\n", + "Pour l'exploration et la création de noeud, nous faisons une exploration des noeuds avec le plus de $1$ en premier. \\\n", + "Lorsque nous arrivons sur un noeud, nous l'enlevons de la pile des noeuds. Nous calculons ensuite s'il est \"trivial\", s'il dépasse déjà l'une des bornes ou s'il est divisible.\n", + "- Dans le cas ou il est \"trivial\", nous calculons sa valeur et remplaçons la meilleure valeur si nécessaire.\n", + "- Dans le cas ou il dépasse une borne, nous l'abandonnons.\n", + "- Et enfin, si il est séparable, nous le séparons en deux en ajoutant la pile des noeuds les deux sous noeuds. Exemple : $[(0,-1,-1)]$ -> $[(0,0,-1) ; (0,1,-1)]$\n", + "\n", + "## Justification de nos résultats\n", + "\n", + "Notre algorithme trouve bien les solutions idéales. Pour le test `uncorrelated_span/knapPI_11_20_1000_1_-1428` Nous avons trouvé la valeur optimale qui est de `1428` en `2131` exploration de noeuds là ou un algorithme classique de parcours de l'ensemble des noeuds aurait exploré `1048576` noeuds, et ce sans dépasser la capacité maximale du sac bien evidemment.\n", + "\n", + "## Analyse\n", + "\n", + "Nous avaons réalisé les calculs de deux manières différentes:\n", + "- La première fois, nous mettions a jour les coordonnées de $x$ dans l'ordre des indices.\n", + "- La deuxième fois, nous avons d'abbord trié les vecteurs `weight` et `price` par leur rapport $prix / poids$. Cela revient à mettre à jour les coordonnées de $x$ dans l'ordre du meilleur rapport $prix / poids$. (NB : le vecteur résultat n'est pas ré-ordonné comme l'était les `price` et `weight` initialement)\n", + "\n", + "Pour le fichier `weakly_correlated/knapPI_2_50_1000_1_-1396`, nous avons comparé les résultats :\n", + "| Valeurs | Sans tri | Avec tri |\n", + "| :------------- | :------- | :------- |\n", + "| Optimal value | 1396 | 1396 |\n", + "| Nodes explored | 3053865 | 846447 |\n", + "| time | 29.372s | 4.480s |\n", + "\n", + "Il serait possible d'optimiser encore plus notre algorithme en ne stockant que les changements de valeurs dans l'abre, au lieu de stocker le vecteurs $x$ en entier à chaque noeud de l'arbre. Cela permetterait d'économiser pas mal de RAM.\n", + "\n", + "## Bonus: Comparaison avec GLPK\n", + "\n", + "La huitième cellule de ce notebook permet de convertir les `.opb` en `.lp`, et permet ainsi de résoudre les problèmes en utilisant GLPK.\n", + "Via les fichiers `test.jl` et `distrib.py` on lance sur toutes les machines de l'enseeiht (la nuit), la résolution de chaque problème via GLPK et via notre implémentation. On observe alors que GLPK obtient la même valeur optimale finale que notre implémentation, cependant le nombre de noeuds explorés par GLPK ainsi que son temps d'exécution sont nettement inférieur à ceux de notre algorithme." + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Julia 1.7.0", + "language": "julia", + "name": "julia-1.7" + }, + "language_info": { + "file_extension": ".jl", + "mimetype": "application/julia", + "name": "julia", + "version": "1.7.0" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/sujet.pdf b/sujet.pdf new file mode 100644 index 0000000..7f6863f Binary files /dev/null and b/sujet.pdf differ diff --git a/test.jl b/test.jl new file mode 100644 index 0000000..73803fe --- /dev/null +++ b/test.jl @@ -0,0 +1,242 @@ +function readKnapInstance(filename) + price = [] + weight = [] + capacity = -1 + open(filename) do f + for i = 1:3 + tok = split(readline(f)) + if (tok[1] == "ListPrices=") + for i = 2:(length(tok)-1) + push!(price, parse(Int64, tok[i])) + end + elseif (tok[1] == "ListWeights=") + for i = 2:(length(tok)-1) + push!(weight, parse(Int64, tok[i])) + end + elseif (tok[1] == "Capacity=") + capacity = parse(Int64, tok[2]) + else + println("Unknown read :", tok) + end + end + end + return price, weight, capacity +end + +function TestsSondabilite_relaxlin(x, price, weight, capacity, BestProfit, Bestsol, affich) + TA, TO, TR = false, false, false + + if (!Constraints(x, weight, capacity)) # Test de faisabilite + TA = true + if affich + println("TA\n") + end + + elseif (Objective(x, price) <= BestProfit) # Test d'optimalite + TO = true + if affich + println("TO\n") + end + + elseif (AllDef(x)) # Test de resolution + TR = true + if affich + println("TR : solution ", " de profit ", Objective(x, price), "\n") + end + if (Objective(x, price) >= BestProfit) # Le profit de la solution trouvée est meilleur que les autres + if affich + println("\t-> Cette solution a un meilleur profit.\n") + end + # On remplace la solution et le profit par les nouvelles valeurs + Bestsol = x + BestProfit = Objective(x, price) + else + if affich + println("\t-> Cette solution est moins bonne.\n") + end + end + + elseif affich + println("non sondable\n") + end + + return TA, TO, TR, Bestsol, BestProfit +end + +function SeparerNoeud_relaxlin(price, listvars, listvals) + # Le noeud est non-sondable. Appliquer le critère de séparation pour le séparer en sous-noeuds + # Cas du noeud le plus à gauche + # On sépare le noeud actuel en 2 sous-noeuds + predX = pop!(listvars) + nextX0 = copy(predX) + nextX1 = copy(predX) + + # On initialise leurs valeurs à zéro + val0 = 0 + val1 = 0 + + # On fixe la nouvelle variable des deux sous-noeuds + n = length(predX) + for i = 1:n + if predX[i] == -1 + # L'un a zéro + nextX0[i] = 0 + # L'autre a un + nextX1[i] = 1 + + # On calcule leurs valeurs + val0 = Objective(nextX0, price) + val1 = Objective(nextX1, price) + break + end + end + + # On ajoute les sous-noeuds a la pile des noeuds a explorer + push!(listvars, nextX0) + push!(listvars, nextX1) + + # On ajoute aussi leurs valeurs + push!(listvals, val0) + push!(listvals, val1) + + return listvars, listvals +end + +function ExplorerAutreNoeud_relaxlin(listvars, listvals) + # Le noeud est sondable, on l'enlève de la pile des noeuds à sonder + + stop = false + if (length(listvars) > 1) + # On passe au noeud suivant + var = pop!(listvars) + val = pop!(listvals) + else + # Il n'y a pas d'autre noeud + stop = true + end + + return listvars, listvals, stop +end + +# Fonction objectif que l'on souhaite maximiser/minimiser (évalué dans le meilleur des cas) +Objective(x, price) = + sum( + if x[i] < 0 + price[i] + else + price[i] * x[i] + end + for i = 1:length(x) + ) + +# Fonction permettant de vérfier toutes les contraintes du modèle (dans le meilleur des cas) +Constraints(x, weight, capacity) = + sum( + if x[i] < 0 + 0 + else + weight[i] * x[i] + end + for i = 1:length(x) + ) <= capacity + +# Fonction qui nous dis si toutes les variables de x sont fixées +function AllDef(x) + for i = 1:length(x) + if x[i] < 0 + return false + end + end + return true +end + +function SolveKnapInstance(filename) + + stop = false + affich = true + + # Extraction des données + price, weight, capacity = readKnapInstance(filename) + + tri = true + if tri + couples = zip(price, weight) + couples = sort(collect(couples), by = x -> x[1] / x[2]) + unzip(a) = map(x->getfield.(a, x), fieldnames(eltype(a))) + price, weight = unzip(couples) + end + + expected = split(split(filename, "-")[2], ".")[1] + nodes_nb = length(price) + nodes_max = 2^nodes_nb + + if affich + println("---", filename, "---") + println("Capacity = ", capacity) + println("Number of objects = ", nodes_nb) + println("Expected optimal value = ", expected) + println("Maximum number of nodes = ", nodes_max) + println() + end + + # Liste des variable pour naviguer de noeuds en noeuds + listvars = [] + listvals = [] + listnodes = [] + + # La meilleur solution et sa valeur + BestProfit = -1 + LastBestProfit = -1 + Bestsol = [] + + # Compter le nombre de noeud explorés + current_node_number = 0 + + # On ajoute le premier noeud à explorer (la racine) + push!(listvars, [-1 for p in price]) + push!(listvals, Objective(last(listvars), price)) + push!(listnodes, 1) + newnodeid = 2 + + while (!stop) + + # Le noeud actuel + x = last(listvars) + + if affich && LastBestProfit != BestProfit + println("Node n°", current_node_number, ": \tBestProfit = ", BestProfit) + LastBestProfit = BestProfit + end + + # Test de sondabilité du noeud actuel + # -> On mets a jour la solution et sa valeur si besoin + TA, TO, TR, Bestsol, BestProfit = TestsSondabilite_relaxlin(x, price, weight, capacity, BestProfit, Bestsol, false) + + is_node_sondable = TA || TO || TR + if (!is_node_sondable) + # Le noeud n'est pas sondable, on le sépare en 2 sous-noeuds + listvars, listvals = SeparerNoeud_relaxlin(price, listvars, listvals) + newnodeid += 2 + else + # Le noeud est sondable, on passe au noeud suivant + listvars, listvals, stop = ExplorerAutreNoeud_relaxlin(listvars, listvals) + end + + current_node_number += 1 + end + + println("\n--------------------------------------------------") + println("Expected optimal value = ", expected) + println("Optimal value = ", BestProfit) + println("Optimal x = ", Bestsol) + println("Maximum number of nodes = ", nodes_max) + println("Nodes explored = ", current_node_number) + if parse(Int64, expected) == BestProfit + return 0 + else + return 1 + end + +end + +SolveKnapInstance(ARGS[1])