{-# LANGUAGE BangPatterns #-} {-# OPTIONS_GHC -O2 #-} module Main where import Control.Monad (replicateM, replicateM_) import Data.ByteString.Char8 qualified as BS import Data.Maybe (fromJust) tuplify2 (x : y : _) = (x, y) tuplify2 _ = undefined -- Input functions with ByteString readInt = fst . fromJust . BS.readInteger readIntTuple = tuplify2 . map readInt . BS.words readIntList = map readInt . BS.words getInt = readInt <$> BS.getLine getIntList = readIntList <$> BS.getLine getIntNList n = map readIntList <$> replicateM (fromIntegral n) BS.getLine getIntMatrix = map readIntList . BS.lines <$> BS.getContents getIntTuple = readIntTuple <$> BS.getLine getIntNTuples n = map readIntTuple <$> replicateM (fromIntegral n) BS.getLine getIntTuples = map readIntTuple . BS.lines <$> BS.getContents -- | -- mwf(n,m,a,b,c,d) = max_{0 <= x < n} a*x + b*floor((c*x + d)/m) について、 -- mwf(n,m,a,b,c,d) <= z かどうかを判定して Bool を返す。 -- -- 前提: -- - z,n,m,a,b,c,d は整数 -- - n > 0, m > 0 -- -- 計算量/メモリ: -- - 時間: O(log m)(ユークリッド互除法的再帰による構造縮約) -- - 追加メモリ: O(1) mwfLeq :: Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Bool mwfLeq z n0 m0 a0 b0 c0 d0 = loop n0 m0 a0 b0 c0 d0 (-z) where loop :: Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Bool loop !n !m !a !b !c !d !s = let (q1, c') = c `divMod` m !a1 = a + b * q1 (q2, d') = d `divMod` m !sum_acc = s + b * q2 in ( sum_acc <= 0 && ( (a1 <= 0 && b <= 0) || ( let !ymax = (c' * (n - 1) + d') `div` m in if ymax == 0 then sum_acc + a1 * (n - 1) <= 0 else if a1 >= 0 then ((sum_acc + a1 * (n - 1) + b * ymax) <= 0) && ((b >= 0) || loop ymax c' b a1 m (m - d' - 1) sum_acc) else loop ymax c' b a1 m (m - d' - 1) (sum_acc + a1 + b) ) ) ) -- | -- max_{l <= x < r} a*x + b*floor((c*x + d)/m) <= z かどうかを返す述語版。 mwfLrLeq :: Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Bool mwfLrLeq z l r m a b c d = let !n = r - l (q, d') = (c * l + d) `divMod` m in mwfLeq (z - a * l - b * q) n m a b c d' -- | -- Δ(D,A,B,x) = floor(x/D) - floor( (floor(x/A) * floor(A*B/D)) / B ) において、 -- u_min(D,A,B,K) = min { u >= 0 | Δ(D,A,B,u*D) > K } を半開区間二分探索 [0, A'BK+2) で求め、 -- x_min(D,A,B,K) = min { x >= 0 | Δ(D,A,B,x) > K } = u_min(D,A,B,K)*D を返す(解なしは -1)。 -- -- 前提: -- * D > 0, A > 0, B > 0, K >= 0(整数) computeXminLeq :: Integer -> Integer -> Integer -> Integer -> Integer computeXminLeq d a b k = let g = gcd d a !d' = d `div` g !a' = a `div` g (m', r') = (a' * b) `divMod` d' !tdelta = b * k in -- 解なし判定 if r' == 0 && d' * k + 1 >= a' then (-1) else let lo0 = if r' == 0 then k + 1 else max ((a' * b * k + r' - (a' - 1) * m') `div` r') (k + 1) hi0 = if r' == 0 then a' else ((a' * b * k + r') `div` r') + 1 -- 二分探索で u_min を求める(F(u) <= T_Δ を満たす最大の u) loop lo hi | lo + 1 >= hi = d * lo | otherwise = let mid = (lo + hi) `div` 2 ok = mwfLrLeq tdelta lo mid a' b (-m') d' 0 in if ok then loop mid hi else loop lo mid in loop lo0 hi0 -- | -- 検算用 Δ(D,A,B,x)。 deltaVal :: Integer -> Integer -> Integer -> Integer -> Integer deltaVal d a b x = let p = x `div` d m = a * b `div` d q = x `div` a * m `div` b in p - q -- | -- 入出力(複数ケース) solve :: IO () solve = do t <- getInt replicateM_ (fromIntegral t :: Int) $ do [d, a, b, k] <- getIntList print (computeXminLeq d a b k) main :: IO () main = solve