ibexopt et variables entières

[bneveu]

J'ai vu dans #168 qu'on pouvait utiliser ibexopt avec des variables entières en utilisant la contrainte floor. les résultats que l'on obtient sont parfois un reel "proche" d'un entier. Est-ce normal et suffit-il d'arrondir à l'entier le plus proche pour avoir la solution ?

exemple exint.bch

variables x in [-10 ,1]; minimize x^2+x; constraints x= floor(x); end

ibexopt exint.bch

optimization successful!

f* in [-2.00000001004e-08,-6.43762076979e-09] (best bound)

x* = (-0.999999993562) (best feasible point)

relative precision on f: 0.678118963128 absolute precision on f: 1.35623793308e-08 [passed] cpu time used: 0.00458300000001s number of cells: 12

[gchabert]

Bertrand, Utilise plutôt directement la contrainte integer(x). Cela te permettra d'utiliser le mode rigoureux et d'avoir un encadrement de l'entier. Je ne l'ai pas encore documenté, je me le note !

[bneveu]

Ok pour cet exemple, mais le mode rigor ne marche pas toujours pour l'exemple suivant (optimum sur une borne) , il met un warning et ne resout pas en 10s.

exint0.bch variables x in [-10 ,-1]; minimize x^2+x; constraints x= floor(x); end

ibexopt --rigor exint0.bch running............

warning: too many active constraints, cannot prove feasibility -> loup lost! time limit 10s. reached

f* in [0,inf] (best bound)

x* = -- (no feasible point found)

---

le mode normal marche

ibexopt exint0.bch f* in [-0,0] (best bound)

x* = (-1) (best feasible point)

[bneveu]

l'exemple précédent marche avec le mode rigor et la contrainte integer(x) à la place de x=floor(x) . il y a juste le warning qui reste.

ibexopt --rigor exint0.bch warning: too many active constraints, cannot prove feasibility -> loup lost! optimization successful!

f* in [-0,0] (best bound)

x* in (<-1, -1>) (best feasible point)

[gchabert]

C'est normal que le mode rigoureux échoue avec la contrainte x=floor(x) qui est singulière sur tout entier x.

Pour le warning avec interger(x), je pense que cela doit apparaître lorsqu'il essaye de traiter la borne du domaine x=1. Il faudrait mettre une trace pour en être sûr. Mais il trouve bien un minimum global certifié, donc tout va bien.

[bneveu]

Il faudrait, je pense enlever ce warning, qui peut remplir des pages entières dans la résolution de certains problèmes minlp

[trombettoni]

Bonjour Gilles,

Bertrand est motivé pour faire du mixte depuis notre réunion de mardi avec les gens du CEDRIC qui était très intéressante (on pourra faire un point par mail sur ce que la collaboration peut apporter, notamment grâce à leur expertise sur les systèmes quadratiques - convexes ou non).

Juste avant cet échange sur le forum de Github, on avait réfléchi à une variante de IbexOpt assez simple où on ajoutait de manière un peu sale/lourde des appels à du "ctcInteger(x) sur toutes les variables" à différentes étapes, par exemple les grosses étapes de la contraction, après HC4, après ACID(HC4), après polytopeHull et HC4 dans le point-fixe, etc.

J'ai l'impression, Gilles, que ton idée est de considérer la cte integer(x) comme les autres ctes numériques, ce qui est effectivement plus élégant. On voit l'avantage sur HC4 où à chaque HC4R d'une cte "normale" qui touche au domaine d'une var entière x, la boucle de propagation remet la cte integer(x) dans la queue...

Mais il y a plein de questions du fait que ce n'est pas une cte très classique (continue, dérivable), etc.

Je pense que ce serait bien qu'on se fasse une visio avec Bertrand et ceux qui ont implanté cette cte integer(x), au moins toi donc, pour bien comprendre comment elle est implantée, ce que ça implique dans PolytopeHull , dans le mode rigoureux et surtout dans LoupFinder (ce n'est pas implanté dans INHC4 si j'ai bien compris, mais quid de InXTaylor ?).

