PARI's commands - Page 204

CRGreathouse · 2011-05-24, 18:49

Faster version:

Code:

pick(n,k)=binomial(if(n<k,n+k-1,n),k)

Basically, your function is binomial(n+k-1, k) if n<k and binomial(n,k) otherwise. I could have written

Code:

pick(n,k)=if(n<k,binomial(n+k-1,k), binomial(n,k))

but I chose to combine the similar parts and use just

Code:

pick(n,k)=binomial(if(n<k,n+k-1,n),k)

When PARI can compute the binomial directly instead of indirectly through the factorial it can save a lot of steps. When you calculate ${9\choose2}=\frac{9!}{7!2!}=\frac{9\cdot8\cdot7\cdot6\cdot5\cdot4\cdot3\cdot2\cdot1}{(7\cdot6\cdot5\cdot4\cdot3\cdot2\cdot1)(2\cdot1)}$
you probably just cross out the 7 * ... * 1 on top and bottom. When you give PARI binomial(9,2) it can do the same; when you give it 9!/(7!*2!) it's forced to compute the large numbers and then divide.

Xitami · 2011-05-24, 20:42

Code:

binom(n,k)={ my(t=i=1);
	if(2*k>n, k=n-k);
	while(i<=k,
		t=(t*n)/i;
		i++; n--
	);
	t
}

-----------------------------

science_man_88 · 2011-05-24, 22:17

Quote:

Originally Posted by CRGreathouse

Faster version:

Code:

pick(n,k)=binomial(if(n<k,n+k-1,n),k)

Basically, your function is binomial(n+k-1, k) if n<k and binomial(n,k) otherwise. I could have written

Code:

pick(n,k)=if(n<k,binomial(n+k-1,k), binomial(n,k))

but I chose to combine the similar parts and use just

Code:

pick(n,k)=binomial(if(n<k,n+k-1,n),k)

When PARI can compute the binomial directly instead of indirectly through the factorial it can save a lot of steps. When you calculate ${9\choose2}=\frac{9!}{7!2!}=\frac{9\cdot8\cdot7\cdot6\cdot5\cdot4\cdot3\cdot2\cdot1}{(7\cdot6\cdot5\cdot4\cdot3\cdot2\cdot1)(2\cdot1)}$
you probably just cross out the 7 * ... * 1 on top and bottom. When you give PARI binomial(9,2) it can do the same; when you give it 9!/(7!*2!) it's forced to compute the large numbers and then divide.

thanks:

Code:

(19:14)>allocatemem(1)
  *** allocatemem: Warning: new stack size = 1024 (0.001 Mbytes).
(19:15)>pick(11,1000000000)
%269 = 275573207396384843474431809413623724923088088349707046408843005955895496034675000001
(19:15)>allocatemem(931987423)
  *** allocatemem: Warning: new stack size = 931987420 (888.812 Mbytes).
(19:15)>pick1(11,1000000000)
  *** _!: the PARI stack overflows !
  current stack size: 931987420 (888.812 Mbytes)
  [hint] you can increase GP stack with allocatemem()

I changed my original to pick1 for this test.

science_man_88 · 2011-05-24, 23:27

found something I find interesting now CRG:

Code:

(20:20)>for(a=1,38,print1(sum(x=1,a,pick(Me(x)%6,Me(x+1)%6)%4)%4","))
0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,

okay it works out this way for all primes it looks like.

science_man_88 · 2011-05-25, 00:48

Quote:

Originally Posted by science_man_88

found something I find interesting now CRG:

Code:

(20:20)>for(a=1,38,print1(sum(x=1,a,pick(Me(x)%6,Me(x+1)%6)%4)%4","))
0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,

okay it works out this way for all primes it looks like.

I made an extended version of MeVec to test this and so far the exponents for x=40 to 47( or guessing) fit this pattern.

CRGreathouse · 2011-05-25, 05:07

What happens if you reverse the order of the last 8?

science_man_88 · 2011-05-25, 11:37

Quote:

Originally Posted by CRGreathouse

What happens if you reverse the order of the last 8?

nothing I bet because late last night I realized why the pattern is set up it's because pick(a,b) with a = 1 or 5, b=1,5 will be 1 mod 4. so the sum is mostly a sum of ones ( which will only cycle mod anything) the 2 that start it 2 and 3 and 3 and 5 are 0 and 1 mod 4 so 0+1 is the start and +1 for each step after. so really all this states is that the first 3 exponents are used to set it up.

science_man_88 · 2011-05-25, 22:57

Quote:

Originally Posted by science_man_88

found something I find interesting now CRG:

Code:

(20:20)>for(a=1,38,print1(sum(x=1,a,pick(Me(x)%6,Me(x+1)%6)%4)%4","))
0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,

okay it works out this way for all primes it looks like.

I've got another pattern to figure out before revealing.

science_man_88 · 2011-05-25, 23:11

Quote:

Originally Posted by science_man_88

I've got another pattern to figure out before revealing.

okay I gave in:

Code:

