解題區/解題區 - UVa

UVa 10498 - Happiness

contents

  1. 1. Problem
  2. 2. Sample Input
  3. 3. Sample Output
  4. 4. Solution

Problem

每個人對於每種食材都有其滿意度，目標滿意度總和不超過上限。

準備一人份食材時，每單位的食材都各自有不同的花費，盡可能花最多的錢。

請問總共要花多少錢。

Sample Input

1
2
3
4
5
3 3
1 0.67 1.67
1 2 1 430
3 0 2 460
1 4 0 420

Sample Output

1
Nasa can spend 1354 taka.

Solution

每個人都花費相同的金額，因此把目標函數最大化後，乘上總人數即可。

現在要把問題轉換成線性規劃的模型，隨後套用 Simplex algorithm (單純形法)，可以參考 李宇骞《浅谈信息学竞赛中的线性规划——简洁高效的单纯形法实现与应用》的簡報，概略說明模型的架構。資料網路上很多，不談論標準的處理程序。

接著閒聊自己的看法。

利用數形轉換的概念，從幾何概念下手，從三個變數來看，單純形法相當於在一個空間上，拿一個物體與平面找接觸，為了使目標函數最大化，會將物體表面的頂點帶入嘗試。

但是要將所有頂點找出來是相當消耗時間，因此爬行法可能會稍微好一點，也就是說，沿著物體表面，藉由一條邊爬行到更好位置上的頂點。由於線性約束，使得物體一定是凸多面體，在更高維度也是一樣的，不斷地爬升到更好的頂點後，最後一定會在最高點停止。

為了使爬行速度加快，每一次找爬升最多的鄰近點，相當於根據某一個約束，盡可能放大某一個維度上的變數 (有可能會將別的變數變小)，使得目標函數增加的值盡可能地多。根據原始的約束寫法，要找到計算仍然相當困難。

定義一個轉軸操作，讓整個問題的空間座標進行變換，方便撰寫。

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
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
#include <stdio.h>
#include <math.h>
#include <string.h>
#include <algorithm>
using namespace std;
const int MAXN = 64;
const int MAXM = 128;
const double eps = 1e-10, INF = 1e+30;
/*
	n: #variable
	m: #inequality constraint
	#inequality constraint format : \sum{a(j, i) xi} <= bj
	for all xi >= 0
 */
class Simplex {
public:
    int n, m;
    int N[MAXN], M[MAXM];
    double a[MAXM][MAXN], b[MAXM], c[MAXN];
    double xval[MAXN], z;
    void init(int n, int m) {
        memset(a, 0, sizeof(a));
        memset(b, 0, sizeof(b));
        memset(c, 0, sizeof(c));
        memset(xval, 0, sizeof(xval));
        z = 0;
        this->n = n + m;
        this->m = m;
        
        for (int i = 0; i < m; i++)
            a[i][n+i] = 1;
        for (int i = 0; i < n; i++)
            N[i] = i;
        for (int i = 0; i < m; i++)
            M[i] = n + i;
    }
    int cmpZero(double v) {
        if (fabs(v) < eps)	return 0;
        return v < 0 ? -1 : 1;
    }
    void compute_xval() {
        for (int i = 0; i < n - m; i++) {
            xval[i] = 0;
            for (int j = 0; j < m; j++) {
                if (M[j] == i)
                    xval[i] = b[j];
            }
        }
    }
    void print() {
        printf("c ");
        for (int i = 0; i < n; i++)
            printf("%6.3lf ", c[i]);
        printf("%6.3lf\n", z);
        
        for (int i = 0; i < m; i++) {
            printf("a ");
            for (int j = 0; j < n; j++)
                printf("%6.3lf ", a[i][j]);
            printf("%6.3lf\n", b[i]);
        }
        puts("--");
        
        compute_xval();
        printf("x ");
        for (int i = 0; i < n - m; i++)
            printf("%6.3lf ", xval[i]);
        puts("");
    }
    void pivot(int row, int col) {
        double e = a[row][col];
        // normalization
        b[row] /= e;
        for (int i = 0; i < n; i++)
            a[row][i] /= e;
        // gaussian elimination
        for (int i = 0; i < m; i++) {
            if (i == row)
                continue;
            double t = a[i][col];
            b[i] = b[i] - t * b[row];
            for (int j = 0; j < n; j++)
                a[i][j] = a[i][j] - t * a[row][j];
        }
        z -= b[row] * c[col];
        
        double t = c[col];
        for (int i = 0; i < n; i++)
            c[i] = c[i] - t * a[row][i];
        swap(N[col], M[row]);
//        print();
    }
    void simplex() {
        double mx;
        int select_col, select_row, mn_row;
        // climb
        do {
            mx = -INF, select_col = select_row = mn_row = -1;
            for (int i = 0; i < n; i++) {
                if (cmpZero(c[i]) > 0) {
                    double mn = INF;
                    // find minmum increase with inequality constraint
                    for (int j = 0; j < m; j++) {
                        if (cmpZero(a[j][i]) > 0) {
                            if (b[j] / a[j][i] < mn)
                                mn = b[j] / a[j][i], mn_row = j;
                        }
                    }
                    // find maxmum climb with axis
                    if (mx < mn * c[i]) {
                        mx = mn * c[i];
                        select_col = i, select_row = mn_row;
                    }
                }
            }
            if (select_col != -1)
                pivot(select_row, select_col);
        } while (select_col != -1);
        
        compute_xval();
    }
} g;
int main() {
    int n, m;
    while (scanf("%d %d", &n, &m) == 2) {
        g.init(n, m);
        
        for (int i = 0; i < n; i++)
            scanf("%lf", &g.c[i]);
        for (int i = 0; i < m; i++) {
            for (int j = 0; j < n; j++)
                scanf("%lf", &g.a[i][j]);
            scanf("%lf", &g.b[i]);
        }
        
        g.simplex();
        
//        for (int i = 0; i < n; i++) {
//            printf("x[%d] = %lf\n", i, g.xval[i]);
//        }
        printf("Nasa can spend %d taka.\n", (int) ceil(-g.z * m));
    }
    return 0;
}
/*
 Max x1 + 14*x2 + 6*x3
 s.t.
 x1 + x2 + x3 	<= 4
 x1              <= 2
 x3              <= 3
 3*x2   + x3     <= 6
 x1 , x2 , x3 , x4 , x5 , x6 , x7 >= 0
 
 3 4
 1 14 6
 1 1 1  4
 1 0 0  2
 0 0 1  3
 0 3 1  6
 
 2 2
 1 1
 2 1 12
 1 2 9
 
 2 2
 -1 -1
 2 1 12
 1 2 9
 
 3 3
 1 0.67 1.67
 1 2 1 430
 3 0 2 460
 1 4 0 420
 
 3 1
 30 10 20
 30 1 3 30
 
 3 3
 30 10 20
 30 1 3 30
 30 1 3 30
 30 1 3 30
 */