{-# LANGUAGE BangPatterns #-} {-# OPTIONS_GHC -O2 #-} module Main where import Data.Maybe (fromJust) import qualified Data.ByteString.Char8 as BS import Control.Monad (replicateM, replicateM_) 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 {-| Max Weighted Floor (mwf) mwf(n,m,a,b,c,d) = max_{0 <= x < n} a*x + b*floor((c*x + d)/m) を返す非再帰(反復)実装。 前提: - n > 0, m > 0 計算量/メモリ: - 時間: O(log m)(ユークリッド互除法的再帰による構造縮約) - 追加メモリ: O(1) -} mwf :: Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer mwf n0 m0 a0 b0 c0 d0 = let !sum0 = 0 !max0 = b0 * (d0 `div` m0) -- x = 0 のとき in go n0 m0 a0 b0 c0 d0 sum0 max0 where go :: Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer go !n !m !a !b !c !d !sumAcc !maxAcc = let (q1, c') = c `divMod` m !a' = a + b * q1 (q2, d') = d `divMod` m !sum' = sumAcc + b * q2 !max' = max maxAcc sum' !ymax = (c' * (n - 1) + d') `div` m in if ymax == 0 then max max' (sum' + a' * (n - 1)) else if a' >= 0 then let !max'' = max max' (sum' + a' * (n - 1) + b * ymax) in go ymax c' b a' m (m - d' - 1) sum' max'' else let !sum'' = sum' + a' + b in go ymax c' b a' m (m - d' - 1) sum'' max' {-| mwf(n,m,a,b,c,d) <= z かどうかを判定して Bool を返す。 (上の mwf と同じ前提・計算量) -} mwfLeq :: Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Bool mwfLeq z n0 m0 a0 b0 c0 d0 = go n0 m0 a0 b0 c0 d0 (-z) where go :: Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Bool go !n !m !a !b !c !d !sumAcc = let (q1, c') = c `divMod` m !a' = a + b * q1 (q2, d') = d `divMod` m !sum' = sumAcc + b * q2 in ((sum' <= 0) && (let !ymax = (c' * (n - 1) + d') `div` m in if ymax == 0 then (sum' + a' * (n - 1)) <= 0 else (sum' + max 0 (a' * (n - 1)) + max 0 (b * ymax) <= 0) || (if a' >= 0 then ((sum' + a' * (n - 1) + b * ymax) <= 0) && go ymax c' b a' m (m - d' - 1) sum' else let !sum'' = sum' + a' + b in go ymax c' b a' m (m - d' - 1) sum''))) {-| max_{l <= x < r} a*x + b*floor((c*x + d)/m) を返す。 既存の mwf(n,m,...)(0 <= x < n)を用いる。 前提: l < r, m > 0 -} mwfLr :: Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer -> Integer mwfLr l r m a b c d = let !n = r - l (q, d') = (c * l + d) `divMod` m in a * l + b * q + mwf n m a b c d' {-| max_{L <= x < R} (...) <= 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 hi0 = a' * b * k + 2 -- [0, hi) の上限 -- 二分探索で u_min を求める(F(u) <= T_Δ を満たす最大の u) go 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 go mid hi else go lo mid in go 0 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