Sieving for Bsmooth values

General Usage Performance
10

My refactor of this would be

function bitlength(x::BigInt)::Int16
    return length(digits(x;base=2))-1
end

const ONE = BigInt(1)

function isSmooth2(num::BigInt, rang::Int, small_primes, logs)
	rn = BigInt[]
	bump = 2^(rang-1)
	stump = 2^rang+ONE
	sta = num - bump - ONE
	cabin = zeros(Int16,stump)
	#sieve using factorbase ignoring the first 100 
	#because those are expensive to mark
	for i in 101:6600
		stu=small_primes[i] - (sta%small_primes[i])
		cabin[stu:small_primes[i]:stump] .+= logs[i]
	end
	#sieve using the extended factorbase
	#for i in np
	#	stu=np[i]-(sta%np[i])
	#	cabin[stu:np[i]:stump] .+=logs2[i]
	#end
	#find bsmooth numbers
	lmax = bitlength(num)-16
	for i in eachindex(cabin)
		if cabin[i]>lmax
			push!(rn,i+sta)
		end
	end
	#expand search for Bsmooth values
	cs=0
	while isempty(rn)
		if cs<0
			cs-=1
		else
			cs+=1
		end
		sta+=stump*cs
		cabin = zeros(Int16,stump)
		for i in 101:6600
			stu=small_primes[i]-(sta%small_primes[i])
			cabin[stu:small_primes[i]:stump] .+= logs[i]
		end
		#for i in np
		#	stu=np[i]-(sta%np[i])
		#	cabin[stu:np[i]:stump] .+=logs2[i]
		#end
		for i in eachindex(cabin)
			if cabin[i]>lmax
				push!(rn,i+sta)
			end
		end
		#now search in the the block of numbers 
		#in the opposite direction
		cs*=-1
	end
	return rn
end

function make_logs(small_primes)
    logs = zeros(Int16,length(small_primes))
    guf = 0
    sqrtof2 = sqrt(2)
    #calculate log base 2 this time rounded....
    for i in 1:6600
        if small_primes[i] > sqrtof2
            guf += 1
            sqrtof2 *= 2.0
        end
        logs[i] = guf
    end
    logs
end

small_primes = primes(66103)
logs = make_logs(small_primes)

for i in 1:10
    x::BigInt= rand(BigInt(2)<<90:BigInt(2)<<100)
    println(x)
    display(isSmooth2(x, 20, small_primes, logs))
end

which takes 0.4 seconds per issmooth call. I’ve also removed TimerOutputs. since IMO it’s usually better to profile code (e.g. using @profview from ProfileView) rather than inserting points to time manually.

4 Likes