How do you integrate `f(x)=(x^3+5x)/((x^2-1)(x-3)(x+9))` using partial fractions?

The denominator factors into distinct linear factors as:

`(x-1)(x+1)(x-3)(x+9)`

So we can look for a partial fraction decomposition of the form:

`(x^3+5x)/((x^2-1)(x-3)(x+9))`

`=A/(x-1)+B/(x+1)+C/(x-3)+D/(x+9)`

`=(A(x+1)(x-3)(x+9)+B(x-1)(x-3)(x+9)+C(x^2-1)(x+9)+D(x^2-1)(x-3))/((x^2-1)(x-3)(x+9))`

`=(A(x^3+7x^2-21x-27)+B(x^3+5x^2-33x+27)+C(x^3+9x^2-x-9)+D(x^3-3x^2-x+3))/((x^2-1)(x-3)(x+9))`

`=((A+B+C+D)x^3+(7A+5B+9C-3D)x^2+(-21A-33B-C-D)x+(-27A+27B-9C+3D))/((x^2-1)(x-3)(x+9))`

Equating coefficients, we get the following system of linear equations:

`{ (A+B+C+D = 1), (7A+5B+9C-3D = 0), (-21A-33B-C-D = 5), (-27A+27B-9C+3D = 0) :}`

This can be expressed in matrix form as:

`((1, 1, 1, 1, 1), (7, 5, 9, -3, 0), (-21, -33, -1, -1, 5), (-27, 27, -9, 3, 0))`

Perform a sequence of row operations on the matrix to get the left hand `4 xx 4` submatrix into the form of an identity matrix. Then the right hand column will contain the values of `A, B, C` and `D`.

Add/subtract multiples of the first row to/from the others to get:

`((1, 1, 1, 1, 1), (0, -2, 2, -10, -7), (0, -12, 20, 20, 26), (0, 54, 18, 30, 27))`

Add/subtract multiples of the second row to/from the third and fourth rows to get:

`((1, 1, 1, 1, 1), (0, -2, 2, -10, -7), (0, 0, 8, 80, 68), (0, 0, 72, -240, -162))`

Subtract `9` times the third row from the fourth row to get:

`((1, 1, 1, 1, 1), (0, -2, 2, -10, -7), (0, 0, 8, 80, 68), (0, 0, 0, -960, -774))`

Divide the `2`nd, `3`rd and `4`th rows by values to make the diagonal all `1`'s:

`((1, 1, 1, 1, 1), (0, 1, -1, 5, 7/2), (0, 0, 1, 10, 17/2), (0, 0, 0, 1, 129/160))`

Subtract row `2` from row `1` to get:

`((1, 0, 2, -4, -5/2), (0, 1, -1, 5, 7/2), (0, 0, 1, 10, 17/2), (0, 0, 0, 1, 129/160))`

Add/subtract multiples of row `3` to/from rows `1` and `2` to get:

`((1, 0, 0, -24, -39/2), (0, 1, 0, 15, 12), (0, 0, 1, 10, 17/2), (0, 0, 0, 1, 129/160))`

Add/subtract multiples of row `4` to/from rows `1`, `2` and `3` to get:

`((1, 0, 0, 0, -3/20), (0, 1, 0, 0, -3/32), (0, 0, 1, 0, 7/16), (0, 0, 0, 1, 129/160))`

Hence:

`{ (A=-3/20), (B=-3/32), (C=7/16), (D=129/160) :}`

So:

`int (x^3+5x)/((x^2-1)(x-3)(x+9)) dx`

`=int A/(x-1)+B/(x+1)+C/(x-3)+D/(x+9) dx`

`=-3/20 ln abs(x-1) -3/32 ln abs(x+1) + 7/16 ln abs(x-3) + 129/160 ln abs(x+9) + "constant"`