Herhangi Matematiksel Şeyler: Euler Projesi 269. Soru

En az bir tam sayı köke sahip polinomlar

P(x)=0 denkleminin bir çözümüne P(x) polinomunun bir kökü ya da sıfırı denir.
Bir n sayısının basamaklarını katsayı kabul eden polinomu P_n ile gösterelim.
Örneğin P₅₇₀₃(x)=5x³+7x²+3.

Şunları gösterebiliriz:

P_n(0), n sayısının son basamağıdır.
P_n(1), n sayısının basamakları toplamıdır.
P_n(10), n sayısıdır.

k sayısından küçük, P_n polinomunun en az bir kökünün olduğu n pozitif tam sayılarının adedini Z(k) ile tanımlayalım.

Z(100.000)=14696 olduğu doğrulanabilir.

Z(10¹⁶) kaçtır?
Cevap: 1311109198529286

C/C++:
#include <stdio.h>
#include <algorithm>

#define L 16
#define BIG 1000000000LL // 1e9

using namespace std;

typedef struct
{
   int s;
   int offset;
}List;

int compare (const void* source, const void* dest)
{
   List* listSource = (List*) source;
   List* listDest = (List*) dest;
   return (listSource->s) - (listDest->s);
}

int main()  {
    
    int h,i,j,d2,n,N,r,p2,s,s2,P2,st,en,iter,add,k,coeff[20];
    long long int tot=0,sum,x,mul,p10,ppow,pr,res;
    List* myList;
    myList=(List*)(malloc)(100000000*sizeof(List));
    
    d2=L/2;
    for(r=-9;r<=-3;r++)  {
        printf("r=%d\n",r);
        myList[0].s=0;
        myList[0].offset=0;
        ppow=1;
        p10=1;
        n=1;
        for(i=1;i<=d2;i++)  {
            N=n;
            for(h=0;h<N;h++)
                for(j=1;j<10;j++)  {
                    myList[n].s=myList[h].s+j*ppow;
                    myList[n].offset=myList[h].offset+j*p10;
                    n++;
               }
            ppow*=r;
            p10*=10;
        }

qsort(myList,n,sizeof(List),compare);
        
        P2=1;
        while(P2<n)  P2<<=1;
        
        ppow=1;
        for(i=1;i<=L-d2;i++)  ppow*=10;
        pr=1;
        for(i=1;i<=d2;i++)    pr*=r;
        
        for(i=0;i<ppow;i++)  {
            sum=0;
            x=i;
            mul=pr;
            for(j=1;j<=L-d2;j++)  sum+=mul*(x%10),x/=10,mul*=r;
            if(sum<BIG&&-BIG<sum){
            
            sum=-sum;
            p2=P2;
            s=-1;
            while(p2)  {
                  s2=s+p2;
                  if(s2<n&&myList[s2].s<sum)  s=s2;
                  p2>>=1;
            }
            st=s+1;
            s=0;
            p2=P2;
            while(p2)  {
                  s2=s+p2;
                  if(s2<n&&myList[s2].s<=sum) s=s2;
                  p2>>=1;
            }
            en=s;
            
            if(st<=en&&myList[st].s==sum)  {
               x=i;
               for(j=d2;j<L;j++)  coeff[j]=x%10,x/=10;
               
               for(iter=st;iter<=en;iter++)  {
               
               x=myList[iter].offset;
               for(j=0;j<d2;j++)  coeff[j]=x%10,x/=10;
               add=1;
               for(j=0;j>r;j--)  {
                   res=0;
                   for(k=L-1;k>=0;k--)  res=j*res+coeff[k];
                   if(res==0)  {
                      add=0;
                      break;
                   }
               }
               tot+=add;
               }
           }}
       }
       printf("tot=%lld\n",tot);
  }
  return 0;
}

