fork download
  1. #include<iostream>
  2. #include<cstdlib>
  3. #include<cstring>
  4. #include<algorithm>
  5. #include<cmath>
  6. #include<stdio.h>
  7. #include<vector>
  8. #include<set>
  9. using namespace std;
  10. const int mod=1e9+7;
  11. typedef long long LL;
  12. int power(int a,int n)
  13. {
  14. if(n==0)
  15. return 1;
  16. int ret=power(a,n/2);
  17. ret=1LL*ret*ret%mod;
  18. if(n%2)
  19. ret=1LL*ret*a%mod;
  20. return ret;
  21. }
  22.  
  23. #define inf 0x3f3f3f3f
  24. #define mem(a,b) memset(a,b,sizeof(a))
  25. #define REP(i,a,b) for(int i=a; i<=b; ++i)
  26. #define PER(i,a,b) for(int i=a; i>=b; --i)
  27. #define MP make_pair
  28. #define PB push_back
  29. #define fi first
  30. #define se second
  31. typedef pair<int,int> pii;
  32. int ex[1101],fac[1101],f[2101],inv[1101];
  33. const double PI = acos(-1.0);
  34. const int Bin = 11;
  35. const int P=1e9+7;
  36. const int maxn = (1<<Bin);
  37. const int md=1e9+7;
  38. const int mmd =(1<<15);
  39. vector<int> vec[12];
  40. struct Comp {
  41. double re, im;
  42. Comp(): re(0), im(0) {}
  43. Comp(const double _re, const double _im) { re = _re; im = _im; }
  44. Comp operator+(const Comp &a) { return Comp(re+a.re, im+a.im); }
  45. Comp operator-(const Comp &a) { return Comp(re-a.re, im-a.im); }
  46. Comp operator*(const Comp &a) { return Comp(re*a.re-im*a.im, re*a.im+im*a.re); }
  47. Comp operator*(const double &a) { return Comp(re*a, im*a); }
  48. Comp operator/(const double &a) { return Comp(re/a, im/a); }
  49. void operator*=(const Comp &a) { (*this) = (*this) * a; }
  50. void operator+=(const Comp &a) { re += a.re; im += a.im; }
  51. void operator*=(const double &a) { re*=a; im*=a; }
  52. void operator/=(const double &a) { re/=a; im/=a; }
  53. Comp conj() { return Comp(re, -im); }
  54. };
  55.  
  56. vector<Comp> ww[Bin + 2];
  57. void Init() {
  58. for(int i = 1, l = 2; l <= maxn; ++i, l <<= 1) {
  59. for(int k = 0; k < l/2; ++k)
  60. ww[i].push_back(Comp(cos(2*PI*k/l),sin(2*PI*k/l)));
  61. }
  62. }
  63.  
  64. void FFT(Comp *a, int len, int type) {
  65. int i, j, k, l, tt;
  66. for(i = 1, j = len>>1; i < len-1; ++i) {
  67. if(i < j) swap(a[i], a[j]);
  68. k = len>>1;
  69. while(j >= k) j -= k, k >>= 1;
  70. j += k;
  71. }
  72. Comp var, u, v;
  73. for(l = 2, tt = 1; l <= len; l <<= 1, ++tt) {
  74. for(k = 0; k < len; k += l) {
  75. for(i = k; i < k+l/2; ++i) {
  76. var = ww[tt][i-k];
  77. if(type == -1) var.im = -var.im;
  78. u = a[i], v = var * a[i+l/2];
  79. a[i] = u+v, a[i+l/2] = u-v;
  80. }
  81. }
  82. }
  83. if(type == -1) for(i = 0; i < len; ++i) a[i] /= len;
  84. }
  85. Comp A[maxn + 5], B[maxn + 5], C[maxn + 5], DD[maxn + 5];
  86. int aa[maxn*2],bb[maxn*2],cc[maxn*2],dd[maxn*2];
  87. void Conv(int *a, int *b, int *c, int len) {
  88. for(int i = 0; i < len; ++i) {
  89. A[i] = Comp(a[i]>>15, b[i]>>15);
  90. B[i] = Comp(a[i]&(mmd-1), b[i]&(mmd-1));
  91. }
  92. FFT(A, len, 1); FFT(B, len, 1);
  93. for(int i = 0; i < len; ++i) {
  94. int j = (len-i) & (len-1);
  95. C[i] = (A[i]*A[i] - (A[j]*A[j]).conj()) * Comp(0, -0.25)+(A[i]*B[i] - (A[j]*B[j]).conj()) * Comp(0.5, 0);
  96. DD[i] = (B[i]*B[i] - (B[j]*B[j]).conj()) * Comp(0, -0.25);
  97. }
  98. FFT(C, len, -1); FFT(DD, len, -1);
  99. for(int i = 0; i < len; ++i) {
  100. c[i] = (((LL(C[i].re+0.5) %mod)<<30)%mod + ((LL(C[i].im+0.5)%mod ) <<15)%mod+ LL(DD[i].re+0.5))%mod ;
  101. }
  102. }
  103. inline int read(){
  104. int f=1,x=0;char ch=getchar();
  105. while (ch<'0'||ch>'9'){if (ch=='-')f=-1;ch=getchar();}
  106. while (ch>='0'&&ch<='9'){x=(x<<3)+(x<<1)+ch-'0';ch=getchar();}
  107. return f*x;
  108. }
  109.  
  110. inline int qpow(int x,int p){int ans=1;
  111. for(;p;p>>=1,x=1ll*x*x%P)if(p&1)ans=1ll*ans*x%P;
  112. return ans;
  113. }
  114. int temp[maxn],a[maxn],b[maxn],c[maxn],d[maxn],tmp[maxn];
  115. inline void GetInv(int *a,int *b,int n){
  116. if (n==1) return void(b[0]=qpow(a[0],P-2));
  117. GetInv(a,b,n>>1);
  118. for (int i=0;i<n;++i)tmp[i]=a[i],tmp[n+i]=0;
  119. int k;for (k=1;k<=n;k<<=1);
  120. for (int i=n;i<k;++i)b[i]=0;
  121. Conv(tmp,b,tmp,k);
  122. for(int i=n;i<k;i++)
  123. tmp[i]=0;
  124. Conv(tmp,b,tmp,k);
  125. for (int i=0;i<n;++i)b[i]=(2*b[i]%mod-tmp[i]+mod)%mod,b[n+i]=0;
  126. }
  127. inline void Getln(int *a,int *b,int *c,int n){
  128. for (int i=0;i<n;++i)b[i]=1ll*(i+1)*a[i+1]%P;
  129. GetInv(a,c,n);
  130. for (int i=0;i<n;++i)tmp[i]=c[i],tmp[n+i]=0;
  131. int k;for (k=1;k<=n;k<<=1);
  132. Conv(tmp,b,tmp,k);
  133. for (int i=1;i<=n;++i)b[i]=1ll*tmp[i-1]*qpow(i,P-2)%P,b[i+n]=0;b[0]=0;
  134. }
  135.  
  136. inline void Getexp(int *a,int *b,int *c,int *d,int n){
  137. if (n==1)return void(b[0]=1);
  138. Getexp(a,b,c,d,n>>1);
  139. Getln(b,c,d,n);
  140. for (int i=0;i<n;++i)tmp[i]=(a[i]-c[i]+P)%P,tmp[n+i]=0;++tmp[0];
  141. int k;for (k=1;k<=n;k<<=1);
  142. Conv(tmp,b,tmp,k);
  143. for (int i=0;i<n;++i)b[i]=tmp[i],b[i+n]=0;
  144. }
  145. void got(int n,int deep)
  146. {
  147. int i;
  148. if(n==1)
  149. {
  150. vec[deep].push_back(1);
  151. vec[deep].push_back(1);
  152. return;
  153. }
  154. got(n>>1,deep+1);
  155. vec[deep].resize(n+1);
  156. int d=n>>1;
  157. for(i=0;i<=n;i++)
  158. {
  159. if(i==0)
  160. ex[i]=1;
  161. else
  162. ex[i]=1LL*ex[i-1]*d%mod;
  163. }
  164. int len=1;
  165. while(len<=d+d)
  166. len<<=1;
  167. for(i=0;i<len;i++)
  168. {
  169. if(i<=(n>>1))
  170. {
  171. aa[i]=1LL*ex[(n>>1)-i]*inv[(n>>1)-i]%mod;
  172. bb[i]=1LL*fac[i]*vec[deep+1][i]%mod;
  173. }
  174. else
  175. bb[i]=aa[i]=0;
  176. }
  177. Conv(aa,bb,cc,len);
  178. for(i=0;i<=d;i++)
  179. cc[i]=1LL*cc[i+(n>>1)]*inv[i]%mod;
  180. for(i=d+1;i<len;i++)
  181. cc[i]=0;
  182. for(i=0;i<len;i++)
  183. {
  184. if(i<=(n>>1))
  185. aa[i]=vec[deep+1][i];
  186. else
  187. aa[i]=0;
  188. }
  189. Conv(aa,cc,bb,len);
  190. if(n&1)
  191. {
  192. for(i=0;i<=n;i++)
  193. {
  194. vec[deep][i]=1LL*bb[i]*n%mod;
  195. vec[deep][i]=(vec[deep][i]+bb[i-1])%mod;
  196. }
  197. }
  198. else
  199. {
  200. for(i=0;i<=n;i++)
  201. vec[deep][i]=bb[i];
  202. }
  203. }
  204. int CC(int n,int m)
  205. {
  206. return 1LL*fac[n]*inv[m]%mod*inv[n-m]%mod;
  207. }
  208. long long S,T;
  209. int n,m;
  210. int main()
  211. {
  212. Init();
  213. scanf("%lld%lld%d%d",&S,&T,&n,&m);
  214. int i,j,s,p,q;
  215. for(i=0;i<=m-n+1;i++)
  216. {
  217. if(i==0)
  218. fac[i]=1;
  219. else
  220. fac[i]=1LL*fac[i-1]*i%mod;
  221. }
  222. inv[m-n+1]=power(fac[m-n+1],mod-2);
  223. for(i=m-n;i>=0;i--)
  224. inv[i]=1LL*inv[i+1]*(i+1)%mod;
  225. int len=1;
  226. while(len<=m-n)
  227. len<<=1;
  228. got(m-n,0);
  229. S=(S-(m-n)+mod)%mod;
  230. T%=mod;
  231. memset(a,0,sizeof(a));
  232. memset(b,0,sizeof(b));
  233. for(i=1;i<=m-n;i++)
  234. b[i]=inv[i+1];
  235. b[0]=1;
  236. GetInv(b,a,len);
  237. for(i=0;i<=m-n+1;i++)
  238. {
  239. if(i==0)
  240. ex[i]=1;
  241. else
  242. ex[i]=1LL*ex[i-1]*(T+1)%mod;
  243. }
  244. for(i=0;i<=m-n;i++)
  245. {
  246. b[i]=1LL*(ex[i+1]-1)*inv[i+1]%mod;
  247. if(b[i]<0)
  248. b[i]+=mod;
  249. }
  250. for(i=m-n+1;i<len<<1;i++)
  251. a[i]=b[i]=0;
  252. Conv(a,b,f,len<<1);
  253. for(i=0;i<=m-n;i++)
  254. {
  255. f[i]=1LL*f[i]*fac[i]%mod;
  256. if(f[i]<0)
  257. f[i]+=mod;
  258. }
  259. memset(a,0,sizeof(a));
  260. memset(b,0,sizeof(b));
  261. memset(c,0,sizeof(c));
  262. int x,y;
  263. for(i=0;i<=m-n;i++)
  264. {
  265. a[i]=1LL*f[i]*inv[i]%mod;
  266. if(i%2)
  267. a[i]=(mod-a[i])%mod;
  268. if(i==0)
  269. {
  270. x=power(a[i],mod-2);
  271. y=power(a[i],n);
  272. }
  273. a[i]=1LL*a[i]*x%mod;
  274. }
  275. Getln(a,b,c,len);
  276. for(i=0;i<=m-n;i++)
  277. b[i]=1LL*b[i]*n%mod;
  278. for(i=m-n+1;i<len;i++)
  279. b[i]=0;
  280. memset(a,0,sizeof(a));
  281. memset(c,0,sizeof(c));
  282. memset(d,0,sizeof(d));
  283. Getexp(b,a,c,d,len);
  284. for(i=0;i<=m-n;i++)
  285. a[i]=1LL*a[i]*y%mod;
  286. int awp=1;
  287. for(i=0;i<=m-n;i++)
  288. {
  289. b[i]=1LL*awp*inv[i]%mod;
  290. awp=1LL*awp*S%mod;
  291. }
  292. for(i=m-n+1;i<len<<1;i++)
  293. a[i]=b[i]=0;
  294. Conv(a,b,c,len<<1);
  295. int ans=0;
  296. for(i=0;i<=m-n;i++)
  297. ans=(ans+1LL*c[i]*vec[0][i]%mod*fac[i])%mod;
  298. ans=1LL*ans*inv[m-n]%mod;
  299. printf("%d\n",ans);
  300. return 0;
  301. }
  302.  
Success #stdin #stdout 0s 15512KB
stdin
1000000000000000000 100000 1000000000 1000001000
stdout
846395340