#include <bits/stdc++.h>
using namespace std;
const int INF = (1 << 30) - 1;
const long double EPS = 1e-9;
const long double PI = fabs(atan2(0.0, -1.0));
const int MOD = 1000000007;
inline int get(int n, int e) { return (n - e) / 2 + 1; }
int memo[4005][4005];
inline int gcd(int a, int b) {
  int &res = memo[a][b];
  if (res != -1) return res;
  return res = (b ? gcd(b, a % b) : a);
}
int cnt[2][2];
void solve() {
  int w, h;
  cin >> w >> h;
  for (int ex = 0; ex < 2; ex++)
    for (int ey = 0; ey < 2; ey++)
      cnt[ex][ey] = (get(w, ex) * 1LL * get(h, ey)) % MOD;
  int res = 0;
  for (int mask = 0; mask < (1 << 6); mask++) {
    int x1 = (((mask) & (1 << (0))) == 0 ? false : true);
    int y1 = (((mask) & (1 << (1))) == 0 ? false : true);
    int x2 = (((mask) & (1 << (2))) == 0 ? false : true);
    int y2 = (((mask) & (1 << (3))) == 0 ? false : true);
    int x3 = (((mask) & (1 << (4))) == 0 ? false : true);
    int y3 = (((mask) & (1 << (5))) == 0 ? false : true);
    int vx1 = x2 - x1, vy1 = y2 - y1;
    int vx2 = x3 - x1, vy2 = y3 - y1;
    int s = abs(vx1 * vy2 - vy1 * vx2) & 1;
    if (s == 0) {
      res +=
          (((cnt[x1][y1] * 1LL * cnt[x2][y2]) % MOD) * 1LL * cnt[x3][y3]) % MOD;
      if (res >= MOD) res -= MOD;
    }
  }
  res -= ((w + 1) * (h + 1)) % MOD;
  if (res < 0) res += MOD;
  res -= ((w + 1) * (h + 1)) * 1LL * ((w + 1) * (h + 1) - 1) * 3 % MOD;
  if (res < 0) res += MOD;
  res -= ((w + 1) * 1LL * (h + 1) * 1LL * h * 1LL * (h - 1)) % MOD;
  if (res < 0) res += MOD;
  res -= ((h + 1) * 1LL * (w + 1) * 1LL * w * 1LL * (w - 1)) % MOD;
  if (res < 0) res += MOD;
  memset(memo, -1, sizeof(memo));
  for (int dx = 1; dx <= w; dx++) {
    for (int dy = 1; dy <= h; dy++) {
      int d = gcd(dx, dy);
      if (d == 0) continue;
      int cnt = ((d - 1) * 12LL * (w - dx + 1) * 1LL * (h - dy + 1)) % MOD;
      res -= cnt;
      if (res < 0) res += MOD;
    }
  }
  cout << res;
}
int main() {
  solve();
  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;
const int INF = (1 << 30) - 1;
const long double EPS = 1e-9;
const long double PI = fabs(atan2(0.0, -1.0));
const int MOD = 1000000007;
inline int get(int n, int e) { return (n - e) / 2 + 1; }
int memo[4005][4005];
inline int gcd(int a, int b) {
  int &res = memo[a][b];
  if (res != -1) return res;
  return res = (b ? gcd(b, a % b) : a);
}
int cnt[2][2];
void solve() {
  int w, h;
  cin >> w >> h;
  for (int ex = 0; ex < 2; ex++)
    for (int ey = 0; ey < 2; ey++)
      cnt[ex][ey] = (get(w, ex) * 1LL * get(h, ey)) % MOD;
  int res = 0;
  for (int mask = 0; mask < (1 << 6); mask++) {
    int x1 = (((mask) & (1 << (0))) == 0 ? false : true);
    int y1 = (((mask) & (1 << (1))) == 0 ? false : true);
    int x2 = (((mask) & (1 << (2))) == 0 ? false : true);
    int y2 = (((mask) & (1 << (3))) == 0 ? false : true);
    int x3 = (((mask) & (1 << (4))) == 0 ? false : true);
    int y3 = (((mask) & (1 << (5))) == 0 ? false : true);
    int vx1 = x2 - x1, vy1 = y2 - y1;
    int vx2 = x3 - x1, vy2 = y3 - y1;
    int s = abs(vx1 * vy2 - vy1 * vx2) & 1;
    if (s == 0) {
      res +=
          (((cnt[x1][y1] * 1LL * cnt[x2][y2]) % MOD) * 1LL * cnt[x3][y3]) % MOD;
      if (res >= MOD) res -= MOD;
    }
  }
  res -= ((w + 1) * (h + 1)) % MOD;
  if (res < 0) res += MOD;
  res -= ((w + 1) * (h + 1)) * 1LL * ((w + 1) * (h + 1) - 1) * 3 % MOD;
  if (res < 0) res += MOD;
  res -= ((w + 1) * 1LL * (h + 1) * 1LL * h * 1LL * (h - 1)) % MOD;
  if (res < 0) res += MOD;
  res -= ((h + 1) * 1LL * (w + 1) * 1LL * w * 1LL * (w - 1)) % MOD;
  if (res < 0) res += MOD;
  memset(memo, -1, sizeof(memo));
  for (int dx = 1; dx <= w; dx++) {
    for (int dy = 1; dy <= h; dy++) {
      int d = gcd(dx, dy);
      if (d == 0) continue;
      int cnt = ((d - 1) * 12LL * (w - dx + 1) * 1LL * (h - dy + 1)) % MOD;
      res -= cnt;
      if (res < 0) res += MOD;
    }
  }
  cout << res;
}
int main() {
  solve();
  return 0;
}

```