結果
| 問題 |
No.2443 特殊線形群の標準表現
|
| コンテスト | |
| ユーザー |
tonegawa
|
| 提出日時 | 2023-08-25 22:02:33 |
| 言語 | C++17 (gcc 13.3.0 + boost 1.87.0) |
| 結果 |
WA
|
| 実行時間 | - |
| コード長 | 39,924 bytes |
| コンパイル時間 | 1,516 ms |
| コンパイル使用メモリ | 151,132 KB |
| 最終ジャッジ日時 | 2025-02-16 14:01:10 |
|
ジャッジサーバーID (参考情報) |
judge3 / judge5 |
(要ログイン)
| ファイルパターン | 結果 |
|---|---|
| other | AC * 6 WA * 15 |
ソースコード
#line 1 ".lib/template.hpp"
#include <iostream>
#include <string>
#include <vector>
#include <array>
#include <tuple>
#include <stack>
#include <queue>
#include <deque>
#include <algorithm>
#include <set>
#include <map>
#include <unordered_set>
#include <unordered_map>
#include <bitset>
#include <cmath>
#include <functional>
#include <cassert>
#include <climits>
#include <iomanip>
#include <numeric>
#include <memory>
#include <random>
#include <thread>
#include <chrono>
#define allof(obj) (obj).begin(), (obj).end()
#define range(i, l, r) for(int i=l;i<r;i++)
#define bit_subset(i, S) for(int i=S, zero_cnt=0;(zero_cnt+=i==S)<2;i=(i-1)&S)
#define bit_kpop(i, n, k) for(int i=(1<<k)-1,x_bit,y_bit;i<(1<<n);x_bit=(i&-i),y_bit=i+x_bit,i=(!i?(1<<n):((i&~y_bit)/x_bit>>1)|y_bit))
#define bit_kth(i, k) ((i >> k)&1)
#define bit_highest(i) (i?63-__builtin_clzll(i):-1)
#define bit_lowest(i) (i?__builtin_ctzll(i):-1)
#define sleepms(t) std::this_thread::sleep_for(std::chrono::milliseconds(t))
using ll = long long;
using ld = long double;
using ul = uint64_t;
using pi = std::pair<int, int>;
using pl = std::pair<ll, ll>;
using namespace std;
template<typename F, typename S>
std::ostream &operator<<(std::ostream &dest, const std::pair<F, S> &p){
dest << p.first << ' ' << p.second;
return dest;
}
template<typename T>
std::ostream &operator<<(std::ostream &dest, const std::vector<std::vector<T>> &v){
int sz = v.size();
if(sz==0) return dest;
for(int i=0;i<sz;i++){
int m = v[i].size();
for(int j=0;j<m;j++) dest << v[i][j] << (i!=sz-1&&j==m-1?'\n':' ');
}
return dest;
}
template<typename T>
std::ostream &operator<<(std::ostream &dest, const std::vector<T> &v){
int sz = v.size();
if(sz==0) return dest;
for(int i=0;i<sz-1;i++) dest << v[i] << ' ';
dest << v[sz-1];
return dest;
}
template<typename T, size_t sz>
std::ostream &operator<<(std::ostream &dest, const std::array<T, sz> &v){
if(sz==0) return dest;
for(int i=0;i<sz-1;i++) dest << v[i] << ' ';
dest << v[sz-1];
return dest;
}
template<typename T>
std::ostream &operator<<(std::ostream &dest, const std::set<T> &v){
for(auto itr=v.begin();itr!=v.end();){
dest << *itr;
itr++;
if(itr!=v.end()) dest << ' ';
}
return dest;
}
template<typename T, typename E>
std::ostream &operator<<(std::ostream &dest, const std::map<T, E> &v){
for(auto itr=v.begin();itr!=v.end();){
dest << '(' << itr->first << ", " << itr->second << ')';
itr++;
if(itr!=v.end()) dest << '\n';
}
return dest;
}
template<typename T>
vector<T> make_vec(size_t sz, T val){return std::vector<T>(sz, val);}
template<typename T, typename... Tail>
auto make_vec(size_t sz, Tail ...tail){
return std::vector<decltype(make_vec<T>(tail...))>(sz, make_vec<T>(tail...));
}
template<typename T>
vector<T> read_vec(size_t sz){
std::vector<T> v(sz);
for(int i=0;i<(int)sz;i++) std::cin >> v[i];
return v;
}
template<typename T, typename... Tail>
auto read_vec(size_t sz, Tail ...tail){
auto v = std::vector<decltype(read_vec<T>(tail...))>(sz);
for(int i=0;i<(int)sz;i++) v[i] = read_vec<T>(tail...);
return v;
}
void io_init(){
std::cin.tie(nullptr);
std::ios::sync_with_stdio(false);
}
#line 1 ".lib/math/mod.hpp"
#line 6 ".lib/math/mod.hpp"
#include <type_traits>
#line 8 ".lib/math/mod.hpp"
#include <ostream>
#line 1 ".lib/math/minior/mod_base.hpp"
#line 4 ".lib/math/minior/mod_base.hpp"
// @param m `1 <= m`
constexpr long long safe_mod(long long x, long long m){
x %= m;
if (x < 0) x += m;
return x;
}
struct barrett{
unsigned int _m;
unsigned long long im;
explicit barrett(unsigned int m) : _m(m), im((unsigned long long)(-1) / m + 1){}
unsigned int umod()const{return _m;}
unsigned int mul(unsigned int a, unsigned int b)const{
unsigned long long z = a;
z *= b;
#ifdef _MSC_VER
unsigned long long x;
_umul128(z, im, &x);
#else
unsigned long long x = (unsigned long long)(((unsigned __int128)(z) * im) >> 64);
#endif
unsigned long long y = x * _m;
return (unsigned int)(z - y + (z < y ? _m : 0));
}
};
// @param n `0 <= n`
// @param m `1 <= m`
constexpr long long pow_mod_constexpr(long long x, long long n, int m){
if(m == 1) return 0;
unsigned int _m = (unsigned int)(m);
unsigned long long r = 1;
unsigned long long y = safe_mod(x, m);
while(n){
if (n & 1) r = (r * y) % _m;
y = (y * y) % _m;
n >>= 1;
}
return r;
}
constexpr bool is_prime_constexpr(int n) {
if (n <= 1) return false;
if (n == 2 || n == 7 || n == 61) return true;
if (n % 2 == 0) return false;
long long d = n - 1;
while (d % 2 == 0) d /= 2;
constexpr long long bases[3] = {2, 7, 61};
for(long long a : bases){
long long t = d;
long long y = pow_mod_constexpr(a, t, n);
while(t != n - 1 && y != 1 && y != n - 1){
y = y * y % n;
t <<= 1;
}
if(y != n - 1 && t % 2 == 0){
return false;
}
}
return true;
}
template<int n>
constexpr bool is_prime = is_prime_constexpr(n);
constexpr int primitive_root_constexpr(int m){
if(m == 2) return 1;
if(m == 167772161) return 3;
if(m == 469762049) return 3;
if(m == 754974721) return 11;
if(m == 998244353) return 3;
int divs[20] = {};
divs[0] = 2;
int cnt = 1;
int x = (m - 1) / 2;
while (x % 2 == 0) x /= 2;
for(int i = 3; (long long)(i)*i <= x; i += 2){
if(x % i == 0){
divs[cnt++] = i;
while(x % i == 0){
x /= i;
}
}
}
if(x > 1) divs[cnt++] = x;
for(int g = 2;; g++){
bool ok = true;
for(int i = 0; i < cnt; i++){
if(pow_mod_constexpr(g, (m - 1) / divs[i], m) == 1){
ok = false;
break;
}
}
if(ok)return g;
}
}
template <int m>
constexpr int primitive_root = primitive_root_constexpr(m);
int ceil_pow2(int n){
int x = 0;
while ((1U << x) < (unsigned int)(n)) x++;
return x;
}
int bsf(unsigned int n){
return __builtin_ctz(n);
}
// @param b `1 <= b`
// @return pair(g, x) s.t. g = gcd(a, b), xa = g (mod b), 0 <= x < b/g
constexpr std::pair<long long, long long> inv_gcd(long long a, long long b){
a = safe_mod(a, b);
if(a == 0) return {b, 0};
long long s = b, t = a;
long long m0 = 0, m1 = 1;
while (t){
long long u = s / t;
s -= t * u;
m0 -= m1 * u;
auto tmp = s;
s = t;
t = tmp;
tmp = m0;
m0 = m1;
m1 = tmp;
}
if(m0 < 0) m0 += b / s;
return {s, m0};
}
#line 13 ".lib/math/mod.hpp"
template<int m>
long long modpow(long long a, long long b){
assert(0 <= b);
assert(0 < m);
a = safe_mod(a, m);
long long ret = 1;
while(b){
if(b & 1) ret = (ret * a) % m;
a = (a * a) % m;
b >>= 1;
}
return ret;
}
// @param 0 <= b, 0 < m
long long modpow(long long a, long long b, int m){
assert(0 <= b);
assert(0 < m);
a = safe_mod(a, m);
long long ret = 1;
while(b){
if(b & 1) ret = (ret * a) % m;
a = (a * a) % m;
b >>= 1;
}
return ret;
}
struct modint_base {};
struct static_modint_base : modint_base {};
template <int m, std::enable_if_t<(1 <= m)>* = nullptr>
struct static_modint : static_modint_base{
using mint = static_modint;
public:
static constexpr int mod(){return m;}
static mint raw(int v) {
mint x;
x._v = v;
return x;
}
static_modint(): _v(0){}
template <class T>
static_modint(T v){
long long x = v % (long long)umod();
if (x < 0) x += umod();
_v = x;
}
unsigned int val()const{return _v;}
mint& operator++(){
_v++;
if (_v == umod()) _v = 0;
return *this;
}
mint& operator--(){
if (_v == 0) _v = umod();
_v--;
return *this;
}
mint operator++(int){
mint result = *this;
++*this;
return result;
}
mint operator--(int){
mint result = *this;
--*this;
return result;
}
mint& operator+=(const mint& rhs){
_v += rhs._v;
if (_v >= umod()) _v -= umod();
return *this;
}
mint& operator-=(const mint& rhs){
_v -= rhs._v;
if (_v >= umod()) _v += umod();
return *this;
}
mint& operator*=(const mint& rhs){
unsigned long long z = _v;
z *= rhs._v;
_v = (unsigned int)(z % umod());
return *this;
}
mint& operator/=(const mint& rhs){return *this = *this * rhs.inv();}
mint operator+()const{return *this;}
mint operator-()const{return mint() - *this;}
mint pow(long long n)const{
assert(0 <= n);
mint x = *this, r = 1;
while(n){
if (n & 1) r *= x;
x *= x;
n >>= 1;
}
return r;
}
mint inv()const{
if(prime){
assert(_v);
return pow(umod() - 2);
}else{
auto eg = inv_gcd(_v, m);
assert(eg.first == 1);
return eg.second;
}
}
friend mint operator+(const mint& lhs, const mint& rhs){return mint(lhs) += rhs;}
friend mint operator-(const mint& lhs, const mint& rhs){return mint(lhs) -= rhs;}
friend mint operator*(const mint& lhs, const mint& rhs){return mint(lhs) *= rhs;}
friend mint operator/(const mint& lhs, const mint& rhs){return mint(lhs) /= rhs;}
friend bool operator==(const mint& lhs, const mint& rhs){return lhs._v == rhs._v;}
friend bool operator!=(const mint& lhs, const mint& rhs){return lhs._v != rhs._v;}
private:
unsigned int _v;
static constexpr unsigned int umod(){return m;}
static constexpr bool prime = is_prime<m>;
};
template<int id>
struct dynamic_modint : modint_base{
using mint = dynamic_modint;
public:
static int mod(){return (int)(bt.umod());}
static void set_mod(int m){
assert(1 <= m);
bt = barrett(m);
}
static mint raw(int v){
mint x;
x._v = v;
return x;
}
dynamic_modint(): _v(0){}
template <class T>
dynamic_modint(T v){
long long x = v % (long long)(mod());
if (x < 0) x += mod();
_v = x;
}
unsigned int val()const{return _v;}
mint& operator++(){
_v++;
if(_v == umod()) _v = 0;
return *this;
}
mint& operator--(){
if (_v == 0) _v = umod();
_v--;
return *this;
}
mint operator++(int){
mint result = *this;
++*this;
return result;
}
mint operator--(int){
mint result = *this;
--*this;
return result;
}
mint& operator+=(const mint& rhs){
_v += rhs._v;
if(_v >= umod()) _v -= umod();
return *this;
}
mint& operator-=(const mint& rhs){
_v += mod() - rhs._v;
if(_v >= umod()) _v -= umod();
return *this;
}
mint& operator*=(const mint& rhs){
_v = bt.mul(_v, rhs._v);
return *this;
}
mint& operator/=(const mint& rhs){return *this = *this * rhs.inv();}
mint operator+()const{return *this;}
mint operator-()const{return mint() - *this;}
mint pow(long long n)const{
assert(0 <= n);
mint x = *this, r = 1;
while(n){
if (n & 1) r *= x;
x *= x;
n >>= 1;
}
return r;
}
mint inv()const{
auto eg = inv_gcd(_v, mod());
assert(eg.first == 1);
return eg.second;
}
friend mint operator+(const mint& lhs, const mint& rhs){return mint(lhs) += rhs;}
friend mint operator-(const mint& lhs, const mint& rhs){return mint(lhs) -= rhs;}
friend mint operator*(const mint& lhs, const mint& rhs){return mint(lhs) *= rhs;}
friend mint operator/(const mint& lhs, const mint& rhs){return mint(lhs) /= rhs;}
friend bool operator==(const mint& lhs, const mint& rhs){return lhs._v == rhs._v;}
friend bool operator!=(const mint& lhs, const mint& rhs){return lhs._v != rhs._v;}
private:
unsigned int _v;
static barrett bt;
static unsigned int umod(){return bt.umod();}
};
template <int id>
barrett dynamic_modint<id>::bt(998244353);
using modint = dynamic_modint<-1>;
using modint998244353 = static_modint<998244353>;
using modint1000000007 = static_modint<1000000007>;
template <class T>
using is_modint = std::is_base_of<modint_base, T>;
template <class T>
using is_modint_t = std::enable_if_t<is_modint<T>::value>;
template <class T>
using is_static_modint = std::is_base_of<static_modint_base, T>;
template <class T>
using is_static_modint_t = std::enable_if_t<is_static_modint<T>::value>;
template <class> struct is_dynamic_modint : public std::false_type {};
template <int id>
struct is_dynamic_modint<dynamic_modint<id>> : public std::true_type {};
template <class T>
using is_dynamic_modint_t = std::enable_if_t<is_dynamic_modint<T>::value>;
template<int m>
std::ostream &operator<<(std::ostream &dest, const static_modint<m> &a){
dest << a.val();
return dest;
}
template<int id>
std::ostream &operator<<(std::ostream &dest, const dynamic_modint<id> &a){
dest << a.val();
return dest;
}
// 0 <= n < m <= int_max
// 前処理 O(n + log(m))
// 各種計算 O(1)
// 変数 <= n
template<typename mint, is_modint<mint>* = nullptr>
struct modcomb{
private:
int n;
std::vector<mint> f, i, fi;
void init(int _n){
assert(0 <= _n && _n < mint::mod());
if(_n < f.size()) return;
n = _n;
f.resize(n + 1), i.resize(n + 1), fi.resize(n + 1);
f[0] = fi[0] = mint(1);
if(n) f[1] = fi[1] = i[1] = mint(1);
for(int j = 2; j <= n; j++) f[j] = f[j - 1] * j;
fi[n] = f[n].inv();
for(int j = n; j >= 2; j--){
fi[j - 1] = fi[j] * j;
i[j] = f[j - 1] * fi[j];
}
}
public:
modcomb(): n(-1){}
modcomb(int _n){
init(_n);
}
void recalc(int _n){
init(std::min(mint::mod() - 1, 1 << ceil_pow2(_n)));
}
mint comb(int a, int b){
if((a < 0) || (b < 0) || (a < b)) return 0;
return f[a] * fi[a - b] * fi[b];
}
mint perm(int a, int b){
if((a < 0) || (b < 0) || (a < b)) return 0;
return f[a] * fi[a - b];
}
mint fac(int x){
assert(0 <= x && x <= n);
return f[x];
}
mint inv(int x){
assert(0 < x && x <= n);
return i[x];
}
mint finv(int x){
assert(0 <= x && x <= n);
return fi[x];
}
};
template<typename mint, is_modint<mint>* = nullptr>
struct modpow_table{
std::vector<mint> v;
// x^maxkまで計算できる
modpow_table(){}
void init(int x, int maxk){
v.resize(maxk + 1);
v[0] = 1;
for(int i = 1; i <= maxk; i++) v[i] = v[i - 1] * x;
}
mint pow(int k){
assert(0 <= k && k < v.size());
return v[k];
}
};
#line 1 ".lib/math/matrix/matrix_mod.hpp"
#line 5 ".lib/math/matrix/matrix_mod.hpp"
template<typename mint>
struct matrix_mod{
int n, m;
using _mint = mint;
private:
using vec = std::vector<mint>;
using matrix = matrix_mod<mint>;
// n × k 行列と k × m 行列の積(n × m行列)
// K == 0だと壊れる
static matrix __mul_mat(const matrix &vl, const matrix &vr){
int N = vl.n, K = vl.m, M = vr.m;
assert(K == vr.n);
assert(K);
if(N == 0) return matrix(0, M);
if(M == 0) return matrix(N, 0);
auto vr_t = vr.t();
matrix ret(N, M, 0);
for(int i = 0; i < N; i++){
for(int j = 0; j < M; j++){
__int128_t S = 0;
for(int k = 0; k < K; k++){
S += (long long)vl.val[i][k].val() * vr_t[j][k].val();
}
ret[i][j] = S % mint::mod();
}
}
return ret;
}
// n × m 行列と n × m 行列の和(n × m行列)
static void __add_mat_inplace(matrix &vl, const matrix &vr){
assert(vl.n == vr.n && vl.m == vr.m);
int N = vl.n, M = vl.m;
for(int i = 0; i < N; i++){
for(int j = 0; j < M; j++){
vl[i][j] += vr[i][j];
}
}
}
// n × m 行列と n × m 行列の差(n × m行列)
static void __sub_mat_inplace(matrix &vl, const matrix &vr){
assert(vl.n == vr.n && vl.m == vr.m);
int N = vl.n, M = vl.m;
for(int i = 0; i < N; i++){
for(int j = 0; j < M; j++){
vl[i][j] -= vr[i][j];
}
}
}
static void __mul_val_inplace(matrix &vl, mint vr){
int N = vl.n, M = vl.m;
for(int i = 0; i < N; i++){
for(int j = 0; j < M; j++){
vl[i][j] *= vr;
}
}
}
static void __add_val_inplace(matrix &vl, mint vr){
int N = vl.n, M = vl.m;
for(int i = 0; i < N; i++){
for(int j = 0; j < M; j++){
vl[i][j] += vr;
}
}
}
static void __sub_val_inplace(matrix &vl, mint vr){
int N = vl.n, M = vl.m;
for(int i = 0; i < N; i++){
for(int j = 0; j < M; j++){
vl[i][j] -= vr;
}
}
}
std::vector<vec> val;
public:
matrix_mod(): n(0), m(0){}
matrix_mod(int _n, int _m, mint x = mint(0)) : n(_n), m(_m), val(_n, vec(_m, x)){}
matrix_mod(const matrix_mod &v) : n(v.n), m(v.m), val(v.val){}
matrix_mod(const vec &v): n(1), m(v.size()), val(1, vec(v.size())){val[0] = v;}
matrix_mod(const std::vector<vec> &v): n(v.size()), m(v[0].size()), val(v){}
matrix_mod operator + (const matrix_mod &vr){matrix_mod tmp(*this); return tmp += vr;}
matrix_mod operator - (const matrix_mod &vr){matrix_mod tmp(*this); return tmp -= vr;}
matrix_mod operator * (const matrix_mod &vr){return __mul_mat(*this, vr);}
matrix_mod operator ^ (const long long vr){return pow(vr);}
matrix_mod operator * (const mint vr){matrix_mod tmp(*this); return tmp *= vr;}
matrix_mod operator += (const matrix_mod &vr){__add_mat_inplace(*this, vr); return *this;}
matrix_mod operator -= (const matrix_mod &vr){__sub_mat_inplace(*this, vr); return *this;}
matrix_mod operator *= (const matrix_mod &vr){return (*this) = __mul_mat(*this, vr);}
matrix_mod operator ^= (const long long vr){return (*this) = pow(vr);}
matrix_mod operator *= (const mint vr){__mul_val_inplace(*this, vr); return *this;}
vec& operator [] (const int i){return val[i];}
// n次の単位行列
static matrix_mod eye(int n){
matrix_mod ret(n, n, 0);
for(int i = 0; i < n; i++) ret[i][i] = mint(1);
return ret;
}
void print()const{
for(int i = 0; i < n; i++){
for(int j = 0; j < m; j++){
std::cout << val[i][j] << (j == m - 1 ? '\n' : ' ');
}
}
}
matrix_mod pow(long long k){
assert(n && m && n == m); // 正方行列でなければならない
matrix_mod ret = eye(n); // k == 0の場合単位行列を返す
matrix_mod m(*this);
while(k){
if(k & 1) ret *= m;
m *= m;
k >>= 1;
}
return ret;
}
// 転置
matrix_mod t()const{
matrix_mod ret(m, n, 0);
for(int i = 0; i < n; i++){
for(int j = 0; j < m; j++){
ret[j][i] = val[i][j];
}
}
return ret;
}
//掃き出し法で上三角行列を作る, {変形後の行列、行のスワップ回数}を返す O(NM^2)
std::pair<matrix_mod, int> gaussian_elimination(){
matrix_mod v(*this);
int row = 0;//確定していない行
int swp = 0;
for(int i = 0; i < m && row < n; i++){
//i列目が0でない行を探す
int r = -1;
for(int j = row; j < n; j++){
if(v[j][i].val()){
r = j;
break;
}
}
if(r == -1) continue;
if(r != row){
swp++;
std::swap(v[r], v[row]);
}
//i列目が0でない行の処理
for(int j = row + 1; j < n; j++){
if(v[j][i].val() == 0) continue;
mint x = v[j][i] / v[row][i];
for(int k = i; k < m; k++){
v[j][k] -= x * v[row][k];
}
}
row++;
}
return {v, swp};
}
//掃き出し法で上三角行列を作る, {変形後の行列、行のスワップ回数}を返す O(NM^2 * log mod)
std::pair<matrix_mod, int> gaussian_elimination_arbitrary_mod(){
matrix_mod v(*this);
int row = 0;//確定していない行
int swp = 0;
for(int i = 0; i < m && row < n; i++){
//i列目が0でない行を探す
int r = -1;
for(int j = row; j < n; j++){
if(v[j][i].val()){
r = j;
break;
}
}
if(r == -1) continue;
if(r != row){
swp++;
std::swap(v[r], v[row]);
}
//i列目が0でない行の処理
for(int j = row + 1; j < n; j++){
while(v[j][i].val() != 0){
if(v[row][i].val() > v[j][i].val()){
swp++;
std::swap(v[row], v[j]);
}
int x = v[j][i].val() / v[row][i].val();
for(int k = i; k < m; k++){
v[j][k] -= x * v[row][k];
}
}
}
row++;
}
return {v, swp};
}
//すでに上三角行列になっていることが前提
int rank(){
int cnt = 0;
for(int i = 0; i < n; i++, cnt++){
bool f = false;
for(int j = i; j < m; j++){
if(val[i][j].val()){
f = true;
break;
}
}
if(!f) break;
}
return cnt;
}
// 行列式 O(N^3)
mint det(){
assert(n == m); // 正方行列のみ
auto [tmp, swp] = gaussian_elimination();
mint res(1);
for(int i = 0; i < n; i++) res *= tmp[i][i];
return swp & 1 ? -res : res;
}
// 行列式 O(N^3 * log mod)
mint det_arbitrary_mod(){
assert(n == m); // 正方行列のみ
auto [tmp, swp] = gaussian_elimination_arbitrary_mod();
mint res(1);
for(int i = 0; i < n; i++) res *= tmp[i][i];
return swp & 1 ? -res : res;
}
// (n, m) + (n, l) -> (n, m + l) 横に結合
matrix_mod concat_horizontal(matrix_mod vr){
assert(n == vr.n);
matrix_mod res(*this);
for(int i = 0; i < n; i++){
res[i].insert(res[i].end(), vr[i].begin(), vr[i].end());
}
res.m += vr.m;
return res;
}
// (n, m) + (l, m) -> (n + l, m) 縦に結合
matrix_mod concat_vertical(matrix_mod vr){
assert(m == vr.m);
matrix_mod res(*this);
for(int i = 0; i < vr.n; i++) res.val.push_back(vr[i]);
res.n += vr.n;
return res;
}
// (n, m) -> (n, k), (n, m - k)
std::pair<matrix_mod, matrix_mod> split_horizontal(int k){
assert(0 <= k && k <= m);
matrix_mod a(n, k), b(n, m - k);
for(int i = 0; i < n; i++){
for(int j = 0; j < k; j++){
a[i][j] = val[i][j];
}
}
for(int i = 0; i < n; i++){
for(int j = 0; j < m - k; j++){
b[i][j] = val[i][j + k];
}
}
return {a, b};
}
// (n, m) -> (k, m), (n - k, m)
std::pair<matrix_mod, matrix_mod> split_vertical(int k){
assert(0 <= k && k <= n);
matrix_mod a(k, m), b(n - k, m);
for(int i = 0; i < k; i++){
for(int j = 0; j < m; j++){
a[i][j] = val[i][j];
}
}
for(int i = 0; i < n - k; i++){
for(int j = 0; j < m; j++){
b[i][j] = val[k + i][j];
}
}
return {a, b};
}
matrix_mod inv(){
assert(n == m);
auto [tmp, swp] = concat_horizontal(eye(n)).gaussian_elimination();
for(int i = 0; i < n; i++){
mint x = tmp[i][i];
if(!x.val()) return matrix_mod{};// 存在しない
x = x.inv();
for(int j = 0; j < 2 * n; j++) tmp[i][j] *= x;
}
for(int i = n - 1; i >= 0; i--){
for(int j = i + 1; j < n; j++){
if(!tmp[i][j].val()) continue;
mint c = tmp[i][j];
for(int k = j; k < 2 * n; k++){
tmp[i][k] -= c * tmp[j][k];
}
}
}
return tmp.split_horizontal(n).second;
}
// https://ja.wikipedia.org/wiki/LU%E5%88%86%E8%A7%A3
// キャッシュのためにuを転置して実装
std::pair<matrix_mod, matrix_mod> lu_decomposition(){
matrix_mod l = eye(n), u(n, n);
for(int i = 0; i < n; i++){
// u[i][i]を決定
u[i][i] = val[i][i];
for(int j = 0; j < i; j++) u[i][i] -= l[i][j] * u[i][j];
if(u[i][i].val() == 0) return {matrix_mod{}, matrix_mod{}}; // 不可能
mint iuii = u[i][i].inv();
// l[0, n)[i]を決定
for(int j = i + 1; j < n; j++){
l[j][i] = val[j][i];
for(int k = 0; k < i; k++) l[j][i] -= l[j][k] * u[i][k];
l[j][i] *= iuii;
}
// u[i][0, n)を決定
for(int j = i + 1; j < n; j++){
u[j][i] = val[i][j];
for(int k = 0; k < i; k++) u[j][i] -= l[i][k] * u[j][k];
}
}
u = u.t();
return {l, u};
}
// Ax = b
// (n, m) * (m, 1) -> (n, 1)
// を満たす連立方程式を解く, 解空間の次元、(rank*変数)の基底を返す
// 解空間の基底は任意のt_iについてA * (v1t1 + v2t2 ...) = 0を満たす
// つまり plus + res[0]t_0 + res[1]t_1 + res[2]t_2...は全て解を満たす
// 解が存在しない場合解空間の次元として-1を返す
std::tuple<int, matrix_mod, vec> system_of_linear_equations(const vec &vr){
assert(vr.size() == n);
matrix_mod tmp = concat_horizontal(matrix_mod(vr).t()).gaussian_elimination().first;
//解空間の次元 = 変数の数 - 階数
int r = tmp.rank();
std::vector<int> fc(r, -1);//各行に初めて非零要素が現れる列
for(int i = 0; i < r; i++){
mint tmp_inv;
bool f = false;
for(int j = i; j < tmp.m; j++){
if(tmp[i][j].val() == 0) continue;
if(j == tmp.m - 1 && !f){
return {-1, matrix_mod{}, vec{}}; // 解なし
}
if(!f){
tmp_inv = tmp[i][j].inv();
fc[i] = j;
f = true;
}
tmp[i][j] = tmp[i][j] * tmp_inv;
}
}
int d = tmp.m - 1 - r, v = tmp.m - 1;
vec plus(v, 0);
for(int i = r - 1; i >= 0; i--){
int idx = fc[i];
assert(idx != -1);
plus[idx] = tmp[i][v];
for(int j = idx + 1; j < v; j++){
plus[idx] -= plus[j] * tmp[i][j];
}
}
matrix_mod res(d, v, 0);
std::vector<bool> not_fc(v, true);
for(int i = 0; i < r; i++) not_fc[fc[i]] = false;
for(int i = 0, j = 0; i < v; i++) if(not_fc[i]) res[j++][i] = 1;
for(int i = r - 1; i >= 0; i--){ //各行に1つまだ確定していない変数が現れる
int col = fc[i];
assert(col != -1);
assert(tmp[i][col].val() == 1);
for(int k = 0; k < d; k++){ // 次元
for(int j = col + 1; j < v; j++){ // すでに確定した要素
res[k][col] -= res[k][j] * tmp[i][j];
}
}
}
return {d, res, plus};
}
};
#line 1 ".lib/data_structure/segment_tree/segment_tree.hpp"
#line 1 ".lib/monoid.hpp"
#include <limits>
#line 8 ".lib/monoid.hpp"
struct point_min_range_min{
template<typename T>
static T id(){
return std::numeric_limits<T>::max();
}
template<typename T>
static T update(T a, T b){
return std::min(a, b);
}
template<typename T>
static T merge(T a, T b){
return std::min(a, b);
}
};
struct point_min_range_second_min{
template<typename T>
static T id(){
return {std::numeric_limits<long long>::max(), std::numeric_limits<long long>::max()};
}
template<typename T>
static T update(T a, T b){
if(a.first <= b.first) return {a.first, std::min(a.second, b.first)};
else return {b.first, std::min(a.first, b.second)};
}
template<typename T>
static T merge(T a, T b){
if(a.first <= b.first) return {a.first, std::min(a.second, b.first)};
else return {b.first, std::min(a.first, b.second)};
}
};
struct point_max_range_max{
template<typename T>
static T id(){
return std::numeric_limits<T>::min();
}
template<typename T>
static T update(T a, T b){
return std::max(a, b);
}
template<typename T>
static T merge(T a, T b){
return std::max(a, b);
}
template<typename T>
static T flip(T a){
return a;
}
};
struct point_max_range_second_max{
template<typename T>
static T id(){
return {std::numeric_limits<long long>::min(), std::numeric_limits<long long>::min()};
}
template<typename T>
static T update(T a, T b){
if(a.first >= b.first) return {a.first, std::min(a.second, b.first)};
else return {b.first, std::min(a.first, b.second)};
}
template<typename T>
static T merge(T a, T b){
if(a.first >= b.first) return {a.first, std::min(a.second, b.first)};
else return {b.first, std::min(a.first, b.second)};
}
};
struct point_mul_range_mul{
template<typename T>
static T id(){
return 1;
}
template<typename T>
static T update(T a, T b){
return a * b;
}
template<typename T>
static T merge(T a, T b){
return a * b;
}
};
struct point_add_range_min{
template<typename T>
static T id(){
return std::numeric_limits<T>::max();
}
template<typename T>
static T update(T a, T b){
if(a == id<T>()) return b;
else if(b == id<T>()) return a;
return a + b;
}
template<typename T>
static T merge(T a, T b){
return std::min(a, b);
}
};
struct point_add_range_max{
template<typename T>
static T id(){
return std::numeric_limits<T>::min();
}
template<typename T>
static T update(T a, T b){
if(a == id<T>()) return b;
else if(b == id<T>()) return a;
return a + b;
}
template<typename T>
static T merge(T a, T b){
return std::max(a, b);
}
};
struct point_add_range_sum{
template<typename T>
static T id(){
return 0;
}
template<typename T>
static T update(T a, T b){
return a + b;
}
template<typename T>
static T merge(T a, T b){
return a + b;
}
template<typename T>
static T flip(T a){
return a;
}
};
struct point_set_range_composite{
static constexpr int mod = 998244353;
template<typename T>
static T id(){
return {1, 0};
}
template<typename T>
static T update(T a, T b){
return b;
}
template<typename T>
static T merge(T a, T b){
int xy = ((long long)a.first * b.first) % mod;
int ab = ((long long)a.second * b.first + b.second) % mod;
return {xy, ab};
}
};
struct point_set_range_composite_flip{
static constexpr int mod = 998244353;
template<typename T>
static T id(){
return {1, 0, 0};
}
template<typename T>
static T update(T a, T b){
return b;
}
template<typename T>
static T flip(T a){
return {a[0], a[2], a[1]};
}
template<typename T>
static T merge(T a, T b){
int xy = ((long long)a[0] * b[0]) % mod;
int ab = ((long long)a[1] * b[0] + b[1]) % mod;
int ba = ((long long)b[2] * a[0] + a[2]) % mod;
return {xy, ab, ba};
}
};
struct point_add_range_gcd{
template<typename Val>
static Val __binary_gcd(Val a, Val b){
if(!a || !b) return !a ? b : a;
if(a < 0) a *= -1;
if(b < 0) b *= -1;
int shift = __builtin_ctzll(a | b); // [1] gcd(2a', 2b') = 2 * gcd(a', b')
a >>= __builtin_ctzll(a);
do{
// if b is odd
// gcd(2a', b) = gcd(a', b), if a = 2a'(a is even)
// gcd(a, b) = gcd(|a - b|, min(a, b)), if a is odd
b >>= __builtin_ctzll(b); // make b odd
if(a > b) std::swap(a, b);
b -= a;
}while(b); // gcd(a, 0) = a
return a << shift; // [1]
}
template<typename Val>
static Val id(){
return 0;
}
template<typename Val>
static Val update(Val a, Val b){
return a + b;
}
template<typename Val>
static Val merge(Val a, Val b){
return __binary_gcd(a, b);
}
};
// 区間set, これまでにsetした物の中ならどれかを取得
struct range_set_get_any{
template<typename Val>
static Val id1(){
return nullptr;
}
template<typename Lazy>
static Lazy id2(){
return nullptr;
}
template<typename Lazy>
static Lazy propagate(Lazy l, Lazy x){
return (x == nullptr ? l : x);
}
template<typename Val, typename Lazy>
static Val apply(Val v, Lazy x, int l, int r){
return (x == nullptr ? v : x);
}
};
struct range_add_range_sum{
template<typename T>
static T id1(){
return T(0);
}
template<typename E>
static E id2(){
return E(0);
}
template<typename T>
static T merge(T a, T b){
return a + b;
}
template<typename T, typename E>
static T apply(T a, E b, int l, int r){
return a + b * (r - l);
}
template<typename E>
static E propagate(E a, E b){
return a + b;
}
template<typename T>
static T flip(T a){
return a;
}
};
struct range_max_range_max{
template<typename T>
static T id1(){
return std::numeric_limits<T>::min();
}
template<typename E>
static E id2(){
return std::numeric_limits<E>::min();
}
template<typename T>
static T merge(T a, T b){
return std::max(a, b);
}
template<typename T, typename E>
static T apply(T a, E b, int l, int r){
return std::max(a, b);
}
template<typename E>
static E propagate(E a, E b){
return std::max(a, b);
}
};
struct range_add_range_max{
template<typename T>
static T id1(){
return std::numeric_limits<T>::min();
}
template<typename E>
static E id2(){
return 0;
}
template<typename T>
static T merge(T a, T b){
return std::max(a, b);
}
template<typename T, typename E>
static T apply(T a, E b, int l, int r){
if(a == id1<T>()) return b * (r - l);
return a + b;
}
template<typename E>
static E propagate(E a, E b){
return a + b;
}
};
struct range_add_range_min{
template<typename T>
static T id1(){
return std::numeric_limits<T>::max();
}
template<typename E>
static E id2(){
return 0;
}
template<typename T>
static T merge(T a, T b){
return std::min(a, b);
}
template<typename T, typename E>
static T apply(T a, E b, int l, int r){
if(a == id1<T>()) return b * (r - l);
return a + b;
}
template<typename E>
static E propagate(E a, E b){
return a + b;
}
};
struct range_add_range_argmin{
template<typename T>
static T id1(){
return {std::numeric_limits<long long>::max(), -1} ;
}
template<typename E>
static E id2(){
return 0;
}
template<typename T>
static T merge(T a, T b){
return std::min(a, b);
}
template<typename T, typename E>
static T apply(T a, E b, int l, int r){
if(a == id1<T>()) return a;
return {a.first + b, a.second};
}
template<typename E>
static E propagate(E a, E b){
return a + b;
}
};
/*
#include <array>
struct range_affine_range_sum{
const static long long MOD = 998244353;
template<typename T>
static T id1(){
return 0;
}
template<typename E>
static E id2(){
return E{1, 0};
}
template<typename T>
static T merge(T a, T b){
return (a + b) % MOD;
}
template<typename T, typename E>
static T apply(T a, E b, int l, int r){
return (a * b[0] + b[1] * (r - l)) % MOD;
}
template<typename E>
static E propagate(E a, E b){
return E{(a[0] * b[0]) % MOD, (a[1] * b[0] + b[1]) % MOD};
}
};
*/
struct range_affine_range_sum{
const static int MOD = 998244353;
template<typename T>
static T id1(){
return 0;
}
template<typename E>
static E id2(){
return E{1, 0};
}
template<typename T>
static T merge(T a, T b){
return (a + b) % MOD;
}
template<typename T, typename E>
static T apply(T a, E b, int l, int r){
return ((long long)a * b.first + (long long)b.second * (r - l)) % MOD;
}
template<typename E>
static E propagate(E a, E b){
return E{((long long)a.first * b.first) % MOD, ((long long)a.second * b.first + b.second) % MOD};
}
};
struct range_add_range_min_count{
static constexpr int INF = std::numeric_limits<int>::max();
template<typename T>
static T id1(){
return {INF, 0};
}
template<typename E>
static E id2(){
return 0;
}
template<typename T>
static T merge(T a, T b){
if(a.first != b.first) return a.first < b.first ? a : b;
return {a.first, a.second + b.second};
}
template<typename T, typename E>
static T apply(T a, E b, int l, int r){
if(a.first == INF) return {b, r - l};
return {a.first + b, a.second};
}
template<typename E>
static E propagate(E a, E b){
return a + b;
}
};
struct range_composite_lct{
static constexpr int MOD = 998244353;
template<typename T>
// 1x + 0, 1x + 0
static T id1(){
return std::array<int, 3>{1, 0, 0};
}
// no use
template<typename E>
static E id2(){
return false;
}
// b(a(x)), a(b(x))
template<typename T>
static T merge(T a, T b){
int ba1 = ((long long)b[0] * a[0]) % MOD;
int ba2 = ((long long)b[0] * a[1] + b[1]) % MOD;
int ab2 = ((long long)a[0] * b[2] + a[2]) % MOD;
return std::array<int, 3>{ba1, ba2, ab2};
}
// no use
template<typename T, typename E>
static T apply(T a, E b, int l, int r){
return a;
}
// no use
template<typename E>
static E propagate(E a, E b){
return false;
}
//
template<typename T>
static T flip(T a){
return std::array<int, 3>{a[0], a[2], a[1]};
}
};
#line 5 ".lib/data_structure/segment_tree/segment_tree.hpp"
#include <cstdint>
#line 10 ".lib/data_structure/segment_tree/segment_tree.hpp"
template<typename monoid, typename Val>
struct segment_tree{
int N, M;
std::vector<Val> sum;
int ceil_pow2(int y){
int x = 0;
while ((1U << x) < (unsigned int)(y)) x++;
return x;
};
segment_tree(){}
segment_tree(int n): N(n), M(1 << ceil_pow2(N)), sum(2 * M - 1, monoid::template id<Val>()){}
segment_tree(const std::vector<Val> v): N(v.size()), M(1 << ceil_pow2(N)), sum(2 * M - 1, monoid::template id<Val>()){
std::copy(v.begin(), v.end(), sum.begin() + M - 1);
for(int i = M - 2; i >= 0; i--){
sum[i] = monoid::template merge<Val>(sum[i * 2 + 1], sum[i * 2 + 2]);
}
}
int size(){
return N;
}
void set(int k, Val x){
assert(0 <= k && k < N);
k += M - 1;
sum[k] = x;
while(k){
k = (k - 1) >> 1;
sum[k] = monoid::template merge<Val>(sum[k * 2 + 1], sum[k * 2 + 2]);
}
}
Val get(int k){
assert(0 <= k && k < N);
return sum[M - 1 + k];
}
void update(int k, Val x){
assert(0 <= k && k < N);
set(k, monoid::template update<Val>(sum[k + M - 1], x));
}
Val query(int l, int r){
l = std::max(l, 0), r = std::min(r, N);
assert(l <= r);
l += M, r += M;
Val L = monoid::template id<Val>(), R = L;
while(l < r){
if(l & 1) L = monoid::template merge<Val>(L, sum[(l++) - 1]);
if(r & 1) R = monoid::template merge<Val>(sum[(--r) - 1], R);
l >>= 1;
r >>= 1;
}
return monoid::template merge<Val>(L, R);
}
Val query_all(){
return sum[0];
}
// f(sum[l, r])がtrueになる最左のr. ない場合は-1
template<typename F>
int bisect_from_left(int l, const F &f){
assert(0 <= l);
assert(!f(monoid::template id<Val>()));
if(l >= N) return -1;
l += M;
Val ret = monoid::template id<Val>();
do{
while(l % 2 == 0) l >>= 1;
if(f(monoid::template merge(ret, sum[l - 1]))){
while(l < M){
l = 2 * l;
if(!f(monoid::template merge(ret, sum[l - 1]))){
ret = monoid::template merge(ret, sum[l - 1]);
l++;
}
}
return l - M;
}
ret = monoid::template merge(ret, sum[l - 1]);
l++;
}while((l & -l) != l);
return -1;
}
// f(sum[l, r])がtrueになる最右のl. ない場合は-1
template<typename F>
int bisect_from_right(int r, const F &f){
assert(0 <= r && r < N);
assert(!f(monoid::template id<Val>()));
r++;
r += M;
Val ret = monoid::template id<Val>();
do{
r--;
while(r > 1 && (r % 2)) r >>= 1;
if(f(monoid::template merge<Val>(sum[r - 1], ret))){
while(r < M){
r = 2 * r + 1;
if(!f(monoid::template merge<Val>(sum[r - 1], ret))){
ret = monoid::template merge<Val>(sum[r - 1], ret);
r--;
}
}
return r - M;
}
}while((r & -r) != r);
return -1;
}
};
#line 5 "a.cpp"
using mint = dynamic_modint<0>;
using M = matrix_mod<mint>;
struct range_matmul{
template<typename T>
static T id(){
return M::eye(2);
}
template<typename T>
static T update(T a, T b){
return b * a;
}
template<typename T>
static T merge(T a, T b){
return b * a;
}
};
int main(){
io_init();
ll n, b, q;
std::cin >> n >> b >> q;
mint::set_mod(b);
std::vector<M> v(n);
range(i, 0, n){
int x, y, z, w;
std::cin >> x >> y;
std::cin >> z >> w;
v[i] = M(2, 2);
v[i][0][0] = x;
v[i][0][1] = y;
v[i][1][0] = z;
v[i][1][1] = w;
}
segment_tree<range_matmul, M> seg(v);
range(i, 0, q){
int l, r, x, y;
std::cin >> l >> r >> x >> y;
l--;
auto tmp = M(2, 1);
tmp[0][0] = x;
tmp[1][0] = y;
auto ans = seg.query(l, r) * tmp;
std::cout << ans[0][0] << " " << ans[1][0] << '\n';
}
}
tonegawa