Integrate $ \int_0^\infty \left( ( x A+I)^{-1} A - \frac{1}{c+x} I \right)\, \mathrm dx $ where $A$ is positive-definite and $c>0$
This might be helpful for you
$(xA+I)^{-1}$ can be presented as an integral using Laplace transform.
Step 1: Since \begin{eqnarray} (A+sI)^{-1}&=&\int_{0}^{\infty}e^{-At}e^{-stI}dt \end{eqnarray} then \begin{eqnarray} (xA+I)^{-1}A&=&\frac{1}{x}\int_{0}^{\infty}e^{-At}Ae^{-\frac{t}{x}I}dt\\ &=&\int_{0}^{\infty}e^{-At}A\frac{e^{-\frac{t}{x}}}{x}dt \end{eqnarray} where $x=\frac{1}{s}$ NOTE: I used the fact that $e^{-\frac{t}{x} I}=Ie^{-\frac{t}{x}}$ since I is the identity matrix.
Step 2: For any positive definite matrix $A$, it is true that:
$$\int_{0}^{\infty}e^{-At}Adt=-e^{-At}|_0^\infty =I$$ then $$\frac{I}{c+x}=\int_{0}^{\infty}\frac{e^{-At}A}{c+x}dt$$
Combine step 1 and step 2 to get
$$(A+sI)^{-1}A-\frac{I}{c+x}=\int_{0}^{\infty}e^{-At}\bigg[\frac{e^{-\frac{t}{x}}}{x}-\frac{1}{c+x}\bigg]Adt$$
Now, integrate both sides with respect to x gives
\begin{eqnarray} \int_{0}^{\infty}\bigg[(A+sI)^{-1}A-\frac{I}{c+x}\bigg]dx&=&\int_{0}^{\infty}e^{-At}\int_{0}^{\infty}\bigg[\frac{e^{-\frac{t}{x}}}{x}-\frac{1}{c+x}\bigg]dxdt\times A \quad...(1)\\ \end{eqnarray}
Solving the scalar integration solves the problem $$I_1(t)=\int_0^\infty\bigg[\frac{e^{-\frac{t}{x}}}{x}-\frac{1}{c+x}\bigg]dx$$
This integral appears in the scalar form as you evaluate it:
\begin{eqnarray} \log(a)+\log(c)&=&\int_{0}^{\infty}e^{-at}\int_{0}^{\infty}\bigg[\frac{e^{-\frac{t}{x}}}{x}-\frac{1}{c+x}\bigg]dxdt\times a\\ &=&\int_{0}^{\infty}e^{-at}I_1(t)dt\times a\\ \end{eqnarray}
Let $a$ be the Laplace varaible, then replacing $a$ by $s$ in the last equations gives \begin{eqnarray} I_1(t)&=& \mathcal{L}^{-1}\left[\frac{\log(s)}{s}+\frac{\log(c)}{s}\right]\\ &=& \frac{1}{t^2}+\log(c)\quad ..(2) \end{eqnarray}
Since A is symmetric it can be diagonalized without any Jordan block, i.e. $A=T^{-1}DT$ where $D=diag\{\lambda_1, ..., \lambda_n\}$. Then \begin{eqnarray} e^{-At}=e^{-T^{-1}DTt}=T^{-1}e^{-Dt}T=T^{-1}\left[\begin{array}{cccc}e^{-\lambda_1t}&0 &... & 0\\ 0 & e^{-\lambda_2t} &...& 0\\ 0 & 0 & ... & e^{-\lambda_nt} \end{array}\right]T \end{eqnarray}
Now, return back to Eq. (1) and substitute Eq. (2) to get:
\begin{eqnarray}
\int_{0}^{\infty}\bigg[(A&+&sI)^{-1}A-\frac{I}{c+x}\bigg]dx=\int_{0}^{\infty}e^{-At} \bigg[\frac{1}{t^2}+\log(c)\bigg]dt A\\
&=&T^{-1}\int_{0}^{\infty}e^{-Dt} \bigg[\frac{1}{t^2}+\log(c)\bigg]dt TA \\
&=&T^{-1}\int_{0}^{\infty}\left[\begin{array}{cccc}e^{-\lambda_1t}&0 &... & 0\\
0 & e^{-\lambda_2t} &...& 0\\
0 & 0 & ... & e^{-\lambda_nt} \end{array}\right]
\bigg[\frac{1}{t^2}+\log(c)\bigg]dt TA \\
&=&T^{-1}\left[\begin{array}{cccc}\int_{0}^{\infty}\big[\frac{1}{t^2}+\log(c)\big]e^{-\lambda_1t}dt&0 &... & 0\\
0 & .... &...& 0\\
0 & 0 & ... & \int_{0}^{\infty}\big[\frac{1}{t^2}+\log(c)\big]e^{-\lambda_nt}dt \end{array}\right] TA \\
\end{eqnarray}
Every diagonal entry of the matrix is the laplace transform with a variable $\lambda_i$, using Eq. (2) they can be written as:
\begin{eqnarray} \int_{0}^{\infty}\big[\frac{1}{t^2}+\log(c)\big]e^{-\lambda_it}&=&\mathcal{L}\{\big[\frac{1}{t^2}+\log(c)\big]\}\\ &=&\frac{\log(\lambda_i)}{\lambda_i}+\frac{\log(c)}{\lambda_i} \end{eqnarray}
Finally, the solution of your problem will be: \begin{eqnarray} \int_{0}^{\infty}\bigg[(A&+&sI)^{-1}A-\frac{I}{c+x}\bigg]dx\\ &=&T^{-1}\left[\begin{array}{cccc}\bigg[\frac{\log(\lambda_1)}{\lambda_1}+\frac{\log(c)}{\lambda_1}\bigg]&0 &... & 0\\ 0 & .... &...& 0\\ 0 & 0 & ... & \bigg[\frac{\log(\lambda_n)}{\lambda_n}+\frac{\log(c)}{\lambda_n} \bigg]\end{array}\right] TA \\ \end{eqnarray}
NOTE:(Laplace inverse derivation in equation (2)) \begin{eqnarray} \mathcal{L}^{-1}\left[\frac{\log(s)}{s}+\frac{\log(c)}{s}\right] =\frac{1}{t^2}+\log(c)\quad ..(2) \end{eqnarray}
It is well known that if $$\mathcal{L}\{g(t)\}=G(s)$$ then
$$\mathcal{L}\{tg(t)\}=\frac{-dG(s)}{ds}$$ and $$\mathcal{L}\{\int_0^t g(\tau)d\tau\}=\frac{G(s)}{s}$$
then $$\mathcal{L}\{t\int_0^t g(\tau)d\tau\}=\frac{-d\bigg[\frac{G(s)}{s}\bigg]}{ds}$$ Substitute $G(s)\left[\frac{\log(s)}{s}+\frac{\log(c)}{s}\right]$ and take the Laplace inverse (which is very easy now) for both sides then divide by t and differentiate with respect to t get g(t).
Regards,
As $A$ is PSD, by a (constant) change of basis, all of the integrands in the family $$ \left\{(xA+I)^{-1}A-\frac1{c+x}I:\ x\ge0\right\} $$ are simultaneously diagonalisable. Hence the integral reduces to the scalar case and it exists if and only if $A$ is positive definite (i.e. $A$ is also nonsingular). Also, when it exists, the integral is equal to $\log A+\log(c)I$.