#include <bits/stdc++.h>
using namespace std;
using ll = long long;
using vi = vector<int>;
using ii = pair<int, int>;
static ll P;
template <typename T>
struct intmod {
  intmod() {}
  constexpr intmod(T t) : x((t + P) % P) {}
  T value() { return x; }
  bool operator!=(const intmod<T> i) { return x != i.x; }
  bool operator==(const intmod<T> i) { return x == i.x; }
  intmod<T> &operator+=(const intmod<T> i) {
    x = (x + i.x) % P;
    return *this;
  }
  intmod<T> &operator-=(const intmod<T> i) {
    x = (x + P - i.x) % P;
    return *this;
  }
  intmod<T> &operator*=(const intmod<T> i) {
    x = ((long long)x * i.x) % P;
    return *this;
  }
  intmod<T> &operator/=(const intmod<T> i) {
    x = ((long long)x * i.inverse().x) % P;
    return *this;
  }
  intmod<T> operator+(const intmod<T> i) const {
    auto j = *this;
    return j += i;
  }
  intmod<T> operator-(const intmod<T> i) const {
    auto j = *this;
    return j -= i;
  }
  intmod<T> operator*(const intmod<T> i) const {
    auto j = *this;
    return j *= i;
  }
  intmod<T> operator/(const intmod<T> i) const {
    auto j = *this;
    return j /= i;
  }
  intmod<T> operator-() const {
    intmod<T> n;
    n.x = P - x;
    return n;
  }
  intmod<T> inverse() const {
    T a = x, b = P;
    T aa = 1, ab = 0;
    T ba = 0, bb = 1;
    while (a) {
      T q = b / a;
      T r = b % a;
      ba -= aa * q;
      bb -= ab * q;
      swap(ba, aa);
      swap(bb, ab);
      b = a;
      a = r;
    }
    return intmod<T>(aa) + intmod<T>(ba);
  }
  static const intmod<T> zero;
  static const intmod<T> unity;

 private:
  T x;
};
using rem = intmod<ll>;
static ll x;
static ll ans;
static ll count(ll first, ll period) {
  ll c = 0;
  if (first > x)
    c = 0;
  else {
    ll over = x - first + 1;
    c = (over + period - 1) / period;
  }
  ((void)0);
  return c;
}
template <typename T>
static T egcd(T a, T b, T &n, T &m) {
  T aa, ab, ba, bb;
  ab = ba = 0;
  aa = bb = 1;
  while (a) {
    T q = b / a;
    T r = b % a;
    ba -= aa * q;
    bb -= ab * q;
    swap(aa, ba);
    swap(ab, bb);
    b = a;
    a = r;
  }
  n = aa + ba;
  m = ab + bb;
  return b;
}
static void add_pair(ll n, ll phi, ll o) {
  ((void)0);
  rem j = (rem(n) - rem(o)) / rem(phi);
  ans += count(j.value() * phi + o, phi * P);
}
int main() {
  ll _a, _b;
  cin >> _a >> _b >> P >> x;
  rem a(_a), b(_b);
  static rem as[1000100];
  static ll index[1000100];
  ll i;
  ll phi = 0;
  fill(index, index + P, -1);
  as[0] = rem(1);
  index[1] = 0;
  for ((i) = (1); (i) < (P); (i)++) {
    phi++;
    as[i] = as[i - 1] * a;
    ((void)0);
    if (as[i] == as[0]) break;
    index[as[i].value()] = i;
  }
  ((void)0);
  assert(as[phi] == as[0]);
  ll quo = x / P;
  ll res = x % P;
  ((void)0);
  ((void)0);
  ll n;
  for ((n) = (1); (n) < (P); (n)++) {
    rem apn = b / rem(n);
    ll pos = index[apn.value()];
    if (pos < 0) continue;
    add_pair(n, phi, pos);
  }
  cout << ans << endl;
  return 0;
}

### Prompt

In cpp, your task is to solve the following problem:

Given an integer x. Your task is to find out how many positive integers n (1 ≤ n ≤ x) satisfy $$$n ⋅ a^n ≡ b   (mod\;p), where a, b, p$$$ are all known constants.

