Formalization

One may assume that \(A_1, \dots, A_N, B_1, \dots, B_N\) are pairwise distinct. (Arbitrary tie-breaking allows general cases.)

Let \((C_1, \dots, C_{2N})\) be the sequence obtained by arranging \(A_1, \dots, A_N, B_1, \dots, B_N\) in ascending order, and define \(a_i\) and \(b_i\) so that \(A_i = C_{a_i}\) and \(B_i = C_{b_i}\).

Let \(g_A(i,k)\) be the probability that the \(k\)-th pancake from the top has the size \(A_i\). Then the polynomial \(P_{a_i} := \sum_{k=1}^{N} g(i,k) x^{k-1}\) equals the product of the following for \(j \in \{1, \dots, N\}\):

• \(\frac{1}{2}\) if \(j = i\)
• \(x\) if \(j \neq i\) and \(A_j, B_j < A_i\)
• \(\frac{1}{2} + \frac{1}{2}x\) if \(j \neq i\) and \(A_j < A_i < B_j\)
• \(1\) if \(j \neq i\) and \(A_i < A_j, B_j\)

Same applies to \(B_i\). Define \(P_{b_i}\) in the same way.

\(\sum_{k=1}^{N} f(k) x^{k-1} = \sum_{i = 1}^{2N} C_{i} P_{i}\), so it suffices to find the sum of these \(2N\) polynomials.

Formalization 2

Out of nowhere, define sequences of rational functions \(L\), \(M\), and \(R\), each of length \(2N\), as the following for \(1 \leq i \leq N\):

• \(L_{a_i} = \frac{1}{2} + \frac{1}{2} x\)
• \(L_{b_i} = x / (\frac{1}{2} + \frac{1}{2} x)\)
• \(R_{a_i} = 1\)
• \(R_{b_i} = 1\)
• \(M_{a_i} = A_i\)
• \(M_{b_i} = B_i / (\frac{1}{2} + \frac{1}{2} x)\)

Then

\[C_{i} P_{i} = \left( \prod_{l < i} L_l \cdot M_i \cdot \prod_{r > i} R_r \right)\]

so the sought polynomial equals

\[\sum_{i = 1}^{2N} C_{i} P_{i} = \sum_{i=1}^{2N} \left( \prod_{l < i} L_l \cdot M_i \cdot \prod_{r > i} R_r \right).\]

Fast computation

For an interval \([l, r)\), correspond a triplet of rational functions \((\mathrm{prodL}, \mathrm{prodR}, \mathrm{sumLMR})\) defined as:

• \(\mathrm{prodL} := \prod_{l \leq i < r} L_i\)
• \(\mathrm{prodR} := \prod_{l \leq i < r} R_i\)
• \(\mathrm{sumLMR} := \sum_{l \leq i < r} \left( \prod_{l \leq l' < i} L_{l'} \cdot M_i \cdot \prod_{i < r' < r} R_{r'} \right)\)

One can obtain the triplet for \([l, r)\) based on the triplets for \([l, m)\) and \([m, r)\).

Multiplying all elements of \(L\), \(M\), and \(R\) by \((\frac{1}{2} + \frac{1}{2} x)\) multiplies the result exactly by \((\frac{1}{2} + \frac{1}{2} x)^{2N}\), and each element is a degree-\(2\) polynomial. Obtaining the triplet for \([l, r)\) based on the triplets for \([l, m)\) and \([m, r)\) costs \(O((r-l) \log (r-l))\) time (because degree-\(n\) polynomials can be multiplied in \(O(n \log n)\) time with Fast Fourier Transform), so the divide-and-conquer algorithm yields the triplet for \([0,N)\), and thus the sought polynomial.

Sample code (C++)

#pragma GCC optimize("unroll-loops")

#include <iostream>
using std::cin;
using std::cout;
using std::cerr;
using std::endl;
#include <vector>
using std::vector;
using std::pair;
#include <array>
using std::array;
#include <cassert>
using std::min;
using std::max;

#include <atcoder/modint.hpp>
using mint = atcoder::modint998244353;
#include <atcoder/convolution.hpp>
using atcoder::convolution;


#ifdef DEBUG
const int debug = 1;
#else
const int debug = 0;
#endif

using ll = int64_t;
using P = pair<ll, ll>;