Tu pourrais avoir une petite heure la semaine prochaine sur le sujet ?

Bon WE à tous en attendant !

[raphaelchenouard]

Bonjour,

De notre côté aussi, on regarde des choses sur le mixte. On a commencé à définir un bisecteur adapté, mais ce n'est pas forcément facile de récupérer les variables entières dans un système, car seule la contrainte integer permet de les identifier si on part de minibex.

Bon week-end, Raphaël.

[gchabert]

Salut les gars,

Gilles, voici quelques explications. La contrainte Minibex "integer(x)" est en fait immédiatement traduite dans Ibex par "saw(x)=0" avec "saw" la fonction en vert sur la figure suivante:

Comme tu peux le voir, la fonction est discontinue entre 2 entiers mais tout à fait continue et dérivable autour de chaque entier. Cela permet donc aussi d'appliquer un Newton avec cette contrainte, d'où la possibilité d'avoir des solutions certifiées! (mode rigoureux). Ce qui a été implémenté avec la fonction saw dans ibex, c'est

• l'évaluation
• le backward (donc HC4 l'accepte)
• le gradient (donc inXTaylor l'accepte).

Par contre, en effet, pas d'implémentation pour inHC4, mais je peux regarder si vous voulez. Ca paraît simple mais il y a le cas des points au milieu des entiers qui est un peu pénible à gérer proprement en terme d'arrondi. J'avais d'ailleurs un peu galéré pour le backward.

Raphaël, effectivement, il n'y a pas de moyen direct de reconnaître une variable entière dans Ibex mais je dirais que c'est presque le but recherché: on traite tout en continu. Je comprends que ça puisse rendre les choses un peu laborieuses si vous souhaitez traiter différemment les variables en fonction de leur type, typiquement pour un bissecteur. Finalement, je trouve ta solution assez adaptée.

On peut prévoir une petite réunion si vous avez encore des questions.

[bneveu]

J'ai essayé quelques problèmes MINLP de http://www.minlplib.org/instances.html

Pour le problème gear4, ibexopt en mode rigueur trouve la solution pour les variables entières, mais n'arrive pas à la prouver avec une précision de 1.e-6. Il y des petites boites non bissectables qu'il n'arrive pas à éliminer et qui du coup empêchent la remontée du loup voici le problème et la trace En mode non rigoureux, on trouve une solution approchée (la contrainte d'intégrité étant relâchée) prouvée

---

variables i1 in [12,60]; i2 in [12,60]; i3 in [12,60]; i4 in [12,60]; x6 in [0,1.e8]; x7 in [0,1.e8];

minimize x6+x7;

constraints integer(i1); integer(i2); integer(i3); integer(i4); -1000000i1i2/(i3*i4) - x6 + x7 + 144279.32477276 =0;

end

---

../../../bin/ibexopt gear4.bch -a 1.e-6 -r 1.e-6 --trace --rigor -t 100

**** setup **** file loaded: gear4.bch rel-eps-f: 1.00001e-06 (relative precision on objective) abs-eps-f: 1.00001e-06 (absolute precision on objective) rigor mode: ON (feasibility of equalities certified) output COV file: gear4.cov timeout: 100s trace: ON

warning: inHC4 disabled (unimplemented operator)

---

running............

                 loup= 99144279.3248

uplo= 0 loup= 61144229.7527 loup= 30144254.5388 loup= 14644266.9318 loup= 13788531.6123 loup= 6466405.46853 loup= 2805342.39665 loup= 974810.860711 loup= 4257.57868696 loup= 3359.77486098 loup= 3359.77486098 loup= 3358.84520982 loup= 3358.84520982 loup= 3269.67863073 loup= 3269.21380516 loup= 1792.46355922 loup= 1551.47845509 loup= 1551.47845509 loup= 1550.54880393 loup= 1550.54880393 loup= 300.542816745 loup= 299.613165594 loup= 245.175072813 loup= 244.245421662 loup= 244.245421662 loup= 137.03993092 loup= 136.110279768 loup= 116.424349663 loup= 49.4851931424 loup= 33.0773314643 loup= 1.643428474 uplo= 0.410856707643 uplo= 0.464825575781 uplo= 0.821713415285 unprocessable tiny box: now uplo<=1.64324968274 uplo=0.821713415285 uplo= 0.929651151562 unprocessable tiny box: now uplo<=1.64319935212 uplo=0.929651151562 uplo= 1.23257012293 uplo= 1.28653899107 uplo= 1.39447672735 uplo= 1.39447672735

time limit 100s. reached

f* in [1.64319935212,1.643428474] (best bound)

x* in (<16, 16> ; <19, 19> ; <43, 43> ; <49, 49> ; <-0, 0> ; [1.64342847393734, 1.643428473995548]) (best feasible point)

relative precision on f: 0.00013943645017 absolute precision on f: 0.000229121884579 cpu time used: 99.9992900001s number of cells: 199256

---

../../../bin/ibexopt gear4.bch -a 1.e-6 -r 1.e-6 --trace -t 100

**** setup **** file loaded: gear4.bch rel-eps-f: 1.00001e-06 (relative precision on objective) abs-eps-f: 1.00001e-06 (absolute precision on objective) output COV file: gear4.cov timeout: 100s trace: ON

warning: inHC4 disabled (unimplemented operator)

---

running............

uplo= 0 uplo= 0.372529029847 uplo= 0.745058059693 uplo= 1.11758708954 uplo= 1.49011611939 loup= 1.64322226451 loup= 1.64321998822 loup= 1.64320100468 loup= 1.64319936037

optimization successful!

f* in [1.64319771717,1.64319936037] (best bound)

x* = (18.9999999901 ; 15.9999999901 ; 43.00000001 ; 49.00000001 ; 2.13763573034e-09 ; 1.64319935823) (best feasible point)

relative precision on f: 9.99999999902e-07 [passed] absolute precision on f: 1.64319771701e-06 cpu time used: 50.8511030001s number of cells: 199588

[trombettoni]

Salut les jeunes,

Merci Chabs pour les explications.

Très rigolo (et élégant) l'utilisation d'une cte integer(x) et d'une implémentation utilisant cette fonction saw(.)... Qui a inventé ça ?? C'est toi, Chabs ? C'est vrai que c'est "so interval" comme approche...

Je pense qu'une réunion pour préciser les conséquences de ce choix sera la bienvenue, mais autant attendre le premier retour sur expérience de Bertrand qui expérimente sur des instances de : http://minlplib.org/instances.html

J'aurai déjà deux remarques qu'on pourra discuter à cette réunion :

a) Pourquoi utiliser integer(x) := saw(x) == 0 plutôt qu'une fonction continue comme integer(x) := sin(\Pi x) == 0 ? (ou sin (2 \Pi * x) ==0)

