#include <bits/stdc++.h>
using namespace std;
const long double pi = 3.14159265359;
template <typename T>
T abs(T x) {
  return x > 0 ? x : -x;
}
template <typename T>
T sqr(T x) {
  return x * x;
}
template <typename T>
void chmin(T &x, T y) {
  x = min(x, y);
}
template <typename T>
void chmax(T &x, T y) {
  x = max(x, y);
}
const long long mod = 1000000007;
long long add(long long x, long long y) { return (x + y) % mod; }
long long sub(long long x, long long y) { return (x - y + mod) % mod; }
long long mul(long long x, long long y) { return (x * y) % mod; }
const int maxn = 4005;
int g[maxn][maxn];
int gcd(int a, int b) { return a ? gcd(b % a, a) : b; }
long long brute(int n, int m) {
  int res = 0;
  for (int i = 0; i <= n; i++) {
    for (int j = 0; j <= m; j++) {
      for (int i2 = 0; i2 <= n; i2++) {
        for (int j2 = 0; j2 <= m; j2++) {
          for (int i3 = 0; i3 <= n; i3++) {
            for (int j3 = 0; j3 <= m; j3++) {
              long long s = (i2 - i) * (j3 - j) - (j2 - j) * (i3 - i);
              if (s && s % 2 == 0) {
                res++;
              }
            }
          }
        }
      }
    }
  }
  return res;
}
int main() {
  srand(time(NULL));
  int n, m;
  cin >> n >> m;
  long long ans = 0;
  for (int i = 0; i < maxn; i++) {
    for (int j = 0; j < maxn; j++) {
      g[i][j] = 1;
    }
  }
  for (int i = 1; i < maxn; i++) {
    for (int j = 1; j < maxn; j++) {
      if (g[i][j] == 1) {
        for (int k = 2; i * k < maxn && j * k < maxn; k++) {
          g[i * k][j * k] = k;
        }
      }
    }
  }
  for (int i = 0; i <= n; i++) {
    int a = 0 - i;
    int b = n - i;
    long long u[2];
    if (a & 1) {
      u[1] = 1 + (b - a) / 2;
      u[0] = b - a + 1 - u[1];
    } else {
      u[0] = 1 + (b - a) / 2;
      u[1] = b - a + 1 - u[0];
    }
    for (int j = 0; j <= m; j++) {
      int c = 0 - j;
      int d = m - j;
      long long v[2];
      if (c & 1) {
        v[1] = 1 + (d - c) / 2;
        v[0] = d - c + 1 - v[1];
      } else {
        v[0] = 1 + (d - c) / 2;
        v[1] = d - c + 1 - v[0];
      }
      long long s = (u[0] + u[1]) * (v[0] + v[1]) % mod;
      long long x = s - u[1] * v[1] + mod;
      x %= mod;
      ans = (ans + x * x) % mod;
      x = u[1] * v[1];
      x %= mod;
      ans = (ans + x * x) % mod;
    }
  }
  for (int i = 1; i <= n; i++) {
    for (int j = 1; j <= m; j++) {
      long long x = g[i][j];
      long long d = (mod - (x - 1)) * 12LL;
      d *= (n + 1 - i) * (m + 1 - j);
      d %= mod;
      ans = (ans + d) % mod;
    }
  }
  ans -= (n + 1) * 1LL * n % mod * (n - 1) % mod * (m + 1) % mod;
  ans -= (m + 1) * 1LL * m % mod * (m - 1) % mod * (n + 1) % mod;
  ans %= mod;
  ans += mod;
  ans %= mod;
  long long x = (n + 1) * 1LL * (m + 1) % mod;
  ans -= x * (x + mod - 1) % mod * 3LL % mod;
  ans %= mod;
  ans += mod;
  ans %= mod;
  ans -= x;
  ans %= mod;
  ans += mod;
  ans %= mod;
  cout << ans << endl;
  return 0;
}

### Prompt

