Problem: Everyone knows what the Fibonacci sequence is. This sequence can be defined by the recurrence relation: 

F1 = 1, F2 = 2, Fi = Fi - 1 + Fi - 2 (i > 2).

We'll define a new number sequence Ai(k) by the formula: 

Ai(k) = Fi × ik (i ≥ 1).

In this problem, your task is to calculate the following sum: A1(k) + A2(k) + ... + An(k). The answer can be very large, so print it modulo 1000000007 (109 + 7).

Input

The first line contains two space-separated integers n, k (1 ≤ n ≤ 1017; 1 ≤ k ≤ 40).

Output

Print a single integer — the sum of the first n elements of the sequence Ai(k) modulo 1000000007 (109 + 7).

Examples

Input

1 1


Output

1


Input

4 1


Output

34


Input

5 2


Output

316


Input

7 4


Output

73825

Solution (in CPP):

Solution:
#include <bits/stdc++.h>
using namespace std;
long long int mod = 1e9 + 7;
long long int exp(long long int x, long long int y) {
  if (y < 0) return 0;
  long long int ret = 1;
  for (; y; y >>= 1, x = (x * x) % mod) {
    if (y & 1) ret = (ret * x) % mod;
  }
  return ret;
}
long long int pexp(long long int x, long long int y) {
  long long int ret = 1;
  for (; y; y >>= 1, x = (x * x)) {
    if (y & 1) ret = (ret * x);
  }
  return ret;
}
long long int gcd(long long int x, long long int y) {
  if (!x || !y) return x + y;
  return x % y == 0 ? y : gcd(y, x % y);
}
long long int lcm(long long int x, long long int y) {
  return x * (y / gcd(x, y));
}
long long int bsum(long long int u, long long int b) {
  long long int ret = 0;
  if (u < b) return u;
  while (u) {
    ret += u % b;
    u /= b;
  }
  return ret;
}
long long int prival(long long int u, long long int p) {
  int cn = 0;
  while (u % p == 0) {
    cn++;
    u = u / p;
  }
  return cn;
}
long long int minv(long long int a, long long int b) {
  return 1 < a ? b - minv(b % a, a) * b / a : 1;
}
long long int extm(long long int a, long long int b) {
  long long int ret = 0;
  while (a != 0) {
    if (a % 2 == 1) {
      ret += b;
      ret %= mod;
    }
    a >>= 1;
    b = (2 * b) % mod;
  }
  return ret;
}
long long int eaphi(long long int x) {
  long long int t = x, ret = x, i;
  for (i = 2; i * i <= x; i++) {
    if (t % i == 0) {
      ret -= ret / i;
      while (t % i == 0) t /= i;
    }
  }
  if (t != 1) ret -= ret / t;
  return ret;
}
bool isprime(long long int x) {
  if (x <= 1) return false;
  for (long long int i = 2; i * i <= x; i++) {
    if (x % i == 0) {
      return false;
    }
  }
  return true;
}
struct mat {
  long long int a[99][99];
};
mat MM, RES;
long long int n, k, ans;
int sz;
long long int vec[999], fac[1111], invfac[1111];
mat mul(mat U, mat V) {
  mat RET;
  int i, j, k;
  for (i = 1; i <= sz; i++) {
    for (j = 1; j <= sz; j++) {
      RET.a[i][j] = 0;
    }
  }
  for (i = 1; i <= sz; i++) {
    for (j = 1; j <= sz; j++) {
      for (k = 1; k <= sz; k++) {
        RET.a[i][j] += (U.a[i][k] * V.a[k][j]) % mod;
      }
      RET.a[i][j] %= mod;
    }
  }
  return RET;
}
mat pow(mat U, long long int po) {
  mat RET;
  int i, j;
  for (i = 1; i <= sz; i++) {
    for (j = 1; j <= sz; j++) {
      if (i == j)
        RET.a[i][j] = 1;
      else
        RET.a[i][j] = 0;
    }
  }
  while (po) {
    if (po % 2 == 0) {
      po = po / 2;
      U = mul(U, U);
    } else if (po % 2 == 1) {
      po--;
      RET = mul(RET, U);
    }
  }
  return RET;
}
long long int ncr(long long int u, long long int v) {
  if (u < 0 || v < 0 || u < v) return 0;
  long long int ret = fac[u];
  ret = (ret * invfac[v]) % mod;
  ret = (ret * invfac[u - v]) % mod;
  return ret;
}
int main(void) {
  long long int i, j;
  cin >> n >> k;
  if (n == 1) {
    printf("1");
    return 0;
  }
  sz = 2 * k + 3;
  fac[0] = 1;
  invfac[0] = 1;
  for (i = 1; i <= 100; i++) {
    fac[i] = (i * fac[i - 1]) % mod;
    invfac[i] = exp(fac[i], mod - 2);
  }
  for (i = 1; i <= k + 1; i++) {
    vec[i] = exp(2, i) + 1;
  }
  for (i = k + 2; i <= 2 * k + 3; i++) {
    vec[i] = 1;
  }
  for (i = 1; i <= k + 1; i++) {
    for (j = 1; j <= i; j++) {
      MM.a[i][j] = ncr(i - 1, j - 1);
    }
    for (j = k + 1 + 1; j <= k + 1 + i; j++) {
      MM.a[i][j] = (exp(2, i + k + 1 - j) * ncr(i - 1, j - k - 2)) % mod;
    }
    MM.a[i][2 * k + 3] = exp(2, i - 1) + 1;
  }
  for (i = k + 2; i <= 2 * k + 2; i++) {
    MM.a[i][i - k - 1] = 1;
  }
  MM.a[2 * k + 3][2 * k + 3] = 1;
  RES = pow(MM, n - 2);
  for (i = 1; i <= 2 * k + 3; i++) {
    ans += (vec[i] * RES.a[k + 1][i]) % mod;
  }
  ans = ans % mod;
  cout << ans;
  return 0;
}