#include <bits/stdc++.h>
#pragma comment(linker, "/STACK:16777216")
using namespace std;
const double pi = 3.14159265358979323846264338327950288419716939937511;
const double eps = 1e-10;
char ch_ch_ch[1 << 20];
inline string gs() {
  scanf("%s", ch_ch_ch);
  return string(ch_ch_ch);
}
inline string gl() {
  gets(ch_ch_ch);
  return string(ch_ch_ch);
}
inline int gi() {
  int x;
  scanf("%d", &x);
  return x;
}
int gcdp[4001][4001];
int gcd(int a, int b) {
  if (gcdp[a][b]) return gcdp[a][b];
  return gcdp[a][b] = (b ? gcd(b, a % b) : a);
}
int w, h;
const int mod = 1000000007;
int precalc[2][2][2][2];
void solution() {
  w = gi();
  h = gi();
  for (int x1 = 0; x1 < (2); ++x1)
    for (int x2 = 0; x2 < (2); ++x2)
      for (int y1 = 0; y1 < (2); ++y1)
        for (int y2 = 0; y2 < (2); ++y2) {
          int cx1 = w / 2 + ((!x1) ? 1 : (w & 1));
          int cx2 = w / 2 + ((!x2) ? 1 : (w & 1));
          int cy1 = h / 2 + ((!y1) ? 1 : (h & 1));
          int cy2 = h / 2 + ((!y2) ? 1 : (h & 1));
          precalc[x1][y1][x2][y2] = (cx1 * 1ll * cy1 * cx2 * cy2) % mod;
        }
  int res = 0;
  for (int x = 0; x < (w + 1); ++x)
    for (int y = 0; y < (h + 1); ++y) {
      int xp = x & 1, yp = y & 1;
      for (int x1 = 0; x1 < (2); ++x1)
        for (int x2 = 0; x2 < (2); ++x2)
          for (int y1 = 0; y1 < (2); ++y1)
            for (int y2 = 0; y2 < (2); ++y2)
              if (!(abs((x1 - xp) * (y2 - yp) + (x2 - xp) * (y1 - yp)) & 1)) {
                res = (res + precalc[x1][y1][x2][y2]);
                if (res >= mod) res -= mod;
              }
    }
  int deg = 0;
  for (int dx = 0; dx < (w + 1); ++dx)
    for (int dy = 0; dy < (h + 1); ++dy) {
      int cnt = (w - dx + 1) * (h - dy + 1);
      int pointsCnt = (gcd(dx, dy) + 1);
      int variants = 0;
      if (pointsCnt == 1)
        variants = 1;
      else
        variants = 6 * (pointsCnt - 1);
      if (dx != 0 && dy != 0) cnt *= 2;
      deg = (deg + cnt * 1ll * variants) % mod;
    }
  cout << (res + mod - deg) % mod << endl;
}
int main(int argc, char** argv) {
  solution();
  return 0;
}

### Prompt

Generate a cpp solution 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>
#pragma comment(linker, "/STACK:16777216")
using namespace std;
const double pi = 3.14159265358979323846264338327950288419716939937511;
const double eps = 1e-10;
char ch_ch_ch[1 << 20];
inline string gs() {
  scanf("%s", ch_ch_ch);
  return string(ch_ch_ch);
}
inline string gl() {
  gets(ch_ch_ch);
  return string(ch_ch_ch);
}
inline int gi() {
  int x;
  scanf("%d", &x);
  return x;
}
int gcdp[4001][4001];
int gcd(int a, int b) {
  if (gcdp[a][b]) return gcdp[a][b];
  return gcdp[a][b] = (b ? gcd(b, a % b) : a);
}
int w, h;
const int mod = 1000000007;
int precalc[2][2][2][2];
void solution() {
  w = gi();
  h = gi();
  for (int x1 = 0; x1 < (2); ++x1)
    for (int x2 = 0; x2 < (2); ++x2)
      for (int y1 = 0; y1 < (2); ++y1)
        for (int y2 = 0; y2 < (2); ++y2) {
          int cx1 = w / 2 + ((!x1) ? 1 : (w & 1));
          int cx2 = w / 2 + ((!x2) ? 1 : (w & 1));
          int cy1 = h / 2 + ((!y1) ? 1 : (h & 1));
          int cy2 = h / 2 + ((!y2) ? 1 : (h & 1));
          precalc[x1][y1][x2][y2] = (cx1 * 1ll * cy1 * cx2 * cy2) % mod;
        }
  int res = 0;
  for (int x = 0; x < (w + 1); ++x)
    for (int y = 0; y < (h + 1); ++y) {
      int xp = x & 1, yp = y & 1;
      for (int x1 = 0; x1 < (2); ++x1)
        for (int x2 = 0; x2 < (2); ++x2)
          for (int y1 = 0; y1 < (2); ++y1)
            for (int y2 = 0; y2 < (2); ++y2)
              if (!(abs((x1 - xp) * (y2 - yp) + (x2 - xp) * (y1 - yp)) & 1)) {
                res = (res + precalc[x1][y1][x2][y2]);
                if (res >= mod) res -= mod;
              }
    }
  int deg = 0;
  for (int dx = 0; dx < (w + 1); ++dx)
    for (int dy = 0; dy < (h + 1); ++dy) {
      int cnt = (w - dx + 1) * (h - dy + 1);
      int pointsCnt = (gcd(dx, dy) + 1);
      int variants = 0;
      if (pointsCnt == 1)
        variants = 1;
      else
        variants = 6 * (pointsCnt - 1);
      if (dx != 0 && dy != 0) cnt *= 2;
      deg = (deg + cnt * 1ll * variants) % mod;
    }
  cout << (res + mod - deg) % mod << endl;
}
int main(int argc, char** argv) {
  solution();
  return 0;
}

```