Python:
def euler269(ndigits,x):
   t=0 
   D={0:1}
   for r in range(0,ndigits,2):
      print r
      D2=defaultdict(int)
      if r==0 and x!=0: startd=1
      else: startd=0
      for d in range(startd,10):
         e=(x**r)*d
         for s,v in D.items():
            D2[s+e]+=v
      D=D2
   Deven=D
   D={0:1}
   for r in range(1,ndigits,2):
      print r
      D2=defaultdict(int)
      for d in range(10):
         e=(x**r)*d
         for s,v in D.items():
            D2[s+e]+=v
      D=D2
   Dodd=D
   for s,n in Deven.items():
      t+=Dodd[-s]*n
   return t

def e269M(x,r,N,n11):
    t=0
   if N==0:
      if n11==0: return 0
      if r!=0: return 0
      return 1
   for d in range(10):
      r2=d+r
      n11b=(d-n11)
      if r2%x==0:
         t+=e269M(x,r2/-x,N-1,n11b)
   return t
e269M=Memoize(e269M)
def euler269c(D,N,n11,startd,tot):
   t=0
   if N==0:
      if n11==0: return 0
      for r in D.values():
         if r==0: return 1
      return 0

for d in range(startd,10):
      D2={}
      for x,r in D.items():
         r2=d+r
         if r2%x==0:
            D2[x]=r2/-x
      n11b=(d-n11)
      if len(D2)==0: continue
      if len(D2)==1:
         t+=e269M(D2.keys()[0],D2.values()[0],N-1,n11b)
      else:
         t+=euler269c(D2,N-1,n11b,0,tot*10+d)
   return t
t=0
N=16
t+=euler269(N,-1)
D={}
for x in range(2,10):
   D[x]=0
t+=euler269c(D,N,0,1,0)
print 10**(N-1)+t

Ruby:
# compute number of positive integers <= 10**digits where the corresponding
# polynomial has integer roots, for specified combination of roots
def count_m(digits, roots)
  # compute minimum and maximum nondivergent residuals
  min = roots.collect do |x|
    case x
    when 0
      0
    when -1
      -1.0/0
    else
      (-8.0 / (x**2 - 1)).ceil
    end
  end
  max = roots.collect do |x|
    case x
    when 0
      9
    when -1
      1.0/0
    else
      (-9*x - 1) / (x**2 - 1)
    end
  end

# res0 maps combinations of residuals to count, start with initial digits
  res0 = Hash.new(0)
  1.upto([9,max.min].min) { |i| res0[[i]*roots.length] = 1 }

# process for each number of digits
  count = 0
  zero = [0]*roots.length
  2.upto(digits) do |d|
    res1 = Hash.new(0)
    res0.each do |r0,cnt|
      # compute new in-range residuals
      min1 = r0.each_index.collect{|i| min[i] - r0[i]*roots[i]}.max
      max1 = r0.each_index.collect{|i| max[i] - r0[i]*roots[i]}.min
      [0,min1].max.upto([9,max1].min) do |d|
        r1 = r0.each_index.collect{|i| r0[i]*roots[i]+d}
        if roots[0] == 0 && d != digits
          count += cnt if r1 == zero
          r1[0] = -1            # collapse counts, avoid false zero
        end
        res1[r1] += cnt
      end
    end
    res0 = res1
    count += res1[zero]
  end

count += 1 if roots.length == 1 && roots[0] == 0 # bug fix
  count
end

def visit(digits, roots)
  count = count_m(digits, roots)
  return 0 if count == 0

(roots.last-1).downto(-9) do |r|
    count -= visit(digits, roots+[r])
  end
  count
end

# count via inclusion-exclusion
def count_all(digits)
  count = 0
  0.downto(-9) do |r|
    count += visit(digits, [r])
  end
  count
end
    
if __FILE__ == $0
  p count_all(16)
end

Java:
import java.util.HashMap;
import java.util.Map;

public class Euler269 {
 static int digits = 16;
 static int z, a, b;
 static Map<Integer, Long>[] memo = new Map[digits + 1];
 
