Let's calculate how many operations we need. The number of pairs of operations is $$$\left\lfloor\frac{a}{x + y}\right\rfloor$$$. Additionally, we add $$$+1$$$ if $$$(a \bmod (x + y)) \ge x$$$, because in this case, the first player will be able to make one more move. Note that the first part is always even, which means the answer depends only on whether $$$(a \bmod (x + y)) \ge x$$$, and in this case, the answer is YES.

#include <bits/stdc++.h>

using namespace std;

void solve() {
    int x, y, a;
    cin >> x >> y >> a;
    if (a % (x + y) < x) {
        cout << "NO\n";
    } else {
        cout << "YES\n";
    }
}

int main() {
    ios_base::sync_with_stdio(0);
    cin.tie(0);

    int t;
    cin >> t;
    while (t--)
        solve();
}

We can transform the problem into: select a row / column and replace the leftmost / topmost $$$0$$$ with $$$1$$$. Then the solution is clear: if, for each $$$1$$$ in the grid, there does not exist any $$$0$$$ on its top, or does not exist any $$$0$$$ on its left, then the answer is YES; otherwise NO.

But how to prove that if the condition is satisfied, we can always construct a valid series of ball-pushing operations? We can consider the problem in reverse order: given the final state, if, for an $$$1$$$ in the grid, there does not exist any $$$0$$$ on its top, or does not exist any $$$0$$$ on its left, then we can turn it into $$$0$$$, and our goal is to turn all $$$1$$$ s into $$$0$$$ s. We can always reach the goal if we choose $$$1$$$ s in descending order of $$$i+j$$$, where $$$(i,j)$$$ is their coordinates in the grid.

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll t,n,m,a[59][59],vis[59][59];
char s[59];
inline ll read(){
	ll s = 0,w = 1;
	char ch = getchar();
	while (ch > '9' || ch < '0'){ if (ch == '-') w = -1; ch = getchar();}
	while (ch <= '9' && ch >= '0') s = (s << 1) + (s << 3) + (ch ^ 48),ch = getchar();
	return s * w;
}
int main(){
	t = read();
	while (t --){
		n = read(),m = read();
		for (ll i = 1;i <= n;i += 1){
			scanf("%s",s + 1);
			for (ll j = 1;j <= m;j += 1){
				a[i][j] = s[j] - '0';
				vis[i][j] = 0;
			}
		}
		for (ll i = 1;i <= n;i += 1){
			for (ll j = 1;j <= m;j += 1){
				if (!a[i][j]) break;
				vis[i][j] = 1;
			}
		}
		for (ll j = 1;j <= m;j += 1){
			for (ll i = 1;i <= n;i += 1){
				if (!a[i][j]) break;
				vis[i][j] = 1;
			}
		}
		bool fl = 1;
		for (ll i = 1;i <= n && fl;i += 1){
			for (ll j = 1;j <= m;j += 1){
				if (a[i][j] && !vis[i][j]){
					fl = 0;
					break;
				}
			}
		}
		puts(fl ? "YES" : "NO");
	}
	return 0;
}

Original idea from myee has $$$4$$$ cases, which is a medium difficulty problem using deletable heap. You may find the original statement in zh-cn on GitHub later.

The Codeforces version by FairyWinx has only $$$2$$$ cases, in order to meet the difficulty as Div2C.

The En/Ru statement changes for times after checking. Sorry for the inconvenience.

Firstly we can get the distance of each cell by BFS or Math Way: $$$d(x,y)=x+y+2[x\bmod3=y\bmod3=2]$$$.

The proof is trivial.

Then we can sort out the first $$$4n$$$ table cell with the lowest distance in $$$O(n)$$$ time: by BFS or enumerating possible $$$d(x,y)=k$$$.

Suppose we record whether a table cell is occupied or not by a bool array. Then we can check whether a cell can be allowed to go by someone $$$o=0/1$$$ in $$$O(1)$$$ time.

Suppose the shortest table cell we can choose is $$$x_0$$$ -th and $$$x_1$$$ -th ones, then $$$x_1\le x_0\le4n$$$. We can keep trying the smallest $$$x_o$$$ until it's okay.

The time complexity is $$$O(n)$$$.

#include <bits/stdc++.h>

using uint = unsigned;
using bol = bool;

const uint T=2000;
const uint R=200000;

uint Ord[R+5];
uint Nxt(uint p,uint o)
{
    uint x=p/T,y=p%T;
    if(o&1)y=y/3*3+3-y%3;
    if(o&2)x=x/3*3+3-x%3;
    return x*T+y;
}

bol G[2][T*T];uint Ans[2];

int main()
{
#ifdef MYEE
    freopen("QAQ.in","r",stdin);
    freopen("QAQ.out","w",stdout);
#endif
    for(uint i=2,tp=0;i<T&&tp<R;i++)if(i%3==2)
    {
        for(uint x=1;x<i&&tp<R;x+=3)Ord[tp++]=x*T+i-x;
    }
    else if(!(i%3))
    {
        for(uint x=1;x<=i-2&&tp<R;x++)
            if(x%3==1)Ord[tp++]=x*T+i-x;
            else if(x%3==2)Ord[tp++]=x*T+i-x-2,Ord[tp++]=x*T+i-x;
        Ord[tp++]=(i-1)*T+1;
    }
    uint t;scanf("%u",&t);
    while(t--)
    {
        uint q;scanf("%u",&q);
        for(uint i=0;i<q*4;i++)G[0][Ord[i]]=G[1][Ord[i]]=0;
        Ans[0]=Ans[1]=0;
        while(q--)
        {
            uint t;scanf("%u",&t);while(G[t][Ord[Ans[t]]])Ans[t]++;
            t=Ord[Ans[t]];printf("%u %u\n",t/T,t%T);
            G[1][t]=1;for(uint o=0;o<4;o++)G[0][Nxt(t,o)]=1;
        }
    }
    return 0;
}

