4
$\begingroup$

I have one case (Nn=1) for 7 Polynomials as following:

p1= a[0]*b[0] + c[0]*d[0] + e[0]*f[0] + g[0]*h[0];
p2= a[1]*c[1] + b[1]*d[1] + e[1]*g[1] + f[1]*h[1]; 
p3= a[2]*d[2] + b[2]*c[2] + e[2]*h[2] + f[2]*g[2]; 
p4= a[3]*e[3] + b[3]*f[3] + c[3]*g[3] + d[3]*h[3]; 
p5= a[4]*f[4] + b[4]*e[4] + c[4]*h[4] + d[4]*g[4]; 
p6= a[5]*g[5] + b[5]*h[5] + c[5]*e[5] + d[5]*f[5]; 
p7= a[6]*h[6] + b[6]*g[6] + c[6]*f[6] + d[6]*e[6];

and I want to pick up 4 terms from the 7*4=28 terms such as a[0]*b[0], c[0]*d[0] and so on, which gives me the output in the form of a[x1]*b[x2]*c[x3]*d[x4]*e[x5]*f[x6]*g[x7]*h[x8] (we can later write it as FF[x1,x2,x3,x4,x5,x6,x7,x8]).

The way I do is as following:

outputstemp={};
outputs={};
FFnCaseList={a[c1_]*b[c2_]*c[c3_]*d[c4_]*e[c5_]*f[c6_]*g[c7_]*h[c8_]->FF[c1,c2,c3,c4,c5,c6,c7,c8], a[_]->0, b[_]->0, c[_]->0,  d[_]->0, e[_]->0, f[_]->0, g[_]->0, h[_]->0,  x_[_]^n_->0};

AllRows = (p1+p2+p3+p4+p5+p6+p7)^4; 
AppendTo[outputstemp, ExpandAll[AllRows]];

Timing[AppendTo[outputs, outputstemp/.FFnCaseList];]

the time is {10.2344, Null}

So is there any quick way to do this? If I have many different cases (for example Nn=6000) of 7 Polynomials, then I can make a loop that run similar things Nn times (6000*10s). That will takes roughly 16.7 hours!!!

Any comments or suggestions are appreciated! Thank you very much!

$\endgroup$
3
$\begingroup$

I would use Coefficient after mopdifying so that terms of the desired form have a unique power in a new variable.

p1 = a[0]*b[0] + c[0]*d[0] + e[0]*f[0] + g[0]*h[0];
p2 = a[1]*c[1] + b[1]*d[1] + e[1]*g[1] + f[1]*h[1];
p3 = a[2]*d[2] + b[2]*c[2] + e[2]*h[2] + f[2]*g[2];
p4 = a[3]*e[3] + b[3]*f[3] + c[3]*g[3] + d[3]*h[3];
p5 = a[4]*f[4] + b[4]*e[4] + c[4]*h[4] + d[4]*g[4];
p6 = a[5]*g[5] + b[5]*h[5] + c[5]*e[5] + d[5]*f[5];
p7 = a[6]*h[6] + b[6]*g[6] + c[6]*f[6] + d[6]*e[6];

Multiply each variable by t^2^j where j is the index of that variable in the list {a,b,c,d,e,f,g,h}.

polys = {p1, p2, p3, p4, p5, p6, p7} /. 
  Thread[{v_a, v_b, v_c, v_d, v_e, v_f, v_g, v_h} -> 
    v*t^(2^Range[0, 7])]

(* Out[13]= {t^3 a[0] b[0] + t^12 c[0] d[0] + t^48 e[0] f[0] + 
  t^192 g[0] h[0], 
 t^5 a[1] c[1] + t^10 b[1] d[1] + t^80 e[1] g[1] + t^160 f[1] h[1], 
 t^6 b[2] c[2] + t^9 a[2] d[2] + t^96 f[2] g[2] + t^144 e[2] h[2], 
 t^17 a[3] e[3] + t^34 b[3] f[3] + t^68 c[3] g[3] + t^136 d[3] h[3], 
 t^18 b[4] e[4] + t^33 a[4] f[4] + t^72 d[4] g[4] + t^132 c[4] h[4], 
 t^20 c[5] e[5] + t^40 d[5] f[5] + t^65 a[5] g[5] + t^130 b[5] h[5], 
     t^24 d[6] e[6] + t^36 c[6] f[6] + t^66 b[6] g[6] + t^129 a[6] h[6]} *)

