4
$\begingroup$

My question is closely related to this one but the accepted answer does not meet my needs.

I have a rational function $f(x)$ in one variable which admits a partial fraction decomposition of the form $$\sum_k\frac{A_k}{(x-\alpha_k)^{d_k}}$$ unknown in advance. I would like to produce from $f(x)$ the list of triples $\left\{ (A_k,\alpha_k,d_k)\right\}$ or three lists $\left\{A_k\right\}$, $\left\{\alpha_k\right\}$, $\left\{d_k\right\}$. Can this be done, say, using the function Apart?

More exactly, my problem with the existing answer is that it only works in pathological examples where the polynomials appearing in the denominator after the application of Apart are monic. For example, it already fails to produce $\left\{\left\{3/2\right\},\left\{1/2\right\}\right\}$ for $f(x)=3/(2x-1)$. Is there a simple way to correct this?

Update: Consider a test function $f(x)=\frac{3}{(2x-1)(4x-3)^2}$ which has a decomposition $$f(x)=\frac{3}{2\left(x-\frac12\right)}+\frac{3}{8\left(x-\frac34\right)^2}-\frac{3}{2\left(x-\frac34\right)}.$$ In this case the code should produce three lists $$A=\left\{\tfrac32,\tfrac38,-\tfrac32\right\},\qquad \alpha=\left\{\tfrac12,\tfrac34,\tfrac34\right\},\qquad d=\left\{1,2,1\right\}.$$ Finally I have managed to produce something that works:

f[k_] := 3/((2*k - 1)*(4*k - 3)^2); 
numer := Numerator[List @@ Apart[f[k], k]]; 
denom := Denominator[List @@ Apart[f[k], k]]; 
poles := First /@ (k /. Solve[#1 == 0, {k}] & ) /@ denom; 
expos := Exponent[denom, k]; 
coefs := Limit[denom/k^expos, k -> Infinity]; 
{numer/coefs, poles, expos}

with the expected result

enter image description here

This is okay for my purposes but maybe there is a nicer solution that would not involve taking the limits? (in principle all information is already contained in the result of application of List@@Apart)

$\endgroup$
2

2 Answers 2

1
$\begingroup$

Adjust the variable k as needed:

pfd = Apart[3/((2*k - 1)*(4*k - 3)^2), k];
pfdData = Replace[List @@ pfd,
   a_.*Power[c_. k + b_., d_] :> {a*c^d, -b/c, -d},
   1]
desiredLists = Transpose@pfdData
(*
{{3/2, 1/2, 1}, {3/8, 3/4, 2}, {-(3/2), 3/4, 1}}
{{3/2, 3/8, -(3/2)}, {1/2, 3/4, 3/4}, {1, 2, 1}}
*)

Alternative, using the undocumented function Integrate`ComplexApart[], which does the decomposition in the form desired by the OP:

pfd = Integrate`ComplexApart[3/((2*k - 1)*(4*k - 3)^2), k];
pfdData = 
 Replace[List @@ pfd, a_.*Power[k + b_., d_] :> {a, -b, -d}, 1]
desiredLists = Transpose@pfdData
(*
{{3/8, 3/4, 2}, {-(3/2), 3/4, 1}, {3/2, 1/2, 1}}
{{3/8, -(3/2), 3/2}, {3/4, 3/4, 1/2}, {2, 1, 1}}
*)
$\endgroup$
2
  • $\begingroup$ Alternative to the alternative: Integrate`ComplexApart[3/((2*k - 1)*(4*k - 3)^2), k] /. {Plus | Times | Power -> List, k -> Nothing} // Map@Flatten // Transpose (assumes the orderless arguments of Plus and Times always occur in the right order). $\endgroup$
    – Michael E2
    Jul 27, 2022 at 22:51
  • $\begingroup$ This is perfect, thank you very much! $\endgroup$ Jul 27, 2022 at 23:01
3
$\begingroup$

To reuse the previous answer

k = 4;

(Format[#[n_]] := Subscript[#, n]) & /@ {A, a, d};

(f[x_] = 
   Total[Array[A, k]/(x - Array[a, k])^Array[d, k]]) // TraditionalForm

enter image description here

dList = Cases[f[x], Power[b_, -exp_] :> exp, Infinity]

enter image description here

f2[x_] = f[x] /. Power[b_, -exp_] :> Power[b, -1]

enter image description here

pf[func_, var_Symbol] := 
 Module[{af = List @@ Apart[func, var]}, {Numerator /@ af, 
   var - (Denominator /@ af)}]

{Alist, aList} = pf[f2[x] // Together, x]

enter image description here

EDIT:

Clear["Global`*"]

kmax = 4;

(Format[#[n_]] := Subscript[#, n]) & /@ {a, b, c, d};

(f[x_] = Sum[b[k]/(c[k]*x - a[k])^d[k], {k, 1, kmax}]) //
  TraditionalForm

enter image description here

({αList, Alist, dList} = 
 Transpose@Cases[List @@ f[x], b_.*(c_.*x - a_)^d_ :> {a/c, b/c^(-d), -d}, 1]) // TraditionalForm

enter image description here

f[x] == Total[Alist/(x - αList)^dList] // 
 Simplify[#, Element[Array[d, kmax], PositiveIntegers]] &

(* True *)
$\endgroup$
3
  • 1
    $\begingroup$ Thank you very much for your answer. However my main problem with the previous code is not the inclusion of higher order poles but the fact that the coefficients of the leading powers of $x$ in the denominators are in general different from $1$. If I understand correctly, this is not corrected in the new code (e.g. for $f(x)=3/(2x-1)$ Denominator/@Apart[3/(2x-1)] gives $2x-1$ and not $x-1/2$ as I would like). $\endgroup$ Jul 27, 2022 at 16:24
  • $\begingroup$ Edit your question to show your actual problem. Your example does not match the form that you said you were working with. $\endgroup$
    – Bob Hanlon
    Jul 27, 2022 at 16:28
  • 1
    $\begingroup$ It does (I never said $A_k$ or $\alpha_k$ are integers) but I will try to make it more clear. $\endgroup$ Jul 27, 2022 at 16:31

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.