mersenneforum.org (https://www.mersenneforum.org/index.php)
-   Puzzles (https://www.mersenneforum.org/forumdisplay.php?f=18)
-   -   Sums of all Squares (https://www.mersenneforum.org/showthread.php?t=13181)

 davar55 2010-03-18 14:00

Sums of all Squares

2^2 + 3^2 + 5^2 + ... + p^2 = 10[sup]m[/sup]K

What is the smallest prime p such that
the sum of squares of all primes up to p
is a multiple of 10 (or 100 or 1000).

 axn 2010-03-18 14:26

[SPOILER]s=0;forprime(p=2,1000,s=s+p^2;if(Mod(s,10)==0, print(p, ":",s)))
907:37464550
967:44505630
977:46403000
991:48351370[/SPOILER]

 davar55 2010-03-18 14:32

Nice and simple and quick reply -- thanks.

 CRGreathouse 2010-03-18 17:11

907, 977, 977, 36643, 1067749, 17777197, 71622461, 2389799983, ...

The next term (if one exists) is more than 4 trillion.

Not yet in the OEIS.

[URL="http://www.research.att.com/%7Enjas/sequences/"]http://www.research.att.com/~njas/sequences/[/URL]

[quote]Search: [B]907,[/B] [spoiler][B] 967, 977, 991[/B] [/spoiler]

I am sorry, but the terms do not match anything in the table.

If your sequence is of general interest, please send it to me using the [URL="http://www.research.att.com/%7Enjas/sequences/Submit.html"]form provided[/URL] and I will (probably) add it to the data base! Include a brief description and if possible enough terms to fill 3 lines on the screen. I need a minimum of 4 terms.[/quote]I think it qualifies. Also, I'm fond of OEIS entries with relatively large initial terms -- especially when the next few terms are so closely spaced as in this one. (Might it set some record in that regard -- highest ratio of initial term to average spacing of next n terms, for n = 3?)

I'd be glad to submit it, but I think it should be one of you guys.

How about generalizing to other bases?

 axn 2010-03-18 19:02

I think CRG's sequence is more "worthy". It is also the solution of OP.

I can think of two ways to generalize: to other bases and other powers (other than squares).

[quote=axn;208810]I think CRG's sequence is more "worthy".[/quote][spoiler]Even with the omitted and repeated (just a typo) terms? :smile:[/spoiler]

What I had in mind was a submission with the best of both your contributions.

[quote]I can think of two ways to generalize: to other bases and other powers (other than squares).[/quote]Oh, wow ... bases 2-16 or so, powers to, say, ninth ==> 135 sequences.

 axn 2010-03-18 19:59

[QUOTE=cheesehead;208818]Even with the omitted and repeated terms?[/QUOTE]
That sequence is the first occurrence of 10^n. 977 repeats (not a typo!) because it ends in 000 and comes before any other 00. So it stands at positions 2 & 3.

EDIT:- Mine is merely the first four occurrences of 10

[quote=axn;208819]That sequence is the first occurrence of 10^n. 977 repeats because it ends in 000 and comes before any other 00. So it stands at positions 2 & 3.[/quote]:doh!:

(Sorry, CRG)

But that doubles the potential number of sequences. 270.

 CRGreathouse 2010-03-19 01:28

[QUOTE=CRGreathouse;208783]The next term (if one exists) is more than 4 trillion.[/QUOTE]

:blush: That's *billion*, not trillion. Now my search limit is 50 billion, giving me

907, 977, 977, 36643, 1067749, 17777197, 71622461, 2389799983, 31252968359, 49460594569, ...

The nth term is very roughly n * log 10 * 10^n, so I was pretty lucky getting that last term. The next one will probably need over 2 trillion. Anyone up to the task? I don't actually have a good segmented sieve coded at the moment...

 bsquared 2010-03-29 16:42

[QUOTE=CRGreathouse;208849] Anyone up to the task? I don't actually have a good segmented sieve coded at the moment...[/QUOTE]

I bit. :smile:

I have [url=http://www.mersenneforum.org/showpost.php?p=209641&postcount=359]recently[/url] spent some time with my sieve, so decided to give this a shot.

I just started a run to 2 trillion. Here is the output so far:
[CODE]found primes in range 0 to 1000000000 in elapsed time = 7.0227
**** 907 is 0 mod 10 ****
**** 977 is 0 mod 100 ****
**** 977 is 0 mod 1000 ****
**** 36643 is 0 mod 10000 ****
**** 1067749 is 0 mod 100000 ****
**** 17777197 is 0 mod 1000000 ****
**** 71622461 is 0 mod 10000000 ****
sum of squares complete in elapsed time = 8.5178, sum is 16352255694497179054764665

found primes in range 1000000000 to 2000000000 in elapsed time = 5.9418
sum of squares complete in elapsed time = 7.9423, sum is 126512354351558021982865866

found primes in range 2000000000 to 3000000000 in elapsed time = 5.9503
**** 2389799983 is 0 mod 100000000 ****
sum of squares complete in elapsed time = 7.7389, sum is 418923904898718760122282892

found primes in range 3000000000 to 4000000000 in elapsed time = 5.8990
sum of squares complete in elapsed time = 7.6150, sum is 979895993641271252685833855

found primes in range 4000000000 to 5000000000 in elapsed time = 5.8293
sum of squares complete in elapsed time = 7.4966, sum is 1894402266333772221759233898
[/CODE]

With these timing trends, should have a result in 7 hours or so.

- ben.

 bsquared 2010-03-29 16:46

update...

[CODE]found primes in range 30000000000 to 31000000000 in elapsed time = 5.4245
sum of squares complete in elapsed time = 7.1620, sum is 416903941002774697723222981803

found primes in range 31000000000 to 32000000000 in elapsed time = 5.5802
**** 31252968359 is 0 mod 1000000000 ****
sum of squares complete in elapsed time = 7.1759, sum is 457955303775896882861615585442

found primes in range 32000000000 to 33000000000 in elapsed time = 5.5252
sum of squares complete in elapsed time = 7.1461, sum is 501598601070515778427418232428
[/CODE]

and another update:
[CODE]found primes in range 47000000000 to 48000000000 in elapsed time = 5.5268
sum of squares complete in elapsed time = 7.0336, sum is 1519756369296424391708040649758

found primes in range 48000000000 to 49000000000 in elapsed time = 5.4239
sum of squares complete in elapsed time = 7.0347, sum is 1615357580573805620690452754303

found primes in range 49000000000 to 50000000000 in elapsed time = 5.4835
**** 49460594569 is 0 mod 1410065408 ****
sum of squares complete in elapsed time = 6.8852, sum is 1714863031171407826702942323341

found primes in range 50000000000 to 51000000000 in elapsed time = 5.4299
[/CODE]

I'm stupidly using %u to print the modulus, but it is stored internally as a 64 bit integer... so it is only a printing error.

 bsquared 2010-03-29 22:57

Right on cue

we have this result:

[CODE]found primes in range 1915000000000 to 1916000000000 in elapsed time = 5.0965
**** 1915014433303 is 0 mod 1215752192 ****
sum of squares complete in elapsed time = 5.0550, sum is 83903230112675776937166385335972895
[/CODE]

So the sum of primes squared up to 1915014433303 is zero mod 100e9. I'm processing about 1 billion numbers per second, so assuming the trend of this sequence holds, to go to the next value at around 20 trillion would take a couple days.

 bsquared 2010-03-30 02:35

[quote=bsquared;210004]...I'm processing about 1 billion numbers per second...[/quote]

That should have been 1e9 every 10 sec, but the time estimate is still about right. It's probably silly, but a run to 20 trillion is ongoing...

Here was the sum of all prime squares up to 1915014433303:
83775363722237720731978600000000000

I'm keeping a file with the sums every 1e9, in case anyone wants to extend the sequence after I get tired of it or for double checks.

 bsquared 2010-03-30 12:57

I see that the death knell of this thread has been sounded :smile:

Maybe the 12th member of the sequence is still of interest...

$$\Sigma_{p = 2}^{4076200167673} p^2 = 786646994677132840800629000000000000$$

 CRGreathouse 2010-03-30 15:07

[QUOTE=bsquared;210051]I see that the death knell of this thread has been sounded :smile:

Maybe the 12th member of the sequence is still of interest...

$$\Sigma_{p = 2}^{4076200167673} p^2 = 786646994677132840800629000000000000$$[/QUOTE]

Are you going to submit this to the OEIS?

 bsquared 2010-03-30 15:40

I have no problem doing so... but I'm not the "discoverer" of this sequence. I'll defer to you or davar55 if you would rather take the credit.

 CRGreathouse 2010-03-30 15:44

I PM'd davar55.

 bsquared 2010-03-30 15:52

Also, [url=http://www.research.att.com/~njas/sequences/A024450]this sequence[/url] could now be greatly extended.

 CRGreathouse 2010-03-30 16:12

[QUOTE=bsquared;210069]Also, [url=http://www.research.att.com/~njas/sequences/A024450]this sequence[/url] could now be greatly extended.[/QUOTE]

Sure. It's quite easy to extend, but tradition limits b-files to 10,000 entries. If you'd like you can extend it to that (or I can), but I wouldn't go beyond. The current b-file has 5000.

As it happens I never computed that sequence for these calculations -- I used pure modular arithmetic. (If I used BCD I could have avoided this while keeping speed high...)

 bsquared 2010-03-30 17:01

[QUOTE=CRGreathouse;210070] Sure. It's quite easy to extend, but tradition limits b-files to 10,000 entries. If you'd like you can extend it to that (or I can), but I wouldn't go beyond. The current b-file has 5000.
[/QUOTE]
I guess I was thinking of a link to tables, or something. But that would require me to generate and host those tables. I'll save that for some day when I'm bored ;)

[QUOTE=CRGreathouse;210070]
As it happens I never computed that sequence for these calculations -- I used pure modular arithmetic. (If I used BCD I could have avoided this while keeping speed high...) [/QUOTE]

Yeah, that's definitely faster, but using pure modular arithmetic would require you to start the sum over for each new modulus, right?

 CRGreathouse 2010-03-30 17:17

[QUOTE=bsquared;210073]Yeah, that's definitely faster, but using pure modular arithmetic would require you to start the sum over for each new modulus, right?[/QUOTE]

Right. You could do two modili at a time without much penalty, though, with appropriate lookup tables and bit operations. Note that you only need to compare (and hence reduce) every 8 primes, each term (other than the first) has index = 5 (mod 8).

 CRGreathouse 2010-03-30 17:27

[QUOTE=bsquared;210073]I guess I was thinking of a link to tables, or something. But that would require me to generate and host those tables. I'll save that for some day when I'm bored ;)[/QUOTE]

Maybe you can submit a sequence of every millionth term or something. Not exactly classy, but it has precedent (A080128, say).

76304519151822049179, 671924965564646162227, 2393465488665494654963, 5889405149040404480379, 11834774513923727795971, 20925456417823033330259, ...

 bsquared 2010-03-30 19:00

[QUOTE=CRGreathouse;210078]Right. You could do two modili at a time without much penalty, though, with appropriate lookup tables and bit operations. Note that you only need to compare (and hence reduce) every 8 primes, each term (other than the first) has index = 5 (mod 8).[/QUOTE]

I just realized my code to do the sum of prime squares routines was hugely inefficient, in that I was using YAFU's built in arbitrary precision functions to do the squaring and summing. In reality we only need fixed precision of, say, 3 64 bit limbs (192 bits should be plenty to represent the sum). Implementing this made my sum of prime squares routine about 35 times faster.

I also store the highest power of 2 dividing the power of 10 modulus, which makes for a very quick pre-test of divisibility by 10 (logical AND followed by a predictable branch) and makes full precision divisions *extremely* rare. Doing things this way is actually faster than using pure modular arithmetic, since we almost never have to perform a division.

As a side benefit, it's easy to build tables of prime sums, or prime square sums, and we also don't have to restart a sum to test for a new power of 10 modulus.

 bsquared 2010-03-30 19:02

[QUOTE=CRGreathouse;210079] Not exactly classy ...[/QUOTE]

Right up my alley, then :smile:

 bsquared 2010-03-30 19:43

I must have some sort of tinkerers disease... can't leave well enough alone.

This disease was causing me to be offended by how long it was taking to compute the primes in a range of 1e9. So I tinkered... and doubled the speed :smile:

before:
[CODE]found 40609038 primes in range 49000000000 to 50000000000 in elapsed time = 5.4835
**** 49460594569 is 0 mod 1410065408 ****
sum of squares complete in elapsed time = 6.8852, sum is 1714863031171407826702942323341
[/CODE]

after:
[CODE]found 40609038 primes in range 49000000000 to 50000000000 in elapsed time = 2.8866
**** 49460594569 is 0 mod 10000000000 ****
sum of squares complete in elapsed time = 0.1639, sum is 1714863031171407826702942323341
[/CODE]

which of course is completely useless, but now I feel better.

 davar55 2010-03-30 20:00

This is fine work by all of you. If you wish to submit the sequence to
oeis, please go ahead. I couldn't do justice to the calculations, which
I'm really impressed by. Joint discovery (attribution) is fine.

 CRGreathouse 2010-03-30 20:14

[QUOTE=davar55;210100]I couldn't do justice to the calculations, which
I'm really impressed by.[/QUOTE]

I'll echo this. I just used a one-line Pari script to discover mine. :smile:

 bsquared 2010-04-02 14:41

Now in OEIS: [URL="http://www.research.att.com/~njas/sequences/A174862"]A174862[/URL]

 davar55 2010-04-06 20:28

With all the work done on the OP, it shouldn't be too hard
to generalize the problem a bit.

I think cubes.

2^3 + 3^3 + 5^3 + ... + p^3 = 10[sup]m[/sup]K

What is the smallest prime p such that
the sum of cubes of all primes up to p
is a multiple of 10 (or 100 or 1000 or 10000 or ...).

I'm also curious about how these (squares and cubes) results
compare to first powers (sum of primes themselves).

Since these series depend on the properties of a number
in base ten, they could be considered recreational --
interesting but not necessarily useful. Still, perhaps the
sequence of sequences can someday be used to derive some
important number theoretic fact. That's one of the purposes
of the oeis.

 petrw1 2010-04-06 22:33

[QUOTE=davar55;210785]I'm also curious about how these (squares and cubes) results compare to first powers (sum of primes themselves).[/QUOTE]

A start on the first powers for the first 6,000,000 primes - p<=104,395,301

5 10
23 100
35677 63731000
106853 515530000
632501 15570900000

 bsquared 2010-04-06 22:36

[QUOTE=petrw1;210796]A start on the first powers for the first 1,000,000 primes - p<=15,485,863

5 10
23 100
35677 63731000
106853 515530000
632501 15570900000[/QUOTE]

 petrw1 2010-04-06 22:44

Crap...I missed the next 2: Bug alert!

 bsquared 2010-04-07 00:34

Crap...I missed the next 2: Bug alert![/quote]

What program are you using? It's entirely possible the bug is in my code.

 axn 2010-04-07 00:48

[QUOTE=bsquared;210804]What program are you using? It's entirely possible the bug is in my code.[/QUOTE]

Pari concurs with your first seven

[CODE]p=2;s=0;ten=10;for(i=2,10000000,s=s+p;if(Mod(s,ten)==0, ten=ten*10;print(i,":",p, ":",s)); p=nextprime(p+1))

4:5:10
10:23:100
3796:35677:63731000
10184:106853:515530000
51532:632501:15570900000
1926966:31190879:29057028000000
3471032:58369153:98078160000000[/CODE]

 axn 2010-04-07 01:01

While we're at it, some results for the cubes:
[CODE]3:5:160
51:233:143309500
1095:8783:167992435025000
2739:24763:9495929161130000
401785:5828099:18803849605481106073200000
616801:9229931:114943299218925309364000000
14335805:262707241:62239590622437034770047320000000
[/CODE]

 petrw1 2010-04-07 03:40

[QUOTE=bsquared;210804]What program are you using? It's entirely possible the bug is in my code.[/QUOTE]

Plain, ordinary, BASIC.

My guess is I lost precision when the sum got too high.

 bsquared 2010-04-07 04:46

[QUOTE=axn;210806]While we're at it, some results for the cubes:
[CODE]3:5:160
51:233:143309500
1095:8783:167992435025000
2739:24763:9495929161130000
401785:5828099:18803849605481106073200000
616801:9229931:114943299218925309364000000
14335805:262707241:62239590622437034770047320000000
[/CODE][/QUOTE]

Doublecheck plus a few more:
[CODE]**** 10 divides prime cube sum up to 5, sum = 160 ****
**** 100 divides prime cube sum up to 233, sum = 143309500 ****
**** 1000 divides prime cube sum up to 8783, sum = 167992435025000 ****
**** 10000 divides prime cube sum up to 24763, sum = 9495929161130000 ****
**** 100000 divides prime cube sum up to 5828099, sum = 18803849605481106073200000 ****
**** 1000000 divides prime cube sum up to 9229931, sum = 114943299218925309364000000 ****
**** 10000000 divides prime cube sum up to 262707241, sum = 62239590622437034770047320000000 ****
**** 100000000 divides prime cube sum up to 7717488553, sum is 39389391603365585735745579849700000000 ****
**** 1000000000 divides prime cube sum up to 34529828929, sum = 14800565770732540706707662233175000000000 ****
**** 10000000000 divides prime cube sum up to 311995561321, sum = 90365528187658782254536155073531290000000000 ****

[/CODE]

By morning the search should be at (all primes below) 10 trillion.

 bsquared 2010-04-07 14:08

[QUOTE=bsquared;210814]
By morning the search should be at (all primes below) 10 trillion.[/QUOTE]

Done. Here's the 11th term:

[CODE]**** 100000000000 divides prime cube sum up to 549120448879, sum is 848814744633978332442418792098769600000000000 ****[/CODE]

 bsquared 2010-04-07 14:17

[QUOTE=petrw1;210813]Plain, ordinary, BASIC.

My guess is I lost precision when the sum got too high.

Yes, that could be the issue. I'm using 3 64 bit words to hold the sum, which may not even be big enough for the cube sum problem eventually. 192 bits can hold a 58 digit number without loss of precision and the cube sum up to 10 trillion is already at 49 digits.

A single 64 bit integer runs out of precision at about 1.8e19, which is enough for the prime sum up to about 9 billion. Since you found the point at (632501 15570900000), you must be using more than 32 bits... but with 64 bits you should be fine for a few more points of the sequence. So maybe its an issue with your prime generation as well. A useful double check is to see if the count of primes you've summed agrees with known prime counts. For instance [URL="http://www.trnicely.net/pi/pix_0000.htm"]here[/URL].

 petrw1 2010-04-07 17:34

YUP- precision.

[QUOTE=petrw1;210813]Plain, ordinary, BASIC.

My guess is I lost precision when the sum got too high.

I can match the first 7 now.

 davar55 2010-04-08 14:28

Another way to extend this problem is to use
a base other than ten. I think binary.

2^2 + 3^2 + 5^2 + ... + p^2 = 2[sup]m[/sup]K

What is the smallest prime p such that
the sum of squares of all primes up to p
is a multiple of 2 (or 4 or 8 or 16 or ...).

This question can also be asked of first powers
or cubes of primes.

Since we basically compute in decimal or binary,
if there's an interesting number theoretic fact here
we may find it in one of these two related
sequences of sequences.

[quote=cheesehead;208818]Oh, wow ... bases 2-16 or so, powers to, say, ninth ==> 135 sequences.[/quote]...

 alexhiggins732 2010-04-10 01:35

[quote=bsquared;210847]Done. Here's the 11th term:

[code]**** 100000000000 divides prime cube sum up to 549120448879, sum is 848814744633978332442418792098769600000000000 ****[/code]I updated my [URL="http://sites.google.com/site/bbuhrow/home/sums-of-prime-squares"]webpage [/URL]too.[/quote]

n Pmax Sum(P^2) from 2 to Pmax; a multiple of 10^n
1 907 37464550[B]
2[/B][B] 977[/B][B] 46403000[/B][B]
3[/B][B] 977[/B][B] 46403000[/B]

You have duplicated 2 and 3 on your site. I also noted that the OEIS sequence didn't have example code. Perhaps you could submit the expression you have entered into YAFU.

 bsquared 2010-04-10 01:48

[quote=alexhiggins732;211243] n Pmax Sum(P^2) from 2 to Pmax; a multiple of 10^n
1 907 37464550
[B]2[/B][B] 977[/B][B] 46403000[/B]
[B]3[/B][B] 977[/B][B] 46403000[/B]
You have duplicated 2 and 3 on your site.
[/quote]
That's not an error... see[URL="http://www.mersenneforum.org/showpost.php?p=208819&postcount=8"][COLOR=#0066cc] post 8[/COLOR][/URL].
[quote=alexhiggins732;211243]
I also noted that the OEIS sequence didn't have example code. Perhaps you could submit the expression you have entered into YAFU. [/quote]

I could submit the PARI code others have used, for instance CRGreathouse:
[CODE]
s=Mod(0,10^4);forprime(p=2,1e9,s+=p^2;if(!s,return(p)))
[/CODE]

or axn, in [URL="http://www.mersenneforum.org/showpost.php?p=210805&postcount=37"]post 37[/URL]:
[CODE]p=2;s=0;ten=10;for(i=2,10000000,s=s+p;if(Mod(s,ten)==0, ten=ten*10;print(i,":",p, ":",s)); p=nextprime(p+1))
[/CODE]

The YAFU code is a little too unwieldy for OEIS... it is a hundred lines or so of C code mixed in with some custom assembly and leveraging thousands of lines of more C code implementing the sieve of eratosthenes. I have no problem making it public, but wasn't aware of any interest in it.

 alexhiggins732 2010-04-10 02:00

[quote=bsquared;211245]That's not an error... see[URL="http://www.mersenneforum.org/showpost.php?p=208819&postcount=8"][COLOR=#0066cc] post 8[/COLOR][/URL].

I could submit the PARI code others have used, for instance CRGreathouse:
[code]
s=Mod(0,10^4);forprime(p=2,1e9,s+=p^2;if(!s,return(p)))
[/code]or axn, in [URL="http://www.mersenneforum.org/showpost.php?p=210805&postcount=37"]post 37[/URL]:
[code]p=2;s=0;ten=10;for(i=2,10000000,s=s+p;if(Mod(s,ten)==0, ten=ten*10;print(i,":",p, ":",s)); p=nextprime(p+1))
[/code]The YAFU code is a little too unwieldy for OEIS... it is a hundred lines or so of C code mixed in with some custom assembly and leveraging thousands of lines of more C code implementing the sieve of eratosthenes. I have no problem making it public, but wasn't aware of any interest in it.[/quote]

The more examples the merrier. Personally, I feel code examples help me to understand the sequences better, as I think can look at code an understand it easier then an explanation using a mathematical explanation.

The best example may be I frankly could not understand the first mathematical notation I read to explain the Lucas-Lehmer test, but looking at the code it understood the notation that was used immediately.

As far as the YAFU code, I was under the impression that you a simple PARI like expression was entered into YAFU.

I am also interested in how all primes up to 10 trillion have been calculated, or have I read it wrong and that is the sum of the primes is over 10 Trillion.

 bsquared 2010-04-10 02:11

[quote=alexhiggins732;211247] The more examples the merrier. Personally, I feel code examples help me to understand the sequences better, as I think can look at code an understand it easier then an explanation using a mathematical explanation.

The best example may be I frankly could not understand the first mathematical notation I read to explain the Lucas-Lehmer test, but looking at the code it understood the notation that was used immediately.

As far as the YAFU code, I was under the impression that you a simple PARI like expression was entered into YAFU. [/quote]

YAFU is nowhere near as powerful an environment as PARI. The only reason I used it at all was because it has a rather fast prime generation routine. But I didn't "use" it by constructing some combination of built-in functions. It probably wasn't clear, but I "used" it by writing new code which leveraged functions elsewhere in the code base.

The way it's written, I doubt the YAFU code would make anything clearer to you :smile:.

[quote=alexhiggins732;211247]
I am also interested in how all primes up to 10 trillion have been calculated, or have I read it wrong and that is the sum of the primes is over 10 Trillion. [/quote]

Going on 100 trillion now, actually :smile:.

The sum is a number with 26 digits, and the square sum is a number with 39 digits. The primes are calculated using the sieve of eratosthenes. That code is public domain... just download YAFU-1.18 and browse around through the soe.c file. Although it's not well documented :blush:.

I compute the primes in batches of 1 billion integers. At a height of 10 trillion, that is 32 million primes or so per batch. It's then pretty straightfoward to sum them and check each for divisiblity by a power of 10 (sweeping many optimizations under the rug, for clarity). On my (admittedly fast) computer, each batch takes about 3 seconds for prime computation, and .2 seconds for the sum/check. This is the only code that's not public yet. Although yafu-1.19 will probably retain it when it comes out.

 alexhiggins732 2010-04-10 02:31

[quote=bsquared;211249]YAFU is nowhere near as powerful an environment as PARI. The only reasons I used it at all was because it has a rather fast prime generation routine. But I didn't "use" it by constructing some combination of built-in functions. It probably wasn't clear, but I "used" it by writing new code which leveraged functions elsewhere in the code base.

The way it's written, I doubt the YAFU code would make anything clearer to you :smile:.[/quote]

I do not seek clarification on this, I understand this. I was making a generalized statement and further it may not be clear immediately reading the description of the sequence. However, a code example helps with that clarification.

Here's .Net Code If you wish to post as well.

VB
[code]dim ten as bigint=10:dim sum as bigint=2:for each prime in primesbelow(10000000):if ten.mod(sum)=0 then : ten*=10 : console.writeline("{0}:{1}:{2}",i,p,s): end if: sum+=prime: next[/code]C#
[code]
bigint ten = 10;bigint sum = 2; foreach (var prime in primesbelow(10000000)) { if (ten.mod(sum) == 0) {ten *= 10; console.writeline("{0}:{1}:{2}", i, p, s); } sum += prime; }}
[/code]

 alexhiggins732 2010-04-10 02:46

[quote=bsquared;211249]YAFU is nowhere near as powerful an environment as PARI. The only reason I used it at all was because it has a rather fast prime generation routine. But I didn't "use" it by constructing some combination of built-in functions. It probably wasn't clear, but I "used" it by writing new code which leveraged functions elsewhere in the code base.

The way it's written, I doubt the YAFU code would make anything clearer to you :smile:.

Going on 100 trillion now, actually :smile:.

The sum is a number with 26 digits, and the square sum is a number with 39 digits. The primes are calculated using the sieve of eratosthenes. That code is public domain... just download YAFU-1.18 and browse around through the soe.c file. Although it's not well documented :blush:.

I compute the primes in batches of 1 billion integers. At a height of 10 trillion, that is 32 million primes or so per batch. It's then pretty straightfoward to sum them and check each for divisiblity by a power of 10 (sweeping many optimizations under the rug, for clarity). On my (admittedly fast) computer, each batch takes about 3 seconds for prime computation, and .2 seconds for the sum/check. This is the only code that's not public yet. Although yafu-1.19 will probably retain it when it comes out.[/quote]

Must be very fast. I use a "similar sliding window" technique if you wish to call it that. However, after each window I need to save the primes to a .bin file in the form of 6k+=1, (eg, skip 2 and 3 and 5,7,11,13,17,19,23,25 =11110000, for each bit and after each iteration start at the beginning of the binary file and calculate the start index of each prime upto $$sqrt(size)$$and sieve the current window, saving results to binary and and continuing.

However, for all primes below 2^32 the file size is 178,956,968 bytes (170 MB) and it doubles with each power of two which led me to ask how how you could sieve up to 2^44 which would require a huge amount of memory even when saved in the most compact possible form (1 bit for each 6K+=1)

Which brings up another point. If the distribution of primes are indeed random, how then is it possible to compress the aforementioned file and receive a compression ratio of 60% (eg 170mb compresses to 102MB)?

 CRGreathouse 2010-04-10 03:40

If you use my code, be sure to change "4" to "n" (and perhaps prefix "a(n)=").

[QUOTE=bsquared;211245]The YAFU code is a little too unwieldy for OEIS... it is a hundred lines or so of C code mixed in with some custom assembly and leveraging thousands of lines of more C code implementing the sieve of eratosthenes. I have no problem making it public, but wasn't aware of any interest in it.[/QUOTE]

I would like to see it, actually. There are plenty of such submissions to the OEIS -- short programs in the %o (or %p or %t, of course) lines, long programs as links to source. But since you have a page set up, why not put it there?

 CRGreathouse 2010-04-10 03:49

[QUOTE=alexhiggins732;211252]However, for all primes below 2^32 the file size is 178,956,968 bytes (170 MB) and it doubles with each power of two which led me to ask how how you could sieve up to 2^44 which would require a huge amount of memory even when saved in the most compact possible form (1 bit for each 6K+=1)[/QUOTE]

Why save the primes? Generate a batch, do what you need to do, then discard them. The only primes you need to keep from iteration to iteration are the sieving primes. If you want to test up to 100 trillion, that's just the first 664,579 (2.5 MB using 32-bit ints). A bit array of a billion integers takes a billion bits, or 119 MB.

 bsquared 2010-04-10 04:35

[quote=CRGreathouse;211255]Why save the primes? Generate a batch, do what you need to do, then discard them. The only primes you need to keep from iteration to iteration are the sieving primes. If you want to test up to 100 trillion, that's just the first 664,579 (2.5 MB using 32-bit ints). A bit array of a billion integers takes a billion bits, or 119 MB.[/quote]

Exactly. Although I don't even keep the sieving primes because that's overhead I don't want to manage and its only a hundredth of a second to compute them. Using -v -v prints some more info in yafu-1.18. For example, here is a yafu command to count the primes in a range from 1e13 to 1.0001e13:

[code]yafu "primes(10000000000000,10001000000000)" -v -v
[/code]

which yeilds this output:
[code]
elapsed time for seed primes = 0.0144
sieving range 9999999999840 to 10001006632800
using 227655 primes, max prime = 3162436
using 8 residue classes
lines have 4194304 bytes and 33554432 flags
lines broken into 128 blocks of size 32768
blocks contain 262144 flags and cover 7864320 primes
bucket sieving 194320 primes > 393216
allocating space for 62341 hits per bucket
using 45026396 bytes for sieving storage
elapsed time = 1.0610
ans = 33405006
[/code]

Counting this range only takes a second. Even with bucket sieving, all told it only uses 45 MB of RAM, because the sieve doesn't store a bit for every number, just the ones not divisible by 2, 3, or 5 (i.e. I line sieve the residue classes mod 30).

Once I know the locations of the primes in the bit arrays, I compute them, but only the 33 odd million of them, of course, not a billion. Stored as 64 bit integers, this takes about 250 MB of RAM and another 2 seconds or so.

 bsquared 2010-04-10 05:33

1 Attachment(s)
[quote=CRGreathouse;211254]
I would like to see it, actually. There are plenty of such submissions to the OEIS -- short programs in the %o (or %p or %t, of course) lines, long programs as links to source. But since you have a page set up, why not put it there?[/quote]

Thanks for your interest! I did have fun with the code for this puzzle, so it's neat to be able to share it. I haven't put it on my site yet because it's not useable as a standalone program due to the calls to internal yafu functions. I don't know if its worth the effort to repackage it as a standalone program, but it will be available within yafu once I release a new version.

To save on the forum database, its attached rather than posted, but keep in mind that this won't compile as is.

- ben.

 bsquared 2010-04-10 05:40

The routine to do prime cubes is very similar, except for the assembly to do the cube/sum:

[code]
[SIZE=2][COLOR=#808080][SIZE=2][COLOR=#808080]ASM_G ([/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"movq %1, %%rcx \n\t" /* store prime */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"mulq %%rcx \n\t" /* square it */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"movq %%rax, %%r8 \n\t" /* save p^2 lo (a) */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"movq %%rdx, %%r9 \n\t" /* save p^2 hi (d) */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"mulq %%rcx \n\t" /* p * a */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"movq %%rax, %%r10 \n\t" /* save p*a lo (apa) */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"movq %%rdx, %%r11 \n\t" /* save p*a hi (apd) */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"movq %%r9, %%rax \n\t" /* p * d */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"mulq %%rcx \n\t" /* lo part in rax (dpa), hi in rdx (dpd) */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"addq %%r10, (%%rbx) \n\t" /* sum0 = sum0 + apa */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"adcq %%rax, 8(%%rbx) \n\t" /* sum1 = sum1 + dpa + carry */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"adcq %%rdx, 16(%%rbx) \n\t" /* sum2 = sum2 + dpd + carry */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"addq %%r11, 8(%%rbx) \n\t" /* sum1 = sum1 + apd */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]"adcq \$0, 16(%%rbx) \n\t" /* sum2 = sum2 + carry */[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]: [/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]: "b"(sum->val), "a"(PRIMES[j])[/COLOR][/SIZE]
[SIZE=2][COLOR=#808080]: "rcx", "rdx", "r8", "r9", "r10", "r11", "memory", "cc"); [/COLOR][/SIZE]
[/COLOR][/SIZE][/code]

 alexhiggins732 2010-04-10 13:49

1 Attachment(s)
[quote=bsquared;211249]YAFU

Just download YAFU-1.18 and browse around through the soe.c file. Although it's not well documented :blush:.

I compute the primes in batches of 1 billion integers. At a height of 10 trillion, that is 32 million primes or so per batch. It's then pretty straightfoward to sum them and check each for divisiblity by a power of 10 (sweeping many optimizations under the rug, for clarity). On my (admittedly fast) computer, each batch takes about 3 seconds for prime computation, and .2 seconds for the sum/check. This is the only code that's not public yet. Although yafu-1.19 will probably retain it when it comes out.[/quote]

Downloaded the source and fired up the Visual Studio Debugger after a bit of .... and still a bit hard to see exactly what is going on because pointers to the arrays don't allow you to see the contents of the array on the memory address. While c/c++ is second language I can't say I fully understand pointer arithmetic entirely, but I am sure if I could figure it out.

So being the crank I am, I did what any self-respecting crank would do, I rewrote it in VB.NET, but have a few problems perhaps you can shine some light on.

I tried, for the most part, to use the exact vb equivalent of the c code
1) The multi-threading is not functional, yet. I understand the logic that you are implementing so I just may rewrite that and deviate from adhering strictly to the original code. I couldn't really figure out what thread_soedata_t->thread_id was being used for. I am not sure what the vb equivalent of CreateThread would be. In any case, I can rewrite the threading. Each thread is assigned a command, waits until it is signal to run then notifies the master thread when done, rinse lather repeat.

2) When you declare a local variable of a structure member a value copy is created on the heap and the original structure is not updated. If this where declared a class, it would work like the c implementation but I didn't what to do that. So I had to rewrite as follows:

[CODE]
Public Sub primes_from_lineflags(ByRef t As thread_soedata_t)
' this is no good. it create value copy on the heap so the original
' thread_soedata_t members do not get updated.
'Dim ddata As soe_dynamicdata_t = t.ddata
'Dim sdata As soe_staticdata_t = t.sdata
' work around is to pass a pointer to the structure and access the
' the members using the pointer like on the next line.
Dim line() As Byte = t.ddata.line
[/CODE]

Problem is, I think that there is somewhere else in the code the value copy of the structure is being accessed and I can't figure out where. Stepping through the yafu source code and the vb source code everything matches up except, the number of primes that are returned. Only thing I can think of is due to pointer arithmetic some values the pointer arrays are different than the vb version.

However, not being able to see that values are in those arrays in the debugger and being a crank that is trying to learn what is going on from the code I don't know if the values in the vb arrays are correct.

Maybe my misinterpretation of how pointers work is the problem.

For example, interpreted flagblock += BLOCKSIZE; as increment the index of the flagblock + blocksize. Since this is illegal in vb, to accommodate I created a new variable fob and instead of flagblock + BLOCKSIZE, its fob+BLOCKSIZE and then when accessing flagblock its flagblock[x+fob]. Is this wrong?

Perhaps the simplest solution is to attach the code I have and hopefully someone with MSVS can debug it and give me a hint whats wrong. Its just a little frustrating after rewriting 2000 lines of code....

Oh well, maybe I'll write it in C# so I can use the pointer arithmetic and just import the DLL into my vb projects :smile:

 bsquared 2010-04-10 15:59

[quote=alexhiggins732;211298]So being the crank I am, I did what any self-respecting crank would do, I rewrote it in VB.NET...
[/quote]

That's a pretty ambitious spur-of-the-moment undertaking. I have no idea if its possible (I don't know much VB.NET so I won't be of any use in helping you with your port).

As you've seen, pointers are used abundantly and I employ all sorts of C hacks like casting void pointers and adding pointer values. Not sure if equivalent code is possible in VB.

It may have been easier to package the existing C code as a .dll library for use in a VB program.

 alexhiggins732 2010-04-10 16:25

[quote=bsquared;211309]That's a pretty ambitious spur-of-the-moment undertaking. I have no idea if its possible (I don't know much VB.NET so I won't be of any use in helping you with your port).

As you've seen, pointers are used abundantly and I employ all sorts of C hacks like casting void pointers and adding pointer values. Not sure if equivalent code is possible in VB.

It may have been easier to package the existing C code as a .dll library for use in a VB program.[/quote]

Yeah tell me about it. That is nothing like my "sliding window" sieve or even the seo. In fact Atkins and Euler have sieve named after them and and made only simple modifications.

[B]YOUR SIEVE IS A BEAST:smile: I NOMINATE THE NAME THE SIEVE OF BUHROW[/B]

COM Interop on windows is quite a drag on performance. I current use a LIBGMP wrapper and just for i=0 to 100000000 takes several minutes.

Actually comparing c/c++/vb performance is quite comparable when dealing with 32 bit words. 16 and 64 may be slightly slower and perhaps the bounds checking performed under the hood would be another performance hit, but nothing compared to the cost of com calls.

Hopefully someone with VS will take a look and give it a shot. I can upload and entire project file. It would simply be a matter of stepping though each program. I mean their almost exact line for line.

FYI, tiny_soe works perfect but then again that is a pretty simple implementation. spSOE fails, first reporting 121 as prime.

 CRGreathouse 2010-04-10 17:19

[QUOTE=alexhiggins732;211311]Yeah tell me about it. That is nothing like my "sliding window" sieve or even the seo. In fact Atkins and Euler have sieve named after them and and made only simple modifications. [/QUOTE]

You think the Sieve of Atkin is a simple modification?

 bsquared 2010-04-10 18:34

[quote=alexhiggins732;211311]Yeah tell me about it. That is nothing like my "sliding window" sieve or even the seo. In fact Atkins and Euler have sieve named after them and and made only simple modifications.

[B]YOUR SIEVE IS A BEAST:smile: I NOMINATE THE NAME THE SIEVE OF BUHROW[/B]

[/quote]

I'm flattered :smile:, but at the end of the day its just another sieve of Eratosthenes. The best that can be said about it is that it uses some optimizations that perhaps haven't been implemented before. Even so, in all but very special cases it isn't even the fastest SOE... ecprime fills those shoes.

 alexhiggins732 2010-04-10 19:43

[quote=CRGreathouse;211315]You think the Sieve of Atkin is a simple modification?[/quote]

While it may be entirely more effecient for larger input and the theory behind it can be quite complex, it is indeed a relatively simple modification compared to the one implement in Yafu.

[code]
Private Shared Sub ParralleSieve(x as Integer)

Dim xx = x * x
For y As Integer = 1 To sqrt
Dim yy = y * y
Dim n = 4 * xx + yy
If n <= max AndAlso (n Mod 12 = 1 OrElse n Mod 12 = 5) Then
isPrime(n) = Not isPrime(n)
End If

n = 3 * xx + yy
If n <= max AndAlso n Mod 12 = 7 Then
isPrime(n) = Not isPrime(n)
End If

n = 3 * xx - yy
If x > y AndAlso n <= max AndAlso n Mod 12 = 11 Then
isPrime(n) = Not isPrime(n)
End If
Next

End Sub
Private Shared Function SieveOfAtkins(ByVal max As Integer) As List(Of Integer)
Dim isPrime = New BitArray(CInt(max) + 1, False)
Dim sqrt = CInt(Math.Sqrt(max))

Dim primes = New List(Of Integer)()
For n As Integer = 5 To sqrt
If isPrime(n) Then
Dim nn As Integer = n * n
Dim k As Integer = nn
While k <= max
isPrime(k) = False
k += nn
End While
End If
Next

For n As Integer = sqrt + 1 To max
If isPrime(n) Then
End If
Next

Return primes
End Function

[/code]

 CRGreathouse 2010-04-10 20:19

[QUOTE=alexhiggins732;211323]While it may be entirely more effecient for larger input and the theory behind it can be quite complex, it is indeed a relatively simple modification compared to the one implement in Yafu.[/QUOTE]

[code]
Private Shared Sub ParralleSieve(x as Integer)

Dim xx = x * x
For y As Integer = 1 To sqrt
Dim yy = y * y
Dim n = 4 * xx + yy
If n <= max AndAlso (n Mod 12 = 1 OrElse n Mod 12 = 5) Then
isPrime(n) = Not isPrime(n)
End If

n = 3 * xx + yy
If n <= max AndAlso n Mod 12 = 7 Then
isPrime(n) = Not isPrime(n)
End If

n = 3 * xx - yy
If x > y AndAlso n <= max AndAlso n Mod 12 = 11 Then
isPrime(n) = Not isPrime(n)
End If
Next

End Sub
Private Shared Function SieveOfAtkins(ByVal max As Integer) As List(Of Integer)
Dim isPrime = New BitArray(CInt(max) + 1, False)
Dim sqrt = CInt(Math.Sqrt(max))

Dim primes = New List(Of Integer)()
For n As Integer = 5 To sqrt
If isPrime(n) Then
Dim nn As Integer = n * n
Dim k As Integer = nn
While k <= max
isPrime(k) = False
k += nn
End While
End If
Next

For n As Integer = sqrt + 1 To max
If isPrime(n) Then
End If
Next

Return primes
End Function

[/code][/QUOTE]

As far as I can tell, this is wrong. But even if it works (I can't get it to give correct answers), it's not more than a pale reflection of the sieve. This takes at least [TEX]n(\log n)^{1+\varepsilon}[/TEX] time, while the 'true' Sieve of Atkin takes only [TEX]O(n/\log\log n)[/TEX] time. Similarly, yours takes about [TEX]n[/TEX] bits of memory, while the real one takes [TEX]n^{0.5+\varepsilon}[/TEX] bits.

The second is not difficult, but the first one takes real magic. By omitting that part you do a great disservice to the inventors of the sieve, Bernstein and Atkin.

But if you think that implementation tricks are what make the algorithm, then you should see [url=http://cr.yp.to/primegen.html]Bernstein's implementation[/url].

Also: the name is Atkin, not Atkins.

 alexhiggins732 2010-04-10 21:07

I concur, it is indeed a great sieve indeed and there are many implementations of it out there. For those interested here is the theory behind it and other sieves.

[url]http://cr.yp.to/papers/primesieves.pdf[/url]

 flouran 2010-04-10 21:38

[QUOTE=CRGreathouse;211325]
Also: the name is Atkin, not Atkins.[/QUOTE]
I guess Atkin is often mistaken for the diet :smile:

 bsquared 2010-04-12 00:55

I've now completed the sum of each sequence (primes, prime squares, and prime cubes) up to 1e14, and found the 12th member of the cube sequence:

[CODE]**** 1000000000000 divides prime cube sum up to 33777547344991, sum is 10532010356822624092227361649102207021134000000000000 ****[/CODE]

[quote=flouran;211328]I guess Atkin is often mistaken for the diet :smile:[/quote]Well, the diet is a sieve, screening out certain types of calories, isn't it? :smile:

 davieddy 2010-04-12 16:38

I thought a calory was 4.2 Joule, but nutritionists use
the term to refer to a kilocalory.

I don't think there are different types for different foods.
Even the "type of energy" it measures is the same for
everything.

David

 alexhiggins732 2010-04-12 17:17

[quote=flouran;211328]I guess Atkin is often mistaken for the diet :smile:[/quote]

[url]http://wwwhomes.uni-bielefeld.de/achim/prime_sieve.html[/url] has it labeled Sieve of Atkins and points the the PDF I referenced above.

Being a crank and all, its practically a given that I'll have many more misnomers than this :smile:.

 Mini-Geek 2010-04-12 17:29

[quote=davieddy;211501]I thought a calory was 4.2 Joule, but nutritionists use
the term to refer to a kilocalory.

I don't think there are different types for different foods.
Even the "type of energy" it measures is the same for
everything.

David[/quote]
In America, the most commonly quoted rule is:
1 Calorie (with a capital C) = 1 kilocalorie = 1000 calories (with a lower-case c) = 4.2 kJ

 davieddy 2010-04-12 17:59

Yes - I deliberated over the capitalization and singular
of calories. I don't think things were thoroughly standardized
pre SI.
I was weaned on the foot poundel, and considered the erg a
newfangled Napoleonic invention:smile:

David

PS is there any concensus about kB, KB, kb or Kb for killerbight?

 xilman 2010-04-12 18:15

[quote=Mini-Geek;211512]In America, the most commonly quoted rule is:
1 Calorie (with a capital C) = 1 kilocalorie = 1000 calories (with a lower-case c) = 4.2 kJ[/quote]True in the UK too.

Though I'd phrase the relationship as 1000 calories is approximately 4.2kJ :smile:

Paul

Atkins sieve diet v1.1.2 bug-fix release

[quote=cheesehead;211487]Well, the diet is a sieve, screening out [/quote]the calories borne by[quote]certain types of[/quote]food [quote], isn't it? :smile:[/quote]

 davieddy 2010-04-12 19:11

Physical Units

I feel a new thread coming on, initiated by this latest
(hijacking) exchange.

David

Science and Technology would be the appropriate
location. Paul?

David

 S485122 2010-04-12 20:25

[QUOTE=davieddy;211516]PS is there any concensus about kB, KB, kb or Kb for killerbight?[/QUOTE]Yes there is a standard. First of all the "k" should be lower case for all units, see the site of the BIPM that defines the SI system, more precisely the page [url=http://www.bipm.org/en/si/prefixes.html]SI prefixes[/url]... Except for 2[sup]10[/sup] (1024) bytes where it is KiB with a capital "K" an "i" to mark the power of two "kilo" of 1024 bytes. 10[sup]3[/sup] (1000) bytes is a kB. See the interesting discussion at [url=http://physics.nist.gov/cuu/Units/binary.html]Prefixes for binary multiples[/URL].

Then bi signifies bits and B signifies bytes.

So depending on the quantity you speak of : it should be 1 KiB for 1024 bytes and 1 kB for 1000 bytes.

1 kb represents 1000 bits.

KB and Kb are no units that I know of.

Jacob

 davieddy 2010-04-12 22:01

THX Jacob.
I will either defer to my previous communique,
or get William Blipp on the job:smile:

David

 only_human 2010-04-15 01:19

[quote=davar55;210978]Another way to extend this problem is to use
a base other than ten. I think binary.

2^2 + 3^2 + 5^2 + ... + p^2 = 2[sup]m[/sup]K

What is the smallest prime p such that
the sum of squares of all primes up to p
is a multiple of 2 (or 4 or 8 or 16 or ...).

This question can also be asked of first powers
or cubes of primes.

Since we basically compute in decimal or binary,
if there's an interesting number theoretic fact here
we may find it in one of these two related
sequences of sequences.[/quote]

Back to this part, I wanted to make a comment or two about sums of odd squares in binary that is too concrete to be nice mathematics but I have mused over recently. Speaking of odd primes (or any odd numbers actually), the squares are all 1 (mod 8), so a sum of these squares will include all these 1's added together as part of the result, therefore:

With regards to the sum of N odd squares (or in your case, squares of odd primes),
2[sup]m[/sup]K == N (mod 8)

I recently noticed that any odd square is one more than 8 times a triangular number:
n[sup]2[/sup] = 8t + 1

Trivial facts:
The nth triangular number is n(n+1)/2
Squares are the sum of N consecutive odd integers
Triangular Numbers are the sum of N consecutive integers
Two consecutive triangular numbers added together is equal to a square

Some aspect of these trivial facts might come in handy; calling a particular triangular number t[sub]i[/sub], and the next triangular number in sequence t[sub]i+1[/sub], we can look at twin primes and consider the sum of their squares:
p[sup]2[/sup] + (p+2)[sup]2[/sup] = 8t[sub]i[/sub] + 1 + 8t[sub]i+1[/sub] + 1
Since t[sub]i[/sub] + t[sub]i+1[/sub] is itself a square that means that the sum of the squares of a the pair twinned primes is two more than a eight times a square. I don't know a use for this except that it is elementary fact accessible to me and it applies to all consecutive odd numbers but of primes, only to twin primes because they are the only consecutive odd numbers that are also prime.

 science_man_88 2010-04-20 14:05

[QUOTE=Mini-Geek;211512]In America, the most commonly quoted rule is:
1 Calorie (with a capital C) = 1 kilocalorie = 1000 calories (with a lower-case c) = 4.2 kJ[/QUOTE]
actually more exact it's 4.184 I think

as to the code if i know what to help with I do have a ASM forum I have gone to maybe they could turn it to straight asm ( apparently ASM is now changed with different OS but if we can get it to a form that will work on most OS'es would that be worthy ?)

 science_man_88 2010-04-20 14:29

I know some assembly but not really enough to help plus to assemble anything I'll have to use my links to masm32 download, and possibly [url]http://win32assembly.online.fr/tut2.html[/url]

 bsquared 2010-04-20 15:53

[QUOTE=science_man_88;212580]

as to the code if i know what to help with I do have a ASM forum I have gone to maybe they could turn it to straight asm ( apparently ASM is now changed with different OS but if we can get it to a form that will work on most OS'es would that be worthy ?)[/QUOTE]

It sounds like you are offering to "help" by asking other people to do some work for you. Also it seems you are somewhat missing the point of this thread; it is a puzzle, not a DC project, and doesn't need any help per se. You are free to do what you like to try to solve the puzzle or one of the variants. If you decide to work on the sum of primes/squares/cubes mod 10^n, then I'd be grateful for any doublechecking you want to do, or you can extend the sequences using the .csv files from my website.

 science_man_88 2010-04-20 16:48

a=0;b=0;c=0;forprime(x=1,200,a=a+x;b=b+x^2;c=c+x^3;print(a","b","c)) is the best pari code I can see except forprime would need a bigger addprimes() to figure this out.

 science_man_88 2010-04-20 16:50

a=0;b=0;c=0;for(x=1,200,if(isprime(x),a=a+x;b=b+x^2;c=c+x^3;print(a","b","c))) is probably better actually can you turn this into assembly if you want an exe ?

 science_man_88 2010-04-20 16:54

to create the exponents without multiplying using counter variables and addition might speed it up.

 bsquared 2010-04-20 16:56

[QUOTE=science_man_88;212596]a=0;b=0;c=0;forprime(x=1,200,a=a+x;b=b+x^2;c=c+x^3;print(a","b","c)) is the best pari code I can see except forprime would need a bigger addprimes() to figure this out.[/QUOTE]

Well, you forgot to check for multiples of 10^n (or a different variant). Also, other people in this thread have already posted PARI code. Its a start though. Now just increase your upper bound by 12 orders of magnitude or so and you'll be up to the state of the art.

 bsquared 2010-04-20 17:00

[QUOTE=science_man_88;212599]to create the exponents without multiplying using counter variables and addition might speed it up.[/QUOTE]

Try this for a prime near 10e12 and say that again.

 science_man_88 2010-04-20 17:02

well log(x)/log(10) could fix the check but[B] I'm not smart enough[/B] to implement them in assembly I've tried to learn it before and to do that I'll probably need to redownload masm.

 science_man_88 2010-04-20 17:05

bsqaured couldn't you use data to have constant powers of 10 then negate to check if it's a multiple of a power of 10 ?

 bsquared 2010-04-20 17:12

[QUOTE=science_man_88;212604]bsqaured couldn't you use data to have constant powers of 10 then negate to check if it's a multiple of a power of 10 ?[/QUOTE]

I'm not sure I understand what you mean. If you're trying to ask if there are shortcuts to doing a full check (for a power of 10) every iteration then yes, there are. The code I posted earlier does a quick check using a single bit operation most of the time, and a full check using a single division occasionally, supported by a few variables.

 science_man_88 2010-04-20 17:18

for a given a could you do this type of thing for log at least:

a=10;b=10;for(counter=0,10,if(a>0,a=a-b));print(counter);

 bsquared 2010-04-20 17:21

[QUOTE=science_man_88;212603]well log(x)/log(10) could fix the check but[B] I'm not smart enough[/B] to implement them in assembly I've tried to learn it before and to do that I'll probably need to redownload masm.[/QUOTE]

Assembly isn't the final answer for speed. Given the quality of optimizing compilers and complex, pipelined, out-of-order cpus nowadays you are almost universally better off using a higher level language to get something done, and done fast enough.

 science_man_88 2010-04-20 17:24

Well I got into programming in vb.net then C then assembly ( some of each). Too bad I haven't ever figured half the stuff out in vb.net ( I'd need MSVS to do any) best I've ever done in vb.net is the basics of a OS deleting virus I never really compiled.

 bsquared 2010-04-20 17:28

[QUOTE=science_man_88;212609] best I've ever done in vb.net is the basics of a OS deleting virus I never really compiled.[/QUOTE]

You're admitting in a public forum that you write viruses? Ballsy move.

 science_man_88 2010-04-20 17:36

NO! I don't write them in fact I haven't used vb.net in a while but I learned enough basics to make a code that could act as one. An array here a for loop there a few commands to do this and that a call to shutdown only to get rid of it and outside backups by the system or something similar take the OS down.

 bsquared 2010-04-20 17:40

Well I can't say as I'm any less confused, but if I misinterpreted I apologize.

Maybe try to do a better job saying what you mean...

 science_man_88 2010-04-20 17:48

best I've ever done in vb.net is the basics of a OS deleting virus [B][I][U]I never really compiled.[/U][/I][/B]
You're admitting in a public forum that you write viruses? Ballsy move.

 science_man_88 2010-04-20 20:20

If you wanted an exe couldn't you compile pari code with the pari compiler ?

 science_man_88 2010-04-20 23:26

"Also available is gp2c, the GP-to-C compiler, which compiles GP scripts into the C language and transparently loads the resulting functions into gp. The advantage of this is that gp2c-compiled scripts will typically run three to four times faster."

not that Wikipedia is reliable but I'm game to try this if you all are

 science_man_88 2010-04-20 23:28

then we just compile to an exe is that's the route you'd like to take.The first PARI code listed is more efficient than mine we could use that if it 's optimised.

 science_man_88 2010-04-20 23:34

or we could compile get the optimised c code turn it to assembly and assemble

All times are UTC. The time now is 05:19.