結果
問題 | No.886 Direct |
ユーザー |
![]() |
提出日時 | 2021-07-10 19:20:55 |
言語 | Crystal (1.14.0) |
結果 |
RE
|
実行時間 | - |
コード長 | 16,490 bytes |
コンパイル時間 | 12,690 ms |
コンパイル使用メモリ | 295,020 KB |
実行使用メモリ | 153,600 KB |
最終ジャッジ日時 | 2024-07-02 02:45:16 |
合計ジャッジ時間 | 14,894 ms |
ジャッジサーバーID (参考情報) |
judge2 / judge4 |
(要ログイン)
ファイルパターン | 結果 |
---|---|
sample | AC * 4 |
other | AC * 21 RE * 11 |
ソースコード
# require "/math/Mint"# require "../atcoder/src/Math"# ac-library.cr by hakatashi https://github.com/google/ac-library.cr## Copyright 2021 Google LLC## Licensed under the Apache License, Version 2.0 (the "License");# you may not use this file except in compliance with the License.# You may obtain a copy of the License at## https://www.apache.org/licenses/LICENSE-2.0## Unless required by applicable law or agreed to in writing, software# distributed under the License is distributed on an "AS IS" BASIS,# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.# See the License for the specific language governing permissions and# limitations under the License.module AtCoder# Implements [ACL's Math library](https://atcoder.github.io/ac-library/master/document_en/math.html)module Mathdef self.extended_gcd(a, b)last_remainder, remainder = a.abs, b.absx, last_x, y, last_y = 0_i64, 1_i64, 1_i64, 0_i64while remainder != 0new_last_remainder = remainderquotient, remainder = last_remainder.divmod(remainder)last_remainder = new_last_remainderx, last_x = last_x - quotient * x, xy, last_y = last_y - quotient * y, yendreturn last_remainder, last_x * (a < 0 ? -1 : 1)end# Implements atcoder::inv_mod(value, modulo).def self.inv_mod(value, modulo)gcd, inv = extended_gcd(value, modulo)if gcd != 1raise ArgumentError.new("#{value} and #{modulo} are not coprime")endinv % moduloend# Simplified AtCoder::Math.pow_mod with support of Int64def self.pow_mod(base, exponent, modulo)if exponent == 0return base.class.zero + 1endif base == 0return baseendb = exponent > 0 ? base : inv_mod(base, modulo)e = exponent.absret = 1_i64while e > 0if e % 2 == 1ret = mul_mod(ret, b, modulo)endb = mul_mod(b, b, modulo)e //= 2endretend# Caluculates a * b % mod without overflow detection@[AlwaysInline]def self.mul_mod(a : Int64, b : Int64, mod : Int64)if mod < Int32::MAXreturn a * b % modend# 31-bit widtha_high = (a >> 32).to_u64# 32-bit widtha_low = (a & 0xFFFFFFFF).to_u64# 31-bit widthb_high = (b >> 32).to_u64# 32-bit widthb_low = (b & 0xFFFFFFFF).to_u64# 31-bit + 32-bit + 1-bit = 64-bitc = a_high * b_low + b_high * a_lowc_high = c >> 32c_low = c & 0xFFFFFFFF# 31-bit + 31-bitres_high = a_high * b_high + c_high# 32-bit + 32-bitres_low = a_low * b_lowres_low_high = res_low >> 32res_low_low = res_low & 0xFFFFFFFF# Overflowif res_low_high + c_low >= 0x100000000res_high += 1endres_low = (((res_low_high + c_low) & 0xFFFFFFFF) << 32) | res_low_low(((res_high.to_i128 << 64) | res_low) % mod).to_i64end@[AlwaysInline]def self.mul_mod(a, b, mod)typeof(mod).new(a.to_i64 * b % mod)end# Implements atcoder::crt(remainders, modulos).def self.crt(remainders, modulos)raise ArgumentError.new unless remainders.size == modulos.sizetotal_modulo = 1_i64answer = 0_i64remainders.zip(modulos).each do |(remainder, modulo)|gcd, p = extended_gcd(total_modulo, modulo)if (remainder - answer) % gcd != 0return 0_i64, 0_i64endtmp = (remainder - answer) // gcd * p % (modulo // gcd)answer += total_modulo * tmptotal_modulo *= modulo // gcdendreturn answer % total_modulo, total_moduloend# Implements atcoder::floor_sum(n, m, a, b).def self.floor_sum(n, m, a, b)n, m, a, b = n.to_i64, m.to_i64, a.to_i64, b.to_i64res = 0_i64if a < 0a2 = a % mres -= n * (n - 1) // 2 * ((a2 - a) // m)a = a2endif b < 0b2 = b % mres -= n * ((b2 - b) // m)b = b2endres + floor_sum_unsigned(n, m, a, b)endprivate def self.floor_sum_unsigned(n, m, a, b)res = 0_i64loop doif a >= mres += n * (n - 1) // 2 * (a // m)a = a % mendif b >= mres += n * (b // m)b = b % mendy_max = a * n + bbreak if y_max < mn = y_max // mb = y_max % mm, a = a, mendresendendendmacro static_modint(name, mod)struct {{name}}MOD = {{mod}}i64def self.zeronewenddef self.raw(value : Int64)result = newresult.value = valueresultendgetter value : Int64def initialize@value = 0i64enddef initialize(value)@value = value.to_i64 % MODenddef initialize(m : self)@value = m.valueendprotected def value=(value : Int64)@value = valueenddef ==(m : self)value == m.valueenddef ==(m)value == menddef + : selfselfenddef - : selfself.class.raw(value != 0 ? MOD &- value : 0i64)enddef +(v)self + self.class.new(v)enddef +(m : self)x = value &+ m.valuex &-= MOD if x >= MODself.class.raw(x)enddef -(v)self - self.class.new(v)enddef -(m : self)x = value &- m.valuex &+= MOD if x < 0self.class.raw(x)enddef *(v)self * self.class.new(v)enddef *(m : self)self.class.new(value &* m.value)enddef /(v)self / self.class.new(v)enddef /(m : self)raise DivisionByZeroError.new if m.value == 0a, b, u, v = m.value, MOD, 1i64, 0i64while b != 0t = a // ba &-= t &* ba, b = b, au &-= t &* vu, v = v, uendself.class.new(value &* u)enddef //(v)self / venddef **(exponent : Int)t, res = self, self.class.raw(1i64)while exponent > 0res *= t if exponent & 1 == 1t *= texponent >>= 1endresend{% for op in %w[< <= > >=] %}def {{op.id}}(other)raise NotImplementedError.new({{op}})end{% end %}def invself.class.raw AtCoder::Math.inv_mod(value, MOD)enddef succself.class.raw(value != MOD &- 1 ? value &+ 1 : 0i64)enddef predself.class.raw(value != 0 ? value &- 1 : MOD &- 1)enddef absselfenddef abs2self * selfenddef to_i64 : Int64valueenddelegate to_s, to: @valuedelegate inspect, to: @valueend{% to = ("to_" + name.stringify.downcase.gsub(/mint|modint/, "m")).id %}struct Int{% for op in %w[+ - * / //] %}def {{op.id}}(value : {{name}}){{to}} {{op.id}} valueend{% end %}{% for op in %w[< <= > >=] %}def {{op.id}}(m : {{name}})raise NotImplementedError.new({{op}})end{% end %}def {{to}} : {{name}}{{name}}.new(self)endendclass Stringdef {{to}} : {{name}}{{name}}.new(self)endendendstatic_modint(Mint, 1000000007)static_modint(Mint2, 998244353)# require "/math/GCD"# require "/atcoder/src/Prime"# ac-library.cr by hakatashi https://github.com/google/ac-library.cr## Copyright 2021 Google LLC## Licensed under the Apache License, Version 2.0 (the "License");# you may not use this file except in compliance with the License.# You may obtain a copy of the License at## https://www.apache.org/licenses/LICENSE-2.0## Unless required by applicable law or agreed to in writing, software# distributed under the License is distributed on an "AS IS" BASIS,# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.# See the License for the specific language governing permissions and# limitations under the License.# require "./Math.cr"module AtCoder# Implements [Ruby's Prime library](https://ruby-doc.com/stdlib/libdoc/prime/rdoc/Prime.html).## ```# AtCoder::Prime.first(7) # => [2, 3, 5, 7, 11, 13, 17]# ```module Primeextend selfinclude Enumerable(Int64)@@primes = [2_i64, 3_i64, 5_i64, 7_i64, 11_i64, 13_i64, 17_i64, 19_i64,23_i64, 29_i64, 31_i64, 37_i64, 41_i64, 43_i64, 47_i64,53_i64, 59_i64, 61_i64, 67_i64, 71_i64, 73_i64, 79_i64,83_i64, 89_i64, 97_i64, 101_i64,]def eachindex = 0loop doyield get_nth_prime(index)index += 1endenddef prime_division(value : Int)raise DivisionByZeroError.new if value == 0int = typeof(value)factors = [] of Tuple(typeof(value), typeof(value))if value < 0value = value.absfactors << {int.new(-1), int.new(1)}enduntil prime?(value) || value == 1factor = valueuntil prime?(factor)factor = find_factor(factor)endcount = 0while value % factor == 0value //= factorcount += 1endfactors << {int.new(factor), int.new(count)}endif value > 1factors << {value, int.new(1)}endfactors.sort_by! { |(factor, _)| factor }endprivate def find_factor(n : Int)# Factor of 4 cannot be discovered by Pollard's Rho with f(x) = x^x+1if n == 4typeof(n).new(2)elsepollard_rho(n).not_nil!endend# Get single factor by Pollard's Rho Algorithmprivate def pollard_rho(n : Int)typeof(n).new(1).upto(n) do |i|x = iy = pollard_random_f(x, n)loop dox = pollard_random_f(x, n)y = pollard_random_f(pollard_random_f(y, n), n)gcd = (x - y).gcd(n)if gcd == nbreakendif gcd != 1return gcdendendendendprivate def pollard_random_f(n : Int, mod : Int)(AtCoder::Math.mul_mod(n, n, mod) + 1) % modendprivate def extract_prime_division_base(prime_divisions_class : Array({T, T}).class) forall TTenddef int_from_prime_division(prime_divisions : Array({Int, Int}))int_class = extract_prime_division_base(prime_divisions.class)prime_divisions.reduce(int_class.new(1)) { |i, (factor, exponent)| i * factor ** exponent }enddef prime?(value : Int)# Obvious patternsreturn false if value < 2return true if value <= 3return false if value.even?return true if value < 9if value < 0xffffreturn false unless typeof(value).new(30).gcd(value % 30) == 17.step(by: 30, to: value) do |base|break if base * base > valueif {0, 4, 6, 10, 12, 16, 22, 24}.any? { |i| value % (base + i) == 0 }return falseendendreturn trueendmiller_rabin(value.to_i64)endprivate def miller_rabin(value)d = value - 1s = 0_i64until d.odd?d >>= 1s += 1endmiller_rabin_bases(value).each do |base|next if base == valuex = AtCoder::Math.pow_mod(base.to_i64, d, value)next if x == 1 || x == value - 1is_composite = s.times.all? dox = AtCoder::Math.mul_mod(x, x, value)x != value - 1endreturn false if is_compositeendtrueend# We can reduce time complexity of Miller-Rabin tests by testing against# predefined bases which is enough to test against primarity in the given range.# https://en.wikipedia.org/wiki/Miller%E2%80%93Rabin_primality_test# ameba:disable Metrics/CyclomaticComplexityprivate def miller_rabin_bases(value)casewhen value < 1_373_653_i64[2, 3]when value < 9_080_191_i64[31, 73]when value < 25_326_001_i64[2, 3, 5]when value < 3_215_031_751_i64[2, 3, 5, 7]when value < 4_759_123_141_i64[2, 7, 61]when value < 1_122_004_669_633_i64[2, 13, 23, 1662803]when value < 2_152_302_898_747_i64[2, 3, 5, 7, 11]when value < 3_474_749_660_383_i64[2, 3, 5, 7, 11, 13]when value < 341_550_071_728_321_i64[2, 3, 5, 7, 11, 13, 17]when value < 3_825_123_056_546_413_051_i64[2, 3, 5, 7, 11, 13, 17, 19, 23]else[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37]endendprivate def get_nth_prime(n)while @@primes.size <= ngenerate_primesend@@primes[n]end# Doubles the size of the cached prime array and performs the# Sieve of Eratosthenes on it.private def generate_primesnew_primes_size = @@primes.size < 1_000_000 ? @@primes.size : 1_000_000new_primes = Array(Int64).new(new_primes_size) { |i| @@primes.last + (i + 1) * 2 }new_primes_max = new_primes.last@@primes.each do |prime|next if prime == 2break if prime * prime > new_primes_max# Here I use the technique of the Sieve of Sundaram. We can# only test against the odd multiple of the given prime.# min_composite is the minimum number that is greater than# the last confirmed prime, and is an odd multiple of# the given prime.min_multiple = ((@@primes.last // prime + 1) // 2 * 2 + 1) * primemin_multiple.step(by: prime * 2, to: new_primes_max) do |multiple|index = new_primes_size - (new_primes_max - multiple) // 2 - 1new_primes[index] = 0_i64endend@@primes.concat(new_primes.reject(0_i64))endprivate struct EachDivisor(T)include Enumerable(T)def initialize(@exponential_factors : Array(Array(T)))enddef eachArray.each_product(@exponential_factors) do |factors|yield factors.reduce { |a, b| a * b }endendend# Returns an enumerator that iterates through the all positive divisors of# the given number. **The order is not guaranteed.**# Not in the original Ruby's Prime library.## ```# AtCoder::Prime.each_divisor(20) do |n|# puts n# end # => Puts 1, 2, 4, 5, 10, and 20## AtCoder::Prime.each_divisor(10).map { |n| 1.0 / n }.to_a # => [1.0, 0.5, 0.2, 0.1]# ```def each_divisor(value : Int)raise ArgumentError.new unless value > 0factors = prime_division(value)if value == 1exponential_factors = [[value]]elseexponential_factors = factors.map do |(factor, count)|cnt = typeof(value).zero + 1Array(typeof(value)).new(count + 1) do |i|cnt_copy = cntif i < countcnt *= factorendcnt_copyendendendEachDivisor(typeof(value)).new(exponential_factors)end# :ditto:def each_divisor(value : T, &block : T ->)each_divisor(value).each(&block)endendendstruct Intdef prime?AtCoder::Prime.prime?(self)endendclass Array(T)# result[i] = Sum_{n | i} a[i]def gcd_zeta!AtCoder::Prime.each do |p|break if p >= sizei, k = size.pred // p, size.pred // p * pwhile k > 0self[i] += self[i * p]i -= 1; k -= pendendselfend# :ditto:def gcd_zetadup.gcd_zeta!end# a[i] = Sum_{n | i} result[i]def gcd_mobius!AtCoder::Prime.each do |p|break if p >= sizei, k = 1, pwhile k < sizeself[i] -= self[k]i += 1; k += pendendselfend# :ditto:def gcd_mobiusdup.gcd_mobius!endendmodule GCDextend self# result[n] = Sum_{gcd(i, j) = n} f[i] * g[j]def convolution(f : Array(T), g : Array(T)) forall Tf.gcd_zeta.zip?(g.gcd_zeta).map { |x, y| (x || T.zero) * (y || T.zero) }.gcd_mobius!endendh, w = read_line.split.map(&.to_i)hh = (0...h).map { |i| (h - i).to_m }ww = (0...w).map { |i| (w - i).to_m }one = h * w.pred + h.pred * wputs GCD.convolution(hh, ww).fetch(1, Mint.zero) * 2 + one