const ll INF = 1'000'000'000'000'000'009LL;

ll n;
vector<ll> a, b;

inline void output (const ll x) {
	cout << x << "\n";
}

vector<mint> operator* (const mint &l, const vector<mint> &r) {
	vector<mint> x(r.size());
	for (ll i = 0; i < (ll)r.size(); i++) x[i] = l * r[i];
	return x;
}
vector<mint> operator* (const vector<mint> &l, const vector<mint> &r) {
	return convolution(l, r);
}
vector<mint> operator+ (const vector<mint> &l, const vector<mint> &r) {
	vector<mint> x(max(l.size(), r.size()), 0);
	for (ll i = 0; i < (ll)l.size(); i++) x[i] += l[i];
	for (ll i = 0; i < (ll)r.size(); i++) x[i] += r[i];
	return x;
}
vector<mint> operator- (const vector<mint> &l, const vector<mint> &r) {
	return l + (mint)(-1) * r;
}

vector<mint> fpsinv (const vector<mint> &x) {
	assert(x[0] != 0);
	vector<mint> y = {x[0].inv()};
	ll k = 1;

	while (k < (ll)x.size()) {
		vector<mint> xx(k*2);
		for (ll i = 0; i < k*2; i++) xx[i] = ((i < x.size()) ? x[i] : 0);

		y = 2 * y - xx * y * y;
		y.resize(k*2);
		k = k*2;
	}
	y.resize(x.size());
	return y;
}


void solve() {
	const mint inv2 = ((mint)2).inv();
	const vector<mint> inv2inv2 = (vector<mint>){inv2, inv2};

	vector<P> c(n*2);
	for (ll i = 0; i < n; i++) {
		c[i*2+0] = {a[i], i*2+0};
		c[i*2+1] = {b[i], i*2+1};
	}
	sort(c.begin(), c.end());

	vector<vector<mint>> cl(n*2), cm(n*2), cr(n*2);
	for (ll i = 0; i < n*2; i++) {
		ll v = c[i].second / 2;
		if (c[i].second % 2 == 0) {
			cl[i] = {inv2 * inv2, inv2, inv2 * inv2};
			cm[i] = {inv2 * a[v] * inv2, inv2 * a[v] * inv2};
			cr[i] = {inv2, inv2};
		} else {
			cl[i] = {0, 1};
			cm[i] = {inv2 * b[v]};
			cr[i] = {inv2, inv2};
		}
	}

	using Fpsdata = struct {
		vector<mint> pl;
		vector<mint> pr;
		vector<mint> sumlmr;
	};
	const auto findlmr = [&](this auto self, const ll l, const ll r) -> Fpsdata {
		if (r-l == 1) {
			Fpsdata ret = {
				pl: cl[l],
				pr: cr[l],
				sumlmr: cm[l]
			};
			return ret;
		}

		ll med = (l + r) / 2;
		Fpsdata fl = self(l, med);
		Fpsdata fr = self(med, r);

		Fpsdata ret = {
			pl: fl.pl * fr.pl,
			pr: fl.pr * fr.pr,
			sumlmr: fl.sumlmr * fr.pr + fl.pl * fr.sumlmr
		};
		return ret;
	};
	Fpsdata wholef = findlmr(0, n*2);

	vector<mint> inv2inv2pow2n(n*2+1);
	mint n2cr = 1;
	for (ll i = 0; i <= n*2; i++) {
		inv2inv2pow2n[i] = inv2.pow(i) * n2cr * inv2.pow(n*2-i);
		n2cr *= (n*2 - i);
		n2cr /= (i + 1);
	}
	// inv2inv2pow2n.resize(n+1);

	vector<mint> xl = wholef.sumlmr;
	xl.resize(n+1);
	vector<mint> xr = fpsinv(inv2inv2pow2n);
	xr.resize(n+1);
	vector<mint> ans = xl * xr;
	const mint pown = ((mint)2).pow(n);
	for (ll i = 0; i < n; i++) {
		cout << (ans[i] * pown).val() << "\n";
	}
}

int main (void) {
	std::cin.tie(nullptr);
	std::ios_base::sync_with_stdio(false);

	cin >> n;
	a.resize(n);
	b.resize(n);
	for (ll i = 0; i < n; i++) {
		cin >> a[i] >> b[i];
	}

	solve();

	
	return 0;
}