Bertrand a déjà expérimenté et ça semble donner des résultats plutôt meilleurs. Je comprends l'intérêt d'avoir une fonction linéaire par morceau, mais le sinus est pas mal foutu non plus (surtout au voisinage des entiers) et a l'avantage d'être continu, déjà implémenté dans Ibex, etc.

b) J'admire l'élégance et la concision de la solution "contrainte integer", mais je ne suis pas entièrement, entièrement convaincu de son efficacité. Pour la "contraction" de integer(x), un arrondi de [x] sur les entiers "fait main" (au dessus pour inf(x) et en dessous pour sup(x)) paraît plus efficace que les calculs sur des fonctions saw ou sin, et éviterait les problèmes de "rigueur" (on peut être rigoureux avec ou sans --rigor), mais l'intérêt de l'implémentation actuelle est surtout le calcul de la dérivée je suppose, dérivée qui offre X-Taylor, In-Xtaylor, voire l'arithmétique affine avec une implémentation par sin(\Pi * x) == 0 (j'imagine pas trop le truc avec saw(x)==0...). Il faudrait peut-être pouvoir avoir le meilleur des deux mondes : des arrondis "faits main" pour la contraction et un passage à saw/sinus pour la dérivée...

Bref, pas mal de trucs rigolos à discuter après un petit retour sur expérience...