(19:34)>for(x=1,47,print1(pick(prime(x),x^2-10)%x","))
0,0,0,3,0,0,0,0,0,4,0,6,0,0,0,0,0,0,0,0,0,0,0,0,0,13,0,20,0,0,0,0,0,0,30,0,0,0,0,0,0,0,0,0,0,0,0,
(19:34)>for(x=1,47,print1(pick(prime(x),x^2-5)%x","))
0,0,2,0,0,2,0,0,0,4,0,9,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,16,0,0,25,0,0,0,13,0,0,28,0,0,0,0,0,
(19:35)>for(x=1,47,print1(pick(MPE[x],x^2-10)%x","))
0,0,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,15,0,0,0,0,0,0,0,0,0,0,
(19:35)>for(x=1,47,print1(pick(MPE[x],x^2-5)%x","))
0,0,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
(19:35)>for(x=1,47,print1(MPE[x]<x^2-5","))
0,0,0,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
(19:47)>for(x=1,47,print1(MPE[x]<x^2-10","))
0,0,0,0,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,

Quote:

MPE[x]!/(k!*(MPE[x]-k)!)
(MPE[x]+k-1)!/((k!)*(MPE[x]-1)!)

is what I consider important

science_man_88 · 2011-05-26, 17:45

Quote:

Originally Posted by science_man_88

okay I gave in:

Code:

(19:34)>for(x=1,47,print1(pick(prime(x),x^2-10)%x","))
0,0,0,3,0,0,0,0,0,4,0,6,0,0,0,0,0,0,0,0,0,0,0,0,0,13,0,20,0,0,0,0,0,0,30,0,0,0,0,0,0,0,0,0,0,0,0,
(19:34)>for(x=1,47,print1(pick(prime(x),x^2-5)%x","))
0,0,2,0,0,2,0,0,0,4,0,9,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,16,0,0,25,0,0,0,13,0,0,28,0,0,0,0,0,
(19:35)>for(x=1,47,print1(pick(MPE[x],x^2-10)%x","))
0,0,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,15,0,0,0,0,0,0,0,0,0,0,
(19:35)>for(x=1,47,print1(pick(MPE[x],x^2-5)%x","))
0,0,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
(19:35)>for(x=1,47,print1(MPE[x]<x^2-5","))
0,0,0,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,
(19:47)>for(x=1,47,print1(MPE[x]<x^2-10","))
0,0,0,0,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,

is what I consider important

the red is for proving not all primes follow the patterns 1 after the other.

the green shows the patterns I think I have ( pick(Me(x),x^2-5) is 0 mod 2 and pick(Me(x),x^2-10) is 0 mod 3).

the orange is to show which path of pick each Me(x) brings us through.

science_man_88 · 2011-06-30, 15:11

not sure how to create this , I know I ignore 2,3,and 5, and I check mod 60 but I have no idea ( the pseudo-code in the Wikipedia article didn't help me as far as I can tell) otherwise.

2011-05-24, 18:49	#2234
CRGreathouse Aug 2006 3·1,993 Posts	Faster version: Code: pick(n,k)=binomial(if(n<k,n+k-1,n),k) Basically, your function is binomial(n+k-1, k) if n<k and binomial(n,k) otherwise. I could have written Code: pick(n,k)=if(n<k,binomial(n+k-1,k), binomial(n,k)) but I chose to combine the similar parts and use just Code: pick(n,k)=binomial(if(n<k,n+k-1,n),k) When PARI can compute the binomial directly instead of indirectly through the factorial it can save a lot of steps. When you calculate ${9\choose2}=\frac{9!}{7!2!}=\frac{9\cdot8\cdot7\cdot6\cdot5\cdot4\cdot3\cdot2\cdot1}{(7\cdot6\cdot5\cdot4\cdot3\cdot2\cdot1)(2\cdot1)}$ you probably just cross out the 7 * ... * 1 on top and bottom. When you give PARI binomial(9,2) it can do the same; when you give it 9!/(7!*2!) it's forced to compute the large numbers and then divide.

2011-05-24, 20:42	#2235
Xitami Apr 2010 2×7 Posts	Code: binom(n,k)={ my(t=i=1); if(2k>n, k=n-k); while(i<=k, t=(tn)/i; i++; n-- ); t } -----------------------------

2011-05-24, 23:27	#2237
science_man_88 "Forget I exist" Jul 2009 Dumbassville 2⁶·131 Posts	found something I find interesting now CRG: Code: (20:20)>for(a=1,38,print1(sum(x=1,a,pick(Me(x)%6,Me(x+1)%6)%4)%4",")) 0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1,2,3,0,1, okay it works out this way for all primes it looks like. Last fiddled with by science_man_88 on 2011-05-24 at 23:29

2011-06-30, 15:11	#2244
science_man_88 "Forget I exist" Jul 2009 Dumbassville 8384₁₀ Posts	sieve of atkins not sure how to create this , I know I ignore 2,3,and 5, and I check mod 60 but I have no idea ( the pseudo-code in the Wikipedia article didn't help me as far as I can tell) otherwise.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Why do I sometimes see all the <> formatting commands when I quote or edit?	cheesehead	Forum Feedback	3	2013-05-25 12:56
Passing commands to PARI on Windows	James Heinrich	Software	2	2012-05-13 19:19
Ubiquity commands	Mini-Geek	Aliquot Sequences	1	2009-09-22 19:33
64-bit Pari?	CRGreathouse	Software	2	2009-03-13 04:22
Are these commands correct?	jasong	Linux	2	2007-10-18 23:40

2011-05-25, 05:07	#2239
CRGreathouse Aug 2006 3·1,993 Posts	What happens if you reverse the order of the last 8?