First, let me share some tragic experiences while writing this article. Initially, I was using Gridea and editing with its editor, but during one writing session, my computer unexpectedly blue-screened. After rebooting, I found that the md file was unreadable, and upon checking, the binary was all 0 emmm (I had saved during the writing process!!). Since then, I switched to hexo emmm.
Symbol Definition#
First, we define some symbols: the number of items is $n$, the volume of the $i$-th item is $v_i$, the value is $w_i$, and the maximum number is $s_i$. The volume of the knapsack is $V$, and the dynamic programming array is $dp$.
Review#
Let's first review the optimization of the multiple knapsack problem (where $x$ is the maximum number of items that can fit in the knapsack when its volume is $j$, i.e., $j / v$). Assume that the volume of the $i$-th item is $v$.
dp[i, j]= max(dp[i - 1, j], dp[i - 1, j - v] + w, dp[i - 1, j - 2v] + 2w, ......,dp[i - 1, j - (x - 1) v] + (x - 1)w, dp[i - 1, j - xv] + xw)
Next, let's consider the relationship between $dp[i,j]$ and $dp[i,j - v]$ in the case of the complete knapsack.
dp[i, j - v]= max(dp[i - 1, j - v], dp[i - 1, j - 2v] + w, dp[i - 1, j - 3v] + 2w, ......, dp[i - 1, j - (x - 1) v] + (x - 2)w, dp[i - 1, j - xv] + (x - 1)w)
(Some students may ask why the last knapsack can still hold $x$ items, but actually, that is not $x$, it is $j-v-(\frac{j-v}{v})v$, and the final state reduces by $xv$.)
Then we are amazed to find that the terms after the first term of $max()$ are so similar to $dp[i,j - v]$, differing only by a $w$, so we can factor out $w$ and replace it with:
dp[i, j]= max(dp[i - 1, j], dp[i][j - v] + w)
Denial#
Now let's look at the multiple knapsack problem; we can also try to optimize it (let's assume that $j$ is much larger than $sv$ at this point).
We first try to see if we can use the complete knapsack method to optimize:
dp[i, j]= max(dp[i - 1, j], dp[i - 1, j - v] + w, dp[i - 1, j - 2v] + 2w, ......,dp[i - 1, j - (s - 1) v] + (s - 1)w, dp[i - 1, j - sv] + sw)
dp[i, j - v]= max(dp[i - 1, j - v], dp[i - 1, j - 2v] + w, dp[i - 1, j - 3v] + 2w, ......, dp[i - 1, j - (s - 1) v] + (s - 2)w, dp[i - 1, j - sv] + (s - 1)w)
It seems to be correct, doesn't it? (However, if we consider the number of terms, we will find that the number of terms in $dp[i, j-v]$ is the same as that in $dp[i,j]$, both having $s$ terms, which does not match the number of terms after the first term of $max()$, thus we cannot use the previous method to optimize.)
New Expansion#
(From this section onward, we will use a one-dimensional array for convenience; if you haven't learned about one-dimensional array representation, go learn it.)
We still consider $dp[j]$ and list the terms.
dp[j]= max(dp[j], dp[j - v] + w, dp[j - 2v] + 2w, ......,dp[j - (s - 1) v] + (s - 1)w, dp[j - sv] + sw)
We observe this term and can easily derive a property: the states that update state $j$ and $j$ are congruent modulo $v$. Thus, we can divide the states from $1$ to $V$ into $v$ classes, namely those congruent to $0$ modulo $v$, those congruent to $1$ modulo $v$, $\cdots$, those congruent to $v-2$ modulo $v$, and those congruent to $v-1$ modulo $v$.
The original dp[1]~dp[V] can be transformed into these $v$ classes (this is a bit hard to imagine, but you must understand it).
dp[0] dp[v] dp[2v] dp[3v] ... dp[kv]
dp[1] dp[1+v] dp[1+2v] dp[1+3v] ... dp[1+kv]
dp[2] dp[2+v] dp[2+2v] dp[2+3v] ... dp[2+kv]
...
dp[v-1] dp[v-1+v] dp[v-1+2v] dp[v-1+3v] ... dp[v-1+kv]
Then we can update separately in these $v$ classes. (Where $... + kv$ is exactly the last state that satisfies, see the code.)
Let’s write such code (first write a two-dimensional version) (if changed to one-dimensional, you need to enumerate states in reverse or use a rolling array).
for(int i = 1;i <= n;i++)
{
for(u = 0;u < v[i];u++)
{
for(int k = 0;u + k * v[i] <= V;k++)
{
int e = u + k * v[i];
for(int t = min(k - s * v[i], u)/*prevent going negative*/;t < k;t += v[i])
{
dp[i][e] = max(dp[i - 1][e], dp[i - 1][t] + (e - t) / v[i] * w[i]);
}
}
}
}
Yes, you read that right, this is a four-layer loop. You might ask, why are there more and more loops, but let’s analyze its time complexity. The first layer loop is $O(n)$, the second layer loop enumerates the starting value of each class, the third layer enumerates each state in that class, and the total time complexity of the second and third layers is $O(V)$, while the last layer enumerates the states that can transition to the current state, which is $O(s)$, so the overall time complexity is still $O(nVs)$.
If the above is not clear, let’s take an example.
Let’s take an example: suppose there are $2$ items, $v_1 = 5$, $w_1 = 6$, $s_1 = 3$, $v_2 = 4$, $w_2 = 8$, $s_2 = 5$, and the knapsack volume is $31$.
When $i = 1$, i.e., enumerating the first item, the classification of the dp states is:
dp[0] dp[5] dp[10] dp[15] dp[20] dp[25] dp[30]
dp[1] dp[6] dp[11] dp[16] dp[21] dp[26] dp[31]
dp[2] dp[7] dp[12] dp[17] dp[22] dp[27]
dp[3] dp[8] dp[13] dp[18] dp[23] dp[28]
dp[4] dp[9] dp[14] dp[19] dp[24] dp[29]
// First (remainder 0) class
dp[1][0] = dp[0][0]
dp[1][5] = max(dp[0][5], dp[0][0] + w[1])
dp[1][10] = max(dp[0][10], dp[0][5] + w[1], dp[0][0] + 2 * w[1])
dp[1][15] = max(dp[0][15], dp[0][10] + w[1], dp[0][5] + 2 * w[1], dp[0][0] + 3 * w[1])
dp[1][20] = max(dp[0][20], dp[0][15] + w[1], dp[0][10] + 2 * w[1], dp[0][5] + 3 * w[1])
// At most select 3, so the maximum volume (second dimension) can reduce at most 3 * v[i]
dp[1][25] = max(dp[0][25], dp[0][20] + w[1], dp[0][15] + 2 * w[1], dp[0][10] + 3 * w[1])
dp[1][30] = max(dp[0][30], dp[0][25] + w[1], dp[0][20] + 2 * w[1], dp[0][15] + 3 * w[1])
// The same applies below
$i = 2$ is the same.
Monotonic Queue Optimization#
You may be wondering, after saying so much, why haven't we seen the shadow of the monotonic queue. In fact, everything we did before was to lay the groundwork for it, and here it comes.
First, what we want to optimize is the fourth layer of the above code: updating the current state from the previous $s$ states (which are not adjacent, with a difference of $v$ between any two).
Let’s think about this intuitively, look at the picture! (qwq
In this picture, $j$ represents the current state, and $j'$ represents the next state, which is $j + v$. The gray boxes represent the states that can update the current state (here $j$ or $j'$).
!!!!!!! Note: The leftmost box in the first row of the image should be $j-sv$, I mistakenly wrote it wrong (too lazy to change it emmm).
We can observe that as the state increases by $v$, the expression dp[j]= max(dp[j], dp[j - v] + w, dp[j - 2v] + 2w, ......,dp[j - (s - 1) v] + (s - 1)w, dp[j - sv] + sw will lose the leftmost state that can update it (the one on the far left in the picture) while also gaining the rightmost state that can update it (the one on the far right in the picture).
Let’s think about it: what data structure can insert at both ends (because $j$ increases by $v$ each time, the max function will gain one more state available to update the current state $j$), delete data (which will also reduce one), and support finding the maximum value in a small ($O(1)$) time complexity (as required by the max function)? Monotonic Queue! (The monotonic queue stores indices, not specific data.)
Some Small Details#
Offset#
dp[0] dp[v] dp[2v] dp[3v] ... dp[kv]
dp[1] dp[1+v] dp[1+2v] dp[1+3v] ... dp[1+kv]
dp[2] dp[2+v] dp[2+2v] dp[2+3v] ... dp[2+kv]
...
dp[v-1] dp[v-1+v] dp[v-1+2v] dp[v-1+3v] ... dp[v-1+kv]
Let’s take any starting point of the class (referring to dp[0] to dp[v-1]), and set it as $u$ to look at.
dp[u] = dp[u]
dp[u+v] = max(dp[u] + w, dp[u + v])
dp[u+2v] = max(dp[u] + 2w, dp[u + v] + w, dp[u + 2v])
dp[u+3v] = max(dp[u] + 3w, dp[u + v] + 2w, dp[u + 2v] + w, dp[u + 3v])
...
We find that the values added to dp[] in the max function are different, which is quite inconvenient for updating values in the monotonic queue. So we need to do something clever.
dp[u] = dp[u]
dp[u+v] = max(dp[u], dp[u + v] - w) + w
dp[u+2v] = max(dp[u], dp[u + v] - w, dp[u + 2v] - 2w) + 2w
dp[u+3v] = max(dp[u], dp[u + v] - w, dp[u + 2v] - 2w, dp[u + 3v] -3w) + 3w
...
This way, it looks much more comfortable. Each time we enter the monotonic queue, the value is dp[u + kv] - k*w. When updating the current state, don’t forget to add $a\times w$, where $a$ is $（\text{current state - front state})/v$. Why? Let’s take a look.
We know that the front of the queue dp[u + kv] - k*w has $u$ fixed, and $k$ can take any maximum value within the range. When we update the current state (let's say $u + bv$), the offset outside the max function is positive $bv$, while dp[u + kv] is less by $kw$. Thus, in total, we need to add $(b - k) * w$, which is $（\text{current state - front state})/v * w$.
The number of elements in the queue should be $s+1$, representing $0$ to $s$ previous states.
Rolling Array#
As mentioned earlier, the multiple knapsack requires either reverse enumeration or rolling (backup) data, but considering the properties of the monotonic queue, reverse enumeration is not very convenient, so we adopt the method of backup arrays.
Code#
Now we can happily write the code.
#include<iostream>
#include<cstdlib>
#include<cstring>// Required header for memcpy function
using namespace std;
const int N = 200005;
int q[N];// Monotonic queue, actually stores indices
int g[N];// Backup array
int dp[N];
int main()
{
int n, V;
cin >> n >> V;
for(int i = 1;i <= n;i++)
{
int v, w, s;// Can choose to read and calculate at the same time, so we don't need to create arrays
memcpy(g, dp, sizeof dp);// Backup
cin >> v >> w >> s;
for(int u = 0;u < v;u++)// Enumerate the starting point of each class
{
int hh = 0, tt = -1;// Because the monotonic queue maintains a closed interval, if tt equals 1
// It indicates that there is one element qwq
for(int k = u;k <= V;k += v)// Here, the meaning of k differs from above
{
while(hh <= tt && (k - q[hh]) > s * v) hh++; // Pop out if exceeding s
while(hh <= tt && g[q[tt]] - (q[tt] - u) / v * w < (g[k] - (k - u) / v * w)) tt--;// Maintain monotonicity of the queue
q[++tt] = k;// Put the current state in
if(hh <= tt) dp[k] = max(dp[k], g[q[hh]] + (k - q[hh]) / v * w);// Update the current state
}
}
}
cout << dp[V];
return 0;
}
The end, let the flowers bloom!