[bneveu]

ces contraintes sont effectivement trop faibles en contraction. En mode rigor, elles n'arrivent pas à éliminer des petites boîtes ne contenant pas de solution entière et ne trouvant pas de loup dans ces boites, ces deviennent non bissectables et empêchent la remontée du uplo. exemple dans gear4.bch la boite suivante [18.99999998999999, 18.99999999000001] ; [15.99999999, 15.99999999000001] ; [49.00000000999999, 49.00000001000001] ; [43.00000000999999, 43.00000001000001] ; [2.220446048083475e-16, 2.220446048083476e-16] ; [1.643199352110968, 1.64319935211097] ; [1.643199352110968, 1.64319935211097]) n'est pas éliminée , alors que les 4 premières variables sont entières. En mode non rigoureux, on trouve des loups dans de telles boites et l'algorithme converge vers une solution non entière.

[trombettoni]

Hum.

Et quid avec un "eps-h" dédié à ces contraintes (en mode non rigoureux) très précis, de l'ordre de 1e-100, soit :

-1e-100 <= saw(x) <= +1e-100

?

[bneveu]

On se rapproche effectivement de l'optimum de la solution entière en baissant eps-h [1.64342102158,1.64342118592] alors que l'optimum serait en 1.643428 pour eps-h = 2.e-10 , mais ensuite on n'arrive plus à trouver de loup dès eps-h = 1.e-10

[bneveu]

J'ai essayé d'ajouter dans optimizer04 l'appel à la contrainte CtcInteger sur les variables. Cela marche bien pour obtenir des solutions entières et certains problèmes sont bien résolus en mode normal et trouvent la solution entière (pour les variables entières) indiquée dans la base minlp. D'autres ne trouvent pas de loup car la solution entière échoue au test is_inner. En effet, elle est souvent sur la frontière de la contrainte , en particulier pour les inégalités ne portant que sur des entiers. Ce test est donc trop fort. Comment pourrait-on lui faire accepter ces points frontières ?

[bneveu]

En implémentant les méthodes is_ineffective et effective_ctrs dans ibex_System, comme suggéré dans les commentaires de ibex_System.cpp et en appelant effective_ctrs dans is_inner à la place de active_ctrs , on arrive à résoudre ces problèmes minlp en acceptant des loup avec des solutions entières. Mais il faut savoir quelles sont les variables entières, il faudrait donc le système connaisse ses variables entières sur lesquelles on appelerait le contracteur CtcInteger (en l'ajoutant avec CtcCompo aux différents contracteurs appelés par l'optimiseur). Pour les essais suivants, j'ai dû pour le moment le faire à la main en l'indiquant pour chaque problème dans optimizer04.cpp ses variables entières.

Par ailleurs, il semble toujours utile d'avoir une relaxation numérique de la contrainte d'intégrité. Pour cela, il semble préférable de remplacer saw(x)=0 par sin(pix=0) : on peut gagner un ordre de grandeur sur certains problèmes. par exemple : ex1252a : ce problème qui a 15 variables réelles et 9 entières (dont 3 booléennes) est résolu: avec CtcInteger en 258s. avec CtcInteger et saw(x_i)=0 en 133s. avec CtcInteger et sin(pix_i)=0 en 14s. pour info, ce problème est résolu par Antigone en 1052s et par Baron en 25s.

J'ai mis quelques problèmes MINLP (les plus petits en nombre de variables de la base http://www.minlplib.org/instances.html comprenant des variables entières et réelles) traduits à la main en minibex dans un nouveau répertoire benchs/optim/benchs-minlp J'ai borné au besoin les variables réelles dans [-1.e8,1.e8] et ajouté des contraintes integer(x) pour indiquer les variables entières.

Je continue à tester ces exemples.