First, let me share some tragic experiences while writing this article. Initially, I was using Gridea and editing with its editor. However, during one writing session, my computer unexpectedly blue-screened. After rebooting, I found that the md file could not be opened, and upon checking, the binary was all zeros emmm (I had saved during the writing process!!). Since then, I switched to hexo emmm.
Symbol Definition#
First, let's define some symbols. The number of items is $n$, the volume of the $i$-th item is $v_i$, its 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. ($x$ is the maximum number of items that can fit in the knapsack with volume $j$, i.e., $j / v$). Let the volume of the $i$-th item be $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)
Now let's think about the relationship between $dp[i,j]$ and $dp[i,j - v]$ in the complete knapsack case.
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, it is not $x$, it is $j-v-(\frac{j-v}{v})v$, and the final state calculated reduces by $xv$.)
Then we are amazed to find that the terms after the first $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. Let's try to see if we can optimize using the complete knapsack method, assuming that $j$ is much larger than $sv$.
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, right? (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 $max()$, thus we cannot use the previous method for optimization.)
New Development#
(From this section onwards, we will use a one-dimensional array for convenience. If you haven't learned about one-dimensional arrays, please study it.)
We still consider $dp[j]$ and list the items.
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 item and easily deduce 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, which are congruent to $0$ modulo $v$, congruent to $1$ modulo $v$, $\cdots$, congruent to $v-2$ modulo $v$, and 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 satisfying the last state, 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)/*to 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 four nested loops. You might ask, why are there more and more loops, but let's analyze its time complexity. The first loop is $O(n)$, the second loop enumerates the starting value of each class, the total time complexity of the second and third loops is $O(V)$, and the last loop enumerates the states that can transition to the current state, which is $O(s)$, so the total time complexity is still $O(nVs)$.
If the above was not clear, let's take an example.
Let's say 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 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 class (remainder 0)
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 by at most 3 times $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])
// Similarly for the second item
Monotonic Queue Optimization#
You might ask, after saying so much, why haven't we seen the shadow of the monotonic queue? Actually, everything we've done before is to lay the groundwork for it, and here it comes.
First, we want to optimize the fourth layer of the above code: updating the current state from the previous $s$ states (which are not adjacent, each differing by $v$).
Let's think about this intuitively, looking at the diagram! (qwq
In this diagram, $j$ represents the current state, and $j'$ represents the next state, i.e., $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 (shown in the leftmost part of the diagram), while also gaining the rightmost state that can update it (shown in the rightmost part of the diagram).
Let's think about what data structure can insert at both ends (since $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 decrease by one), and support a small ($O(1)$) time complexity to find the maximum value (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 a random starting point of the class (referring to dp[0]
to dp[v-1]
), and set it as $u$.
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 the max function after dp[]
are different, which is quite inconvenient for updating values in the monotonic queue. So let's 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 dp[u + kv] - k*w
is fixed at $u$, 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]
has less reduced $kw$, and finally, the max function outside needs to add back $bv$, so the overall result is $(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, for the multiple knapsack, we either enumerate in reverse or roll (backup) data, but considering the properties of the monotonic queue, enumerating in reverse is not very convenient, so we use the backup array method.
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 simultaneously, 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;// Since the monotonic queue maintains a closed interval, if tt equals 1
// it means there is one element qwq
for(int k = u;k <= V;k += v)// Here, k has a different meaning than 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 in 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;
}
End of the article!