/****************************************************************************** * * * COMPUTE MOBIUS FUNCTION (Deleglise & Rivat's algorithm) * * 10/30/15 (dkc) * * * ******************************************************************************/ #include <math.h> int newmobl(unsigned long long a, unsigned long long b, long long *out, unsigned int *table, unsigned int tsize) { unsigned int i,p,count; unsigned long long beg,ps,rem; long long t; count=(unsigned int)(b-a); for (i=0; i<count; i++) out[i]=1; for (i=0; i<tsize; i++) { p=table[i]; ps=(unsigned long long)p*(unsigned long long)p; if (ps>b) goto askip; rem=a-(a/ps)*ps; if (rem!=0) beg=ps-rem; else beg=0; while (beg<count) { out[beg]=0; beg=beg+ps; } rem=a-(a/p)*p; if (rem!=0) beg=p-rem; else beg=0; while (beg<count) { out[beg]=-out[beg]*p; beg=beg+p; } } return(p); askip: for (i=0; i<count; i++) { if (out[i]==0) continue; t=out[i]; if (t<0) t=-t; if ((unsigned long long)t<(i+a)) out[i]=-out[i]; } for (i=0; i<count; i++) { if (out[i]==0) continue; if (out[i]>0) out[i]=1; else out[i]=-1; } return(1); }