O(1)GCD

#include <cstdio>
#include <cstring>
#define SIZE 1000
#define U unsigned int
using namespace std;
  
  int ss[1000001],cnt,p[1000001];
  short g[1001][1001],b[1000001];
  int f[1000001][3];
  int a[2001],bb[2001];

  int pregcd(int x,int y){
      if (y==0) {g[x][y]=x;return(x);}
      if (g[x][y]!=-1) return(g[x][y]);
      g[x][y]=pregcd(y,x%y);
      return(g[x][y]);
  }
  
  void euler(){
      p[1]=1;
      for (int i=2;i<=SIZE*SIZE;i++){
        if (!b[i]){
          p[i]=i;ss[++cnt]=i;
      }    
      for (int j=1;ss[j]*i<=SIZE*SIZE&&j<=cnt;j++){
          b[ss[j]*i]=1;
          p[ss[j]*i]=ss[j];
          if (i%ss[j]==0) break;
      }
    }
  }
  
  int gcd(int x,int y){
      if (x<=SIZE&&y<=SIZE) return(g[x][y]);
      int d=1;
      for (int i=0;i<=2;i++)
        if (p[f[x][i]]==f[x][i]){
          if (y%f[x][i]==0){
            d*=f[x][i];y/=f[x][i];    
        }
      }else{
          int t=g[f[x][i]][y%f[x][i]];
          d*=t;y/=t;
      }
    return(d);  
  }
  
  int main(){
      freopen("a.in","r",stdin);
      
      for (int i=0;i<=SIZE;i++)
        for (int j=0;j<=SIZE;j++)
          g[i][j]=-1;
    for (int i=0;i<=SIZE;i++)
      for (int j=0;j<=SIZE;j++)
        pregcd(i,j);//预处理SIZE以内gcd
        
    euler();
    
    f[1][0]=f[1][1]=f[1][2]=1;
    for (int i=2;i<=SIZE*SIZE;i++){
      memcpy(f[i],f[i/p[i]],sizeof(f[i]));    
      if (f[i][0]<=SIZE/p[i]) f[i][0]*=p[i];else
      if (f[i][1]<=SIZE/p[i]) f[i][1]*=p[i];else
      f[i][2]*=p[i];
    }//预处理所有数的分解
    
    int T,n,m;
    scanf("%d",&T);
    while (T--){
      scanf("%d%d",&n,&m);
      U ans=0;
      for (int i=0;i<n;i++) scanf("%d",&a[i]);
      for (int i=0;i<m;i++) scanf("%d",&bb[i]);
      for (int i=0;i<n;i++)
        for (int j=0;j<m;j++)
          ans+=gcd(a[i],bb[j])^i^j;
      printf("%u
",ans);
    }
  }

BZOJ4454