Input

The only line contains four integers a,b,p,x (2 ≤ p ≤ 10^6+3, 1 ≤ a,b < p, 1 ≤ x ≤ 10^{12}). 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.

### Response

```cpp
#include <bits/stdc++.h>
using namespace std;
using ll = long long;
using vi = vector<int>;
using ii = pair<int, int>;
static ll P;
template <typename T>
struct intmod {
  intmod() {}
  constexpr intmod(T t) : x((t + P) % P) {}
  T value() { return x; }
  bool operator!=(const intmod<T> i) { return x != i.x; }
  bool operator==(const intmod<T> i) { return x == i.x; }
  intmod<T> &operator+=(const intmod<T> i) {
    x = (x + i.x) % P;
    return *this;
  }
  intmod<T> &operator-=(const intmod<T> i) {
    x = (x + P - i.x) % P;
    return *this;
  }
  intmod<T> &operator*=(const intmod<T> i) {
    x = ((long long)x * i.x) % P;
    return *this;
  }
  intmod<T> &operator/=(const intmod<T> i) {
    x = ((long long)x * i.inverse().x) % P;
    return *this;
  }
  intmod<T> operator+(const intmod<T> i) const {
    auto j = *this;
    return j += i;
  }
  intmod<T> operator-(const intmod<T> i) const {
    auto j = *this;
    return j -= i;
  }
  intmod<T> operator*(const intmod<T> i) const {
    auto j = *this;
    return j *= i;
  }
  intmod<T> operator/(const intmod<T> i) const {
    auto j = *this;
    return j /= i;
  }
  intmod<T> operator-() const {
    intmod<T> n;
    n.x = P - x;
    return n;
  }
  intmod<T> inverse() const {
    T a = x, b = P;
    T aa = 1, ab = 0;
    T ba = 0, bb = 1;
    while (a) {
      T q = b / a;
      T r = b % a;
      ba -= aa * q;
      bb -= ab * q;
      swap(ba, aa);
      swap(bb, ab);
      b = a;
      a = r;
    }
    return intmod<T>(aa) + intmod<T>(ba);
  }
  static const intmod<T> zero;
  static const intmod<T> unity;

 private:
  T x;
};
using rem = intmod<ll>;
static ll x;
static ll ans;
static ll count(ll first, ll period) {
  ll c = 0;
  if (first > x)
    c = 0;
  else {
    ll over = x - first + 1;
    c = (over + period - 1) / period;
  }
  ((void)0);
  return c;
}
template <typename T>
static T egcd(T a, T b, T &n, T &m) {
  T aa, ab, ba, bb;
  ab = ba = 0;
  aa = bb = 1;
  while (a) {
    T q = b / a;
    T r = b % a;
    ba -= aa * q;
    bb -= ab * q;
    swap(aa, ba);
    swap(ab, bb);
    b = a;
    a = r;
  }
  n = aa + ba;
  m = ab + bb;
  return b;
}
static void add_pair(ll n, ll phi, ll o) {
  ((void)0);
  rem j = (rem(n) - rem(o)) / rem(phi);
  ans += count(j.value() * phi + o, phi * P);
}
int main() {
  ll _a, _b;
  cin >> _a >> _b >> P >> x;
  rem a(_a), b(_b);
  static rem as[1000100];
  static ll index[1000100];
  ll i;
  ll phi = 0;
  fill(index, index + P, -1);
  as[0] = rem(1);
  index[1] = 0;
  for ((i) = (1); (i) < (P); (i)++) {
    phi++;
    as[i] = as[i - 1] * a;
    ((void)0);
    if (as[i] == as[0]) break;
    index[as[i].value()] = i;
  }
  ((void)0);
  assert(as[phi] == as[0]);
  ll quo = x / P;
  ll res = x % P;
  ((void)0);
  ((void)0);
  ll n;
  for ((n) = (1); (n) < (P); (n)++) {
    rem apn = b / rem(n);
    ll pos = index[apn.value()];
    if (pos < 0) continue;
    add_pair(n, phi, pos);
  }
  cout << ans << endl;
  return 0;
}

```