前置知识：(1) 矩阵树定理。(2) 多项式四则运算（加、减、乘、求逆（牛顿迭代））。
网上的介绍很多。这里就不细讲了。
首先要把这个(gcd)给拆掉，才能用我们熟悉的“矩阵树定理”等一系列方法来解题。
众所周知，(sum_{d|x}varphi(d)=x)。本题中，我们就可以利用这个(varphi)来拆掉(gcd)。具体来说，如果我们用(T)来表示枚举的一棵树，(e_1,e_2,dots ,e_{n-1})来表示(T)包含的边的编号，那么首先，答案可以表示为：
我们把(gcd(w_{e_1},w_{e_2},dots,w_{e_{n-1}}))用上述的公式带入，得到：
然后就可以把这个(d)，放到前面来，得到：
其中，(W)表示最大的权值。后面括号里的部分，相当于【所有【(w_{e_i})是(d)的倍数的边】组成的子图】的所有生成树的边权和。
普通的矩阵树定理，可以用来求所有生成树的边权积之和。其中最为大家所熟知的应用，就是当所有边权都为(1)时，就相当于是求生成树的数量，也就是生成树计数问题。它的实现，就是求矩阵行列式，因为需要用到高斯消元，所以时间复杂度是(O(n^3))的。
但是本题中，要求的是边权和，而不是“边权积之和”。一种朴素的想法是考虑每一条边的贡献，也就是求【包含这条边的生成树数量】。它又等于【原图的生成树数量】减去【原图去掉这条边后的生成树数量】。对每条边都求一次【原图去掉这条边后的生成树数量】，时间复杂度(O(mn^3))。因为外层还要枚举(d)，所以总时间复杂度(O(Wcdot mn^3))，无法通过本题。
求所有生成树的边权和，其实有更好的方法。我们把每条边的边权，看做一个多项式：((1+w_{e_i}x))。其中，(x)不是任何具体的数，它只是一个符号，表示多项式的“一次项”。那么，一个生成树，它的边权之积的“一次项”前的系数，就是这颗生成树的边权和（这可以根据多项式乘法的定义来理解）。
如此以来，求所有生成树的边权和的问题，当我们把边权换成这样一个多项式后，就转化为了求所有生成树的边权积之和的问题。只不过，现在新的边权，不再是一个数，而是一个多项式。所以我们只需要定义出多项式的四则运算，就可以直接用矩阵树定理求解了。另外，我们只关心多项式的“一次项”系数，所以更高的项可以舍去。换句话说，所有的多项式运算，可以在(mod x^2)意义下进行。那么把多项式定义为一个( exttt{pair})就可以了。
多项式的加、减、乘运算都比较简单。不了解的可以看如下代码：
typedef pair<int,int> pii;
#define fi first
#define se second
//为了便于理解，这里就不取模了
pii operator+(const pii& a,const pii& b){
return mk(a.fi+b.fi,a.se+b.se);
}
pii operator-(const pii& a,const pii& b){
return mk(a.fi-b.fi,a.se-b.se);
}
pii operator*(const pii& a,const pii& b){
return mk(a.fi*b.fi,a.se*b.fi+b.se*a.fi);
}
多项式的除法，因为在(mod x^2)意义下进行，所以只需要求出多项式(mod x^2)意义下的逆即可。用经典的倍增法。设原多项式为(A(x))，要求的、它的逆为(F(x)=A^{-1}(x)pmod{x^n})。如果已经求出了(F_0(x)=A^{-1}(x)pmod{x^{lceilfrac{n}{2} ceil}})，那么：
对本题来说，我们就不需要递归了。先求出常数项，也就是求逆元。然后把常数项作为(F_0(x))带入上式即可。具体的代码：
pii inv(const pii& a){
int t=pow_mod(a.fi,MOD-2);
return mk(2*t,0)-mk(t,0)*mk(t,0)*a;
}
pii operator/(const pii& a,const pii& b){
return a*inv(b);
}
还有一个小细节是，高斯消元求行列式的时候，如果当前行的（要作为被除数）的这个多项式，常数项为(0)，上述的求逆的方法就不管用了，因为(0)没有逆元。这种情况下，我们首先考虑，在它下面，找一个常数项不为(0)的行，与它交换（根据行列式的性质，两行交换后行列式要变号，也就是(res:=-res)）。如果它下面所有行常数项都为(0)，那对这一列，我们就不需要管常数项了，直接拿一次项消就行（就把每一项都看成普通的数而不是多项式就行）。
至此，所有四则运算都可以(O(1))实现，所以做一次高斯消元，求出行列式的复杂度就是(O(n^3))的。外层还要枚举(d)，所以总时间复杂度(O(Wn^3))。这个复杂度仍然无法通过本题，我们还需要再对(W)这部分做一些优化。
发现对于一个(d)。如果“是它的倍数”的边，数量小于(n-1)条，那必不可能有任何生成树。特判这一情况后，我们惊喜地发现，复杂度降为：(O(frac{sumsigma_0(w_{e_i})}{n-1}cdot n^3)=O(n^2sumsigma_0(w_{e_i})))。其中(sigma_0(x))表示(x)的约数数量。这个复杂度很好理解，因为每个(d)，要作为约数出现(n-1)次才被计算，所以被计算到的(d)最多只有(frac{sumsigma_0(w_{e_i})}{n-1})个。再看这个复杂度，一个小于等于(W)的数，约数个数是(O(sqrt{W}))级别的，又因为边数是(O(n^2))级别的，所以总复杂度就是(O(n^4sqrt{W}))。事实上，这个(sqrt{W})只是一个理论上限，在本题的数据范围下，打表发现，(w)的约数个数最多为(144)。足以通过本题。
参考代码（在LOJ查看）：
//problem:LOJ3304
#include <bits/stdc++.h>
using namespace std;
#define pb push_back
#define mk make_pair
#define lob lower_bound
#define upb upper_bound
#define fi first
#define se second
#define SZ(x) ((int)(x).size())
typedef unsigned int uint;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int,int> pii;
const int MOD=998244353;
inline int mod1(int x){return x<MOD?x:x-MOD;}
inline int mod2(int x){return x<0?x+MOD:x;}
inline void add(int& x,int y){x=mod1(x+y);}
inline void sub(int& x,int y){x=mod2(x-y);}
inline int pow_mod(int x,int i){int y=1;while(i){if(i&1)y=(ll)y*x%MOD;x=(ll)x*x%MOD;i>>=1;}return y;}
pii operator+(const pii& lhs,const pii& rhs){
return mk(mod1(lhs.fi+rhs.fi),mod1(lhs.se+rhs.se));
}
pii operator-(const pii& lhs,const pii& rhs){
return mk(mod2(lhs.fi-rhs.fi),mod2(lhs.se-rhs.se));
}
pii operator*(const pii& lhs,const pii& rhs){
return mk((ll)lhs.fi*rhs.fi%MOD,((ll)lhs.se*rhs.fi+(ll)rhs.se*lhs.fi)%MOD);
}
pii inv(const pii& a){
//assert(a.fi!=0);
int t=pow_mod(a.fi,MOD-2);
return mk(2LL*t%MOD,0)-mk(t,0)*mk(t,0)*a;//牛顿迭代
}
pii operator/(const pii& lhs,const pii& rhs){
return lhs*inv(rhs);
}
pii operator-(const pii& rhs){
return mk(mod2(-rhs.fi),mod2(-rhs.se));
}//负号
pii& operator+=(pii& lhs,const pii& rhs){
lhs=lhs+rhs;
return lhs;
}
pii& operator-=(pii& lhs,const pii& rhs){
lhs=lhs-rhs;
return lhs;
}
pii& operator*=(pii& lhs,const pii& rhs){
lhs=lhs*rhs;
return lhs;
}
pii& operator/=(pii& lhs,const pii& rhs){
lhs=lhs/rhs;
return lhs;
}
const int MAXW=152501;
int p[MAXW+5],cnt,phi[MAXW+5];
bool v[MAXW+5];
void sieve(int lim){
phi[1]=1;
for(int i=2;i<=lim;++i){
if(!v[i]){
phi[i]=i-1;
p[++cnt]=i;
}
for(int j=1;j<=cnt && p[j]*i<=lim;++j){
v[i*p[j]]=1;
if(i%p[j]==0){
phi[i*p[j]]=phi[i]*p[j];
break;
}
phi[i*p[j]]=phi[i]*phi[p[j]];
}
}
}//线性筛,求phi
const int MAXN=30;
const int MAXM=MAXN*(MAXN-1)/2;
int n,m;
struct Edge{
int u,v,w;
}e[MAXM+5];
struct Matrix{
int n;
pii a[MAXN+5][MAXN+5];
Matrix(){
n=0;
memset(a,0,sizeof(a));
}
};
pii get_det2(Matrix A,int i){
//常数项全部为0
pii res=mk(1,0);
int line=0;
for(int j=i;j<=A.n;++j){
if(A.a[j][i].se!=0){
line=j;
break;
}
}
if(!line)return mk(0,0);
if(line!=i){
for(int j=1;j<=A.n;++j)swap(A.a[i][j],A.a[line][j]);
res=-res;
}
int _inv=pow_mod(A.a[i][i].se,MOD-2);
for(int j=i+1;j<=A.n;++j){
int t=(ll)A.a[j][i].se*_inv%MOD;
for(int k=1;k<=A.n;++k)
sub(A.a[j][k].se,(ll)A.a[i][k].se*t%MOD);
}
res*=A.a[i][i];
return res;
}
pii get_det(Matrix A){
pii res=mk(1,0);
for(int i=1;i<=A.n;++i){
int line=0;
for(int j=i;j<=A.n;++j){
if(A.a[j][i].fi!=0){
line=j;
break;
}
}
if(line==0){
res*=get_det2(A,i);
continue;
}
if(line!=i){
for(int j=1;j<=A.n;++j)swap(A.a[i][j],A.a[line][j]);
res=-res;
}
pii _inv=inv(A.a[i][i]);
for(int j=i+1;j<=A.n;++j){
pii t=A.a[j][i]*_inv;
for(int k=1;k<=A.n;++k)
A.a[j][k]-=A.a[i][k]*t;
}
res*=A.a[i][i];
}
return res;
}
int main() {
cin>>n>>m;
int max_w=0;
for(int i=1;i<=m;++i){
cin>>e[i].u>>e[i].v>>e[i].w;
max_w=max(max_w,e[i].w);
}
sieve(max_w);
int ans=0;
for(int i=1;i<=max_w;++i){
int cnt_e=0;
for(int j=1;j<=m;++j)if(e[j].w%i==0)cnt_e++;
if(cnt_e<n-1)continue;
Matrix mat;
for(int j=1;j<=m;++j)if(e[j].w%i==0){
mat.a[e[j].u][e[j].u]+=mk(1,e[j].w);
mat.a[e[j].v][e[j].v]+=mk(1,e[j].w);
mat.a[e[j].u][e[j].v]-=mk(1,e[j].w);
mat.a[e[j].v][e[j].u]-=mk(1,e[j].w);
}
mat.n=n-1;
pii det=get_det(mat);
ans=((ll)ans+(ll)phi[i]*det.se)%MOD;
}
cout<<ans<<endl;
return 0;
}