# Deleting powers of a common prime except the highest one from a list with further only primes

In a list with factors of which most are prime, some are integer powers of some common prime as base.

All the powers from 1 to imax of that base are present in the list. The value of imax, however, is unknown. I have to delete from the list all powers of that base except the highest.

Equivalently we can e.g.:

• change all the values of the powers of the known base into the base value and multiply them,
• delete all the lower powers of the known common base from the list.

An example with base = 5 and imax = 3:

list1 = {2, 5, 13, 25, 125, 3001, 230686501}
should become {2, 13, 125, 3001, 230686501}


or one with base = 13 and imax = 2:

list2 = {13, 169, 233, 741469, 159607993}
should become {169, 233, 741469, 159607993}


I did find Gather and GatherBy, is it possible to GatherBy common powers of the same (beforehand unknown common) base?

In[1]:= li3 = {2, 5, 11, 5, 7, 5};

In[2]:= Times @@@ Gather[li3]

Out[2]= {2, 125, 11, 7}


Any help is appreciated. Thanks.

-
You may consider changing your nickname for something more "human" – belisarius has settled Dec 23 '13 at 21:25
Should 1 be deleted under any base? – Silvia Dec 24 '13 at 0:17
@Silvia Your example does not fit in the Question. – Rombert Dec 24 '13 at 1:40
In case I misunderstand your quesiton, does this what you want: $\{p_1^{n_{11}},\dots,p_1^{n_{1N}}, p_2^{n_{21}},\dots,p_2^{n_{2M}},\dots\}\xrightarrow{\text{for base }p_2} \{p_1^{n_{11}},\dots,p_1^{n_{1N}}, p_2^{\max(n_{2k}|k=1,\dots,M)},\dots\}$? – Silvia Dec 24 '13 at 2:00
@Silvia: Not precisely. The non-primes in the list are all the powers of some common prime up to a certain max. – Rombert Dec 24 '13 at 2:09

Another one:

f = # ~Complement~ Flatten[Most /@ Divisors[#]] &

f @ {2, 5, 13, 25, 125, 3001, 230686501}
(* {2, 13, 125, 3001, 230686501} *)

f @ {13, 169, 233, 741469, 159607993}
(* {169, 233, 741469, 159607993} *)


In more familiar notation:

f[x_] := Complement[x, Flatten[Most /@ Divisors[x]]]


A faster method

If the list is long, the numbers are large, and performance is important, it may be better to search the list backwards until a non-prime is found. This must be base^imax (assuming the list is sorted). It is then straightforward to construct the set of elements which should be removed from the list.

