Codeforces 919E - Congruence Equation【数论-欧拉降幂】

E. Congruence Equation
time limit per test 3 seconds memory limit per test 256 megabytes

Problem Description

After trying hard for many years, Little Q has finally received an astronaut license. To celebrate the fact, he intends to buy himself a spaceship and make an interstellar travel.

Given an integer x. Your task is to find out how many positive integers n (1 ≤ n ≤ x) satisfy $n*a$^$n$$\equiv$$b(mod$ $p$)

where a, b, p are all known constants.

Input

The only line contains four integers a, b, p, x (2 ≤ p ≤ 106 + 3, 1 ≤ a, b < p, 1 ≤ x ≤ 1012). It is guaranteed that p is a prime.

Output

Print a single integer: the number of possible answers n.

Examples

Input
2 3 5 8

Output
2

Input
4 6 7 13

Output
1

Input
233 233 10007 1

Output
1

Note

In the first sample, we can see that n = 2 and n = 8 are possible answers.

【题意】

求出【1，x】中满足$n*a$^$n$$\equiv$$b(mod$ $p$)的$n$的个数。

【思路】

很容易知道$n$ mod $p$的循环节为$p$

再根据费马小定理：若$c$为质数，$a$^{$p-1$$\equiv$$1$(mod $p$),可以得到$a$$n$的循环节为$p-1$}

那么我们可以得到$n*a$^$n$的循环节即为$P$=$p*(p-1)$

因为$a$^{$n$的循环节为$p-1$，我们可以令 $n=i*(p-1)+j$ $①$}

然后我们根据欧拉降幂公式：若$c$为质数，$a$^{$b$ mod $c$$\equiv$ $a$b%$($c$-$1$)$ (mod $c$)}

那么$n*a$^$n$$\equiv$$n*a$^{$n$%($p$-$1$)}$\equiv$$b$(mod $p$)

所以$n*a$^$j$$\equiv$$b$(mod $p$) $③$

$i*(p-1)+j$$\equiv$$b*a$^$-j$(mod $p$)

$i=j-b*a$^$-j$ $②$

于是我们就可以通过在【$0$，$p-2$】范围内枚举$j$，然后通过上面的②再算出$i$代入①得到一个n，然后再判断+$P$、+$2P$…..后是否在x范围内即可统计个数。

#include <cstdio>
#include <ctime>
#include <cstdlib>
#include <queue>
#include <cstring>
#include <algorithm>
using namespace std;
#define mst(a,b) memset((a),(b),sizeof(a))
#define rush() int T;scanf("%d",&T);while(T--)

typedef long long ll;
const int maxn = 200005;
const ll mod = 1e9+7;
const int INF = 0x3f3f3f3f;
const double eps = 1e-9;

ll a,b,p,x;

ll fast_mod(ll a,ll n,ll Mod)
{
    ll ans=1;
    a%=Mod;
    while(n)
    {
        if(n&1) ans=(ans*a)%Mod;
        a=(a*a)%Mod;
        n>>=1;
    }
    return ans;
}

int main()
{
    scanf("%lld%lld%lld%lld",&a,&b,&p,&x);
    ll ans=0;
    ll Mod=p*(p-1);
    for(ll j=0;j<p-1;j++)
    {
        ll inv=fast_mod(fast_mod(a,j,p),p-2,p);  //a^j的逆元
        ll tmp=inv*b%p;
        ll i=(j-tmp+p)%p;
        ll n=(i*(p-1)+j)%Mod;
        if(n==0) n=Mod;
        ans+=x/Mod+(x%Mod>=n);
    }
    printf("%lld\n",ans);
}