Second derivative of $\sum_i \left(\frac{1}{1+\|Xa_i\|^2}-b_i\right)^2$ Let $A$ be a $n\times m$ matrix, $X$ be a $d \times m$ matrix and $b_i\in \mathbb R$, $i=1\dots n$.
$$
f(x) = 4 XA^\top (I+\Delta_X)^{-3}(B+B\Delta_X-I)A
$$
where $B=\sum_ie_ib_ie_i^\top$ and $\Delta_X=\sum_i e_ie_i^\top AX^\top XA^\top e_ie_i^\top$ are diagonal matrices.
I have found that $f(x)$ is the first derivative of
$$
\sum_i \left(\frac{1}{1+\|Xa_i\|^2}-b_i\right)^2
$$
What is the derivative $\partial f / \partial X$ ?
I'm trying to get it from this:
$$\tag{1}\label{eq1}
\frac{\partial^2}{\partial x^2}\left(\frac{1}{1+g_i(x)}-b_i\right)^2 = 2\left(\frac{3-b_i-b_ig_i(x)}{(1+g_i(x))^4}g'(x)^2 + \frac{b_i+b_ig_i(x)-1}{(1+g_i(x))^3}g''(x)\right)
$$
where $g_i(x)=\|xa_i\|^2$.
But I'm facing some difficulties with $g'(x)^2$, considering $g'(x)$ is a matrix. Besides, I'm not really sure \eqref{eq1} can help.
 A: Before proceeding to the second derivative, let's check on that first derivative.
Define the diagonal matrix $C$ using the Hadamard product as follows: 
$$ \eqalign {
 C &= I\circ (AX^TXA^T) + I \cr
} $$
And let $u$ denote the vector of all ones. 
Then your objective function can be reformulated as 
$$ \eqalign {
 F(X) &= \|C^{-1}u- b\|^2 \cr
} $$
which is easier to work with.  For the first derivative, I obtained
$$ \eqalign {
 f(X) &= \frac {\partial F} {\partial X} \cr
   &= 4 XA^T \bigg[I\circ(C^{-1}bu^TC^{-1} - C^{-2}uu^TC^{-1})\bigg] A \cr
} $$
I'm not sure if this expression is equal to yours or not.
Update with details
Defining $M\!=\!(C^{-1}u\!-\!b)$ for later convenience, let's begin by stating a couple of well-known differentials
$$\eqalign{
 dC^{-1} &= -C^{-1}\,dC\,C^{-1} \cr
 d\,\|M\|^2 &= d\,(M:M) = 2\,M:dM  \cr
} $$
and employ them to find $dF$.  
$$\eqalign{
 dM &= dC^{-1}u  \cr
 dF &= 2\,M:dM  \cr
    &= 2\,Mu^T:dC^{-1}  \cr
    &= -2\,Mu^T : C^{-1}\,dC\,C^{-1}  \cr
    &= -2\,C^{-1}Mu^TC^{-1} : dC  \cr
} $$
The next thing we need is the differential of $C$: 
$$\eqalign{
  C &= I\circ\big(AX^TXA^T\big)+I \cr
  dC &= I\circ\big(A\,\,d(X^TX)\,\,A^T\big) \cr
     &= I\circ\big(A\,\,(dX^TX + X^TdX)\,\,A^T\big) \cr
     &= 2\,I\circ\big(A\,\,{\rm sym}(X^TdX)\,\,A^T\big) \cr
} $$
I should pause and define the Frobenius product and the sym() operator
$$\eqalign{
  A:B &= {\rm tr}(A^TB) \cr
 {\rm sym}(A) &= \frac {1} {2} (A+A^T) \cr
} $$
and a few of their more useful properties, especially with respect to the Hadamard product
$$\eqalign{
 (A\circ B):C &= A:(B\circ C) \cr
  A:{\rm sym}(B) &= {\rm sym}(A):B \cr
  {\rm sym}(I\circ A) &= I\circ A = A\circ I \cr
  A:BC &= B^TA:C \cr
       &= AC^T:B \cr
} $$
Now substitute $dC$ back into the expansion of $dF$.
$$\eqalign{
 dF &= -2\,C^{-1}Mu^TC^{-1} : \big(2\,I\circ(A\,\,{\rm sym}(X^TdX)\,\,A^T)\big)  \cr
    &= -4(\,C^{-1}Mu^TC^{-1})\circ I : \big(A\,\,{\rm sym}(X^TdX)\,\,A^T\big)  \cr
    &= -4\,A^T((\,C^{-1}Mu^TC^{-1})\circ I)A : {\rm sym}(X^TdX)  \cr
    &= -4\,A^T((\,C^{-1}Mu^TC^{-1})\circ I)A : (X^TdX)  \cr
    &= -4\,XA^T((\,C^{-1}Mu^TC^{-1})\circ I)A : dX  \cr
    &= \bigg(\frac {\partial F} {\partial X}\bigg) : dX  \cr
} $$
Update #2
It occurs to me that since $C$ is diagonal and $u$ is a vector of ones, we can get rid of the Hadamard product in the sub-expression 
$$ \eqalign {
  I\circ(C^{-1}Mu^TC^{-1}) &= C^{-1}(I\circ(Mu^T))C^{-1} \cr
       &= C^{-1}{\rm diag}(M)C^{-1} \cr
       &= C^{-1}(C^{-1}-B)C^{-1} \cr
       &= C^{-3}(I - BC) \cr
} $$
where $B = {\rm diag}(b)$, and all of the matrices commute since they're diagonal.  
Plugging this back into the derivative expression
$$ \eqalign {
  \frac {\partial F} {\partial X} &= 4\,XA^TC^{-3}(BC-I)A \cr
} $$
I now see that setting $C=(I+\Delta_X)$  recovers your original expression for the derivative.
$$
$$
Update #3:  The Second Derivative
What we have so far is 
$$ \eqalign {
  f &= 4\,XA^TC^{-3}(BC-I)A \cr
 dC &= I\circ(A(dX^TX+X^TdX)A^T) \cr
} $$
Since $f$ is a matrix, I prefer to denote it by $Y$.
Now the idea is to take the differential
$$ \eqalign {
 dY &= 4(dX)A^TC^{-3}(BC-I)A + 4XA^TC^{-3}B(dC)A \cr
  &+\, 4XA^T\bigg(dC^{-3}\bigg)(BC-I)A \cr
 \cr 
  &= 4(dX)A^TC^{-3}(BC-I)A + 4XA^TC^{-3}B(dC)A \cr
  &+\, 4XA^T\bigg(C^{-3}(dC)C^{-1} + C^{-2}(dC)C^{-2} + C^{-1}(dC)C^{-3}\bigg)(I-BC)A \cr
} $$
And apply the ${\rm vec}()$ operation
$$ \eqalign {
  {\rm vec}(AXB) &= (B^T\otimes A){\rm vec}(X) \cr
} $$
to isolate all of the $dC, dX$ terms as vectors.
$$ \eqalign {
 {\rm vec}(dY) 
  &= 4((A^TC^{-3}(BC-I)A)^T\otimes I){\rm vec}(dX) \cr
  &+ 4\big(A^T(I-BC)C^{-1}\otimes C^{-3}AX^T\big){\rm vec}(dC) \cr
  &+ 4\big(A^T(I-BC)C^{-2}\otimes C^{-2}AX^T\big){\rm vec}(dC) \cr
  &+ 4\big(A^T(I-BC)C^{-3}\otimes C^{-1}AX^T\big){\rm vec}(dC) \cr
  &+ 4(A^T\otimes XA^TC^{-3}B){\rm vec}(dC) \cr
} $$
For notational convenience, collect all the matrix coefficients and re-write this mess as a simple vector equation
$$ \eqalign {
 dy &= G\,dx + H\,dc \cr
} $$
But we're still not done, we need to express $dc$ in terms of $dx$.
However, we must be careful with notation since  ${\rm vec}(dX) = dx$, but ${\rm vec}(dX^T) = P\,dx \neq dx^T$. 

The $P$ matrix is a permutation which depends on the dimension of $dx$.
$$ \eqalign {
  dc &= {\rm vec}(I)\circ{\rm vec}(A\,d\!X^TXA^T+AX^Td\!XA^T) \cr 
     &= {\rm vec}(I)\circ(AX^T\otimes A)Pdx + {\rm vec}(I)\circ(A\otimes AX^T)\,dx \cr 
     &= {\rm diag}({\rm vec}(I))\,\bigg((AX^T\otimes A)P + (A\otimes AX^T)\bigg)\,dx \cr 
     &= K\,dx \cr 
} $$
Finally, we're ready to write down an expression for the second derivative
$$ \eqalign {
 \frac {\partial y} {\partial x} &= G + HK \cr
 &= \frac {\partial\,{\rm vec}(Y)} {\partial\,{\rm vec}(X)} \cr
} $$
