PARI's commands

[science_man_88]

Originally Posted by vittoire

hi, im new to pari/gp and i need some guidance on how to go about storing my factor value into a variable. i have for example 192384 and i need to get 1 of the prime factor of this value.

i tried storing it to a variable "a=factor(192384)" but it stores all the prime value to it. can anyone guide me how i could get 1 of the prime factor please.

where is it in the factorization matrix originally ? what column and what row ?

[CRGreathouse]

You could do factor(192384)[1,1] to get the smallest prime factor, factor(192384)[2,1] to get the second (if it exists), and so forth.

[science_man_88]

I just want people to know that if some of the codes I post elsewhere don't seem correct I now have both 2.5.3 and 2.7.0 and probably another. I love some of these new add-on's to existing functions etc.

[science_man_88]

yes I made another try at aligen that I would hope is at least more inclusive than previous versions:

aligentry(n)={c=[];
forpart(x=n-1,
          if(#x==#vecsort(x,,8) && numdiv(x[1]*x[#x])-2==#x && setminus(divisors(x[1]*x[#x]),Vec(x))==[1,x[1]*x[#x]],
             c=concat(c,x[1]*x[#x])
            )
          ,a=[2,n-1])
      ;c
}

[CRGreathouse]

science_man, I don't understand what that function is trying to do. But looking back through this thread using the search term "aligen" I found this exchange, and maybe a better answer to it will help you:

Originally Posted by science_man_88

is there an easy way to create all the distinct partitions of a number n without the number 1 involved ( though I also have a way with 1 involved I think) ? because I'm trying to speed up my aligen function and it's based on partitions.

Originally Posted by CRGreathouse

Here's some naive code:

[code]partNo1(n)={
select(v->vecmin(v)>1,partitions(n))
};

With newer versions of gp, there's now a faster way to do this:

partNo1(n)={
  partitions(n, [2, n])
};

If you wanted those partitions with only distinct parts, you can do

distinctPartNo1(n)={
  select(v->#v==#Set(v), partitions(n, [2, n]))
};

or better yet

distinctPartNo1(n)={
  my(d=sqrtint(2*n+2));
  if(d*(d+1)/2 > 2*n+2, d--);
  select(v->#v==#Set(v), partitions(n, [2, n], d))
};

which is faster and requires less memory.

[science_man_88]

Originally Posted by CRGreathouse

science_man, I don't understand what that function is trying to do. But looking back through this thread using the search term "aligen" I found this exchange, and maybe a better answer to it will help you:





With newer versions of gp, there's now a faster way to do this:

partNo1(n)={
  partitions(n, [2, n])
};

If you wanted those partitions with only distinct parts, you can do

distinctPartNo1(n)={
  select(v->#v==#Set(v), partitions(n, [2, n]))
};

or better yet

distinctPartNo1(n)={
  my(d=sqrtint(2*n+2));
  if(d*(d+1)/2 > 2*n+2, d--);
  select(v->#v==#Set(v), partitions(n, [2, n], d))
};

which is faster and requires less memory.

then I'm confused why does forpart exists ? it goes over the partitions of it's first value, it then checks if that partition has all distinct members, then for those with the proper number of elements to be the proper divisors without 1 of it's elements at opposite ends, then and only then does it check it those are the true divisors (because that could take a long time I wanted to cut them down first). the "a=" part is the same as the end of your partitions call.

[CRGreathouse]

Originally Posted by science_man_88

then I'm confused why does forpart exists ? it goes over the partitions of it's first value, it then checks if that partition has all distinct members, then for those with the proper number of elements to be the proper divisors without 1 of it's elements at opposite ends, then and only then does it check it those are the true divisors (because that could take a long time I wanted to cut them down first). the "a=" part is the same as the end of your partitions call.

You can certainly use forpart -- and it should be easy to adapt my code to that use.

[science_man_88]

Originally Posted by CRGreathouse

You can certainly use forpart -- and it should be easy to adapt my code to that use.

adding:

floor(solve(y=1,n-1,sum(z=1,floor(y),z)-(n-1)))

into my code cut the time in half and looks to still be as accurate to the older one. is that the type of fix you were talking about ?

[CRGreathouse]

Something like that, yes. Actually I'm a bit surprised solve() doesn't choke on that, wouldn't it be better to use something like

solve(y=1,n-1, y*(y+1)/2-n+1)\1

to remove the flat parts. I guess it doesn't really matter -- this is outside the hot path so speed is nearly irrelevant.

[science_man_88]

Originally Posted by CRGreathouse

Something like that, yes. Actually I'm a bit surprised solve() doesn't choke on that, wouldn't it be better to use something like

solve(y=1,n-1, y*(y+1)/2-n+1)\1

to remove the flat parts. I guess it doesn't really matter -- this is outside the hot path so speed is nearly irrelevant.

yeah yours becomes faster by 19,715 ms by aligentry(130) I guess a quicker way is to check if every member of x is divisible by the primes in x. if not it can't be a list of divisors. parvector and select may be able to do this eliminating some of the checks.

[CRGreathouse]

Your original code:

aligentry0(n)={c=[];
forpart(x=n-1,
          if(#x==#vecsort(x,,8) && numdiv(x[1]*x[#x])-2==#x && setminus(divisors(x[1]*x[#x]),Vec(x))==[1,x[1]*x[#x]],
             c=concat(c,x[1]*x[#x])
            )
          ,a=[2,n-1])
      ;c
};

Cleaned up with indentation, lexical scoping (always use "my" or "local" for variables), Set, a list instead of concatenation (much faster), etc.:

aligentry1(n)={
	my(c=List());
	forpart(x=n-1,
		if(#x==#Set(x) && numdiv(x[1]*x[#x])-2==#x && setminus(divisors(x[1]*x[#x]),Vec(x))==[1,x[1]*x[#x]],
			listput(c, x[1]*x[#x])
		)
	,
		[2,n-2]
	,
		solve(y=1,n-1, y*(y+1)/2-n+1)\1
	);
	Vec(c)
};

Additional thought: you're looking for a partition x1 > ... > xn equal to the proper divisors (> 1) of x1 * xk. So a faster method would be to loop through the primes p which might be the least member, look for partitions of n-p with all parts at least q = nextprime(p+1) and at most q*(n-q)\nextprime(q+1), and then add on p 'by hand'. This trims off a lot of partitions to check, I think.