Pull out terms that have t to the power 1+2+4+...+128, or 255. In order to speed the conversion to the FF[...] notation I replace sums and products by lists (this speeds otherwise slow pattern match). In a more general setting one might have to canonicalize the ordering using Sort.

Timing[
 allvarsmonoms = 
  Expand[Coefficient[(Apply[Plus, polys])^4, t^Total[2^Range[0, 7]]]];
 res0 = Apply[List, allvarsmonoms];
 res1 = res0 /. {Plus -> List, Times -> List};
 res = Apply[Plus, 
   res1 /. {num_, a[v1_], b[v2_], c[v3_], d[v4_], e[v5_], f[v6_], 
      g[v7_], h[v8_]} -> num*FF[v1, v2, v3, v4, v5, v6, v7, v8]];]
Length[res]

(* Out[71]= {0.0625, Null}

Out[72]= 105 *)
$\endgroup$
3
  • $\begingroup$ Thank you very much! If variable in the list is {a,b,c,d,e,f}, should I use v*t^(2^Range[0, 5])? $\endgroup$
    – Xuemei
    Apr 21 '20 at 17:25
  • 1
    $\begingroup$ Yes. For cleaner code, I should have defined the variable symbols in a list as vars and then used Length[vars] for that exponent. $\endgroup$ Apr 21 '20 at 18:04
  • $\begingroup$ Thank you! that speed up a lot. Very interesting! $\endgroup$
    – Xuemei
    Apr 21 '20 at 18:17
5
$\begingroup$

I am sure there is a less awkward method than the following proposal:

n = 4;
expr = p1 + p2 + p3 + p4 + p5 + p6 + p7;

Nest[Map[If[Signature[Head /@ (List @@ #)] == 0, 0, #] &, Expand[expr * #]] &, 1, n] /. 
(a[c1_] b[c2_] c[c3_] d[c4_] e[c5_] f[c6_] g[c7_] h[c8_]) :>
FF[c1, c2, c3, c4, c5, c6, c7, c8]

which yields 24 a[0] b[0] c[0] d[0] e[0] f[0] g[0] h[0] + ... 24 a[6] b[6] c[6] d[6] e[6] f[6] g[6] h[6] in about a hundredth of the time taken for the OP's original code.

$\endgroup$
5
  • $\begingroup$ wow, that's amazing. Well, i didn't saw hundredth of time from my computer but factor of 20 is quite good! Thank you! $\endgroup$
    – Xuemei
    Apr 20 '20 at 18:44
  • 1
    $\begingroup$ changing the replacement rule to {Times -> (# FF[##2] &), Alternatives[a, b, c, d, e, f, g, h] -> Identity} seems to double the speed. (+1) $\endgroup$
    – kglr
    Apr 20 '20 at 19:36
  • $\begingroup$ @kglr, yes, that's cool! Thank you very much! $\endgroup$
    – Xuemei
    Apr 20 '20 at 19:53
  • $\begingroup$ @Xuemei, your "additional" question should surely be a separate question. $\endgroup$ Apr 21 '20 at 0:50
  • 1
    $\begingroup$ @J.M., thank you very much and I have posted it as a separate question here. Could you have a look? would be extremely helpful and thank you a lot! $\endgroup$
    – Xuemei
    Apr 21 '20 at 7:02
3
$\begingroup$
expr  = p1 + p2 + p3 + p4 + p5 + p6 + p7;

res = Total[Multinomial[1, 1, 1, 1] FF @@@ 
       (ReplaceAll[ _[a_] :> a] @ Select[DuplicateFreeQ @* Map[Head]]@
         Apply[Sort@*Join, Subsets[expr /. Plus|Times -> List, {4}],  {-4}])]; // 
           RepeatedTiming // First
 .16

this gives the same result as J.M.'s approach:

res2 = Nest[Map[If[Signature[Head /@ (List @@ #)] == 0, 0, #] &, 
        Expand[expr*#]] &, 1, n] /. 
         (a[c1_] b[c2_] c[c3_] d[c4_] e[c5_] f[c6_] g[c7_] h[c8_]) :> 
           FF[c1, c2, c3, c4, c5, c6, c7, c8]; //  RepeatedTiming // First
0.54
res == res2
True
$\endgroup$
1
  • $\begingroup$ thank you very much and I have posted it as a separate question here. Could you have a look? would be extremely helpful and thank you a lot! $\endgroup$
    – Xuemei
    Apr 21 '20 at 7:05

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.