\\ ---- PARI/GP program for A352739 ----- \\ Peter Munn, April 2022 \\ ----------Abbreviations used---------- \\ AFG: see comment at point of use \\ No = Number \\ PnD = Primitive nonDeficient (see A006039 for definition) \\ RFactor = Replacement Factor \\ Seq = Sequence A352739uptovalue(UpperLimit) = { GetPrime = 1; GetExponent = 2; \\ Fixed values, giving names to the factor list columns. my (SeqIndex); forfactored (FactoredCandidate = 1, UpperLimit, my (FactorList, PrimeDivisor, Prime_sExponent, ComplementaryDivisor, maybeStronglyPnD, AFG, CriticalReplacementPrime, ReplacementPrime, RFactorNo, Candidate_sAbundancy, ComparisonAbundancy); Candidate_sAbundancy = sigma(FactoredCandidate)/FactoredCandidate[1]; if (Candidate_sAbundancy >= 2, maybeStronglyPnD = 1; \\ 1 = True; the candidate is still under consideration FactorList = FactoredCandidate[2]; forstep (FactorNo=omega(FactoredCandidate), 1, -1, \\ Run through all the factor list PrimeDivisor = FactorList[FactorNo,GetPrime]; ComplementaryDivisor = FactoredCandidate[1]/PrimeDivisor; if (sigma(ComplementaryDivisor) >= 2*ComplementaryDivisor, \\ This proper divisor is nondeficient, so the candidate is not primitive maybeStronglyPnD = 0; break , \\ otherwise: apply the replacement prime part of the stronger test Prime_sExponent = FactorList[FactorNo,GetExponent]; AFG = PrimeDivisor * (PrimeDivisor^Prime_sExponent - 1) / (PrimeDivisor - 1); \\ AFG stands for AbundancyFactorGenerator. PrimeDivisor contributes (AFG+1)/AFG to the \\ abundancy of the candidate; the rest comes from the abundancy of ComplementaryDivisor. CriticalReplacementPrime = nextprime(AFG+1); while (FactoredCandidate[1] % CriticalReplacementPrime == 0, CriticalReplacementPrime = nextprime(CriticalReplacementPrime+1) ); \\ We need only test using the CriticalReplacementPrime and the candidate's prime divisors. ReplacementPrime = CriticalReplacementPrime; RFactorNo = omega(FactoredCandidate)+1; while (1, ComparisonAbundancy = sigma(ComplementaryDivisor*ReplacementPrime) / (ComplementaryDivisor*ReplacementPrime); if (Candidate_sAbundancy > ComparisonAbundancy && ComparisonAbundancy >= 2, maybeStronglyPnD=0; break(2) \\ candidate fails the replacement prime test ); if (RFactorNo == 1, break) ; ReplacementPrime = FactorList[RFactorNo--,GetPrime]; ) \\ end of while, running through replacement primes ) \\ endif ); \\ endfor, running through the candidate's factor list (of primes to be replaced) \\ if the candidate remains under consideration, it has passed all the tests if (maybeStronglyPnD, print(SeqIndex++, " ", FactoredCandidate[1])) ) ) \\ end of forfactored, running through the candidates }