 // fill with (x + b)
 static long w(int num_d, int have_for_last) {
  if(num_d < 2) return 0 < have_for_last && have_for_last < 10 ? 1 : 0;
  if(memo[num_d].containsKey(have_for_last)) return memo[num_d].get(have_for_last);
  long count = 0;
  int min_n = (b * 1001 - 1 - have_for_last) / b - 1000;
  int max_n = (b * 1000 + 9 - have_for_last) / b - 1000;
  for(int n = min_n; n <= max_n; ++n) count += w(num_d - 1, n);
  memo[num_d].put(have_for_last, count);
  return count;
 }
 
 // fill with (x + A) (x + B) == x^2 + a x + b 
 static long w2(int num_d, int have_for_last, int have_for_sec_last) {
  if(num_d < 3) return 0 <= have_for_last && have_for_last < 10 && 0 < have_for_sec_last && have_for_sec_last < 10 ? 1 : 0;
  int hk = have_for_last + 0x10000 * have_for_sec_last; 
  if(memo[num_d].containsKey(hk)) return memo[num_d].get(hk);
  long count = 0;
  int min_n = (b * 1001 - 1 - have_for_last) / b - 1000;
  int max_n = (b * 1000 + 9 - have_for_last) / b - 1000;
  for(int n = min_n; n <= max_n; ++n) count += w2(num_d - 1, have_for_sec_last + n * a, n);
  memo[num_d].put(hk, count);
  return count;
 }
 
 // fill with (x + A) (x + B) (x + C) == x^3 + z x^2 + a x + b 
 static long w3(int num_d, int have_for_last, int have_for_sec_last, int have_for_third_last) {
  if(num_d < 4) return 0 <= have_for_last && have_for_last < 10 && 0 <= have_for_sec_last && have_for_sec_last < 10  && 0 < have_for_third_last && have_for_third_last < 10 ? 1 : 0;
  int hk = have_for_last + 0x400 * have_for_sec_last + 0x100000 * have_for_third_last; 
  if(memo[num_d].containsKey(hk)) return memo[num_d].get(hk);
  long count = 0;
  int min_n = (b * 1001 - 1 - have_for_last) / b - 1000;
  int max_n = (b * 1000 + 9 - have_for_last) / b - 1000;
  for(int n = min_n; n <= max_n; ++n) count += w3(num_d - 1, have_for_sec_last + n * a, have_for_third_last + n * z, n);
  memo[num_d].put(hk, count);
  return count;
 }
 
 static void initB(int B) {
  b = B;
  for(int d = 0; d < memo.length; ++d) memo[d] = new HashMap<Integer, Long>();
 }
 
 static void initAB(int A, int B) {
  a = A + B;
  b = A * B;
  for(int d = 0; d < memo.length; ++d) memo[d] = new HashMap<Integer, Long>();
 }
 
 static void initABC(int A, int B, int C) {
  z = A + B + C;
  a = A * B + B * C + A * C;
  b = A * B * C;
  for(int d = 0; d < memo.length; ++d) memo[d] = new HashMap<Integer, Long>();
 }
 
 public static void main(String[] args) {
  long time = System.currentTimeMillis();
  long total = 1;
  for(int d = 1; d < digits; ++d) total *= 10;
  for(int B = 1; B < 10; ++B) {
   initB(B);
   long w = w(digits, 0);
   total += w;
   //System.out.println(B + " : " + w);
  }
  for(int A = 1; A < 10; ++A) {
   for(int B = A + 1; B < 10; ++B) {
    if(A * B >= 10) break;
    initAB(A, B);
    long w = w2(digits, 0, 0);
    total -= w;
    //System.out.println(A + " & " + B + " : " + w);
   }
  }
  for(int A = 1; A < 10; ++A) {
   for(int B = A + 1; B < 10; ++B) {
    for(int C = B + 1; C < 10; ++C) {
     if(A * B * C >= 10) break;
     initABC(A, B, C);
     long w = w3(digits, 0, 0, 0);
     total += w;
     //System.out.println(A + " & " + B + " & " + C + " : " + w);
    }
   }
  }
  time = System.currentTimeMillis() - time;
  System.out.println(total);
  System.out.println("time = " + time / 1000.0 + "s");
  
 }

}

