(in-package :cl-user) (eval-when (:compile-toplevel :load-toplevel :execute) (sb-int:defconstant-eqx opt #+swank '(optimize (speed 3) (safety 2)) #-swank '(optimize (speed 3) (safety 0) (debug 0)) #'equal) #+swank (ql:quickload '(:cl-debug-print :fiveam) :silent t) #-swank (set-dispatch-macro-character #\# #\> (lambda (s c p) (declare (ignore c p)) `(values ,(read s nil nil t))))) #+swank (cl-syntax:use-syntax cl-debug-print:debug-print-syntax) (defmacro define-int-types (&rest bits) `(progn ,@(mapcar (lambda (b) `(deftype ,(intern (format nil "UINT~A" b)) () '(unsigned-byte ,b))) bits) ,@(mapcar (lambda (b) `(deftype ,(intern (format nil "INT~A" b)) () '(signed-byte ,b))) bits))) (define-int-types 2 4 7 8 15 16 31 32 62 63 64) (defconstant +mod+ 998244353) (defmacro dbg (&rest forms) #+swank (if (= (length forms) 1) `(format *error-output* "~A => ~A~%" ',(car forms) ,(car forms)) `(format *error-output* "~A => ~A~%" ',forms `(,,@forms))) #-swank (declare (ignore forms))) (declaim (inline println)) (defun println (obj &optional (stream *standard-output*)) (let ((*read-default-float-format* 'double-float)) (prog1 (princ obj stream) (terpri stream)))) ;; BEGIN_INSERTED_CONTENTS ;;; ;;; Binomial coefficient with mod ;;; build: O(n) ;;; query: O(1) ;;; (defpackage :cp/binomial-coefficient-mod (:use :cl) (:export #:binom #:perm #:multinomial #:stirling2 #:catalan #:+binom-mod+ #:*fact* #:*fact-inv* #:*inv*)) (in-package :cp/binomial-coefficient-mod) ;; TODO: non-global handling (defconstant +binom-size+ 510000) (defconstant +binom-mod+ (if (boundp 'cl-user::+mod+) (symbol-value 'cl-user::+mod+) #.(+ (expt 10 9) 7))) (defparameter *fact* (make-array +binom-size+ :element-type '(unsigned-byte 31)) "table of factorials") (defparameter *fact-inv* (make-array +binom-size+ :element-type '(unsigned-byte 31)) "table of inverses of factorials") (defparameter *inv* (make-array +binom-size+ :element-type '(unsigned-byte 31)) "table of inverses of non-negative integers") (declaim ((simple-array (unsigned-byte 31) (*)) *fact* *fact-inv* *inv*)) (defun initialize-binom () (declare (optimize (speed 3) (safety 0))) (setf (aref *fact* 0) 1 (aref *fact* 1) 1 (aref *fact-inv* 0) 1 (aref *fact-inv* 1) 1 (aref *inv* 1) 1) (loop for i from 2 below +binom-size+ do (setf (aref *fact* i) (mod (* i (aref *fact* (- i 1))) +binom-mod+) (aref *inv* i) (- +binom-mod+ (mod (* (aref *inv* (rem +binom-mod+ i)) (floor +binom-mod+ i)) +binom-mod+)) (aref *fact-inv* i) (mod (* (aref *inv* i) (aref *fact-inv* (- i 1))) +binom-mod+)))) (initialize-binom) (declaim (inline binom)) (defun binom (n k) "Returns nCk." (if (or (< n k) (< n 0) (< k 0)) 0 (mod (* (aref *fact* n) (mod (* (aref *fact-inv* k) (aref *fact-inv* (- n k))) +binom-mod+)) +binom-mod+))) (declaim (inline perm)) (defun perm (n k) "Returns nPk." (if (or (< n k) (< n 0) (< k 0)) 0 (mod (* (aref *fact* n) (aref *fact-inv* (- n k))) +binom-mod+))) ;; TODO: compiler macro or source-transform (declaim (inline multinomial)) (defun multinomial (&rest ks) "Returns the multinomial coefficient K!/k_1!k_2!...k_n! for K = k_1 + k_2 + ... + k_n. K must be equal to or smaller than MOST-POSITIVE-FIXNUM. (multinomial) returns 1." (let ((sum 0) (result 1)) (declare ((integer 0 #.most-positive-fixnum) result sum)) (dolist (k ks) (incf sum k) (setq result (mod (* result (aref *fact-inv* k)) +binom-mod+))) (mod (* result (aref *fact* sum)) +binom-mod+))) (declaim (inline stirling2)) (defun stirling2 (n k) "Returns the stirling number of the second kind S2(n, k). Time complexity is O(klog(n))." (declare ((integer 0 #.most-positive-fixnum) n k)) (labels ((mod-power (base exp) (declare ((integer 0 #.most-positive-fixnum) base exp)) (loop with res of-type (integer 0 #.most-positive-fixnum) = 1 while (> exp 0) when (oddp exp) do (setq res (mod (* res base) +binom-mod+)) do (setq base (mod (* base base) +binom-mod+) exp (ash exp -1)) finally (return res)))) (loop with result of-type fixnum = 0 for i from 0 to k for delta = (mod (* (binom k i) (mod-power i n)) +binom-mod+) when (evenp (- k i)) do (incf result delta) (when (>= result +binom-mod+) (decf result +binom-mod+)) else do (decf result delta) (when (< result 0) (incf result +binom-mod+)) finally (return (mod (* result (aref *fact-inv* k)) +binom-mod+))))) (declaim (inline catalan)) (defun catalan (n) "Returns the N-th Catalan number." (declare ((integer 0 #.most-positive-fixnum) n)) (mod (* (aref *fact* (* 2 n)) (mod (* (aref *fact-inv* (+ n 1)) (aref *fact-inv* n)) +binom-mod+)) +binom-mod+)) (defpackage :cp/mod-power (:use :cl) (:export #:mod-power)) (in-package :cp/mod-power) (declaim (inline mod-power)) (defun mod-power (base power modulus) "Returns BASE^POWER mod MODULUS. Note: 0^0 = 1. BASE := integer POWER, MODULUS := non-negative fixnum" (declare ((integer 0 #.most-positive-fixnum) modulus power) (integer base)) (let ((base (mod base modulus)) (res (mod 1 modulus))) (declare ((integer 0 #.most-positive-fixnum) base res)) (loop while (> power 0) when (oddp power) do (setq res (mod (* res base) modulus)) do (setq base (mod (* base base) modulus) power (ash power -1))) res)) ;;; ;;; Arithmetic operations with static modulus ;;; (defpackage :cp/mod-operations (:use :cl) (:export #:define-mod-operations)) (in-package :cp/mod-operations) ;; NOTE: Currently MOD* and MOD+ doesn't apply MOD when the number of ;; parameters is one. For simplicity I won't fix it for now. (defmacro define-mod-operations (divisor) (let ((mod* (intern "MOD*")) (mod+ (intern "MOD+")) (incfmod (intern "INCFMOD")) (decfmod (intern "DECFMOD")) (mulfmod (intern "MULFMOD"))) `(progn (defun ,mod* (&rest args) (reduce (lambda (x y) (mod (* x y) ,divisor)) args)) (defun ,mod+ (&rest args) (reduce (lambda (x y) (mod (+ x y) ,divisor)) args)) #+sbcl (eval-when (:compile-toplevel :load-toplevel :execute) (locally (declare (sb-ext:muffle-conditions warning)) (sb-c:define-source-transform ,mod* (&rest args) (if (null args) 1 (reduce (lambda (x y) `(mod (* ,x ,y) ,',divisor)) args))) (sb-c:define-source-transform ,mod+ (&rest args) (if (null args) 0 (reduce (lambda (x y) `(mod (+ ,x ,y) ,',divisor)) args))))) (define-modify-macro ,incfmod (delta) (lambda (x y) (mod (+ x y) ,divisor))) (define-modify-macro ,decfmod (delta) (lambda (x y) (mod (- x y) ,divisor))) (define-modify-macro ,mulfmod (multiplier) (lambda (x y) (mod (* x y) ,divisor)))))) ;; BEGIN_USE_PACKAGE (eval-when (:compile-toplevel :load-toplevel :execute) (use-package :cp/mod-operations :cl-user)) (eval-when (:compile-toplevel :load-toplevel :execute) (use-package :cp/mod-power :cl-user)) (eval-when (:compile-toplevel :load-toplevel :execute) (use-package :cp/binomial-coefficient-mod :cl-user)) (in-package :cl-user) ;;; ;;; Body ;;; (define-mod-operations +mod+) (defun solve-small (n m a b) (let ((res 0)) (loop for x from 1 to n do (loop for y from 1 to n for term = 0 when (<= (+ x y) n) do (loop for d from a to b do (dbg x y d) (incfmod term #>(mod* (- m d) (mod-power (- d 1) (- n x y) +mod+)))) (incfmod res (mod* (binom n x) (binom (- n x) y) term)))) res)) (defun solve (n m a b) (let ((res 0) (cumuls (make-array (+ n m) :element-type 'uint31 :initial-element 0))) (dotimes (i (- (length cumuls) 1)) (setf (aref cumuls (+ i 1)) (mod+ (aref cumuls i) (binom (+ i (- (* a (- n 1))) n -2) (- n 2))))) (loop for base from 1 to m for l = (min (+ base a) m) for r = (min (+ base b) m) for width1 = (- l base) for width2 = (- r base) when (>= width1 a) do (let ((value (- (aref cumuls (+ width2 1)) (aref cumuls width1)))) (incfmod res value))) (mod* res (aref *fact* n)))) (defun main () (let* ((n (read)) (m (read)) (a (read)) (b (read))) (println (solve n m a b)))) #-swank (main) ;;; ;;; Test and benchmark ;;; #+swank (defun get-clipbrd () (with-output-to-string (out) #+os-windows (run-program "powershell.exe" '("-Command" "Get-Clipboard") :output out :search t) #+os-unix (run-program "xsel" '("-b" "-o") :output out :search t))) #+swank (defparameter *this-pathname* (uiop:current-lisp-file-pathname)) #+swank (defparameter *dat-pathname* (uiop:merge-pathnames* "test.dat" *this-pathname*)) #+swank (defun run (&optional thing (out *standard-output*)) "THING := null | string | symbol | pathname null: run #'MAIN using the text on clipboard as input. string: run #'MAIN using the string as input. symbol: alias of FIVEAM:RUN!. pathname: run #'MAIN using the text file as input." (let* ((*standard-output* (or out (make-string-output-stream))) (res (etypecase thing (null (with-input-from-string (*standard-input* (delete #\Return (get-clipbrd))) (main))) (string (with-input-from-string (*standard-input* (delete #\Return thing)) (main))) (symbol (5am:run! thing)) (pathname (with-open-file (*standard-input* thing) (main)))))) (if out res (get-output-stream-string *standard-output*)))) #+swank (defun gen-dat () (uiop:with-output-file (out *dat-pathname* :if-exists :supersede) (format out ""))) #+swank (defun bench (&optional (out (make-broadcast-stream))) (time (run *dat-pathname* out))) ;; To run: (5am:run! :sample) #+swank (it.bese.fiveam:test :sample (it.bese.fiveam:is (equal "4 1 2 3 4 " (run "6 5 3 1 4 1 2 1 3 1 4 2 5 3 6 " nil))) (it.bese.fiveam:is (equal "-1 " (run "4 3 4 1 4 1 2 2 3 3 4 " nil))) (it.bese.fiveam:is (equal "6 1 2 4 7 9 10 " (run "10 9 5 1 10 7 2 4 9 5 3 1 9 2 9 8 7 8 5 6 4 10 2 " nil))))