the usual enforcement made here is that A'x + B' / x = 0..............Eqn. (1)
so that the derivative is written as
y' = A - B / x^2
y'' = A' - B' / x^2 + 2 B / x^3
Putting y and its derivatives into the inhomogeneous DE:
x^2 y'' + xy' - y = 4x ln(x)
x^2(A' - B' / x^2 + 2 B / x^3) + x(A - B / x^2) - Ax - B/x = 4x ln(x)
A' x^2 - B' + 2B/x + Ax - B/x - Ax - B/x = 4x ln(x)
A' x^2 - B' = 4x ln(x)
A' - B'/x^2 = 4 ln(x) / x...........Eqn.(2)
Thus we have two equations in two unknowns (A' and B')
A'x + B' / x...= 0........................Eqn.(1)
A' - B' / x^2 = 4 ln(x) / x..........Eqn.(2)
Recast it in matrix form:
[x............1/x...][A'] = [0]
[1.........-1/x^2][B'] = [4 ln(x) / x]
(notice the reason I divided through by x^2 is now apparent, the left-hand matrix now is the Wronskian matrix, it has the functions x and 1/x, and the subsequent row is their derivatives repsectively).
define the above matrices with the labels:
WX = b
Solve this with Cramer's rule (note: det( ) means determinant):
A' = det(b|W) / det(W)
where
det(W) = det(x..............1/x...) = -1/x - 1/x = -2/x
.....................(1...........-1/x^…
det(b|W) = det(0.........................1/x...) = -4 ln(x) / x^2
.........................(4 ln(x) / x.........-1/x^2)
so
A' = det(b|W) / det(W) = (- 4 ln(x) / x^2 ) / (-2/x) = 2 ln(x) / x
and
B' = det(W|b) / det(W)
where
det(W|b) = det(x........................ 0...) = 4 ln(x)
.........................(1............… ln(x) / x)
so
B' = det(W|b) / det(W) = 4 ln(x) / (-2/x) = -2 x ln(x)
in Summary
A' = 2 ln(x) / x --> integrate to get A = ln^2(x)
B' = -2 x ln(x) --> integrate to get B = (1/2)x^2(1 - 2 ln(x)), note this integration is done by parts twice
Thus,
y = Ax + B/x
y = ln^2(x) x + (1/2)x^2(1 - 2 ln(x)) / x
y = ln^2(x) x + (1/2) x (1 - 2 ln(x))
y = ln^2(x) x + x/2 - x ln(x)
y = x ln(x) (ln(x) - 1) + x/2.............[Ans.]