Lemma (Bertrand's postulate): For each positive integer $$$x$$$, there is a prime $$$p$$$ inside the interval $$$[x,2x]$$$.

Find a prime number $$$p$$$ between $$$\lfloor \frac{n}{3} \rfloor$$$ and $$$\lceil \frac{2n}{3} \rceil$$$, and construct the permutation in an alternating way: $$$p, p-1, p+1, p-2, p+2, \dots$$$. Place the remaining numbers arbitrarily.

In this construction, we have $$$c_1 = c_3 = c_5 = \dots = c_{2 \lfloor \frac{n}{3} \rfloor - 1} = p$$$, which meets the requirement with at least $$$\lfloor \frac{n}{3} \rfloor$$$ numbers being prime.

#include<bits/stdc++.h>
using namespace std;

bool isPrime(int x){
	if(x <= 1)
		return 0;
	for(int i = 2 ; i * i <= x ; ++i)
		if(x % i == 0)
			return 0;
	return 1;
}

vector<int> generateSol(int n, int p){
	vector<int> ans;
	ans.push_back(p);
	for(int i = 1 ; i <= n ; ++i){
		if(p - i > 0)
			ans.push_back(p - i);
		if(p + i <= n)
			ans.push_back(p + i);
	}
	return ans;
}

int main(){
	int t;
	cin >> t;
	for(int i = 1 ; i <= t ; ++i){
		int n;
		cin >> n;
		vector < int > ans;
		for(int x = 0 ; ; ++x){
			if(isPrime(n / 2 - x)){
				ans = generateSol(n , n / 2 - x);
				break;
			}
			if(isPrime(n / 2 + x)){
				ans = generateSol(n , n / 2 + x);
				break;
			}
		}
		for(int i = 0 ; i < n ; ++i)
			cout << ans[i] << " \n"[i == n - 1];
	}
	return 0;
}

For each $$$ a_i $$$, calculate the minimum number of cyclic right shifts required to reduce it to zero, denoted as $$$ d_i $$$. The answer is the maximum $$$ d_i $$$ across all $$$ i $$$.

For cyclic problems, we can consider breaking the cycle into a chain (i.e., let $$$ a_{i+n} = a_i $$$, $$$ b_{i+n} = b_i $$$).

Consider constructing a parenthesis sequence $$$ s $$$ by concatenating $$$ a_0 $$$ left parentheses, $$$ b_0 $$$ right parentheses, $$$ a_1 $$$ left parentheses, $$$ b_1 $$$ right parentheses, ... ,$$$ a_{2n-1} $$$ left parentheses and $$$ b_{2n-1} $$$ right parentheses. The operation in the original problem corresponds to the process of matching parentheses in $$$ s $$$:

• "Set the smaller of $$$(a_i, b_i)$$$ to zero and the larger to their difference" $$$\Rightarrow$$$ "Match the parentheses corresponding to $$$(a_i, b_i)$$$".
• "Cyclically right-shift $$$ a $$$" $$$\Rightarrow$$$ "Unmatched left parentheses continue to find matching right parentheses".

The first time $$$ a_i $$$ is reduced to zero depends on which $$$ b_j $$$ ’s right parenthesis matches the leftmost parenthesis corresponding to $$$ a_i $$$.

We can use prefix sums to find matchings of the parenthesis: Let $$$ c_i = \sum\limits_{j=0}^i (a_j - b_j) $$$. Define $$$ p_i $$$ as the smallest index satisfying $$$ p_i \gt i $$$ and $$$ c_{p_i} \leq c_i $$$, then $$$ d_i = p_i - i $$$. We can compute all $$$ p_i $$$ with a monotonic stack.

It should be pointed out that matching parentheses is just an intuitive way of understanding the solution. You can, of course, also draw the above conclusion by direct observation.

The time complexity is $$$ O(\sum n) $$$.

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll t,n,k,a[400009],stk[400009];
inline ll read(){
	ll s = 0,w = 1;
	char ch = getchar();
	while (ch > '9' || ch < '0'){ if (ch == '-') w = -1; ch = getchar();}
	while (ch <= '9' && ch >= '0') s = (s << 1) + (s << 3) + (ch ^ 48),ch = getchar();
	return s * w;
}
int main(){
	t = read();
	while (t --){
		n = read(),k = read();
		for (ll i = 1;i <= n;i += 1) a[i] = read();
		for (ll i = 1;i <= n;i += 1) a[i] = a[i + n] = a[i] - read();
		for (ll i = 1;i <= 2 * n;i += 1) a[i] += a[i - 1];
		ll tp = 0,ans = 0;
		for (ll i = 2 * n;i >= 1;i -= 1){
			while (tp && a[stk[tp]] > a[i]) tp -= 1;
			if (i <= n) ans = max(ans,stk[tp] - i);
			stk[++ tp] = i;
		}
		printf("%lld\n",ans);
	}
	return 0;
}

We can use binary search on the answer $$$ x $$$, and calculate the minimum number of changes needed to make $$$a$$$ zero within $$$ x $$$ rounds.

However, direct computation on the cyclic chain of $$$ 2n $$$ elements is troublesome (due to overlapping contributions). Instead, we can cyclically shift $$$ a $$$ and $$$ b $$$ such that $$$ c_{n-1} = \min\limits_{i=0}^{n-1} c_i$$$, which ensures that no $$$ a_{n-1} \gt 0 $$$ will be cyclically shifted to $$$a_0$$$, enabling us to only consider the first $$$n$$$ elements of the chain.

For the adjusted $$$ a $$$ and $$$ b $$$, we can greedily process $$$ i $$$ from $$$ n - x - 1 $$$ to $$$ -1 $$$: If $$$ \min\limits_{j=i+1}^{i+x}c_j \gt c_i $$$, then we should apply $$$ \min\limits_{j=i+1}^{i+x}c_j - c_i $$$ changes to $$$ a_{i+1} $$$ (assume that $$$c_{-1}=0$$$). We can also use a monotonic stack to maintain the whole process.

Note that if $$$\sum\limits_{i=0}^{n-1} a_i=k$$$, we should output $$$0$$$.

The time complexity is $$$ O(\sum n\log k) $$$.

Bonus: Find an $$$O(\sum n)$$$ solution.

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll t,n,k,pos,a[400009],stk[400009];
inline ll read(){
	ll s = 0,w = 1;
	char ch = getchar();
	while (ch > '9' || ch < '0'){ if (ch == '-') w = -1; ch = getchar();}
	while (ch <= '9' && ch >= '0') s = (s << 1) + (s << 3) + (ch ^ 48),ch = getchar();
	return s * w;
}
bool judge(ll x){
	ll l = 1,r = 0,res = 0;
	for (ll i = pos;i >= pos - n;i -= 1){
		if (l <= r && stk[l] - i > x) l += 1;
		if (pos - i >= x) res += max(0ll,a[stk[l]] - a[i]);
		if (res > k) return 0;
		while (l <= r && a[stk[r]] > a[i]) r -= 1;
		stk[++ r] = i;
	}
	return 1;
}
int main(){
	t = read();
	while (t --){
		n = read(),k = read();
		ll suma = 0;
		for (ll i = 1;i <= n;i += 1) a[i] = read(),suma += a[i];
		for (ll i = 1;i <= n;i += 1) a[i] -= read(),a[i + n] = a[i];
		if (suma == k){puts("0"); continue;}
		for (ll i = 1;i <= 2 * n;i += 1) a[i] += a[i - 1];
		pos = n;
		for (ll i = n + 1;i <= 2 * n;i += 1) if (a[i] < a[pos]) pos = i;
		ll l = 1,r = n,ans = n;
		while (l <= r){
			ll mid = (l + r) >> 1;
			if (judge(mid)) ans = mid,r = mid - 1;
			else l = mid + 1;
		}
		printf("%lld\n",ans);
	}
	return 0;
}

Background

This problem is inspired by a minigame from a recently released party video game, but the story was modified to comply with the title. We tried to isolate the artistic part and the theoretical part, and I apologize if you still feel like working on reading comprehension exercises.

Solution

Noticing that each member's strategy is not maximizing the total outcome (which is the expected number of successful matches between keys and locks), but only maximizing the outcome in her current turn, and equally likely to adopt each of optimal (pure) strategies, it will be ideal if there is some simple form of the strategy that only depends on $$$L$$$ and $$$K$$$.

Analyzing the Strategy

At the very beginning, when there is no information about the relationship between keys and locks, the first member in her first turn must randomly choose a key and a lock.

What about the second member, if the first one is not lucky enough? Since she will not choose the same pair of keys and locks, there are three different types of (pure) strategies:

• Choose a different key to open the same lock;

• Choose the same key to open a different lock;

• Choose a different key and a different lock.

For the first type of strategy, since there are $$$(L+K-1)$$$ keys remaining (after the very first attempt), and among them is only one correct key, the probability of opening the same lock is $$$1/(L+K-1)$$$.

For the second type, it can be proven that the probability is also $$$1/(L+K-1)$$$. Another way to think of this is, add $$$K$$$ imaginary locks and match each of them with a counterfeit key. Then, it is apparent that a new lock matches the selected key with the given probability.

For the third type of strategy, its probability can be computed from either the first type or the second type, as $$$(L+K-1)$$$ different outcomes are equally likely to happen:

\begin{equation} \frac{1-1/(L+K-1)}{L+K-1}=\frac{L+K-2}{(L+K-1)^2}<\frac{1}{L+K-1}. \end{equation}

In conclusion, the second member will either choose the same key or the same lock. Because there are $$$(L - 1)$$$ pairs of the same key and other locks, and $$$(L + K - 1)$$$ pairs of other keys and the same lock, the second member will use the same key with probability $$$(L - 1)/(2L + K - 2)$$$, or check the same lock with probability $$$(L + K - 1)/(2L + K - 2)$$$.

As for all the following attempts, a similar Proof shows that, everyone will imitate the strategy adopted by the second member. That is to say, if the second member decides to continue with the same key, everyone will keep using this same key, until a lock is opened or all the locks are checked; if the second member picks a different key for the same lock, everyone will challenge this same lock, until it is eventually settled (and obviously it will).

Furthermore, whenever a key is identified as counterfeit, or a lock is opened, the information 'collapses', which turns the original problem into a subproblem with either one less counterfeit key (corresponding to $$$K' = K - 1$$$), or both one less valid key and one less lock ($$$L' = L - 1$$$), respectively.

To sum up, the process of the game is:

• Whenever a lock is opened, the original problem is reduced to a subproblem with $$$L' = L - 1$$$;
• Whenever a key is proven counterfeit, the original problem is reduced to a subproblem where $$$K' = K - 1$$$; -The first member randomly chooses a key and a lock to check if they match;
• If the first member fails, the second member adopts the same key with probability $$$(L - 1)/(2L + K - 2)$$$, or investigates the same lock with probability $$$(L + K - 1)/(2L + K - 2)$$$;
• If the second member fails as well, all the following attempts will replicate the choice of whether to employ the same key or the same lock.

Design the Algorithm with Strategies Solved

Knowing the concrete strategies, the task now becomes how to compute the expected values efficiently. As you might have imagined, it can be achieved using dynamic programming due to the existence of overlapping subproblems. But let's examine another approach which might be easier to understand: compute $$$p_i (l, k)$$$, the total probability that the $$$i$$$-th member begins a subproblem with $$$l$$$ locks and $$$k$$$ counterfeit keys, and accumulate $$$e_i$$$'s during the computation.

Transitions for $$$p$$$ itself include (and those from $$$p$$$ to $$$e$$$ will be similar to):

1. The first member manages to open the lock without any prior information;
2. The second member makes a successful guess;
3. Any of the succeeding attempts succeeds;
4. Everyone uses the same key but it fails on all the locks.

For each $$$p_i(l, k)$$$, all the transition can be easily processed in $$$\Theta(1)$$$ time, except the third one. A verbatim implementation requires $$$O(L + K)$$$ time per state, which is apparently not efficient enough. Since $$$N \ll L$$$, a possible work-around is to compute the transition coefficient for each member, instead of each attempt. Intuitively, all attempts should share the same probability of success, regardless of the order. A simplified model is to draw all balls sequentially without replacement from a box of $$$1$$$ red ball and $$$(x-1)$$$ white balls, in which the red ball appears at the $$$i$$$-th attempt with an equal probability of $$$1/x$$$. Of course, the observation for our problem can also be proven mathematically, but let's omit it for simplicity.

This important observation reduces the time complexity of the third type of transition to $$$O(N)$$$, as we can first compute the maximum number of attempts for each member, and then multiply it with the common probability of a single attempt. However, this is still not enough for passing the problem (and we adjusted the constraints to prevent so).

The final trick to pass all testcases requires revisiting the transition coefficients we just computed. Notice that, if all the members have $$$a$$$ attempts in total (but for the third type of transition only, which excludes the first two), then:

• If $$$a$$$ is a multiple of $$$N$$$, then every member can have at most $$$a/N$$$ attempts;

• Otherwise, the first $$$a \bmod N$$$ members (right after the second member, who decides to advance with the same key or the same lock) can have at most $$$\lceil a/N\rceil$$$ attempts, while all the others can have at most $$$\lfloor a/N\rfloor$$$ attempts.

This implies that the transition coefficients consists of at most $$$3$$$ maximal contiguous subsequences of the same coefficients (in fact, if the coefficients are treated as a circular sequence, then there will be exactly $$$1$$$ or $$$2$$$ for the two cases respectively). Using techniques such as prefix sums helps further decrease the complexity of the computation to $$$\Theta(1)$$$ per state.

The total time complexity is $$$O(NLK)$$$.

#include <cstdio>
#include <cstring>
using namespace std;
#define MOD 1000000007
#define MAXN 108
#define MAXL 5004
#define MAXK 27
int e[MAXN], p[2][MAXK][MAXN], inv[2 * MAXL + MAXK];  
inline void _update(int * const a, const int lbd, const int rbd, const int base, const int diff) {
	(a[0] += base) %= MOD;
	(a[lbd] += diff) %= MOD;
	(a[rbd] += MOD - diff) %= MOD;
	return ;
}
#define Hill(_a, _lbd, _rbd, _base, _diff) _update(_a, _lbd, _rbd, _base, _diff)
#define Valley(_a, _lbd, _rbd, _base, _diff) _update(_a, _rbd, _lbd, ((_base) + (_diff)) % MOD, MOD - (_diff))
#define Update(_a, _lbd, _rbd, _base, _diff) { if (_lbd < _rbd) Hill(_a, _lbd, _rbd, _base, _diff); else Valley(_a, _lbd, _rbd, _base, _diff); }
int like() {
	int n, l, k, i, j, a, b, dh, now, nxt, totk, ntry, prob, pd, pb, lbd, rbd;
	scanf("%d %d %d", &n, &l, &k);
	if (n == 1) {
		printf("%d\n", l);
		return 0;
	}
	inv[1] = 1;
	j = 2 * l + k;
	for (i = 2; i <= j; ++i) inv[i] = 1ll * (MOD - MOD / i) * inv[MOD % i] % MOD;
	
	memset(p, 0, sizeof(p));
	memset(e, 0, sizeof(e));
	p[0][0][0] = 1;
	p[0][0][1] = -1;
	now = 0; nxt = 1;
	for (i = 0; i < l; ++i) {
		for (j = 0; j <= k; ++j) {
			totk = l + k - i - j;
			for (a = 1; a < n; ++a) (p[now][j][a] += p[now][j][a - 1]) %= MOD;
			for (a = 0; a < n; ++a) {
				// Match at the first try
				prob = 1ll * p[now][j][a] * inv[totk] % MOD;
				(e[a] += prob) %= MOD;
				if (a != n - 1) (e[a + 1] += MOD - prob) %= MOD;
				(p[nxt][j][(a + 1) % n] += prob) %= MOD;
				if (a != n - 2) (p[nxt][j][(a + 2) % n] += MOD - prob) %= MOD;

				prob = 1ll * prob * inv[totk + l - i - 2] % MOD;
				// Same lock
				pd = 1ll * prob * (ntry = totk - 1) % MOD;
				pb = 1ll * (ntry / n) * pd % MOD;
				if (dh = ntry % n) {
					lbd = (a + 1) % n;
					rbd = (a + dh) % n + 1; 
					Update(e, lbd, rbd, pb, pd);
					(++lbd) %= n;
					++(rbd  %= n); 
					Update(p[nxt][j], lbd, rbd, pb, pd);
				}
				else {
					(e[0] += pb) %= MOD;
					(p[nxt][j][0] += pb) %= MOD;
				}

				// Same key
				pd = 1ll * prob * (ntry = l - i - 1) % MOD;
				if (j < k) {
					pb = 1ll * pd * (k - j) % MOD; 
					b = (a + l - i) % n;
					(p[now][j + 1][b] += pb) %= MOD;
					if (b != n - 1) (p[now][j + 1][b + 1] += MOD - pb) %= MOD;
				}
				pb = 1ll * (ntry / n) * pd % MOD;
				if (dh = ntry % n) {
					lbd = (a + 1) % n;
					rbd = (a + dh) % n + 1; 
					Update(e, lbd, rbd, pb, pd);
					(++lbd) %= n;
					++(rbd  %= n); 
					Update(p[nxt][j], lbd, rbd, pb, pd);
				}
				else {
					(e[0] += pb) %= MOD;
					(p[nxt][j][0] += pb) %= MOD;
				}
			}
		}

		now ^= 1;
		nxt ^= 1;
		memset(p[nxt], 0, sizeof(int) * MAXK * MAXN);
	}
	for (i = 1; i < n; ++i) (e[i] += e[i - 1]) %= MOD;

	for (i = 0; i < n; ++i) printf("%d%c", e[i], " \n"[i == n - 1]);
	return 0;
}
int main() {
	int t;
	scanf("%d", &t);
	while (t--) {
		like();
	}
	return 0;
}

Each operation transform three adjacent characters into one. If the string starts with 11, the answer should be yes since we can first transform the remaining part into a single 0 or 1 and then the value of 11? is always $$$1$$$.

If the string ends with 1, only 101 has no solution, otherwise it starts with 11 or 0, or we can transform 10? in the front of the string into a 0. Then transform the rest part of the string into a single character to make the last operation 0?(anything):1 makes 1.

If an operation transformed a string ending with 0 into one ending with 1, it must be 1?1:0 makes 1.

Consider the case where the strings contains 11 as a substring. After transforming the remaining part, ?11?? or ??11? will be derived. The cases are 01100, 00110, 01110, 10110. All other cases start with 11 or end with 1.

All cases of ??11? have a solution:

(0?0:1)?1:0 makes 1.

(0?1:1)?1:0 makes 1.

1?(0?1:1):0 makes 1.

However, 01100 has no solution.

Since 0?0:1 equals 1, before the last 1, any two adjacent 0's with no 1 in between can be eliminated. Thus, if there are two 1's with an even number of 0's in between, and there are an even number of characters in front the former 1, i.e., (..)*1(00)*1(..)*0 in regular expression, the answer should be yes, since after eliminating the 0's in between, the 1's become adjacent. The parity of the characters before the former 1 guarantees that it will not become 01100.

Another case is that the string ends with 0 and every pair of adjacent 1's has an odd number of 0's in between. Since 0(1?0:1)0 makes 000, 1(0?1:0)1 makes 101, 10(1?0:0)01 makes 10001, there is no way to change the parity of the consecutive 0's to make 11 to change the ending 0 into 1.

In the remaining case, there are even 0's between some pair of adjacent 1's, but the number of characters before the first 1 is odd, there must be only one pair. Suppose the string is 1-indexed. It is because the indices of such pair of 1's must be [odd, even]; while the indices in the previous case must be [even, odd]. Thus, if there are two such pairs of this case [odd, even], in the parity sequence [odd, even, ..., odd, even] there must be an [even, odd] as a substring.

Therefore, any other pair of adjacent 1's must have odd 0's in between, i.e., $$$0(00)^*\color{blue}{(10(00)^*)^*}\color{red}1\color{black}(00)^*\color{red}1\color{blue}{((00)^*01)^*}\color{black}(00)^+$$$ in regular expression. Any operation remains the string in this case, except that the string is $$$0(00)^*\color{blue}{(10(00)^*)^*}\color{red}1\color{red}1\color{black}(00)^+$$$ and the operation is $$$\color{red}1\color{black}?\color{red}1\color{black}:0$$$ makes 1. Since it must ends with 00, it becomes a string ending with 0 and every pair of adjacent 1's has odd 0's in between, which is a case mentioned before that has no solution.

In conclusion, it has a solution if and only of at least one of the two following conditions is met:

• It ends with 1 but it is not 101.
• It has two adjacent 1's such that there are even 0's in between and the number of characters before the former 1 is even.

#include<bits/stdc++.h>
char s[333333];
char val(char x,char y,char z){
	return x=='1'?y:z;
}
std::pair<std::string, char> fold(int l,int r){
	std::string str;
	char v=s[r];
	for(int i=l;i<r;i+=2){
		str+=s[i];
		str+='?';
		str+=s[i+1];
		str+=':';
	}
	str+=s[r];
	for(int i=r-2;i>=l;i-=2){
		v=val(s[i],s[i+1],v);
	}
	return std::make_pair(str,v);
}
void solve(){
	int n;
	scanf("%d",&n);
	n=2*n+1;
	scanf("%s",s+1);
	int last=0;
	for(int i=1;i<=n;++i){
		if(s[i]=='1'&&last>0){
			if(i%2==0&&last&1){
				puts("Yes");
				std::string s4=fold(last+1,i).first;
				auto tmp=fold(i+1,n);
				std::string s5=tmp.first;
				char v5=tmp.second;
				if(last==1){
					printf("1?(%s):(%s)\n",s4.c_str(),s5.c_str());
					return;
				}
				tmp=fold(1,last-2);
				std::string s1=tmp.first;
				char v1=tmp.second;
				if(v1=='1'){
					printf("(%s)?(%c?1:(%s)):(%s)\n",s1.c_str(),s[last-1],s4.c_str(),s5.c_str());
				}
				else{
					printf("((%s)?%c:1)?(%s):(%s)\n",s1.c_str(),s[last-1],s4.c_str(),s5.c_str());
				}
				return;
			}
		}
		if(s[i]=='1'){
			last=i;
		}
	}
	if(s[n]=='1'&&std::string(s+1)!="101"){
		puts("Yes");
		auto left=s[1]=='0'?fold(1,1):fold(1,3);
		auto mid=s[1]=='0'?fold(2,n-1):fold(4,n-1);
		printf("(%s)?(%s):1\n",left.first.c_str(),mid.first.c_str());
		return;
	}
	puts("No");
	return;
}
int main(){
	int T;
	scanf("%d",&T);
	while(T--){
		solve();
	}
	return 0;
}

Firstly, consider a classical dp to calculate $$$dp_{u,k}$$$ as the number of delete $$$k$$$ vertices in subtree $$$u$$$ without considering deleting vertex $$$n-i+1$$$ everyday. This is trivial.

Then, due to the existence of an “automatic deletion” process, we cannot view the selected points within the subtree of $$$u$$$ as operations with only relative order, because they may be illegal. However, we consider the final operation sequence, putting the vertex chosen on the $$$i$$$ th day in the position $$$n-i+1$$$. That is, if we firstly choose $$$4$$$, then $$$3$$$, then $$$1$$$, it will be $$$0431$$$. Then, one vertex $$$u$$$ must be put in position $$$u$$$ or later.

Then, we find that for the subtree of vertex $$$u$$$ with DFS order interval $$$[l,r]$$$, all vertices can't be put in position before $$$l$$$, and it is always legal to put in position after $$$r$$$ if parent-child relations are satisfied.

So we consider a dp $$$f_{u,i,j,k}$$$ which represents that the subtree of $$$u$$$ leaves position $$$i$$$ for ancestors (which means nothing can be placed at position $$$i$$$ or before, and if $$$i=0$$$, then nothing is left for ancestors), there are $$$j$$$ empty spaces (not including the one left for ancestors), and $$$k$$$ vertices need to be put into positions after the subtree of $$$u$$$. In the $$$i=0$$$ case, we merge subtrees, and the subtree with smaller DFS order can put some backward operations into the space of the subtree with bigger DFS order.

Specifically, consider the sons of $$$u$$$ enumerated in reverse order of DFS, and the son is $$$v$$$. Then enumerate $$$w$$$ as how many backward operations of $$$v$$$ will fill in space of $$$u$$$. The transition is $$$f_{u,0,j+j1-w,k+k1-w}\leftarrow f_{u,0,j+j1-w,k+k1-w}+f_{u,0,j,k}\times f_{v,0,j1,k1}\times \binom{j}{w}\times\binom{k+k1-w}{k}$$$. That is, after putting $$$w$$$ operations from backward operations of $$$v$$$ into spaces of $$$u$$$, we merge the spaces of $$$u$$$ and $$$v$$$, then the backward operations of $$$u$$$ and $$$v$$$.

For $$$f_{v,i,j,k}(i\neq 0)$$$, it won't interfere with the transition with its following siblings. And for $$$f_{u,i,j,k}$$$,any operations from $$$v$$$ can't put into the subtree of $$$u$$$ and $$$v$$$, so we only need to merge backward operations of $$$u$$$ and $$$v$$$.

After finishing the transition of subtrees, we consider where to put $$$u$$$.

• If we don't put $$$u$$$, there will be one more space, and we can consider whether to leave it for ancestors if $$$i=0$$$.That is $$$f_{u,i,j,k}\rightarrow f_{u,i,j+1,k},f_{u,0,j,k}\rightarrow f_{u,u,j,k}$$$.
• We put $$$u$$$ into the position left for ancestors, that is $$$f_{u,i,j,k}\rightarrow f_{u,0,j+1,k}$$$. And we can leave another position for ancestors of $$$u$$$, that is $$$f_{u,i,j,k}\rightarrow f_{u,i1,j+1,k}(i1 \lt i)$$$.
• We put $$$u$$$ into position $$$u$$$, in that case, we can't leave position for ancestors of $$$u$$$, that is $$$f_{u,0,j,k}\rightarrow f_{u,0,j,k}$$$.
• We put $$$u$$$ backward, then we can't put anything in the subtree of $$$u$$$, but we can leave position for ancestors of $$$u$$$, that is $$$f_{u,i,siz_u-2,k}\rightarrow f_{u,i,siz_u-1,k+1},f_{u,0,siz_u-1,k}\rightarrow f_{u,0,siz_u,k+1}$$$.

The answer for $$$i$$$ operations is just we leave position $$$n-i+1$$$ for $$$1$$$, and there are $$$n-i$$$ spaces (therefore the spaces are before $$$n-i+1$$$) and no backward operations.

The time complexity is $$$O(\sum n^5)$$$.

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
int read() {
	int x=0,f=1;
	char c=getchar();
	while(c<'0'||c>'9') {
		if(c=='-')f=-1;
		c=getchar();
	}
	while(c>='0'&&c<='9')x=x*10+c-'0',c=getchar();
	return x*f;
}
namespace tokido_saya {
	const int maxn=81,mod=998244353;
	typedef vector<int>::iterator iter;
	int n,lp[maxn],rp[maxn],siz[maxn],t;
	vector<int> v[maxn];
	ll f[maxn][maxn][maxn][maxn],dp[maxn][maxn],fac[maxn],inv[maxn];
	int zc[maxn][maxn][maxn][maxn];
	ll qpow(ll x,int y)
	{
		ll w=1;
		while(y)
		{
			if(y&1)w=w*x%mod;
			x=x*x%mod,y>>=1;
		}
		return w;
	}
	ll C(int x,int y)
	{
		if(y<0||y>x)return 0;
		return fac[x]*inv[x-y]%mod*inv[y]%mod;
	}
	void dfs1(int u)
	{
		lp[u]=rp[u]=u;
		for(iter it=v[u].begin();it!=v[u].end();it++)dfs1(*it),rp[u]=max(rp[u],rp[*it]);
	}
	void dfs2(int u)
	{
		f[u][0][0][0]=1,dp[u][0]=1;
		for(iter it=v[u].begin();it!=v[u].end();it++)
		{
			int v=*it;
			dfs2(v);
			for(int i=siz[u];i>=0;i--)
				for(int j=siz[v];j;j--)dp[u][i+j]=(dp[u][i+j]+dp[u][i]*dp[v][j]%mod*C(i+j,j))%mod;
			for(int j1=0;j1<=siz[v];j1++)
				for(int k1=0;k1<=min(j1+1,siz[v]);k1++)
				{
					memset(zc[j1][k1],0,sizeof(zc[j1][k1]));
					for(int j=0;j<=siz[u];j++)
						for(int k=0;k<=j+1;k++)
							for(int w=0;w<=min(k1,j);w++)zc[j1][k1][j+j1-w][k+k1-w]=(zc[j1][k1][j+j1-w][k+k1-w]+f[u][0][j][k]*C(j,w)%mod*C(k+k1-w,k))%mod;
				}	
			memset(f[u][0],0,sizeof(f[u][0]));
			for(int j=0;j<=siz[v];j++)
				for(int k=0;k<=j+1;k++)
					for(int j1=0;j1<=siz[u]+siz[v];j1++)
						for(int k1=0;k1<=min(siz[u]+siz[v],j1+1);k1++)
							if(zc[j][k][j1][k1])
							{
								const int w=zc[j][k][j1][k1];
								f[u][0][j1][k1]=(f[u][0][j1][k1]+f[v][0][j][k]*w)%mod;
								for(int i=lp[v];i<=rp[v];i++)f[u][i][j1][k1]=(f[u][i][j1][k1]+f[v][i][j][k]*w)%mod;
							}
			for(int i=rp[v]+1;i<=rp[u];i++)
				for(int j=siz[u]-1;j>=i-rp[v]-1;j--)
					for(int k=j+1;k>=0;k--)
					{
						const int val=f[u][i][j][k];
						f[u][i][j][k]=0;
						for(int p=0;p<=siz[v];p++)
							for(int w=0;w<=min(p,j-(i-rp[v]-1));w++)f[u][i][j+siz[v]-w][k+p-w]=(f[u][i][j+siz[v]-w][k+p-w]+dp[v][p]*val%mod*C(j-(i-rp[v]-1),w)%mod*C(k+p-w,k))%mod;
					}
			siz[u]+=siz[v];
		}
		for(int i=siz[u];i>=0;i--)dp[u][i+1]=(dp[u][i+1]+dp[u][i])%mod;
		if(u!=1)for(int j=siz[u];j>=0;j--)
			for(int k=min(j+1,siz[u]);k>=0;k--)
			{
				ll sum=0;
				for(int i=rp[u];i>=lp[u];i--)
				{
					const int w=f[u][i][j][k];
					f[u][i][j+1][k]=(f[u][i][j+1][k]+w)%mod;
					if(j==siz[u]-1)f[u][i][j+1][k+1]=(f[u][i][j+1][k+1]+w)%mod;
					f[u][i][j][k]=sum,sum=(sum+w)%mod;
				}
				const int w=f[u][0][j][k];
				f[u][lp[u]][j][k]=(f[u][lp[u]][j][k]+w)%mod,f[u][0][j+1][k]=(f[u][0][j+1][k]+w+sum)%mod;
				if(j==siz[u])f[u][0][j+1][k+1]=(f[u][0][j+1][k+1]+w)%mod,f[u][lp[u]][j][k+1]=(f[u][lp[u]][j][k+1]+w)%mod;
			}
		siz[u]++;
	}
	int main() {
		int x,y;
		t=read(),fac[0]=1;
		for(int i=1;i<=80;i++)fac[i]=fac[i-1]*i%mod;
		inv[80]=qpow(fac[80],mod-2);
		for(int i=80;i;i--)inv[i-1]=inv[i]*i%mod;
		while(t--)
		{
			n=read(),memset(f,0,sizeof(f)),memset(dp,0,sizeof(dp)),memset(siz,0,sizeof(siz));
			for(int i=1;i<=n;i++)v[i].clear();
			for(int i=1;i<n;i++)
			{
				x=read(),y=read();
				if(x>y)swap(x,y);
				v[x].push_back(y);
			}
			for(int i=1;i<n;i++)sort(v[i].begin(),v[i].end()),reverse(v[i].begin(),v[i].end());
			dfs1(1),dfs2(1);
			for(int i=n;i>1;i--)printf("%lld ",f[1][i][i-2][0]);
			puts("1");
		}
		return 0;
	}
}
int main() {
	return tokido_saya::main();
}

Do we really have to consider the problem on the whole matrix?

Let $$$r$$$ be the number of rows where the bitwise XOR of all the numbers in it is $$$1$$$, and $$$c$$$ be the number of columns where the bitwise XOR of all the numbers in it is $$$1$$$. Our goal is to make both $$$r$$$ and $$$c$$$ zero.

When we change an element in the matrix, both $$$r$$$ and $$$c$$$ change by exactly one. Thus, it's not hard to see that the answer is $$$\max (r,c)$$$.

The time complexity is $$$O(\sum nm)$$$.

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll t,n,m,r,c;
char s[1009];
bool a[1009][1009];
inline ll read(){
	ll s = 0,w = 1;
	char ch = getchar();
	while (ch > '9' || ch < '0'){ if (ch == '-') w = -1; ch = getchar();}
	while (ch <= '9' && ch >= '0') s = (s << 1) + (s << 3) + (ch ^ 48),ch = getchar();
	return s * w;
}
int main(){
	t = read();
	while (t --){
		n = read(),m = read(),r = c = 0;
		for (ll i = 1;i <= n;i += 1){
			scanf("%s",s + 1);
			for (ll j = 1;j <= m;j += 1) a[i][j] = s[j] - '0';
		}
		for (ll i = 1;i <= n;i += 1){
			bool sum = 0;
			for (ll j = 1;j <= m;j += 1) sum ^= a[i][j];
			if (sum) r += 1;
		}
		for (ll j = 1;j <= m;j += 1){
			bool sum = 0;
			for (ll i = 1;i <= n;i += 1) sum ^= a[i][j];
			if (sum) c += 1;
		}
		printf("%lld\n",max(r,c));
	}
	return 0;
}

Consider $$$x$$$ in binary form.

What are all possible values of $$$x$$$ after $$$n+m$$$ operations?

Let's consider how to find the maximum value first.

Consider $$$x$$$ in binary form. Let the number formed by the last $$$n+m$$$ bits of $$$x$$$ be denoted as $$$r$$$, and the remaining higher bits form the number $$$l$$$ (for instance, if $$$x=12=(1100)_2$$$, $$$n=1$$$, $$$m=2$$$, then $$$l=1=(1)_2$$$, $$$r=4=(100)_2$$$), then the final value of $$$x$$$ after $$$n+m$$$ operations will be either $$$l$$$ or $$$l+1$$$, depending on whether a carry-over occurs in $$$r$$$ during the last operation.

If there are no $$$1$$$s in the higher $$$m$$$ bits of $$$r$$$, then the carry-over in $$$r$$$ during the last operation can never occur (operating on some small tests may give you a better understanding). Otherwise, we can choose to perform $$$n$$$ floor operations first followed by $$$m$$$ ceil operations, which can guarantee a carry-over in $$$r$$$ during the last operation.

This demonstrates that performing $$$n$$$ floor operations first followed by $$$m$$$ ceil operations will always yield the maximum value. Similarly, performing $$$m$$$ ceil operations first followed by $$$n$$$ floor operations will always produce the minimum value. We can simply simulate these operations because after $$$O(\log x)$$$ same operations, $$$x$$$ will remain unchanged.

The time complexity is $$$O(T\log x)$$$.

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll t,x,n,m;
inline ll read(){
	ll s = 0,w = 1;
	char ch = getchar();
	while (ch > '9' || ch < '0'){ if (ch == '-') w = -1; ch = getchar();}
	while (ch <= '9' && ch >= '0') s = (s << 1) + (s << 3) + (ch ^ 48),ch = getchar();
	return s * w;
}
ll F(ll x,ll n){
	while (n --){
		if (!x) return x;
		x = (x >> 1);
	}
	return x;
}
ll C(ll x,ll n){
	while (n --){
		if (x <= 1) return x;
		x = ((x + 1) >> 1);
	}
	return x;
}
int main(){
	t = read();
	while (t --){
		x = read(),n = read(),m = read();
		printf("%lld %lld\n",F(C(x,m),n),C(F(x,n),m));
	}
	return 0;
}

Consider $$$x$$$ in binary form.

What are all possible numbers of operations to make $$$x$$$ equal to $$$1$$$?

Consider using dynamic programming to calculate the answer.

Similar to the previous problem, the possible number of operations can only be $$$n-1$$$ or $$$n$$$, depending on whether a carry-over occurs in $$$x$$$ during the $$$(n-1)$$$ -th operation.

We can use dynamic programming to compute the probability of a carry-over occurring during the $$$(n-1)$$$ -th operation.

Let $$$f_i$$$ denote the probability of a carry-over occurring at the $$$i$$$ -th least significant bit. The recurrence relation is:

The final answer is given by $$$n-1+f_{n-1}$$$.

The time complexity is $$$O(\sum n)$$$.

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll inv2 = 5e8 + 4,mod = 1e9 + 7;
ll t,n;
char s[1000009];
inline ll read(){
	ll s = 0,w = 1;
	char ch = getchar();
	while (ch > '9' || ch < '0'){ if (ch == '-') w = -1; ch = getchar();}
	while (ch <= '9' && ch >= '0') s = (s << 1) + (s << 3) + (ch ^ 48),ch = getchar();
	return s * w;
}
int main(){
	t = read();
	while (t --){
		n = read(),scanf("%s",s + 1);
		ll ans = 0;
		for (ll i = n;i > 1;i -= 1) ans = (ans + (s[i] == '1')) * inv2 % mod;
		printf("%lld\n",(n - 1 + ans) % mod);
	}
	return 0;
}

Let's only consider the comparison between adjacent pairs first. Call a pair $$$(a_i,a_{i+1})$$$ an inversion pair if $$$a_i \gt a_{i+1}$$$.

In each operation, how can the number of inversion pairs change?

Under what circumstances can we decrease the maximum number of inversion pairs in every operation?

Let's only consider the comparison between adjacent pairs first. Call a pair $$$(a_i,a_{i+1})$$$ an inversion pair if $$$a_i \gt a_{i+1}$$$. Our goal is to decrease the number of inversion pairs to zero.

If we choose to replace $$$a_l,a_{l+1},\ldots,a_r$$$, only the comparison between $$$(a_{l-1},a_l)$$$ and $$$(a_{r},a_{r+1})$$$ can change, which means we can eliminate at most two inversion pairs in one operation. Let $$$s$$$ be the number of inversion pairs in $$$a_1,a_2,\ldots,a_n$$$, then a lower bound of the answer is $$$l=\left\lceil\dfrac{s}{2}\right\rceil$$$.

Case 1: $$$s \gt 0$$$ and $$$2\mid s$$$

In order to obtain the lower bound $$$l$$$, each operation must eliminate exactly two inversion pairs. This necessitates that we cannot alter any numbers before the first element $$$a_{p_1}$$$ of the first inversion pair or after the second element $$$a_{p_2}$$$ of the last inversion pair. (Why?)

An intuitive idea is that the difference between $$$a_{p_1}$$$ and $$$a_{p_2}$$$ should be as large as possible (more precisely, $$$a_{p_2} - a_{p_1} \geq p_2 - p_1$$$), so that a strictly increasing sequence can be "accommodated" between them. It is not difficult to constructively prove by recursion that this condition is both necessary and sufficient.

In other words, when $$$a_{p_2} - a_{p_1} \geq p_2 - p_1$$$, the answer is $$$l$$$; otherwise, an additional operation is required, making the answer $$$l+1$$$.

Case 2: $$$s=0$$$ or $$$2\not\mid s$$$

It's not difficult to prove that the answer in this case is $$$l$$$. (Why?)

The time complexity is $$$O(\sum n)$$$.

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll t,n,a[200009];
inline ll read(){
	ll s = 0,w = 1;
	char ch = getchar();
	while (ch > '9' || ch < '0'){ if (ch == '-') w = -1; ch = getchar();}
	while (ch <= '9' && ch >= '0') s = (s << 1) + (s << 3) + (ch ^ 48),ch = getchar();
	return s * w;
}
int main(){
	t = read();
	while (t --){
		n = read();
		for (ll i = 1;i <= n;i += 1) a[i] = read();
		ll ans = 0,pos1 = 0,pos2 = 0;
		for (ll i = 1;i < n;i += 1) if (a[i] > a[i + 1]){
			ans += 1,pos2 = i + 1;
			if (!pos1) pos1 = i;
		}
		if ((ans & 1) || (pos1 && a[pos2] - a[pos1] < pos2 - pos1)) printf("%lld\n",(ans >> 1) + 1);
		else printf("%lld\n",ans >> 1);
	}
	return 0;
}

Do we really have to consider the problem on the whole matrix?

Use greedy algorithm and guess some possible conclusions. Can you prove them (or construct some counterexamples) ?

The problem can be rephrased as follows: Let the XOR sum of each row be $$$ r_1, r_2, \ldots, r_n $$$, and the XOR sum of each column be $$$ c_1, c_2, \ldots, c_m $$$. In each operation, we can choose any $$$ 1 \le i \le n $$$, $$$ 1 \le j \le m $$$, and $$$ 0 \le x \le 3 $$$, then let $$$ r_i \leftarrow r_i \oplus x $$$ and $$$ c_j \leftarrow c_j \oplus x $$$. The goal is to determine the minimum number of operations required to make all $$$ r_1, r_2, \ldots, r_n $$$ and $$$ c_1, c_2, \ldots, c_m $$$ zero.

Let $$$ R_x\ ( 0 \le x \le 3 )$$$ denote the count of $$$ x $$$ in $$$ r_1, r_2, \ldots, r_n $$$, and $$$ C_x $$$ denote the count of $$$ x $$$ in $$$ c_1, c_2, \ldots, c_m $$$. A trivial upper bound for the answer is $$$ s = R_1 + R_2 + R_3 + C_1 + C_2 + C_3 $$$, which corresponds to the scenario where each operation only zeros out one of $$$ r_i $$$ and $$$ c_j $$$.

To reduce the number of operations, we can consider grouping the non-zero elements in $$$ r_1, \ldots, r_n, c_1, \ldots, c_m $$$ into disjoint groups where:

1. The XOR sum of each group is zero.
2. Each group contains at least one element from $$$ r $$$ and one from $$$ c $$$.

For such a group, we can reduce the number of operations by one, as there always exists an operation which can zero out both $$$ r_i $$$ and $$$ c_j $$$. Thus, the goal is to maximize the number of such groups, and the minimum number of operations will be $$$s - \text{(number of groups formed)}$$$.

Define the following group types:

• $$$ P_2 $$$: Groups of the form $$$(R_1,C_1)$$$, $$$(R_2,C_2)$$$ or $$$(R_3,C_3)$$$.
• $$$ P_3 $$$: Groups of the form $$$(R_1,R_2,C_3)$$$, $$$(R_1,C_2,C_3)$$$, $$$(R_1,C_2,R_3)$$$, $$$(C_1,C_2,R_3)$$$, $$$(C_1,R_2,R_3)$$$ or $$$(C_1,R_2,C_3)$$$.
• $$$ P_4 $$$: Groups of the form $$$(R_1,R_1,C_2,C_2)$$$, $$$(R_1,R_1,C_3,C_3)$$$, $$$(R_2,R_2,C_1,C_1)$$$, $$$(R_2,R_2,C_3,C_3)$$$, $$$(R_3,R_3,C_1,C_1)$$$ or $$$(R_3,R_3,C_2,C_2)$$$.

Conclusion 1: There exists an optimal solution where all groups are of the form $$$ P_2 $$$, $$$ P_3 $$$, or $$$ P_4 $$$.

Conclusion 2: There exists an optimal solution where all $$$ P_3 $$$ groups share the same structure, and so do all $$$ P_4 $$$ groups.

Conclusion 3: There exists an optimal solution where the elements in $$$ P_4 $$$ groups are a subset of those in $$$ P_3 $$$ groups (e.g., if $$$ P_3 $$$ is $$$ (R_1, R_2, C_3) $$$, then $$$ P_4 $$$ can only be $$$ (R_1, R_1, C_3, C_3) $$$ or $$$ (R_2, R_2, C_3, C_3) $$$). The proof is analogous to Conclusion 2 and is left as an exercise.

Based on these conclusions, the optimal strategy is: first form as many $$$P_2$$$ groups as possible; then form as many $$$P_3$$$ groups as possible from the remaining elements; then form as many $$$P_4$$$ groups as possible from the remaining elements.

After grouping, each group can be zeroed out with one less operation. Note that the remaining ungrouped elements require individual operations.

The time complexity is $$$O(\sum nm)$$$.

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll t,n,m,a[1009][1009];
char s[1009];
vector <ll> row[4],col[4];
inline ll read(){
	ll s = 0,w = 1;
	char ch = getchar();
	while (ch > '9' || ch < '0'){ if (ch == '-') w = -1; ch = getchar();}
	while (ch <= '9' && ch >= '0') s = (s << 1) + (s << 3) + (ch ^ 48),ch = getchar();
	return s * w;
}
int main(){
	t = read();
	while (t --){
		n = read(),m = read();
		for (ll i = 1;i <= 3;i += 1) row[i].clear(),col[i].clear();
		for (ll i = 1;i <= n;i += 1){
			scanf("%s",s + 1);
			for (ll j = 1;j <= m;j += 1) a[i][j] = s[j] - '0';
		}
		for (ll i = 1;i <= n;i += 1){
			ll sum = 0;
			for (ll j = 1;j <= m;j += 1) sum ^= a[i][j];
			if (sum) row[sum].emplace_back(i);
		}
		for (ll j = 1;j <= m;j += 1){
			ll sum = 0;
			for (ll i = 1;i <= n;i += 1) sum ^= a[i][j];
			if (sum) col[sum].emplace_back(j);
		}
		ll ans = 0;
		for (ll i = 1;i <= 3;i += 1){
			while (row[i].size() && col[i].size()){
				ans += 1;
				a[row[i].back()][col[i].back()] ^= i;
				row[i].pop_back(),col[i].pop_back();
			}
		}
		for (ll i = 1;i <= 3;i += 1){
			ll j = i % 3 + 1,k = j % 3 + 1;
			while (row[i].size() && col[j].size() && col[k].size()){
				ans += 2;
				a[row[i].back()][col[j].back()] ^= j;
				a[row[i].back()][col[k].back()] ^= k;
				row[i].pop_back(),col[j].pop_back(),col[k].pop_back();
			}
			while (col[i].size() && row[j].size() && row[k].size()){
				ans += 2;
				a[row[j].back()][col[i].back()] ^= j;
				a[row[k].back()][col[i].back()] ^= k;
				col[i].pop_back(),row[j].pop_back(),row[k].pop_back();
			}
		}
		for (ll i = 1;i <= 3;i += 1) for (ll j = 1;j <= 3;j += 1) if (i != j){
			while (row[i].size() >= 2 && col[j].size() >= 2){
				ll r1 = row[i].back(); row[i].pop_back();
				ll r2 = row[i].back(); row[i].pop_back();
				ll c1 = col[j].back(); col[j].pop_back();
				ll c2 = col[j].back(); col[j].pop_back();
				ans += 3;
				a[r1][c1] ^= i;
				a[r2][c1] ^= (i ^ j);
				a[r2][c2] ^= j;
			}
		}
		for (ll i = 1;i <= 3;i += 1){
			while (row[i].size()) ans += 1,a[row[i].back()][1] ^= i,row[i].pop_back();
			while (col[i].size()) ans += 1,a[1][col[i].back()] ^= i,col[i].pop_back();
		}
		printf("%lld\n",ans);
		for (ll i = 1;i <= n;i += 1,puts("")) for (ll j = 1;j <= m;j += 1) printf("%lld",a[i][j]);
	}
	return 0;
}

There are only $$$O(n\log n)$$$ edges that are truly useful in the whole graph.

First, for vertices with the same weight, they can naturally be treated as a single vertex.

We observe that in the interval of weights $$$(kx,(k+1)x)\ (k\ge 1)$$$, only the smallest $$$y$$$ is useful. It is because, if $$$x \lt y \lt z \lt 2x$$$, connecting edges $$$(x,y)$$$ and $$$(y,z)$$$ will always be better than connecting $$$(x,y)$$$ and $$$(x,z)$$$.

Therefore, for a vertex with weight $$$x$$$, we enumerate $$$k = x , 2 x , \ldots , \lfloor \frac{n}{x}\rfloor\cdot x$$$. For each $$$k$$$, we connect $$$x$$$ to the smallest $$$y$$$ that is greater than or equal to $$$k$$$.

It can be proven that this approach guarantees a connected graph (since connecting only the first edge larger than $$$x$$$ already forms a connected structure).

This method generates $$$O(n\log n)$$$ edges (due to the harmonic series property). Directly applying Kruskal's algorithm would result in a time complexity of $$$O(n\log^2 n)$$$. With radix sort, this can be optimized to $$$O(n\log n)$$$.

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<string.h>
#include<vector>
#include<queue>
#include<map>
#include<ctime>
#include<bitset>
#include<set>
#include<math.h>
//#include<unordered_map>
#define fi first
#define se second
#define mp make_pair
#define pii pair<int,int>
#define pb push_back
#define pil pair<int,long long>
#define pll pair<long long,long long>
#define vi vector<int>
#define vl vector<long long>
#define ci ios::sync_with_stdio(false)
//#define int long long
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
int read(){
	char c=getchar();
	ll x=1,s=0;
	while(c<'0'||c>'9'){
	   if(c=='-')x=-1;c=getchar();
	}
	while(c>='0'&&c<='9'){
	   s=s*10+c-'0';c=getchar();
	}
	return s*x;
}
const int N=5e5+55;
int MAXN=0;
int n,tot,ans,maxx,f[N],siz[N],p[N],bj[N],pre[N],suf[N],pos[N];
struct node{
	int from,to,dis;
}a[N*20];
int cmp(node x,node y){
	return x.dis<y.dis;
}
int find(int x){
	if(f[x]!=x)return f[x]=find(f[x]);
	else return x;
}
void merge(int x,int y){
	x=find(x);y=find(y);
	if(siz[x]<siz[y])swap(x,y);
	siz[x]+=siz[y];
	f[y]=x;
}
int T;
int main(){
	T=read();
	while(T--)
	{
		tot=ans=maxx=0;
		for(int i=1;i<=MAXN;++i)
		{
			f[i]=siz[i]=p[i]=bj[i]=pre[i]=suf[i]=pos[i]=0;
		}
		for(int i=0;i<=MAXN*19;++i)
		a[i].from=a[i].to=a[i].dis=0;
		
		
		
		n=read();
		MAXN=n;
		for(int i=1;i<=n;i++)p[i]=read(),bj[p[i]]=1,maxx=max(maxx,p[i]),MAXN=max(MAXN,p[i]);
		sort(p+1,p+n+1);
		for(int i=1;i<=maxx;i++){
			if(bj[i])pre[i]=i;
			else pre[i]=pre[i-1];
		}
		for(int i=maxx;i>=1;i--){
			if(bj[i])suf[i]=i;
			else suf[i]=suf[i+1];
		}
		for(int i=1;i<=n;i++)pos[p[i]]=i;
		for(int i=1;i<=n;i++){
			if(p[i]==p[i-1]){
				a[++tot]={i,i-1,0};
				continue;
			}
			for(int k=1;k*p[i]<=maxx;k++){
				int tmy=0;
				if(bj[k*p[i]] and k!=1){
					a[++tot]={pos[k*p[i]],i,0};
					continue;
				}
				if(suf[k*p[i]+1]<=(k+1)*p[i])tmy=suf[k*p[i]+1];
				if(!tmy)continue;
				int cost=abs(k*p[i]-tmy);
				a[++tot]={pos[tmy],i,cost};
			}
		}
		sort(a+1,a+tot+1,cmp);
		for(int i=1;i<=n;i++)f[i]=i,siz[i]=1;
		for(int i=1;i<=tot;i++){
			if(find(a[i].from)!=find(a[i].to)){
				merge(a[i].from,a[i].to);
				ans+=a[i].dis;
			}
		}
		cout<<ans<<endl;
	}
	return 0;
}

Move each chip as deep as possible in descending order of labels. Can we simplify the operations after that?

Now we only need to move the chips upward and swap adjacent chips with the same color. What can we do after that?

We can perform dynamic programming (DP) to calculate the numbers. What information do we need to record?

We need to record the color and the length of the top monochromatic segment. How can we perform the DP in $$$\mathcal{O}(m^2)$$$ time complexity?

First, in descending order of labels, move each chip to the deepest possible position successively. It can be shown that this allows each chip to reach the theoretically deepest position it can attain. Consequently, every final state can be obtained from this configuration by performing only upward moves of chips and swaps of adjacent chips with the same color.

Next comes the process of moving chips upward while performing dynamic programming (DP). Let $$$e_u$$$ denote the edge from node $$$u$$$ to its parent. Assume we have decided the order of all the chips in subtree $$$u$$$ (except the ones on $$$e_u$$$), and we need to merge them with the original chips on $$$e_u$$$. If the bottom chip on $$$e_u$$$ shares color with the subtree's top chip, we need to know the length of the subtree's topmost maximal monochromatic contiguous segment to calculate the number of ways to merge the two parts, i.e., if the length of the bottommost segment on $$$e_u$$$ is $$$x$$$ and the length of the topmost segment of the subtree is $$$y$$$, we have $$$\dbinom{x+y}{x}$$$ ways to combine them. So we need to record the color and the length of the top segment.

Let $$$f_{u,0/1,i}$$$ denote the number of permutation schemes on $$$e_u$$$ after moving all the chips in subtree $$$u$$$ to $$$e_u$$$, where:

• The topmost chip on $$$e_u$$$ is black/white (0/1).
• The length of the topmost maximal monochromatic contiguous segment on $$$e_u$$$ is $$$i$$$.

When calculating $$$f_{u,0/1,i}$$$ for all $$$i$$$, we first merge the DP states of all the subtrees of node $$$u$$$, and then merge the original chips on $$$e_u$$$ into the resulting DP states.

Subtree merging: Assume we merge two subtrees $$$v$$$ and $$$w$$$. Let $$$S$$$ and $$$T$$$ denote the number of chips on $$$e_v$$$ and $$$e_w$$$ respectively, and let $$$g_{0/1,i}$$$ be a temporary DP array initialized to all zeros.

Let's consider two cases.

Case 1: Different colors at tops of $$$e_v$$$ and $$$e_w$$$

Suppose the top chip comes from $$$e_v$$$ after merging. Then the maximal monochromatic segment length from $$$e_w$$$ becomes irrelevant, and therefore let $$$h_c=\sum_{i=1}^T f_{w,c,i}$$$. If the top chip of $$$e_w$$$ would split the top segment of $$$e_v$$$, enumerate split positions. The split position should be strictly inside the top segment of $$$e_v$$$.

Otherwise, the position of the top chip of $$$e_w$$$ becomes irrelevant.

With suffix sum optimization, this achieves $$$\mathcal{O}(m^2)$$$ time complexity.

Case 2: Same color at tops

If there are only chips of this color on both edges:

Otherwise, assume the topmost different-colored chip is from subtree $$$w$$$ after merging. Enumerate insertion positions where this chip splits the top segment of the other (not necessarily strictly inside the segment).

(Similar transition equation for $$$v$$$.)

With suffix sum optimization, this forms a tree knapsack DP with $$$\mathcal{O}(m^2)$$$ time complexity.

After computing $$$g_{0/1,i}$$$, we can treat $$$g$$$ as the DP states of a subtree, and continue to perform the process above.

Merging original chips on $$$e_u$$$ into the results: We first account for internal permutations of original monochromatic segments. After that, if the bottom chip on $$$e_u$$$ shares color with merged subtrees' top chip (with segment lengths $$$x$$$ and $$$y$$$), multiply the DP state by $$$\dbinom{x+y}{x}$$$.

Assuming $$$n$$$ and $$$m$$$ are of the same order, the overall time complexity is $$$\mathcal{O}(m^2)$$$.

#include <bits/stdc++.h>
#define eb emplace_back
using namespace std;
typedef long long ll;
const int p=998244353;
int T,n,m,fa[5005],col[5005],d[5005],stk[5005],top,ct[5005][2];
int C[5005][5005],fac[5005],siz[5005],f[5005][2][5005],g[2][5005],ans;
vector<int> to[5005],vec[5005];
inline void diff(int u,int v){
	int S=siz[u],T=siz[v];
	for(int x=0,A,B;x<2;++x){
		A=B=0;
		for(int j=1;j<=T;++j)B=(B+f[v][x^1][j])%p;
		for(int i=S;i>=1;--i){
			g[x][i]=(g[x][i]+((ll)A*C[S+T-i-1][T-1]+(ll)f[u][x][i]*C[S+T-i][T])%p*B)%p;
			A=(A+f[u][x][i])%p;
		}
	}
}
inline void same1(int u,int v){
	int S=siz[u],T=siz[v];
	for(int x=0,A,B;x<2;++x){
		for(int i=1;i<S;++i){
			if(!(A=f[u][x][i]))continue;
			B=0;
			for(int j=T;j>=0;--j){
				B=(B+f[v][x][j])%p;
				g[x][i+j]=(g[x][i+j]+(ll)A*B%p*C[i+j][j]%p*C[S-i-1+T-j][T-j])%p;
			}
		}
	}
}
inline void same(int u,int v){
	same1(u,v),same1(v,u);
	int S=siz[u],T=siz[v];
	for(int x=0;x<2;++x)
		g[x][S+T]=(g[x][S+T]+(ll)f[u][x][S]*f[v][x][T]%p*C[S+T][T])%p;
}
inline void merge(int u,int v){
	for(int i=1;i<=siz[u]+siz[v];++i)g[0][i]=g[1][i]=0;
	diff(u,v),diff(v,u),same(u,v);
	siz[u]+=siz[v];
	for(int i=1;i<=siz[u];++i)
		f[u][0][i]=g[0][i],f[u][1][i]=g[1][i];
}
void dfs(int u){
	siz[u]=0;
	for(auto v:to[u]){
		dfs(v);
		if(!siz[v])continue;
		if(siz[u])merge(u,v);
		else{
			for(int i=1;i<=siz[v];++i)
				f[u][0][i]=f[v][0][i],f[u][1][i]=f[v][1][i];
			siz[u]=siz[v];
		}
	}
	if(vec[u].empty())return;
	int prd=1,L=vec[u].size(),c=vec[u][0],fir=0,lst=0,res=0;
	for(int i=0,j=0;i<L;i=j){
		while(j<L && vec[u][i]==vec[u][j])++j;
		if(i==0)fir=j-i;
		if(j==L)lst=j-i;
		prd=(ll)prd*fac[j-i]%p;
	}
	for(int i=siz[u]+1;i<=siz[u]+L;++i)f[u][0][i]=f[u][1][i]=0;
	if(!siz[u])f[u][vec[u][L-1]][lst]=prd;
	else if(fir==L){
		for(int i=siz[u];i>=1;--i){
			f[u][c][i+L]=(ll)f[u][c][i]*C[i+L][i]%p*fac[L]%p;
			res=(res+f[u][c^1][i])%p;
			f[u][c^1][i]=0;
		}
		for(int i=1;i<L;++i)f[u][c][i]=0;
		f[u][c][L]=(ll)res*fac[L]%p;
	}
	else{
		for(int i=1;i<=siz[u];++i){
			res=(res+(ll)f[u][c][i]*C[i+fir][i]+f[u][c^1][i])%p;
			f[u][0][i]=f[u][1][i]=0;
		}
		f[u][vec[u][L-1]][lst]=(ll)res*prd%p;
	}
	siz[u]+=L;
}
int main(){
	scanf("%d",&T);
	while(T--){
		scanf("%d%d",&n,&m);
		for(int i=0;i<=n;++i){
			to[i].clear(),vec[i].clear();
			ct[i][0]=ct[i][1]=0;
		}
		for(int i=1;i<=n;++i)
			scanf("%d",fa+i),to[fa[i]].eb(i);
		for(int i=1;i<=m;++i)scanf("%d",col+i);
		for(int i=1;i<=m;++i)scanf("%d",d+i);
		for(int i=m,u;i>=1;--i){
			top=0,u=d[i];
			while(u)stk[++top]=u,u=fa[u];
			for(int j=top;j>=1;--j){
				if(j>1 && !ct[stk[j]][col[i]^1])continue;
				vec[stk[j]].eb(col[i]),++ct[stk[j]][col[i]];
				break;
			}
		}
		C[0][0]=fac[0]=1;
		for(int i=1;i<=m;++i){
			for(int j=1;j<=i;++j)
				C[i][j]=(C[i-1][j-1]+C[i-1][j])%p;
			C[i][0]=1;
			fac[i]=(ll)fac[i-1]*i%p;
		}
		dfs(1),ans=0;
		for(int i=1;i<=m;++i)
			ans=((ll)ans+f[1][0][i]+f[1][1][i])%p;
		printf("%d\n",ans);
	}
	return 0;
}

First, it is not difficult to observe that when $$$ n $$$ is odd, a solution exists if and only if $$$ n = 1 $$$.

When $$$ n $$$ is even, a solution always exists.

Consider the following construction:

Let $$$ n = 2m $$$. Use $$$(x, y)$$$ to denote the cell in the $$$ x $$$-th row and $$$ y $$$-th column of the matrix, and assume the top-left corner of the matrix is $$$(1, 1)$$$, while the bottom-right corner is $$$(n, n)$$$.

For $$$(a, b)$$$ and $$$(c, d)$$$, let $$$ k_1 = c - a $$$ and $$$ k_2 = d - b $$$. When $$$ |k_1| = |k_2| $$$, denote $$$(a, b) \sim (c, d)$$$ as:

• $$$(a, b), (a+1, b+1), \dots, (a+k_1, b+k_2)$$$ if $$$ k_1, k_2 \geq 0 $$$.
• $$$(a, b), (a+1, b-1), \dots, (a+k_1, b+k_2)$$$ if $$$ k_1 \geq 0, k_2 \leq 0 $$$.
• $$$(a, b), (a-1, b+1), \dots, (a+k_1, b+k_2)$$$ if $$$ k_1 \leq 0, k_2 \geq 0 $$$.
• $$$(a, b), (a-1, b-1), \dots, (a+k_1, b+k_2)$$$ if $$$ k_1, k_2 \leq 0 $$$.

For any $$$ 1 \leq i \leq m $$$, we alternately fill the numbers $$$ 2i-1 $$$ and $$$ 2i-2 $$$ in the cells along the paths $$$(1, 2i) \sim (n-2i+1, n)$$$, $$$(n-2i+2, n) \sim (n, n-2i+2)$$$, $$$(n, n-2i+1) \sim (2i, 1)$$$, and $$$(2i-1, 1) \sim (1, 2i-1)$$$.

An illustration for $$$ n = 10 $$$:

Next, we prove that this construction satisfies all the requirements of the problem:

• Each row and column of $$$ A $$$ is a permutation of $$$ 0 \sim n-1 $$$.

• $$$ a_{i,j} + a_{i,n-j+1} = n-1 $$$ and $$$ a_{i,j} + a_{n-i+1,j} = n-1 $$$.

• All ordered pairs $$$ \langle a_{i,j}, a_{i,j+1} \rangle $$$ (for $$$ 1 \leq i \leq n, 1 \leq j \lt n $$$) are distinct, and all ordered pairs $$$ \langle a_{i,j}, a_{i+1,j} \rangle $$$ (for $$$ 1 \leq i \lt n, 1 \leq j \leq n $$$) are distinct.

In conclusion, we have provided a construction for even $$$ n $$$ and proven that it satisfies all the problem's requirements.

The time complexity is $$$O(\sum n^2)$$$.

#include<bits/stdc++.h>
using namespace std;
int Test_num,n;
int a[3002][3002];
void solve()
{
	scanf("%d",&n);
	if(n>=3 && (n&1))return (void)(puts("NO"));
	puts("YES");
	if(n==1)return (void)(puts("0"));
	for(int i=0;i<n;i+=2)
	{
		for(int j=1;j<=i+1;++j)a[j][i+2-j]=i+((j&1)^1);
		for(int j=i+2;j<=n;++j)a[j][j-i-1]=i+((j&1)^1);
		for(int j=1;j<=n-i-1;++j)a[j][i+1+j]=i+(j&1);
		for(int j=n-i;j<=n;++j)a[j][2*n-i-j]=i+(j&1);
	}
	for(int i=1;i<=n;++i)for(int j=1;j<=n;++j)printf("%d%c",a[i][j],j==n? '\n':' ');
}
int main()
{
	for(scanf("%d",&Test_num);Test_num--;)solve();
	return 0;
}

If $$$g'(np)\neq g'(n)$$$ , then $$$p\le n$$$ unless $$$n=2,p=3$$$ or $$$n=1,p=2$$$ .

Think about the method of $$$\min 25$$$ sieve.

We define $$$\displaystyle g(n)=\frac n{\varphi(n)}=\prod_{p|n}\frac p{p-1},g'(n)=\lfloor g(n)\rfloor$$$.

It's easy to see that $$$res=\frac{N(N+1)}2-\sum_{i=1}^Ng'(i)\varphi(i)$$$.

Let's consider how to calculate $$$\sum_{i=1}^Ng'(i)\varphi(i)$$$.

$$$\texttt{Key observation 1}$$$：if $$$g'(np)\neq g'(n)$$$, then $$$p\le n$$$ unless $$$n=2,p=3$$$ or $$$n=1,p=2$$$.

Therefore, all the prime factors $$$p$$$ that $$$\gt \sqrt N$$$ wouldn't change the number of $$$g'(n)$$$.

Use a linear sieve to find all the prime numbers $$$\le \sqrt N$$$, following the approach for handling composite numbers in the min25 sieve：

Let's define:

When $$$p\leq\frac{n+(\varphi(n)-n\mod\varphi(n))}{\varphi(n)-n\mod\varphi(n)}$$$, there is $$$g'(np)=g'(n)+1$$$, otherwise $$$g'(np)=g'(n)$$$.

Pre-calculate the prefix of $$$\varphi(p)$$$ among prime numbers, we can calculate $$$s(n,k)$$$ like this:

The answer is $$$\frac{N(N+1)}{2}-1-S(1,1)$$$.

With time complexity $$$\mathcal O\left(\frac{n^{0.75}}{\log n}\right)$$$, it could pass the test cases of $$$N\le 10^{11}$$$.

#include<cstdio>
#include<cmath>
typedef unsigned long long ull;
typedef unsigned uint;
const uint M=1e6+5;
ull N;const uint ANS[11]={0,0,0,1,1,2,2,3,3,6,8};
uint B,n4,top,pri[M],F0[M],F1[M],G0[M],G1[M],SF[M],SG[M];
inline ull Div(const ull&n,const ull&m){
	return double(n)/m;
}
inline ull min(const ull&a,const ull&b){
	return a>b?b:a;
}
inline uint DFS(const ull&n,uint k,const ull&phi){
    const ull&T=Div(N,n);
    const uint&g=int(n*1./phi),&R=min(Div((g+1)*phi,(g+1)*phi-n),min(ull(B),T));
	uint ans=0;
	if(k>top)ans=T<=B?0:g*phi*(SF[n]-SG[B]);
	else if(pri[k]>T)ans=0;
	else if(pri[k]<=R)ans=g*phi*((n<=B?SF[n]:SG[T])-SG[R])+(g+1)*phi*(SG[R]-SG[pri[k-1]]);
	else ans=g*phi*((n<=B?SF[n]:SG[T])-SG[pri[k-1]]);
	for(;k<=top;++k){
		const uint&p=pri[k],&w=p<=R?g+1:g;if(1ull*p*p>T)break;
		for(ull x=n,tp=phi*(p-1);x*p<=N;tp*=p)ans+=DFS(x*=p,k+1,tp)+w*tp;
		ans-=w*phi*(p-1);
	}
	return ans;
}
signed main(){
	scanf("%llu",&N);
	if(N<=10)return!printf("%d",ANS[N]);B=sqrt(N);n4=sqrt(B);
	for(uint i=1;i<=B;++i){
		const ull&w=Div(N,i);
		F0[i]=w-1;F1[i]=w*(w+1)/2-1;
		G0[i]=i-1;G1[i]=i*(i+1ull)/2-1;
	}
	for(uint p=2;p<=n4;++p)if(G0[p]!=G0[p-1]){
		const uint&lim=Div(B,p),&S0=G0[p-1],&S1=G1[p-1];pri[++top]=p;
		for(uint k=1;k<=lim;++k){
			F0[k]-=F0[k*p]-S0;F1[k]-=p*(F1[k*p]-S1);
		}
		for(uint k=lim+1;k<=B;++k){
			const uint&x=Div(N,k*p);F0[k]-=G0[x]-S0;F1[k]-=p*(G1[x]-S1);
		}
		for(uint k=B,e=lim;e>=p;--e)for(uint x=e*p,V0=G0[e]-S0,V1=p*(G1[e]-S1);k>=x;--k){
			G0[k]-=V0;G1[k]-=V1;
		}
	}
	for(uint p=n4+1;p<=B;++p)if(G0[p]!=G0[p-1]){
		const uint&lim=Div(B,p),&T=Div(N,1ll*p*p),&S0=G0[p-1],&S1=G1[p-1];pri[++top]=p;
		for(uint k=1;k<=lim;++k){
			F0[k]-=F0[k*p]-S0;F1[k]-=p*(F1[k*p]-S1);
		}
		for(uint k=lim+1;k<=T;++k){
			const uint&x=Div(N,k*p);F0[k]-=G0[x]-S0;F1[k]-=p*(G1[x]-S1);
		}
	}
	for(uint i=1;i<=B;++i)SF[i]=F1[i]-F0[i],SG[i]=G1[i]-G0[i];
	printf("%u",N*(N+1)/2-1-DFS(1,1,1));
}

Consider the DFS structure of $$$\min 25$$$ sieve.

What kind of nodes are important?

How many of nodes are important?

Consider using Du Sieve to get some important functions.

For further optimization, let's consider the DFS tree of $$$\min25$$$ sieve: $$$N$$$ vertices, rooted on $$$1$$$, the father of a node $$$d$$$ is $$$\frac{d}{\text{maxp}(d)}$$$ where $$$\text{maxp}(d)$$$ denotes the maximum prime factor of $$$d$$$.

Consider all the positions $$$d$$$ that $$$g'(d)\neq g'(fa_d)$$$. Let's divide them into two cases:

1. $$$d$$$ is a leaf, which means $$$d\le N\lt d\times \operatorname{maxp}(d)$$$.
2. $$$d$$$ is not a leaf, which means $$$d\times \operatorname{maxp}(d)\le N$$$.

The answer could be calculated like this:

where $$$G(d,p)$$$ denotes the sum of $$$\varphi$$$ in the subtree of node $$$d$$$, and $$$G(d,p)$$$ could be calculated like this:

We could use DFS to search all the nodes $$$d$$$ which satisfy case 2, and we can calculate its contribution to the answer during our DFS.

For a non-leaf node $$$u$$$, let's consider its child nodes $$$up$$$ which satisfy case 1.

It's easy to see that $$$p$$$ is in an interval, which allows $$$O(1)$$$ calculation of the contributions of these nodes.

As a result, we get the number of nodes $$$d$$$.

//pi[n] means the number of prime number <=n
//Sp[k] The prefix sum of phi(p) among the first k prime numbers
//g[k]=(pri[k]-1.0)/pri[k]
inline ull div(ull n,ull m){
    return n*1./m;
}
inline bool check(double P,ull n,uint k){
    ull x(1);while(P>1&&k<=m&&x*pri[k]<=n)P*=g[k],x*=pri[k++];return P>1;
}
bool DFS1(ull n,uint k,ull phi){
    const ull T=div(N,n);const uint w=int(1.*n/phi)+1; 
    if(check(w*phi*1./n,T,k))return true;
	//if the g' would not change in the subtree, return directly
    const uint R=pi[min(div(w*phi,w*phi-n),min(ull(B),T))];
    //g(n)*p/(p-1)>=int(g(n))+1
    for(;k<=m;++k){
        const uint&p=pri[k];if(p*p>T)break;
        if(k<=R)d[k].push_back(nb(n,phi)),n*p<=K&&(lim=max(lim,k));
        for(ull x=n,tp=phi*(p-1);x*p<=N;tp*=p)if(DFS1(x*=p,k+1,tp))break;
		//The subtree of np^e is isomorphic to the subtree of np^(e+1) , so if we didn't find in the previous one, we do not need to dfs into the latter one
		//This will not effect the efficiency significantly since dfs is not the efficiency bottleneck
    }
    if(k<=R)sum+=phi*(Sp[R]-Sp[k-1]);
	//The contribution of its child nodes which is leaf
    return false;
}

Let's consider how to calculate $$$G(d,p)$$$. Define:

where $$$\text{minp}(i)$$$ denote the minimum prime number of $$$i$$$.

Then we can see that:

It's easy to see that $$$F_0(x)=\frac{\zeta(x-1)}{\zeta(x)}$$$, we can use Du Sieve to get all the $$$f_0(\lfloor \frac Nd\rfloor)$$$ in $$$O(N^{\frac 23})$$$.

According to $$$\texttt{Key observation 1}$$$, $$$\forall d\notin{2,6},\mathrm{maxp}(d)\leq\sqrt d$$$.

Let $$$k=\max_{g'(d)\neq g'(fa_d)\And d\leq B}(\mathrm{maxp}(d))$$$, categorize these nodes $$$d$$$ again:

• (i) $$$\mathrm{maxp}(d)\leq k$$$, Then $$$\mathrm{maxp}(d)\leq \sqrt B$$$.
• (ii) $$$\mathrm{maxp}(d) \gt k$$$, then $$$B\leq\mathrm{maxp^2(d)}\leq d$$$, $$$\lfloor\frac{N}{d}\rfloor \lt \frac NB$$$.

For those $$$d$$$ satisfy condition (i), let's begin with $$$f_0$$$, delete the prime numbers from small to large, until $$$f_k$$$ and calculate the contributions of these $$$d$$$.

We could use Dirichlet convolution to calculate $$$f_{p+1}$$$ from $$$f_p$$$.

Step 1:

Step 2:

These could solve in $$$O\left(\frac{\sqrt N \times \sqrt B}{\log n}\right)$$$.

For those $$$d$$$ satisfy condition (ii), we calculate $$$G(d,p)$$$ in this method:

Because $$$p\geq \sqrt B$$$ and $$$d\geq B$$$, we only need to consider $$$k\le \frac{\log N}{\log B}$$$, which is similar to $$$O(1)$$$.

In this part, we only need $$$f(x)$$$ for $$$x\le B$$$, so to calculate these $$$f$$$ is actually a two-dimensional partial order problem which can be solved using a BIT, time complexity $$$O(\frac NB\log N)$$$.

Since for each $$$d$$$, we get $$$G(d,p)$$$ in $$$O(1)$$$ in case (i), and $$$O(\log N)$$$ in case (ii).

The total time complexity is $$$O\left(N^{\frac 23} + \frac{\sqrt N\times \sqrt B}{\log N}+ \frac NB\log N+D\log N\right)$$$.

Let $$$B=N^{\frac 13}\log^{\frac 43}N$$$, the total time complexity is $$$O\left(N^{\frac 23}+ \frac{N^{\frac 23}}{\log^{\frac 13}N}+D\log N\right)$$$.

Further optimization to $$$n=10^{13}$$$

This part is included in the offline contest.

$$$\texttt{Key observation 2}$$$ : The number of $$$i$$$ in range $$$[1,n]$$$ that satisfy $$$\mathrm{minp}\geq n^{\frac{1}{4}}$$$ is of the order of $$$\frac{n}{\ln n}$$$.

Since we only need to calculate $$$G(\frac{d}{\mathrm{maxp}(d)},\mathrm{maxp}(d))$$$ where $$$\text{maxp}(d)\gt \sqrt B$$$ in case (ii).

We actually only need to consider those $$$\text{minp}(i)\gt \sqrt B$$$, and $$$i\le \frac{d}{\text{maxp}(d)}\le\frac N{\sqrt B}$$$.

Then if $$$B\gt N^{\frac 13}$$$, the actual size to append into the BIT is only $$$O(\frac{N}{B\log N})$$$, so the time complexity of case ii is $$$O(\frac NB)$$$.

Then we could adjust $$$B$$$ to $$$N^{\frac 13}\log^{\frac 23}N$$$, the time complexity could reduce to $$$O\left(\frac {N^{\frac 23}}{\log^{\frac 23}N}+D\log N+N^{\frac 23}\right)$$$.

To further optimize, let's introduce a method that could significantly reduce the time consumption of Du Sieve.

The formula of Du Sieve is $$$\sum_{xy\leq N}f(x)g(y)=S_{f*g}(N)$$$. We can see that $$$\min(x,y)\leq\sqrt{N}$$$, Let $$$B=\lfloor\sqrt{N}\rfloor$$$, the formula could be rewritten like this:

Traditionally, we need $$$[6B,8B]$$$ times of division, but using this method need a maximum of $$$2B$$$ times of division.

As for what we need to calculate in this problem, we can get the formula:

Without using things like $$$\text{std::map}$$$, and use dynamic programming instead of DFS, we can further reduce the time consumption.

Then consider how to optimize Du Sieve to $$$O\left(\frac{N^{\frac{2}{3}}}{\log N}\right)$$$.

Definition: if $$$f$$$ is a multiplicative function, then we define a multiplicative function $$$f_x(p^k)=[p\geq x]f(p^k)$$$.

We can calculate $$$\varphi_{N^{\frac{1}{6}}}(\lfloor\frac {N}{i} \rfloor)$$$ for each $$$i$$$ and then recurrence $$$\varphi$$$ from $$$\varphi_{N^{\frac 16}}$$$. Let $$$n6=N^{\frac{1}{6}}$$$, we can change the formula into:

Where $$$\text{id}(x)=x$$$ and $$$1(x)=1$$$, which can be calculated through $$$\text{min25}$$$ Seive, but we only need to solve $$$p\leq n6$$$ instead of $$$p\leq\sqrt N$$$.

According to $$$\texttt{Key observation 2}$$$, we know that the number of $$$i\le \sqrt N$$$ that satisfy $$$\mathrm{minp}(i)\geq n6$$$ is about $$$\frac{\sqrt N}{\log N}$$$, find them and iterate these numbers in the Du Sieve. The time complexity of this part is $$$O\left(\frac{N^{\frac{2}{3}}}{\log N}\right)$$$.

At last we get a solution of $$$O\left(\frac{N^{\frac{2}{3}}}{\log N}+\frac{\sqrt{NK}}{\log N}+\frac{N}{K}\right)$$$.

#include<bitset>
#include<cstdio>
#include<vector>
#include<cmath>
#include<ctime>
using std::min;
using std::max;
typedef unsigned uint;
typedef unsigned long long ull;
namespace Sol{
    const uint M=4162277+5,Pi=1001258+5,M2=500000000+5;
    uint sum;ull N;
    uint F[M],G[M];
    uint top,Sp[Pi],pri[Pi],pi[M],phi[M2];double g[Pi];
    uint n4,K,B,m,lim,BIT[M];
    double st,ed;
    struct nb{
        ull n;uint phi;
        nb(const ull&n=1,const uint&phi=1):n(n),phi(phi){}
    };std::vector<nb>d[Pi];
    inline ull div(const ull&n,const ull&m){
        return n*1./m;
    }
    inline void Add(uint n,const uint&V){
        while(n<=B)BIT[n]+=V,n+=n&-n;
    }
    inline uint Qry(uint n){
        uint ans(0);while(n)ans+=BIT[n],n-=n&-n;return ans;
    }
    inline bool check(double P,ull n,uint k){
        ull x(1);while(P>1&&k<=m&&x*pri[k]<=n)P*=g[k],x*=pri[k++];return P>1;
    }
    inline bool DFS1(const ull&n,uint k,const ull&phi){
        const ull&T=div(N,n);const uint&w=int(n/phi)+1; 
        if(check(w*phi*1./n,T,k))return true;
        //g(n)*p/(p-1)>=int(g(n))+1
        const uint&R=pi[min(div(w*phi,w*phi-n),min(ull(B),T))];
        for(;k<=m;++k){
            const uint&p=pri[k];if(p*p>T)break;
            if(k<=R)d[k].push_back(nb(n,phi)),n*p<=B&&(lim=max(lim,k));
            for(ull x=n,tp=phi*(p-1);x*p<=N;tp*=p)if(DFS1(x*=p,k+1,tp))break;
        }
        if(k<=R)sum+=phi*(Sp[R]-Sp[k-1]);
        return false;
    }
    inline void DFS2(const uint&n,uint k,const uint&phi){
        const uint&T=div(B,n);Add(n,phi);
        for(;k<=m;++k){
            const uint&p=pri[k],&T2=div(B,p);if(p>T)return;
            for(uint x=n,tp=phi*(p-1);x<=T2;tp*=p)DFS2(x*=p,k+1,tp);
        }
    }
    void sieve(const uint&n){
        phi[1]=1;
        for(uint i=2;i<=n;++i){
            if(!phi[i]){
            	phi[i]=i-1;if(i<=B)pri[++top]=i,g[top]=1.*(i-1)/i,Sp[top]=Sp[top-1]+phi[i],pi[i]=1;
			}
            for(uint x,j=1;j<=top&&(x=i*pri[j])<=n;++j){
			    phi[x]=phi[i]*(pri[j]-1);if(div(i,pri[j])*pri[j]==i){phi[x]+=phi[i];break;}
			}
            if(i<=B)pi[i]+=pi[i-1];
        }
    }
    void Getphi(const ull&n){
    	const uint&B2=div(N,K);
		for(uint i=1;i<=B;++i)G[i]=G[i-1]+phi[i];for(uint i=div(N,B);i>=1;--i)F[B]+=phi[i];
		for(uint i=B-1;i>B2;--i){
			const uint&L=div(N,i+1)+1,&R=div(N,i);F[i]=F[i+1];for(uint k=L;k<=R;++k)F[i]+=phi[k];
		}
    	for(uint i=B2;i>=1;--i){
    		const ull&T=div(N,i);const uint&lim=sqrt(T);F[i]=lim*G[lim]+T*(T-1)/2;
    		for(uint k=2;k<=lim;++k){
    			const uint&x=div(T,k);F[i]-=(i*k<=B?F[i*k]:G[x])+phi[k]*x;
			}
		}
    }
    uint Solve(const ull&n){
    	st=clock();
        N=n;B=sqrt(N);n4=sqrt(B);
        sieve(K=uint(pow(N,2./3)));m=pi[B];
        fprintf(stderr,"sieve(%d) done %lf\n",K,((ed=clock())-st)/1000);st=ed;
        DFS1(1,1,1);
        fprintf(stderr,"DFS done %lf\n",((ed=clock())-st)/1000);st=ed;
        Getphi(N);sum+=F[1];
        fprintf(stderr,"Getphi done %lf\n",((ed=clock())-st)/1000);st=ed;
        for(uint i=1;i<=lim;++i){
            const uint&p=pri[i],&T=B/p;const ull&k=div(N,p);
            for(nb&x:d[i])sum+=x.phi*(x.n<=B?F[x.n]:G[div(N,x.n)]);
            for(uint i=1;i<=T;++i)F[i]-=p*F[i*p];
            for(uint i=T+1;i<=B;++i)F[i]-=p*G[div(k,i)];
            for(uint i=B,j=div(i,p);j;--j)for(uint e=j*p,V=p*G[j];i>=e;--i)G[i]-=V;
            for(uint i=p,j=1;i<=B;++j)for(uint e=min(i+p,B+1),V=G[j];i<e;++i)G[i]+=V;
            for(uint i=B;i>T;--i)F[i]+=G[div(k,i)];
            for(uint i=T;i>=1;--i)F[i]+=F[i*p];
            for(nb&x:d[i])sum-=x.phi*(x.n<=B?F[x.n]:G[div(N,x.n)]);
        }
        fprintf(stderr,"case1 done %lf\n",((ed=clock())-st)/1000);st=ed;
        Add(1,1);
        for(uint i=m;i>lim;--i){
            const uint&p=pri[i];
            for(nb&x:d[i]){
                ull n=N/x.n/p,phi=x.phi*(p-1);
                while(n)sum+=phi*Qry(n),n/=p,phi*=p;
            }
            for(ull x=p,phi=p-1;x<=B;x*=p,phi*=p)DFS2(x,i+1,phi);
        }
        fprintf(stderr,"case2 done %lf\n",((ed=clock())-st)/1000);st=ed;
        return n*(n+1)/2-sum;
    }
}
namespace sol{
	uint top,pri[10005],pos[10005],phi[10005];
	inline uint solve(const uint&N){
		uint sum(0);phi[1]=1;
		for(uint i=2;i<=N;++i){
			if(!pos[i])pri[pos[i]=++top]=i,phi[i]=i-1;
			for(uint x,j=1;j<=pos[i]&&(x=i*pri[j])<=N;++j)phi[x]=(pos[x]=j)==pos[i]?phi[i]*pri[j]:phi[i]*(pri[j]-1);
			sum+=i%phi[i];
		}
		return sum;
	} 
}
signed main(){
    ull n;scanf("%llu",&n);if(n<=10000)return printf("%u\n",sol::solve(n)),0;printf("%u\n",Sol::Solve(n));
}