f2[x_] := Module[{base, imax},
{base, imax} = FactorInteger[Select[Reverse@x, ! PrimeQ[#] &, 1]][[1, 1]];
Complement[x, base^Range[imax - 1]]]

-
That's bloody beautiful and I feel rather foolish for what I wrote now, but I'll leave my answer so all can see how you trounced me so completely. – Mr.Wizard Dec 27 '13 at 22:27
@Simon Thank you for your answer. I have studied your method and understands how it works, up to the end. I am not so impressed as Mr.Wizard appears to be, but that is because I could never have come up with this solution and in this form. In the last step you use infix notation which is completely new to me. I understand the infix part in Mr.Wizard's new Method (standing now above) but yr part is more complicated because of the #'s before and after the infix part. Would you be so kind to explain it to me? – Rombert Dec 28 '13 at 14:57
@Simon Using for instance: #~Complement~Flatten[Most /@ Divisors[#]] &@{5, 25} .Thanking you in advance ... – Rombert Dec 28 '13 at 15:09
@Mr.Wizard, thank you - I think "trounced" is an exaggeration though! - as you've pointed out the code may be concise but the algorithm is inefficient for large imax. – Simon Woods Dec 28 '13 at 16:47
@Rombert, I've added a version using more familiar notation. – Simon Woods Dec 28 '13 at 16:50
base = 5 ;
list = Sort@{2, 5, 13, 25, 125, 3001, 230686501};

sel = list /. x_ /; Mod[x, base] == 0 :> False /. x_?NumericQ :> True /.
{a___, False, b : True ..} :> {a, True, b}
Pick[list, sel]
(*
{2, 13, 125, 3001, 230686501}
*)


Edit

Or:

base = 5 ;
list = {2, 5, 13, 25, 125, 3001, 230686501};
list1 = Prepend[Sort@list, 1];
Rest@Join[#[[1]], #[[2, -1 ;; -1]]] &@ GatherBy[list1, Mod[#, base] == 0 &]
(*
{2, 13, 3001, 230686501, 125}
*)


Or:

base = 5 ;
list = Sort@{2, 5, 13, 25, 125, 3001, 230686501};
DeleteCases[list, Alternatives @@ Most@Cases[list, x_ /; Mod[x, base] == 0]]
(*
{2, 13, 3001, 230686501, 125}
*)

-
I certainly have to study this! Thank you! – Rombert Dec 23 '13 at 20:09
Nice answers, thank you again! OK, you introduce the value of the base. Is it possible, you think, to discover that value without visual inspection, so automatically? – Rombert Dec 23 '13 at 20:49
@user8778 Unless you have more restrictions on the problem, there could be multiple "bases" – belisarius has settled Dec 23 '13 at 20:58
You can automatically detect that "bases"? – Rombert Dec 23 '13 at 21:05

I am answering my interpretation of question. Apologies if I have misunderstood. If the aim is to retain the largest power of a prime (and not the integer with largest power prime) the following will achieve this:

rem[b_Integer, u_Integer] := Module[{lg, mx, rl},
lg = Log[b, #] & /@ u;
mx = Max@Cases[lg, _Integer];
rl = Or[Not[IntegerQ[#]], # == mx] & /@ lg;
Pick[u, rl]]


Testing:

list = Sort@{2, 5, 13, 25, 125, 3001, 230686501, 150, 250};
rem[5,list]


yields:

{2, 13, 125, 150, 250, 3001, 230686501}

Note belisarius code which is based on Mod will select largest multiple of base and will remove all powers in above test case and retain 250. Other counter examples can also show this, e.g remove 250 in above and 150 will be retained and all powers removed.

If the aim is to find integer with largest base prime exponent that is a different question.

-
In my lists are only two ingridients: A: simple primes occurring once and B: powers of a common base prime but they occur also only once and are all present up to the power imax. – Rombert Dec 24 '13 at 8:25
thank you for the clarification. In that case it seems to me Modulus based approach will return a multiple of base that is not a power of base and reject smaller but pure powers of base... – ubpdqn Dec 24 '13 at 9:01
I studied your answer and like it, however had a little problem in getting it to work (I changed u_Integer into u_ in the arguments of your Module). I agree with your remark about the solution of belisarius using Mod. However, in the lists that I am using only Primes and integer powers of primes will occur. I will also investgate the answers of belisarius and Mr.Wizard , so I did not yet make my choice, which will be rather difficult. – Rombert Dec 27 '13 at 12:28
@Rombert ok. My code determines the maximum power of base (b) in the list and then picks from the list either elements which are not integer powers of the base or the maximum ower of the base, where you choose the base, Your case of interest is a prime integer. – ubpdqn Dec 27 '13 at 13:17

I believe the following should work based on your comment:

in my lists are only two ingridients: simple primes occurring onceand powers of a common base prime

## New Methods

After a bit of reflection I think there is a much simpler approach:

strip2[a : {__Integer}] :=
Power @@@ Last /@ GatherBy[Join @@ FactorInteger @ Sort @ a, First] // Sort


Test:

strip2 @ {2, 5, 13, 25, 125, 3001, 230686501}
strip2 @ {13, 169, 233, 741469, 159607993}

{2, 125, 13, 3001, 230686501}

{169, 233, 741469, 159607993}


Here is another method based on Simon's answer. It is only slightly longer than his code yet is more efficient for large $imax$ because it eliminates a large number of factors before sorting:

strip3[a_] := Sort @ GatherBy[Divisors @ a, #[[2]] &][[All, -1, -1]]


### Timings

To justify the existence of this answer after Simon Woods posted a wonderfully clean method I shall illustrate that there can be a performance advantage to my method, specifically in the case of large $imax$ values.

Building a long list b matching the problem description:

SeedRandom[1]
a = Union@RandomPrime[2000, 300];
base = RandomChoice[Take[a, 5]];
imax = 200;
b = Union[a, base^Range[imax]];


Timings:

Do[strip2[b], {500}]  // Timing // First
Do[strip3[b], {500}]  // Timing // First
Do[f[b], {500}]       // Timing // First

0.218

1.17

3.495


This is quite a pathological example as my test list has values up to 1.6*10^60 which is probably quite unlike your application. With an $imax$ of 25 the timings are equivalent.

Also one may find value in working with the gathered factors list if other manipulation is desired.

## Old method

Here is a fall-back in case the new method above does not work well:

strip[a : {__Integer}] :=
Module[{p, c, idx},
{p, c} = GatherBy[Sort@a, PrimeQ];
idx = Position[Outer[Log, p, c], {__Integer}];
p ~Delete~ idx ~Append~ Last[c] // Sort
]

-
I have figured out how your methods (new and old) are working. They both work fine! Btw thank you. I just received a new answer from Simon Woods. Something "compeletely" different. Now I have to figure out the precise working of Complement in this case! How long can I wait with my choosing the "winner"? What is the best way to publish my motivation? I often told my students: the best way to learn something new, is to ask a question (;) – Rombert Dec 27 '13 at 17:10
@Rombert You're welcome. I just saw Simon's answer now and it looks like the clear winner. You can wait as long as you like to Accept an answer, but know that you can also change your Accept at any time, and if someone posts a better answer you can (and arguably should) change to Accept to it. I don't know what you meant by "What is the best way to publish my motivation?" Could you clarify? – Mr.Wizard Dec 27 '13 at 22:24
With publish my motivation, I mean give the argumentation for Accepting that one. Or is it just make yr choice. You know, we are used to say why! Your methods have some advantages in checking whether all the powers of the common base prime are present in the given list or not. So I believe. – Rombert Dec 27 '13 at 22:44
@Rombert You're surely welcome to explain why you chose a particular answer but it's not necessary. I really like clean, concise code so I wish I had thought of Simon's solution. Nevertheless I agree that there are aspects of my (new) method that may be advantageous. In certain circumstances it is more than an order of magnitude faster. I shall update my answer with an example. – Mr.Wizard Dec 27 '13 at 22:58