How to calculate the square root of matrix $A+B$ perturbatively? $A=diag\{\lambda_1,...,\lambda_n\}$ and $\lambda_i>0$, $B$ is a positive definite symmetric matrix and  $max\{B_{ij} \}\ll min\{\lambda_i\}$
Note that the perturbative calculation of square root of $I+B$ is very easy, where $B$ is a small matrix.  
How to calculate the square root of $A+B$ perturbatively?
 A: EDIT. Finally user34669 is right. We assume that $A$ is fixed and $B$ tends to $0$. The following $3$ lines are not correct because, in general, $\sqrt{XY}\not=\sqrt{X}\sqrt{Y}$.
"$\Delta=\sqrt{A+B}-\sqrt{A}\approx (A^{-1/2}B)/2$. We may also write $\Delta\approx (BA^{-1/2})/2$ or in a symmetric form $\Delta\approx (A^{-1/2}B+BA^{-1/2})/4$ and the formula is valid for $A$ symmetric $>0$ and $B$ a small symmetric matrix." 
In fact, the previous approximations give an error in $O(||[A^{1/2},B]||)$, that is not interesting except if $A,B$ commute (we would like an error in $O(||B||^2)$ or at least in $o(||B||)$).
On the other hand, the other part is correct.
The function $S:X\rightarrow \sqrt{X}$ is defined and derivable on the set of SPD matrices. Let $K=DS_A(H)$ be the derivative of $S$ in $A$, where $H$ is a variable SYMMETRIC matrix. Here $SS=I$ implies $K\sqrt{A}+\sqrt{A}K=H$, a Sylvester equation in the unknown $K$. 
i) A closed form. It is known that  $K=\int_0^{\infty}e^{-t\sqrt{A}}He^{-t\sqrt{A}}dt$. That implies that if $B$ is a small symmetric matrix, then $\Delta\approx \int_0^{\infty}e^{-t\sqrt{A}}Be^{-t\sqrt{A}}dt$, a symmetric matrix.
ii) Numerically, it is easier to solve directly the Sylvester equation diagonalizing $A$.
A: See Theorem 3.25 on page 122 (as given by the page numbers) of http://www.renyi.hu/~petz/pdf/matrixPD.pdf
