#include #include #include #include #include #include #include #include #include static const int MOD = 1000000007; using ll = long long; using u32 = unsigned; using u64 = unsigned long long; using namespace std; template constexpr T INF = ::numeric_limits::max()/32*15+208; template struct modint { u32 val; public: static modint raw(int v) { modint x; x.val = v; return x; } modint() : val(0) {} template modint(T v) { ll x = (ll)(v%(ll)(M)); if (x < 0) x += M; val = u32(x); } modint(bool v) { val = ((unsigned int)(v) % M); } modint& operator++() { val++; if (val == M) val = 0; return *this; } modint& operator--() { if (val == 0) val = M; val--; return *this; } modint operator++(int) { modint result = *this; ++*this; return result; } modint operator--(int) { modint result = *this; --*this; return result; } modint& operator+=(const modint& b) { val += b.val; if (val >= M) val -= M; return *this; } modint& operator-=(const modint& b) { val -= b.val; if (val >= M) val += M; return *this; } modint& operator*=(const modint& b) { u64 z = val; z *= b.val; val = (u32)(z % M); return *this; } modint& operator/=(const modint& b) { return *this = *this * b.inv(); } modint operator+() const { return *this; } modint operator-() const { return modint() - *this; } modint pow(long long n) const { modint x = *this, r = 1; while (n) { if (n & 1) r *= x; x *= x; n >>= 1; } return r; } modint inv() const { return pow(M-2); } friend modint operator+(const modint& a, const modint& b) { return modint(a) += b; } friend modint operator-(const modint& a, const modint& b) { return modint(a) -= b; } friend modint operator*(const modint& a, const modint& b) { return modint(a) *= b; } friend modint operator/(const modint& a, const modint& b) { return modint(a) /= b; } friend bool operator==(const modint& a, const modint& b) { return a.val == b.val; } friend bool operator!=(const modint& a, const modint& b) { return a.val != b.val; } }; using mint = modint; class Factorial { vector facts, factinv; public: explicit Factorial(int n) : facts(n+1), factinv(n+1) { facts[0] = 1; for (int i = 1; i < n+1; ++i) facts[i] = facts[i-1] * mint(i); factinv[n] = facts[n].inv(); for (int i = n-1; i >= 0; --i) factinv[i] = factinv[i+1] * mint(i+1); } mint fact(int k) const { if(k >= 0) return facts[k]; else return factinv[-k]; } mint operator[](const int &k) const { if(k >= 0) return facts[k]; else return factinv[-k]; } mint C(int p, int q) const { if(q < 0 || p < q) return 0; return facts[p] * factinv[q] * factinv[p-q]; } mint P(int p, int q) const { if(q < 0 || p < q) return 0; return facts[p] * factinv[p-q]; } mint H(int p, int q) const { if(p < 0 || q < 0) return 0; return q == 0 ? 1 : C(p+q-1, q); } }; struct Prime { // Wheel factorization static constexpr int wheel[] = {4, 2, 4, 2, 4, 6, 2, 6}, wheel2[] = {7, 11, 13, 17, 19, 23, 29, 31}, wheel_sum[] = {0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 3, 3, 3, 3, 4, 4, 5, 5, 5, 5, 6, 6, 6, 6, 6, 6, 7, 7}; static inline int f(int n){ return (n-1)/30*8 + wheel_sum[(n-1)%30]; } static inline int g(int n){ return ((n-1) >> 3)*30 + wheel2[(n-1)&7]; } vector primes; Prime(int M) { primes = {2, 3, 5}; if(M < 7){ while(!primes.empty() && M < primes.back()) primes.pop_back(); return; } int n = f(M), m = g(n), k = f(int(floor(sqrt(M)))); primes.reserve(3+max(0, (int)(n/(log(n)-1.12)))); vector sieve(n+1, true); for (int i = 1; i <= k; ++i) { if(sieve[i]){ ll p = g(i), q = p*p; int j = (i-1)&7; while(q <= m){ sieve[f(q)] = false; q += wheel[j] * p; j = (j+1)&7; } } } for (int i = 1; i <= n; ++i) { if(sieve[i]) primes.emplace_back(g(i)); } } }; constexpr int Prime::wheel[], Prime::wheel2[], Prime::wheel_sum[]; Prime p(200000); template void div_transform(vector &a){ int n = a.size(); for (auto &&i : p.primes) { if(i >= n) break; for (int k = (n-1)/i; k > 0; --k) { a[k] += a[k*i]; } } for (int i = 1; i < n; ++i) a[i] += a[0]; } template void div_itransform(vector &a){ int n = a.size(); for (int i = 1; i < n; ++i) a[i] -= a[0]; for (auto &&i : p.primes) { if(i >= n) break; for (int k = 1; k*i < n; ++k) { a[k] -= a[k*i]; } } } int main() { int h, w; cin >> h >> w; auto f = [&](mint x){ return x*(x-1)*(x-2)/6; }; mint ans = f(mint(h)*w) - mint(h)*f(w) - mint(w)*f(h); int sz = max(h, w); vector a(sz), b(sz); for (int i = 1; i < h; ++i) a[i] = h-i; for (int i = 1; i < w; ++i) b[i] = w-i; div_transform(a); div_transform(b); for (int i = 0; i < sz; ++i) a[i] *= b[i]; div_itransform(a); for (int i = 1; i < sz; ++i) ans -= a[i]*(i-1)*2; cout << ans.val << "\n"; return 0; }