2
$\begingroup$

Consider the following polynomial:

expr = (1 + a[x] + c[x, aa] + c[x, aa] da[x, aa] + 
 da[x, aa]*da[x, aa] + 
 dc[x, aa, bb]*da[x, cc]*da[x, bb]*eps[aa, bb, cc, dd])^3 // 

Expand

Among all the terms, I want to leave only those which contain exactly three a[x] or da[x,yy_] and one c or dc[x,yy_,zz_]; it is okay if they have eps[aa,bb,cc,dd]. At least one of the terms from expr matching these requirements is 3 a[x]^2 c[x, aa] da[x, aa].

Could you please tell me if this is possible in a more or less efficient way in Mathematica (in my realistic example, I deal with thousands of terms)?

Edit

This is my ugly realization. It uses the following answer:

ruleFixedPower[expr_, field_, dfield_, n_] := (Normal@
     Series[expr /. {factor : 
          field[x] | field[x, a_] | dfield[x, a_] | 
           dfield[x, a_, b_] :> factor*$T}, {$T, 0, n}] /. {$T -> 
      1}) - (Normal@
     Series[expr /. {factor : 
          field[x] | field[x, a_] | dfield[x, a_] | 
           dfield[x, a_, b_] :> factor*$T}, {$T, 0, 
       Max[0, n - 1]}] /. {$T -> 1})

Then,

expr1 = ruleFixedPower[expr, a, da, 3];
expr2 = ruleFixedPower[expr1, c, dc, 1]

6 a(x) da(x,bb) da(x,cc) dc(x,aa,bb) eps(aa,bb,cc,dd)+6 a(x) c(x,aa) da(x,aa)^2+3 a(x)^2 c(x,aa) da(x,aa)+6 c(x,aa) da(x,aa)^3

Unfortunately, it is slow as hell for my realistic tasks with many thousands of terms and more than two different structures, the powers of which have to be fixed. I may pre-select to have only the terms with the total degree equal to n, where n is the total degree of the desired output, but still, it is slow.

The derivative expansion would work much faster than this construction with Series, but it obviously introduces wrong expansion coefficients.

To illustrate the slowness, consider the following expansion:

expr = (1 + a[x] + c[x, aa] + a[x]*b[aa] + c[x, aa]* da[x, aa] + 
      dc[x, bb] dx[x, dd] + da[x, aa]*da[x, aa] + 
      dc[x, aa, bb]*da[x, cc]*da[x, bb]*eps[aa, bb, cc, dd])^16 // 
   Expand;
expr // Length

167637

(*Selecting terms not with the common power not higher than n*)
ruleCommonPower[expr_, 
  n_] := (Normal@
    Series[expr /. {factor : 
         a[x] | ca[x, a_] | da[x, a_] | dc[x, a_, b_] :> 
        factor*$T}, {$T, 0, 14}] /. {$T -> 1})
expr0 = ruleCommonPower[expr, 14]; // AbsoluteTiming
expr1 = ruleFixedPower[expr0, a, da, 10]; // AbsoluteTiming
expr2 = ruleFixedPower[expr1, c, dc, 4] // AbsoluteTiming

This code gets stuck.

$\endgroup$
6
  • 1
    $\begingroup$ Please edit your question to include the expression for the expected result. $\endgroup$
    – Bob Hanlon
    Commented Aug 10, 2023 at 20:57
  • $\begingroup$ @BobHanlon : I have provided the desired output together with my (very slow) solution for this problem. $\endgroup$ Commented Aug 10, 2023 at 21:23
  • $\begingroup$ It is not clear how the terms of the shown desired output satisfies the requirement of "contain exactly three a[x] or da[x,yy_] and one dc[x,yy_,zz_]" $\endgroup$
    – Bob Hanlon
    Commented Aug 10, 2023 at 22:50
  • 1
    $\begingroup$ For ruleFixedPower, you ought to be able to use SeriesCoefficient and avoid subtracting two series. $\endgroup$
    – Michael E2
    Commented Aug 10, 2023 at 23:28
  • $\begingroup$ This is very confusing. In the description there is mention of a[x], da[x,yy_] , dc[x,yy_,zz_] but the putative term then given, 3 a[x]^2 c[x, aa] da[x, aa], has no dc[...], contrary to the requirement that there be exactly one (I assume this counting is by multiplicity). Before tackling this I'd want to see the actual objective in unambiguous terms. $\endgroup$ Commented Aug 11, 2023 at 0:14

