#include <bits/stdc++.h>
using namespace std;
template <class T>
inline T abs(T a) {
  return a > 0 ? a : -a;
}
int n;
int m;
const long long mod = 1e9 + 7;
long long ans = 0;
int G[4001][4001];
int gcd(int a, int b) {
  if (G[a][b] >= 0) return G[a][b];
  return (G[a][b] = (b ? gcd(b, a % b) : a));
}
long long brute(long long w, long long h) {
  long long ans = 0;
  for (int x1 = 0; x1 < (w + 1); x1++)
    for (int hdsfjkhgjlsdfhgsdf = 0; hdsfjkhgjlsdfhgsdf < (h + 1);
         hdsfjkhgjlsdfhgsdf++)
      for (int x2 = 0; x2 < (w + 1); x2++)
        for (int y2 = 0; y2 < (h + 1); y2++)
          for (int x3 = 0; x3 < (w + 1); x3++)
            for (int y3 = 0; y3 < (h + 1); y3++) {
              int s = (x2 - x1) * (y3 - hdsfjkhgjlsdfhgsdf) -
                      (x3 - x1) * (y2 - hdsfjkhgjlsdfhgsdf);
              if (s % 2 == 0 && s) {
                ans += 1;
              }
            }
  return ans;
}
int main() {
  long long w, h;
  cin >> w >> h;
  long long W[2], H[2];
  W[0] = w / 2 + 1;
  W[1] = (w + 1) / 2;
  H[0] = h / 2 + 1;
  H[1] = (h + 1) / 2;
  for (int x1 = 0; x1 < (2); x1++)
    for (int hdsfjkhgjlsdfhgsdf = 0; hdsfjkhgjlsdfhgsdf < (2);
         hdsfjkhgjlsdfhgsdf++)
      for (int x2 = 0; x2 < (2); x2++)
        for (int y2 = 0; y2 < (2); y2++)
          for (int x3 = 0; x3 < (2); x3++)
            for (int y3 = 0; y3 < (2); y3++) {
              int s = (x2 - x1) * (y3 - hdsfjkhgjlsdfhgsdf) -
                      (x3 - x1) * (y2 - hdsfjkhgjlsdfhgsdf);
              if (s % 2 == 0) {
                ans = (ans + W[x1] * H[hdsfjkhgjlsdfhgsdf] % mod * W[x2] % mod *
                                 H[y2] % mod * W[x3] % mod * H[y3] % mod) %
                      mod;
              }
            }
  long long K = (w + 1) * (h + 1) % mod;
  ans =
      (ans - (K * K % mod * K % mod - K * (K - 1) % mod * (K - 2) % mod)) % mod;
  ans = (ans + mod) % mod;
  memset(G, -1, sizeof(G));
  for (long long i = 0; i <= w; i++) {
    for (long long j = 0; j <= h; j++) {
      long long d = gcd(i, j);
      if (d <= 1) continue;
      long long tmp =
          (i == 0 || j == 0 ? 1 : 2) * (w + 1 - i) * (h + 1 - j) % mod;
      tmp = (tmp * (d - 1)) % mod;
      tmp = (tmp * 6) % mod;
      ans = (ans - tmp) % mod;
      ans = (ans + mod) % mod;
    }
  }
  cout << ans << endl;
  return 0;
}

### Prompt

Please create a solution in cpp to the following problem:

The Little Elephant is playing with the Cartesian coordinates' system. Most of all he likes playing with integer points. The Little Elephant defines an integer point as a pair of integers (x; y), such that 0 ≤ x ≤ w and 0 ≤ y ≤ h. Thus, the Little Elephant knows only (w + 1)·(h + 1) distinct integer points.

The Little Elephant wants to paint a triangle with vertexes at integer points, the triangle's area must be a positive integer. For that, he needs to find the number of groups of three points that form such triangle. At that, the order of points in a group matters, that is, the group of three points (0;0), (0;2), (2;2) isn't equal to the group (0;2), (0;0), (2;2).

Help the Little Elephant to find the number of groups of three integer points that form a nondegenerate triangle with integer area.

Input

A single line contains two integers w and h (1 ≤ w, h ≤ 4000).

Output

In a single output line print an integer — the remainder of dividing the answer to the problem by 1000000007 (109 + 7).

Examples

Input

2 1


Output

36


Input

2 2


Output

240

### Response

```cpp
#include <bits/stdc++.h>
using namespace std;
template <class T>
inline T abs(T a) {
  return a > 0 ? a : -a;
}
int n;
int m;
const long long mod = 1e9 + 7;
long long ans = 0;
int G[4001][4001];
int gcd(int a, int b) {
  if (G[a][b] >= 0) return G[a][b];
  return (G[a][b] = (b ? gcd(b, a % b) : a));
}
long long brute(long long w, long long h) {
  long long ans = 0;
  for (int x1 = 0; x1 < (w + 1); x1++)
    for (int hdsfjkhgjlsdfhgsdf = 0; hdsfjkhgjlsdfhgsdf < (h + 1);
         hdsfjkhgjlsdfhgsdf++)
      for (int x2 = 0; x2 < (w + 1); x2++)
        for (int y2 = 0; y2 < (h + 1); y2++)
          for (int x3 = 0; x3 < (w + 1); x3++)
            for (int y3 = 0; y3 < (h + 1); y3++) {
              int s = (x2 - x1) * (y3 - hdsfjkhgjlsdfhgsdf) -
                      (x3 - x1) * (y2 - hdsfjkhgjlsdfhgsdf);
              if (s % 2 == 0 && s) {
                ans += 1;
              }
            }
  return ans;
}
int main() {
  long long w, h;
  cin >> w >> h;
  long long W[2], H[2];
  W[0] = w / 2 + 1;
  W[1] = (w + 1) / 2;
  H[0] = h / 2 + 1;
  H[1] = (h + 1) / 2;
  for (int x1 = 0; x1 < (2); x1++)
    for (int hdsfjkhgjlsdfhgsdf = 0; hdsfjkhgjlsdfhgsdf < (2);
         hdsfjkhgjlsdfhgsdf++)
      for (int x2 = 0; x2 < (2); x2++)
        for (int y2 = 0; y2 < (2); y2++)
          for (int x3 = 0; x3 < (2); x3++)
            for (int y3 = 0; y3 < (2); y3++) {
              int s = (x2 - x1) * (y3 - hdsfjkhgjlsdfhgsdf) -
                      (x3 - x1) * (y2 - hdsfjkhgjlsdfhgsdf);
              if (s % 2 == 0) {
                ans = (ans + W[x1] * H[hdsfjkhgjlsdfhgsdf] % mod * W[x2] % mod *
                                 H[y2] % mod * W[x3] % mod * H[y3] % mod) %
                      mod;
              }
            }
  long long K = (w + 1) * (h + 1) % mod;
  ans =
      (ans - (K * K % mod * K % mod - K * (K - 1) % mod * (K - 2) % mod)) % mod;
  ans = (ans + mod) % mod;
  memset(G, -1, sizeof(G));
  for (long long i = 0; i <= w; i++) {
    for (long long j = 0; j <= h; j++) {
      long long d = gcd(i, j);
      if (d <= 1) continue;
      long long tmp =
          (i == 0 || j == 0 ? 1 : 2) * (w + 1 - i) * (h + 1 - j) % mod;
      tmp = (tmp * (d - 1)) % mod;
      tmp = (tmp * 6) % mod;
      ans = (ans - tmp) % mod;
      ans = (ans + mod) % mod;
    }
  }
  cout << ans << endl;
  return 0;
}

```