Decision Optimization

Decision Optimization

Delivers prescriptive analytics capabilities and decision intelligence to improve decision-making.

 View Only

June 2019 IBM Ponder This challenge

  • 1.  June 2019 IBM Ponder This challenge

    Posted Thu July 04, 2019 04:08 AM

    Hi,

     

    the June challenge was http://www.research.ibm.com/haifa/ponderthis/solutions/June2019.html

     

    and an additional challenge is to use OPL CPLEX.

     

     

     

    Let me do that in 2 steps:

    1) Use MIP to find a solution that is less than 1e-5

    2) And then use CPO to turn that into a good solution for 1e-15

     

     

    1)

    {int} p={2,3,5,7};
    int n=1;
    int n2=2;
    int nbRootsMax=4;
    int maxy=10;


    range r=0..n;
    range r2=0..n2;

    {int} optionsy={ si*ftoi(pow(i,a))| i in p,a in r2,si in -1..1};


    sorted {int} smooths={ ftoi(pow(2,a)*pow(3,b)*pow(5,c)*pow(7,d)) | a,b,c,d in r};


    {float} sqr={ sqrt(i) | i in smooths} diff {1.0};

     

     

    dvar int y[sqr] in -maxy..maxy;

    dexpr float delta=sum(i in sqr) y[i]*i;


    subject to
    {
    1<=sum(i in sqr) (y[i]!=0)<=nbRootsMax;
    -1e-6<=delta<=1e-6;


    forall(i in sqr) (sum(j in optionsy) (y[i]==j))==1;
    }

     

    execute
    {

     

    for(var i in sqr) if (y[i]!=0) writeln("sqrt(",i*i,") * ",y[i]);
    writeln();

    }

    {int} goodsquares={sgn(y[i])*y[i]*y[i]*ftoi(round(i*i)) | i in sqr : y[i]!=0};

    execute
    {
    writeln(goodsquares);
    }

    // But let 's use power 3 in order to get the difference less than 10 power -15

    int nbCubeGoodSquares=ftoi(pow(card(goodsquares),3));
    float cubeGoodSquares[i in goodsquares][j in goodsquares][k in goodsquares]=
       i*j*k;
       
    execute
    {
    var output=new IloOplOutputFile("cubeGoodSquares.dat");
    output.write("nbCubeGoodSquares=");
    output.write(nbCubeGoodSquares);
    output.writeln(";");
    output.write("cubeGoodSquares=[");
    for(var i in goodsquares) for(var j in goodsquares) for(k in goodsquares)
      output.writeln(cubeGoodSquares[i][j][k],",");
      output.writeln("];");
      output.close();

    }
       

    which generates cubegoodsquares.dat and

    sqrt(5) * 1
    sqrt(7) * -9
    sqrt(21) * -2
    sqrt(105) * 3

     {5 -567 -84 945}

    2) Use CPO:

    using CP;

    execute
    {
    cp.param.timelimit=120;
    }

    {int} p={2,3,5,7};
    int n=1;

    int maxy=10;

    range r=0..n;

    sorted {int} smooths={ ftoi(pow(2,a)*pow(3,b)*pow(5,c)*pow(7,d)) | a,b,c,d in r};
    {int} sqr=smooths; // diff {1};


    int nbCubeGoodSquares=...;

    range goodsquares=1..nbCubeGoodSquares;
    int cubeGoodSquares[goodsquares]=...;

    dvar int y[goodsquares] in 1..1000000;
    dvar int squaresmooth[goodsquares];
    dvar int s[goodsquares] in -1..1;
    dvar int multiplier[smooths];

    range r4=1..4;

    dvar int power[smooths][r4][p] in 0..8;
    dvar boolean needed[smooths][r4]; // do we need that term ?
    dvar int term[smooths][r4];

    minimize sum(i in smooths,j in r4) needed[i][j];
    subject to
    {
    forall(i in goodsquares)
    {
    squaresmooth[i] in sqr;
    cubeGoodSquares[i]==y[i]*y[i]*s[i]*squaresmooth[i];
    }

    forall(i in smooths) multiplier[i]==sum(j in goodsquares) (squaresmooth[j]==i)*y[j]*s[j];

    // and now decompose the multipliers as sum of smooth numbers

    forall(i in smooths)
    {
    ftoi(abs(multiplier[i]))==sum(j in r4) needed[i][j]*term[i][j];
    forall(j in r4)term[i][j]==prod(k in p)ftoi(pow(k,power[i][j][k]));
    }

    }

    execute
    {
    writeln("the cube");
    for(var i in smooths) if (multiplier[i]!=0) writeln(multiplier[i]," * sqr of ",i);
    }

     

    {int} list1={  term[i][j]*term[i][j]*i     | i in smooths,j in r4 : (multiplier[i]>=1) && (1==needed[i][j])};
    {int} list2={  term[i][j]*term[i][j]*i     | i in smooths,j in r4 : (multiplier[i]<=-1) && (1==needed[i][j])};

    execute
    {
    writeln();
    writeln("So the two lists are ");

    writeln(list1);
    writeln(list2);
    }

    with that .dat and then we get the solution:

     

    So the two lists are
     {2657205 51200000 3732480 8573040 1620304560 6827950080 47840625}
     {14288400 96446700 326592 7533176175 1701 1801088541}

    3) We can even double check with python high precision that we're right:

     

    from decimal import *
    getcontext().prec = 40


    l1=[2657205, 51200000, 3732480, 8573040, 1620304560, 6827950080, 47840625]
    l2=[14288400 ,96446700, 326592 ,7533176175 ,1701, 1801088541]

    s1=sum([Decimal.sqrt(Decimal(i)) for i in l1])
    s2=sum([Decimal.sqrt(Decimal(i)) for i in l2])

    print(s1)
    print(s2)

    print(Decimal(s2-s1));

    which gives

    143446.5588138200170046814527613232168467
    143446.5588138200170043256720431189535309
    -3.557807182042633158E-16

    So the solution is fine.

     

    Again MIP and CPO within OPL CPLEX were great

     

    regards

    https://www.linkedin.com/pulse/puzzles-having-fun-useful-mathematics-alex-fleischer/

     


    #DecisionOptimization
    #OPLusingCPLEXOptimizer