3 Answers 3

1
$\begingroup$

One can create artificial powers in a new variable and then extract only those of interest from the expanded form. For this example we want to have exactly 3 of one type of variable and exactly 1 of another. We can take use E in the exponents of the first type, and integers for the second. Then we are seeking powers of 1+3*E.

expr = (1 + a[x] + c[x, aa] + c[x, aa] da[x, aa] + 
     da[x, aa]*da[x, aa] + 
     dc[x, aa, bb]*da[x, cc]*da[x, bb]*eps[aa, bb, cc, dd])^3;

powers = {E, E, 1, 1};
reps = Thread[{a[x], da[x, yy_], dc[x, yy_, zz_], c[x, yy_]} -> 
    t^powers*{a[x], da[x, yy], dc[x, yy, zz], c[x, yy]}];

e2 = Expand[expr] /. reps;
terms = Cases[e2, t^(1 + 3*E)*aa_ :> aa]

(* Out[213]= {3 a[x]^2 c[x, aa] da[x, aa], 6 a[x] c[x, aa] da[x, aa]^2, 
 6 c[x, aa] da[x, aa]^3, 
 6 a[x] da[x, bb] da[x, cc] dc[x, aa, bb] eps[aa, bb, cc, dd]} *)
$\endgroup$
2
  • $\begingroup$ Thanks! Looks very elegant, and it may be extended to an arbitrary sets of variables. $\endgroup$ Commented Aug 12, 2023 at 12:28
  • $\begingroup$ Thanks again! It works fantastically fast. $\endgroup$ Commented Aug 15, 2023 at 19:16
2
$\begingroup$

I am not sure, if the following selection meets your demand, but its the way to implement such sieves

   Select[expr, ( Count[#, _a | _da] == 3 && Count[#, _d | _dc] == 1 &)]   

     6 a[x] da[x, bb] da[x, cc] dc[x, aa, bb] eps[aa, bb, cc, dd] + 
         6 a[x] c[x, aa] da[x, bb] da[x, cc] dc[x, aa, bb] 
           eps[aa, bb, cc,dd] +
         6 c[x, aa] da[x, aa] da[x, bb] da[x, cc] dc[x, aa, bb] 
       eps[aa, bb, cc, dd] + 
         6 c[x, aa]^2 da[x, aa] da[x, bb] da[x, cc] dc[x, aa, bb] eps[aa, bb, cc, dd] + 
         6 a[x] da[x, aa]^2 da[x, bb] da[x, cc] dc[x, aa, bb]  bb, cc, eps[aa,dd]

Eventually one has to check for the count of additional features as integer multiples, integer powers or

    Count[#,pattern,depth] 

or specifying the variables pattern inside the given functions a by

 Count[#, a[x]|a[y],Infinity]
$\endgroup$
2
$\begingroup$

You could try with pattern matching. The pattern will be complicated and needs to be protected against evaluation by "HoldPattern".

Then we will first search for terms with 3 a or da by the pattern: "___ (a[x] | da[x, _]) ___ (a[x] | da[x, _]) ___ (a[x] | da[x, _]) ___":

tmp = Cases[expr, 
  HoldPattern[___   (a[x] | da[x, _])  ___  (a[x] | 
      da[x, _])  ___  (a[x] | da[x, _])  ___]];

Then with search the results for c or dc by the pattern: "___ (c | dc[x, _, _]) ___":

Cases[tmp, HoldPattern[___ (c | dc[x, _, _]) ___]]

{6 a[x] da[x, bb] da[x, cc] dc[x, aa, bb] eps[aa, bb, cc, dd], 
 6 a[x] c[x, aa] da[x, bb] da[x, cc] dc[x, aa, bb] eps[aa, bb, cc, 
   dd], 6 c[x, aa] da[x, aa] da[x, bb] da[x, cc] dc[x, aa, bb] eps[aa,
    bb, cc, dd], 
 6 a[x] c[x, aa] da[x, aa] da[x, bb] da[x, cc] dc[x, aa, bb] eps[aa, 
   bb, cc, dd], 
 6 c[x, aa]^2 da[x, aa] da[x, bb] da[x, cc] dc[x, aa, bb] eps[aa, bb, 
   cc, dd], 
 6 a[x] da[x, aa]^2 da[x, bb] da[x, cc] dc[x, aa, bb] eps[aa, bb, cc, 
   dd]}
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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