paz_math/murty.cpp at main · pazvirtual/paz_math · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
#include "PAZ_Math"

void paz::murty(const MatRef& costMat, std::size_t numBest, std::vector<std::
    vector<std::size_t>>& rowSols, std::vector<double>& costs)
{
    rowSols.clear();
    costs.clear();

    if(!numBest)
    {
        return;
    }

    std::vector<std::size_t> bestSol;
    const double bestCost = jv(costMat, bestSol);

    if(numBest == 1 || bestCost == paz::inf())
    {
        rowSols = {bestSol};
        costs = {bestCost};
        return;
    }

    const std::size_t rows = costMat.rows();
    const std::size_t cols = costMat.cols();
    if(rows > cols)
    {
        throw std::runtime_error("Matrix is too tall.");
    }

    std::vector<Mat> costMatsList = {costMat};
    std::vector<std::vector<std::size_t>> solsList = {bestSol};
    std::vector<double> costsList = {bestCost};

    while(true)
    {
        const std::size_t idx = std::distance(costsList.begin(), std::
            min_element(costsList.begin(), costsList.end()));
        costs.push_back(costsList[idx]);
        rowSols.push_back(solsList[idx]);

        if(rowSols.size() == numBest)
        {
            break;
        }

        Mat curCostMat = costMatsList[idx];
        const auto curSol = solsList[idx];

        std::swap(costMatsList[idx], costMatsList.back());
        costMatsList.pop_back();
        std::swap(solsList[idx], solsList.back());
        solsList.pop_back();
        std::swap(costsList[idx], costsList.back());
        costsList.pop_back();

        for(std::size_t i = 0; i < curSol.size(); ++i)
        {
            if(curSol[i] != None)
            {
                if(curCostMat(i, curSol[i]) == inf())
                {
                    throw std::logic_error("Current solution should have been s"
                        "kipped due to infinite cost.");
                }

                Mat tempCostMat = curCostMat;
                tempCostMat(i, curSol[i]) = inf();

                std::vector<std::size_t> tempSol;
                const double tempCost = jv(tempCostMat, tempSol);
                if(tempCost < inf() && std::find(tempSol.begin(), tempSol.
                    end(), None) == tempSol.end())
                {
                    costMatsList.push_back(tempCostMat);
                    solsList.push_back(tempSol);
                    costsList.push_back(tempCost);
                }

                const double temp = curCostMat(i, curSol[i]);
                for(std::size_t j = 0; j < cols; ++j)
                {
                    curCostMat(i, j) = inf();
                }
                for(std::size_t j = 0; j < rows; ++j)
                {
                    curCostMat(j, curSol[i]) = inf();
                }
                curCostMat(i, curSol[i]) = temp;
            }
        }

        if(costsList.empty())
        {
            break;
        }
    }
}