1+ #line 1 "cp-algo/min-bundled/cp-algo/math/laurent.hpp"
2+ #line 1 "cp-algo/math/laurent.hpp"
3+ #line 1 "cp-algo/util/big_alloc.hpp"
4+ #include < vector>
5+ #include < cstddef>
6+ #include < iostream>
7+ #if defined(__linux__) || defined(__unix__) || (defined(__APPLE__) && defined(__MACH__))
8+ # define CP_ALGO_USE_MMAP 1
9+ # include < sys/mman.h>
10+ #else
11+ # define CP_ALGO_USE_MMAP 0
12+ #endif
13+ namespace cp_algo {template <typename T,std::size_t Align=32 >
14+ class big_alloc {static_assert (Align>=alignof (void *)," Align must be at least pointer-size" );
15+ static_assert (std::popcount(Align)==1 ," Align must be a power of two" );
16+ public: using value_type=T;
17+ template <class U >struct rebind {using other=big_alloc<U,Align>;};
18+ constexpr bool operator ==(const big_alloc&)const =default ;
19+ constexpr bool operator !=(const big_alloc&)const =default ;
20+ big_alloc ()noexcept =default ;
21+ template <typename U,std::size_t A>
22+ big_alloc (const big_alloc<U,A>&)noexcept {}
23+ [[nodiscard]]T*allocate (std::size_t n){std::size_t padded=round_up (n*sizeof (T));
24+ std::size_t align=std::max<std::size_t >(alignof (T),Align);
25+ #if CP_ALGO_USE_MMAP
26+ if (padded>=MEGABYTE){void *raw=mmap (nullptr ,padded,
27+ PROT_READ|PROT_WRITE,
28+ MAP_PRIVATE|MAP_ANONYMOUS,-1 ,0 );
29+ madvise (raw,padded,MADV_HUGEPAGE);
30+ madvise (raw,padded,MADV_POPULATE_WRITE);
31+ return static_cast <T*>(raw);}
32+ #endif
33+ return static_cast <T*>(::operator new (padded,std::align_val_t (align)));}
34+ void deallocate (T*p,std::size_t n)noexcept {if (!p)return ;
35+ std::size_t padded=round_up (n*sizeof (T));
36+ std::size_t align=std::max<std::size_t >(alignof (T),Align);
37+ #if CP_ALGO_USE_MMAP
38+ if (padded>=MEGABYTE){munmap (p,padded);return ;}
39+ #endif ::operator delete (p,padded,std::align_val_t (align));}
40+ private: static constexpr std::size_t MEGABYTE=1 <<20 ;
41+ static constexpr std::size_t round_up (std::size_t x)noexcept {return (x+Align-1 )/Align*Align;}};
42+ template <typename T>
43+ using big_vector=std::vector<T,big_alloc<T>>;}
44+ #line 4 "cp-algo/math/laurent.hpp"
45+ #include < memory>
46+ #include < optional>
47+ #line 7 "cp-algo/math/laurent.hpp"
48+ #include < cassert>
49+ #include < algorithm>
50+ namespace cp_algo ::math{template <typename T>
51+ struct provider {mutable big_vector<T>cache;
52+ mutable int cache_offset=0 ;
53+ mutable bool initialized=false ;
54+ mutable bool all_cached=false ;
55+ virtual ~provider ()=default ;
56+ virtual std::optional<int >valuation ()const {return std::nullopt ;}
57+ virtual std::optional<int >degree ()const {return std::nullopt ;}
58+ virtual T coeff_lazy (int k)const =0;
59+ virtual void double_up ()const {int old_size=cache.size ();
60+ int new_size=old_size==0 ?1 :2 *old_size;
61+ if (auto deg=degree ()){int max_idx=*deg-cache_offset;
62+ if (max_idx<new_size){new_size=max_idx+1 ;
63+ all_cached=true ;}}
64+ cache.resize (new_size);
65+ for (int i=old_size;i<new_size;i++){cache[i]=coeff_lazy (cache_offset+i);}}
66+ virtual T coeff (int k)const {if (!initialized){auto v=valuation ();
67+ cache_offset=v?*v:0 ;
68+ initialized=true ;}
69+ int idx=k-cache_offset;
70+ if (idx<0 ){return T (0 );}
71+ if (all_cached&&idx>=(int )cache.size ()){return T (0 );}
72+ while (idx>=(int )cache.size ()&&!all_cached){double_up ();}
73+ if (idx<(int )cache.size ()){return cache[idx];}
74+ return T (0 );}
75+ T get (int k)const {return coeff (k);}};
76+ template <typename T>
77+ struct constant_provider :provider<T>{T value;
78+ int offset;
79+ constant_provider (T value,int offset=0 ):value(value),offset(offset){}
80+ std::optional<int >valuation ()const override {return value!=T (0 )?std::optional<int >(offset):std::nullopt ;}
81+ std::optional<int >degree ()const override {return value!=T (0 )?std::optional<int >(offset):std::nullopt ;}
82+ T coeff_lazy (int k)const override {return k==offset?value:T (0 );}
83+ T coeff (int k)const override {return coeff_lazy (k);}};
84+ template <typename T>
85+ struct polynomial_provider :provider<T>{polynomial_provider(big_vector<T>coeffs,int offset=0 ){auto non_zero=[](const T&x){return x!=T (0 );};
86+ auto first=std::ranges::find_if (coeffs,non_zero);
87+ auto last=std::ranges::find_if (coeffs|std::views::reverse,non_zero);
88+ if (first!=coeffs.end ()){int leading_zeros=first-coeffs.begin ();
89+ int trailing_zeros=last-coeffs.rbegin ();
90+ coeffs=big_vector<T>(first,coeffs.end ()-trailing_zeros);
91+ offset+=leading_zeros;}else {coeffs.clear ();}
92+ this ->cache =std::move (coeffs);
93+ this ->cache_offset =offset;
94+ this ->initialized =true ;
95+ this ->all_cached =true ;}
96+ std::optional<int >valuation ()const override {if (this ->cache .empty ())return std::nullopt ;
97+ return this ->cache_offset ;}
98+ std::optional<int >degree ()const override {if (this ->cache .empty ())return std::nullopt ;
99+ return this ->cache_offset +(int )this ->cache .size ()-1 ;}
100+ T coeff_lazy (int k)const override {int idx=k-this ->cache_offset ;
101+ if (idx<0 ||idx>=(int )this ->cache .size ()){return T (0 );}
102+ return this ->cache [idx];}
103+ T coeff (int k)const override {return coeff_lazy (k);}};
104+ template <typename T>
105+ struct unary_provider :provider<T>{std::shared_ptr<provider<T>>operand;
106+ unary_provider (std::shared_ptr<provider<T>>operand):operand(std::move(operand)){}
107+ virtual T transform (T const &a)const =0;
108+ std::optional<int >valuation ()const override {return operand->valuation ();}
109+ std::optional<int >degree ()const override {return operand->degree ();}
110+ T coeff_lazy (int k)const override {return transform (operand->coeff_lazy (k));}
111+ T coeff (int k)const {return transform (operand->coeff (k));}};
112+ template <typename T>
113+ struct binary_provider :provider<T>{enum class bound_type {valuation,degree};
114+ std::shared_ptr<provider<T>>lhs,rhs;
115+ binary_provider (std::shared_ptr<provider<T>>lhs,std::shared_ptr<provider<T>>rhs):lhs(std::move(lhs)),rhs(std::move(rhs)){}
116+ virtual T combine (T const &a,T const &b)const =0;
117+ std::optional<int >valuation ()const override {auto lv=lhs->valuation ();
118+ auto rv=rhs->valuation ();
119+ if (!lv||!rv)return std::nullopt ;
120+ if (*lv!=*rv){return std::min (*lv,*rv);}
121+ T val=combine (lhs->coeff_lazy (*lv),rhs->coeff_lazy (*lv));
122+ if (val!=T (0 )){return *lv;}
123+ return std::nullopt ;}
124+ std::optional<int >degree ()const override {auto ld=lhs->degree ();
125+ auto rd=rhs->degree ();
126+ if (!ld||!rd)return std::nullopt ;
127+ if (*ld!=*rd){return std::max (*ld,*rd);}
128+ T val=combine (lhs->coeff_lazy (*ld),rhs->coeff_lazy (*ld));
129+ if (val!=T (0 )){return *ld;}
130+ return std::nullopt ;}
131+ T coeff_lazy (int k)const override {return combine (lhs->coeff_lazy (k),rhs->coeff_lazy (k));}
132+ T coeff (int k)const {return combine (lhs->coeff (k),rhs->coeff (k));}};
133+ template <typename T>
134+ struct add_provider :binary_provider<T>{using binary_provider<T>::binary_provider;
135+ T combine (T const &a,T const &b)const override {return a+b;}};
136+ template <typename T>
137+ struct subtract_provider :binary_provider<T>{using binary_provider<T>::binary_provider;
138+ T combine (T const &a,T const &b)const override {return a-b;}};
139+ template <typename T>
140+ struct negate_provider :unary_provider<T>{using unary_provider<T>::unary_provider;
141+ T transform (T const &a)const override {return -a;}};
142+ template <typename T>
143+ struct scale_provider :unary_provider<T>{T scalar;
144+ scale_provider (std::shared_ptr<provider<T>>operand,T scalar):unary_provider<T>(std::move(operand)),scalar(scalar){}
145+ T transform (T const &a)const override {return a*scalar;}};
146+ template <typename T>
147+ struct multiply_provider :provider<T>{std::shared_ptr<provider<T>>lhs,rhs;
148+ multiply_provider (std::shared_ptr<provider<T>>lhs,std::shared_ptr<provider<T>>rhs):lhs(std::move(lhs)),rhs(std::move(rhs)){}
149+ std::optional<int >valuation ()const override {auto lv=lhs->valuation ();
150+ auto rv=rhs->valuation ();
151+ if (lv&&rv)return *lv+*rv;
152+ return std::nullopt ;}
153+ std::optional<int >degree ()const override {auto ld=lhs->degree ();
154+ auto rd=rhs->degree ();
155+ if (ld&&rd)return *ld+*rd;
156+ return std::nullopt ;}
157+ T coeff_lazy (int k)const override {T result=T (0 );
158+ for (int j=0 ;j<=k;j++){result+=lhs->coeff (j)*rhs->coeff (k-j);}
159+ return result;}};
160+ template <typename T>
161+ struct laurent {std::shared_ptr<provider<T>>impl;
162+ laurent ():impl(std::make_shared<constant_provider<T>>(T(0 ),0 )){}
163+ laurent (T value,int offset=0 ):impl(std::make_shared<constant_provider<T>>(value,offset)){}
164+ laurent (big_vector<T>coeffs,int offset=0 ):impl(std::make_shared<polynomial_provider<T>>(std::move(coeffs),offset)){}
165+ laurent (std::shared_ptr<provider<T>>impl):impl(std::move(impl)){}
166+ T operator [](int k)const {return impl->get (k);}
167+ laurent operator -()const {return std::make_shared<negate_provider<T>>(impl);}
168+ laurent operator +(const laurent&other)const {return std::make_shared<add_provider<T>>(impl,other.impl );}
169+ laurent operator -(const laurent&other)const {return std::make_shared<subtract_provider<T>>(impl,other.impl );}
170+ laurent operator *(const laurent&other)const {return std::make_shared<multiply_provider<T>>(impl,other.impl );}
171+ laurent&operator +=(const laurent&other){return *this =*this +other;}
172+ laurent&operator -=(const laurent&other){return *this =*this -other;}
173+ laurent&operator *=(const laurent&other){return *this =*this *other;}
174+ laurent operator *(T const &scalar)const {return std::make_shared<scale_provider<T>>(impl,scalar);}
175+ laurent&operator *=(T const &scalar){return *this =*this *scalar;}};
176+ template <typename T>
177+ laurent<T>operator *(T const &scalar,laurent<T>const &series){return series*scalar;}}
0 commit comments