Haskell:
import Data.Ratio ((%))
import Data.Map (toList, fromListWith)

next :: Integral a => [a] -> [a] -> [a]
next [q] [b] = [ceiling $ -b % q .. floor $ (9-b) % q]
next [q1, q0] [b1, b0] = [ceiling low .. floor upp]
    where low = - (b0 + (q1 + q0) * b1) % (q1 * q0)
          upp = low + 9 % (q1 * q0)
next [q2, q1, q0] [b2, b1, b0] = [ceiling low .. floor upp]
        where low = - (b0 + (q0 + q1 + q2) * b1 + (q0*q1 + q1*q2 + q2*q0) * b2) % (q0 * q1 * q2)
              upp = low + 9 % (q0 * q1 * q2)

feasible :: Int -> Bool
feasible x = 0 <= x && x <= 9
verify :: [Int] -> [Int] -> Bool
verify _ [b] = feasible b
verify [q1, q0] [b1, b0] = all feasible [b1, b0 + (q1+q0)*b1]
verify [q2, q1, q0] [b2, b1, b0] = all feasible [b2, b1 + s*b2, b0 + s*b1 + (q0*q1 + q1*q2 + q2*q0)*b2]
        where s = q0 + q1 + q2

count :: Int -> [Int] -> Integer
count (n+1) (0:qs) = count n qs
count n [] = 10^n
count n qs = pred.sum.map snd.filter (verify qs.fst).(!!(n-k)).iterate grow $ [(replicate k 0, 1)]
    where k = length qs
          grow ss = toList.fromListWith (+).concatMap nextS $ ss
          nextS (bs, c) = [(b : init bs, c) | b <- next qs bs]

roots :: Int -> Int -> [[Int]]
roots 0 q = [[]]
roots (n+1) q = concat [map (q':).roots n $ q'+1 | q' <- [q..9]]

p269 :: Int -> Integer
p269 n = sum (map (count n) $ roots 1 0) - sum (map (count n) $ roots 2 0)
         + sum (map (count n) $ roots 3 0) - sum (map (count n.(0:))$ roots 3 1)

main :: IO ()
main = print $ p269 16

Maple:
> maxfind:=(puis,sol)->9*add(sol^(2*i),i=0..iquo(puis,2)):
> minfind:=(puis,sol)->9*add(sol^(2*i+1),i=0..iquo(puis-1,2)):
> finder:=proc(zero,puis,sol)
> option remember;
> if puis=-1 then
>   if zero=0 then return 1
>   else return 0 fi;
> fi;
> if zero+minfind(puis,sol)>0 or zero+maxfind(puis,sol)<0 then return 0 fi;
> return add(finder(zero+sol^puis*i,puis-1,sol),i=0..9)
> end:
> 
> sfinder:=proc(tzero,puis,tsol)
> option remember;
> if puis=-1 then
>   if andmap(x->is(x=0),tzero) then return 1
>   else return 0 fi
> fi;
> if ormap(x->is(finder(x[1],puis,x[2])=0),[seq([tzero[i],tsol[i]],i=1..nops(tzero))]) then return 0 fi;
> return add(sfinder([seq(tzero[j]+tsol[j]^puis*i,j=1..nops(tzero))],puis-1,tsol),i=0..9)
> end:
> 
> with(combinat): t:=time():
> add(-(-1)^i*add(sfinder([0$i],15,j),j in choose([$-9..0],i)),i=1..9)-1;
> time()-t;

1311109198529286
                                52.219

Herhangi Matematiksel Şeyler

Sayfalar

21 Temmuz 2019

Euler Projesi 269. Soru