不妨认为 \(n,q\) 同阶。
考虑根号重构。如果没有第 3 种操作的话，我们每 \(\mathcal O(\sqrt n)\) 操作整体更新一次，每个询问只需要考虑块内的修改所在置换环上有多少 \([l,r]\) 内的数。这个是容易 \(\mathcal O(n\sqrt n)\) 做的。
然后考虑置换环会变怎么办。注意到一个块内真正会变的 \(p\) 值只有 \(\mathcal O(\sqrt n)\) 个，标记这些会变的位置，剩下的位置会形成若干条链，它们本身不会发生改变之类，所以我们直接把整条链缩到一起。这样一个操作块内就只剩 \(\mathcal O(\sqrt n)\) 个点了。
考虑如何统计答案。我们需要的是一个置换环上 \([l,r]\) 内的点的个数。把询问拆成 \([1,r]\) 减去 \([1,l-1]\)，挂到每个缩起来的点上。顺序加入一遍 \(1\sim n\)，每次加入把其对应的缩点 +1。其实就是一个扫描线状物。
时间复杂度 \(\mathcal O(n\sqrt n)\)。我实现的不好，需要卡常。我使用了简易火车头。注意到常数大在缩起来的点数还是有点多，其它地方顺序访问常数不大，所以适当取小块长即可。还有交换数组两维使其顺序访问的优化，我还没加上就过了。
#pragma GCC optimize("Ofast")
#pragma GCC optimize("unroll-loops")
#pragma GCC target("sse,sse2,sse3,ssse3,sse4,popcnt,abm,mmx,avx,avx2,tune=native")
#include <bits/stdc++.h>
#define pb emplace_back
#define fir first
#define sec second
using i64 = long long;
using pii = std::pair<int, int>;
constexpr int maxn = 2e5 + 5, S = 1000 + 5;
int n, p[maxn], bs, q, idx[maxn], cnt, c[maxn], tr[maxn], tot, ed[maxn], coef[S][S], buc[S];
i64 sum[maxn], ans[maxn], res[maxn], a[maxn];
std::vector<std::tuple<int, int, int>> opt, qry;
void solve() {
cnt = tot = 0;
std::fill(idx + 1, idx + 1 + n, 0);
std::fill(ed + 1, ed + 1 + n, 0);
for (auto& [op, x, y] : opt) {
if (op >= 2 && !idx[x]) idx[x] = ++cnt, ed[cnt] = x;
if (op == 3 && !idx[y]) idx[y] = ++cnt, ed[cnt] = y;
}
std::vector<int> st;
for (auto& [op, x, y] : opt) {
if (op >= 2) st.pb(p[x]);
if (op == 3) st.pb(p[y]);
}
for (auto& x : st) {
if (idx[x]) continue ;
++cnt; int y = x;
for (; !idx[y]; y = p[y])
idx[y] = cnt, ed[cnt] = y;
}
qry.clear();
for (auto& [op, x, y] : opt) {
if (op == 1) {
ans[++tot] = sum[y] - sum[x - 1];
qry.pb(x - 1, -1, tot);
qry.pb(y, 1, tot);
}
}
std::sort(qry.begin(), qry.end()); int j = 1;
std::fill(buc + 1, buc + 1 + S, 0);
std::fill(res + 1, res + 1 + S, 0);
for (int i = 1; i <= cnt; ++i)
for (int j = 1; j <= tot; ++j)
coef[i][j] = 0;
for (auto& [x, v, id] : qry) {
while (j <= x) {
++buc[idx[j]];
++j;
}
for (int i = 1; i <= cnt; ++i)
coef[i][id] += v * buc[i];
}
tot = 0;
for (auto& [op, x, y] : opt) {
if (op == 1) {
++tot;
for (int i = 1; i <= cnt; ++i)
ans[tot] += (i64)coef[i][tot] * res[i];
} else if (op == 2) {
int z = idx[x], c = 0;
while (true) {
c += z == idx[x];
if (c > 1) break ;
res[z] += y;
z = idx[p[ed[z]]];
}
} else {
std::swap(p[x], p[y]);
}
}
for (int i = 1; i <= tot; ++i)
printf("%lld\n", ans[i]);
for (int i = 1; i <= n; ++i)
a[i] += res[idx[i]], sum[i] = sum[i - 1] + a[i];
return ;
}
int main() {
scanf("%d", &n);
for (int i = 1; i <= n; ++i)
scanf("%lld", &a[i]), sum[i] = sum[i - 1] + a[i];
for (int i = 1; i <= n; ++i)
scanf("%d", &p[i]);
scanf("%d", &q); bs = sqrt(q); bs = std::max(bs / 2, 1);
for (int i = 1; i <= q; ++i) {
int op, l, r; scanf("%d %d %d", &op, &l, &r);
opt.pb(op, l, r);
if (i % bs == 0 || i == q) solve(), opt.clear();
}
return 0;
}