Please provide a CPP coded solution to the problem described below:

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;
const long double pi = 3.14159265359;
template <typename T>
T abs(T x) {
  return x > 0 ? x : -x;
}
template <typename T>
T sqr(T x) {
  return x * x;
}
template <typename T>
void chmin(T &x, T y) {
  x = min(x, y);
}
template <typename T>
void chmax(T &x, T y) {
  x = max(x, y);
}
const long long mod = 1000000007;
long long add(long long x, long long y) { return (x + y) % mod; }
long long sub(long long x, long long y) { return (x - y + mod) % mod; }
long long mul(long long x, long long y) { return (x * y) % mod; }
const int maxn = 4005;
int g[maxn][maxn];
int gcd(int a, int b) { return a ? gcd(b % a, a) : b; }
long long brute(int n, int m) {
  int res = 0;
  for (int i = 0; i <= n; i++) {
    for (int j = 0; j <= m; j++) {
      for (int i2 = 0; i2 <= n; i2++) {
        for (int j2 = 0; j2 <= m; j2++) {
          for (int i3 = 0; i3 <= n; i3++) {
            for (int j3 = 0; j3 <= m; j3++) {
              long long s = (i2 - i) * (j3 - j) - (j2 - j) * (i3 - i);
              if (s && s % 2 == 0) {
                res++;
              }
            }
          }
        }
      }
    }
  }
  return res;
}
int main() {
  srand(time(NULL));
  int n, m;
  cin >> n >> m;
  long long ans = 0;
  for (int i = 0; i < maxn; i++) {
    for (int j = 0; j < maxn; j++) {
      g[i][j] = 1;
    }
  }
  for (int i = 1; i < maxn; i++) {
    for (int j = 1; j < maxn; j++) {
      if (g[i][j] == 1) {
        for (int k = 2; i * k < maxn && j * k < maxn; k++) {
          g[i * k][j * k] = k;
        }
      }
    }
  }
  for (int i = 0; i <= n; i++) {
    int a = 0 - i;
    int b = n - i;
    long long u[2];
    if (a & 1) {
      u[1] = 1 + (b - a) / 2;
      u[0] = b - a + 1 - u[1];
    } else {
      u[0] = 1 + (b - a) / 2;
      u[1] = b - a + 1 - u[0];
    }
    for (int j = 0; j <= m; j++) {
      int c = 0 - j;
      int d = m - j;
      long long v[2];
      if (c & 1) {
        v[1] = 1 + (d - c) / 2;
        v[0] = d - c + 1 - v[1];
      } else {
        v[0] = 1 + (d - c) / 2;
        v[1] = d - c + 1 - v[0];
      }
      long long s = (u[0] + u[1]) * (v[0] + v[1]) % mod;
      long long x = s - u[1] * v[1] + mod;
      x %= mod;
      ans = (ans + x * x) % mod;
      x = u[1] * v[1];
      x %= mod;
      ans = (ans + x * x) % mod;
    }
  }
  for (int i = 1; i <= n; i++) {
    for (int j = 1; j <= m; j++) {
      long long x = g[i][j];
      long long d = (mod - (x - 1)) * 12LL;
      d *= (n + 1 - i) * (m + 1 - j);
      d %= mod;
      ans = (ans + d) % mod;
    }
  }
  ans -= (n + 1) * 1LL * n % mod * (n - 1) % mod * (m + 1) % mod;
  ans -= (m + 1) * 1LL * m % mod * (m - 1) % mod * (n + 1) % mod;
  ans %= mod;
  ans += mod;
  ans %= mod;
  long long x = (n + 1) * 1LL * (m + 1) % mod;
  ans -= x * (x + mod - 1) % mod * 3LL % mod;
  ans %= mod;
  ans += mod;
  ans %= mod;
  ans -= x;
  ans %= mod;
  ans += mod;
  ans %= mod;
  cout << ans << endl;
  return 0;
}

```