Sums of all Squares - Page 13

science_man_88 · 2010-12-14, 12:46

Quote:

Originally Posted by NBtarheel_33

Isn't this the sum of P^2 for p prime, NOT p^p? This thread has me confused...

If so I've got them beat!

@CRG so in other words prime(x) is almost like:

Code:

for(y=1,#primelist,if(x==y,return(primelist[x])))

science_man_88 · 2010-12-14, 12:49

Code:

(08:47)>a=0;for(m=1,100,for(x=1,#v,a=a+v[x];b=10^(m);if(a%b==0,print(prime(x));break()));a=0)
11
661
4397

where

Code:

(08:47)>v=vector(1000,n,(prime(n)^prime(n))%10)

science_man_88 · 2010-12-14, 13:08

Quote:

Originally Posted by NBtarheel_33

Isn't this the sum of P^2 for p prime, NOT p^p? This thread has me confused...

no changing v to d where it replaces prime(x)^prime(x) with prime(x)^2 gave 907 as the first term.

CRGreathouse · 2010-12-14, 15:39

Quote:

Originally Posted by science_man_88

@CRG so in other words prime(x) is almost like:

Code:

for(y=1,#primelist,if(x==y,return(primelist[x])))

The basic idea: Pari stores the differences between primes up to primelimit, so it can quickly go from one prime to the next in as long as it knows where it is on the list. (It's can't use that for nextprime() because it doesn't know what position the given prime is at.)

So while it cannot jump straight to, say, the 500th prime it can quickly go to 'the next prime' 500 times. But that's still 500 steps.

Now the situation is somewhat complicated by a speedup in the code. It can jump straight to the 1000-th, 2000-th, 3000-th, 4000-th, 5000-th, 6000-th, 10,000-th, 20,000-th, 30,000-th, or 40,000-th prime. So prime(1000) takes only about 60 nanoseconds, while prime(999) takes 2350 nanoseconds (2.35 microseconds).

But as you go beyond 40,000, prime() can take a really long time. prime(10^7) takes about 29 milliseconds, ten thousand times longer than prime(999).

Now, all of this probably sounds fairly abstract, since we're talking about really short times. But when you actually do a real loop you see the difference right away. for(n=1,10^5,prime(n)) takes 5.5 seconds, compared to the 4.5 milliseconds that forprime(p=2,prime(10^5),p) takes. Increase those to 10^6 and forprime takes 44 milliseconds compared to over 20 minutes (I killed the script) with for(n=1,10^6,prime(n)).

science_man_88 · 2010-12-14, 15:50

Quote:

Originally Posted by CRGreathouse

The basic idea: Pari stores the differences between primes up to primelimit, so it can quickly go from one prime to the next in as long as it knows where it is on the list. (It's can't use that for nextprime() because it doesn't know what position the given prime is at.)

So while it cannot jump straight to, say, the 500th prime it can quickly go to 'the next prime' 500 times. But that's still 500 steps.

Now the situation is somewhat complicated by a speedup in the code. It can jump straight to the 1000-th, 2000-th, 3000-th, 4000-th, 5000-th, 6000-th, 10,000-th, 20,000-th, 30,000-th, or 40,000-th prime. So prime(1000) takes only about 60 nanoseconds, while prime(999) takes 2350 nanoseconds (2.35 microseconds).

But as you go beyond 40,000, prime() can take a really long time. prime(10^7) takes about 29 milliseconds, ten thousand times longer than prime(999).

Now, all of this probably sounds fairly abstract, since we're talking about really short times. But when you actually do a real loop you see the difference right away. for(n=1,10^5,prime(n)) takes 5.5 seconds, compared to the 4.5 milliseconds that forprime(p=2,prime(10^5),p) takes. Increase those to 10^6 and forprime takes 44 milliseconds compared to over 20 minutes (I killed the script) with for(n=1,10^6,prime(n)).

yeah the problem with my script is the for loop you are complaining about goes through the indexes of an vector so a forprime would only check prime indexes.

CRGreathouse · 2010-12-14, 16:03

Quote:

Originally Posted by science_man_88

yeah the problem with my script is the for loop you are complaining about goes through the indexes of an vector so a forprime would only check prime indexes.

huh? I was talking about your code in #126:

Code:

a=0;for(m=1,2,for(n=1,1000,a=a+prime(n)^prime(n);if(a%(10^m)==0,print(n":"m)));a=0)

which can be improved like so:

Code:

for(m=1,2,n=0;a=0;forprime(p=2,7919,a+=p^p;n++;if(a%(10^m)==0,print(n":"m))))

Even better, if you're willing to accept a different format for output, you can make it

Code:

for(m=1,2,a=0;forprime(p=2,7919,a+=p^p;if(a%(10^m)==0,print(p" "m))))

Of course none of this matters much in this case, since essentially all of the time is spent creating and storing large numbers. But the general principle is important -- see my post above where you can see that even getting a million primes with prime() can take over 20 minutes.

science_man_88 · 2010-12-14, 16:06

Code:

(09:07)>v=vector(1000,n,(prime(n)^prime(n))%10)
%208 = [4, 7, 5, 3, 1, 3, 7, 9, 7, 9, 1, 7, 1, 7, 3, 3, 9, 1, 3, 1, 3, 9, 7, 9, 7, 1, 7, 3, 9
(09:07)>a=0;for(m=1,100,for(x=1,#v,a=a+v[x];b=10^(m);if(a%b==0,print(prime(x));break()));a=0)
11
661
4397

is my latest code

science_man_88 · 2010-12-14, 16:08

Quote:

Originally Posted by CRGreathouse

huh? I was talking about your code in #126:

Code:

a=0;for(m=1,2,for(n=1,1000,a=a+prime(n)^prime(n);if(a%(10^m)==0,print(n":"m)));a=0)

which can be improved like so:

Code:

for(m=1,2,n=0;a=0;forprime(p=2,7919,a+=p^p;n++;if(a%(10^m)==0,print(n":"m))))

Even better, if you're willing to accept a different format for output, you can make it

Code:

for(m=1,2,a=0;forprime(p=2,7919,a+=p^p;if(a%(10^m)==0,print(p" "m))))

Of course none of this matters much in this case, since essentially all of the time is spent creating and storing large numbers. But the general principle is important -- see my post above where you can see that even getting a million primes with prime() can take over 20 minutes.

Code:

(11:52)>for(m=1,2,a=0;forprime(p=2,7919,a+=p^p;if(a%(10^m)==0,print(p" "m))
11 1
17 1
67 1
83 1
137 1
191 1
197 1
211 1
227 1
257 1
331 1
347 1
401 1
431 1
563 1
571 1
587 1
653 1
661 1
677 1
727 1
751 1
773 1
797 1
811 1
829 1
883 1
977 1
1039 1
1051 1
1087 1
1129 1
1201 1
1217 1
1303 1
1367 1
1453 1
1471 1
1489 1
1523 1
1553 1
1607 1
1777 1
1787 1
1801 1
1873 1
1999 1
2221 1
2281 1
2351 1
2447 1
2521 1
2539 1
2647 1
2689 1
2789 1
2857 1
2879 1
2897 1
3067 1
3169 1
3329 1
3359 1
3457 1
3469 1
3499 1
3529 1
3613 1
3691 1
3803 1
4099 1
4217 1
4241 1
4253 1
4261 1
4289 1
4327 1
4397 1
4421 1
4507 1
4517 1
4591 1
4723 1
4903 1
5009 1
5153 1
5303 1
5393 1
5441 1
5483 1
5519 1
5623 1
5641 1
5717 1
5783 1
6113 1
6163 1
6197 1
6269 1
6299 1
6481 1
6763 1
6781 1
6829 1
6967 1
7079 1
7127 1
7151 1
7393 1
7507 1
7537 1
7583 1
7591 1
7607 1
7639 1
7753 1
751 2
797 2
811 2
883 2
1129 2
1201 2
1367 2
1453 2
1553 2
3529 2
3613 2
3803 2
4507 2
5393 2

just as crappy as it prints too much.

lavalamp · 2010-12-14, 17:21

I checked both sequences (p^2 and p^p) as far as I know how to in Pari (without writing my own prime generating code for primes beyond 4294965247 anyway).

Code for p^2:

Code:

default(primelimit,4294965247);

t=0;	\\ This is the running total for the sum of prime squares.
c=1;	\\ This is the number of zeroes we want at the end for the next output.

forprime(p=2,4294965247,t+=p^2;while(t%10^c==0,write("C:/s.txt",c,", ",p);c++));

Which gives this output:

Code:

1, 907
2, 977
3, 977
4, 36643
5, 1067749
6, 17777197
7, 71622461
8, 2389799983

Then for p^p I modified the first code slightly:

Code:

default(primelimit,4294965247);

t=Mod(0,10^15);	\\ This is the running total for the sum of prime squares.
c=1;		\\ This is the number of zeroes we want at the end for the next output.

forprime(p=2,4294965247,t+=Mod(p,10^15)^p;while(t%10^c==0,write("C:/p.txt",c,", ",p);c++));

That gives this output:

Code:

1, 11
2, 751
3, 1129
4, 361649
5, 361649
6, 12462809
7, 12462809
8, 1273183931
9, 1273183931

The first code took about 3 minutes to run, the second code unsurprisingly took somewhat longer, about 24 minutes.

All results match what everyone else got.

CRGreathouse · 2010-12-14, 17:25

Quote:

Originally Posted by science_man_88

just as crappy as it prints too much.

I was trying to make my code act the same as yours. If (say) you just wanted one answer for each size, do

Code:

find(m)={
  my(s=Mod(0,m));
  forprime(p=2,default(primelimit),
    s+=Mod(p,m)^p;
    if(!s, return(p))
  )
};
{for(n=1,6,
  p=find(10^n);
  print("2^2 + 3^3 + ... + "p"^"p" ends in "n" zero(s).")
)}

CRGreathouse · 2010-12-14, 17:36

Quote:

Originally Posted by lavalamp

All results match what everyone else got.

Glad to hear it.

Quote:

Originally Posted by lavalamp

I checked both sequences (p^2 and p^p) as far as I know how to in Pari (without writing my own prime generating code for primes beyond 4294965247 anyway).

I think I checked to 10^10 without finding any further terms. 64-bit Pari doesn't have the 4294965247 limitation.

I really need to stop being lazy and finish my Pari (not GP) code for finding primes above primelimit (tentatively "forbigprime"); I feel that this would be useful for a lot of people, and I'm sure Karim Belabas would be happy to accept the code if I submitted it.

2010-12-14, 12:49	#134
science_man_88 "Forget I exist" Jul 2009 Dumbassville 8384₁₀ Posts	Code: (08:47)>a=0;for(m=1,100,for(x=1,#v,a=a+v[x];b=10^(m);if(a%b==0,print(prime(x));break()));a=0) 11 661 4397 where Code: (08:47)>v=vector(1000,n,(prime(n)^prime(n))%10)

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Regarding Squares	a1call	Miscellaneous Math	42	2017-02-03 01:29
Basic Number Theory 12: sums of two squares	Nick	Number Theory Discussion Group	0	2016-12-11 11:30
Integers = sums of 2s and 3s.	3.14159	Miscellaneous Math	12	2010-07-21 11:47
Sums of three squares	CRGreathouse	Math	6	2009-11-06 19:20
squares or not squares	m_f_h	Puzzles	45	2007-06-15 17:46