Proof of Descartes' theorem I came across the use of Descartes' theorem while solving a question.I searched it but I  could only find the theorem but not any proof.Even Wikipedia also, just states the theorem!!I want to know the procedure to find the radius of the Soddy Circle??
I apologize if its duplicate and to mention it is not a homework.
 A: Part I - Proof of Soddy-Gosset theorem (generalization of Descartes theorem).
For any integer $d \ge 2$, consider the problem of placing $n = d + 2$ hyper-spheres touching each other in $\mathbb{R}^d$. Let $\vec{x}_i \in \mathbb{R}^d$ and $R_i \in \mathbb{R}$ be the center and radius for the $i^{th}$ sphere. The condition for these spheres touching each other can be expressed as:
$$|\vec{x}_i - \vec{x}_j| = | R_i + R_j | \quad\text{ for }\quad 1 \le i < j \le n$$
or equivalently
$$|\vec{x}_i - \vec{x}_j|^2 = (R_i + R_j)^2 - 4R_iR_j\delta_{ij}\quad\text{ for }\quad 1 \le i, j \le n\tag{*1}$$
where $\;\delta_{ij} = \begin{cases}1,&i = j\\0,& i \ne j\end{cases}\;$ is the 
Kronecker delta.
Since the $n = d+2$ points $\vec{x}_i$ live in $\mathbb{R}^d$, the $d+1$ vectors
$\vec{x}_2 - \vec{x}_1, \vec{x}_3 - \vec{x}_1, \ldots, \vec{x}_n - \vec{x}_1$ are linear
depedent, this means we can find $n-1$ numbers $\beta_2, \beta_3, \ldots, \beta_n$ not all zero such that
$$\sum_{k=2}^n \beta_k (\vec{x}_k - \vec{x}_1 ) = \vec{0}$$
Let $\beta_1 = -(\beta_2 + \ldots + \beta_n)$, we can rewrite this relation in a more symmetric form:
$$
\sum_{k=1}^n \beta_k = 0
\quad\text{ and }\quad
\sum_{k=1}^n \beta_k \vec{x}_k = \vec{0}
\quad\text{ subject to some }\;\; \beta_k \ne 0
$$
If we fix $j$ in $(*1)$, multiple the $i^{th}$ term by $\beta_i$ and then sum over $i$, we get
$$
\sum_{i=1}^n\beta_i |\vec{x}_i|^2
= 
\sum_{i=1}^n\beta_i R_i^2
+ 2 \left( \sum_{i=1}^n \beta_i R_i \right) R_j
- 4R_j^2 \beta_j
$$
This leads to
$$4R_j^2 \beta_j =
2 A R_j  + B
\quad\text{ where }\quad\
\begin{cases}
A &= \sum\limits_{i=1}^n \beta_i R_i\\
B &= \sum\limits_{i=1}^n\beta_i ( R_i^2 - |\vec{x}_i|^2 )
\end{cases}
\tag{*2}
$$
Divide $(*2)$ by $R_j$ and sum over $j$, we get 
$$4A = 2nA + B\sum_{j=1}^n\frac{1}{R_j}\quad\iff\quad A = -\frac{B}{2d}\sum_{j=1}^n\frac{1}{R_j}\tag{*3}$$
A consequence of this is $B$ cannot vanish. Otherwise $B = 0 \implies A = 0$ and $(*2)$ implies all $\beta_j = 0$ which is clearly isn't the case.
Divide $(*2)$ by $R_j^2$ and sum over $j$, we get
$$0 = 4\sum_{j=1}^n \beta_j = 2A\sum_{j=1}^n\frac{1}{R_j} + B\sum_{j=1}^n\frac{1}{R_j^2}$$
Combine with $(*3)$, the RHS becomes
$$B \left( \sum_{j=1}^n \frac{1}{R_j^2} - \frac{1}{d}\left( \sum_{j=1}^n\frac{1}{R_j}\right)^2\right) = 0
\quad\iff\quad
\left( \sum_{j=1}^n\frac{1}{R_j}\right)^2 = d\sum_{j=1}^n \frac{1}{R_j^2}\tag{*4}
$$
The RHS of $(*4)$ is sometimes called Soddy-Gosset theorem. When $d = 2$, it reduces to the Descartes four circle theorem, the theorem we wish to prove:
$$\left( \frac{1}{R_1} + \frac{1}{R_2} + \frac{1}{R_3} + \frac{1}{R_4} \right)^2
= 
2 \left( \frac{1}{R_1^2} + \frac{1}{R_2^2} + \frac{1}{R_3^2} + \frac{1}{R_4^2} \right)
$$
Part II - Construction of inner/outer Soddy hyper-spheres
There is an interesting side-product of the proof in Part I. $\beta_k$ is determined up to an overall scaling factor. If we normalize $\beta_k$ such that $B = 4$, $(*2)$ and $(*3)$ together allows us to derive an explicit expression for $\beta_j$
$$\beta_j = \frac{1}{R_j^2} - \left(\frac{1}{d} \sum_{k=1}^n \frac{1}{R_k}\right)\frac{1}{R_j}\tag{*5}$$
We can use this relation to construct the inner and outer Soddy hyper-spheres.
Assume we already have $n-1 = d+1$ hyper-spheres touching among themselves.
The inner Soddy hyper-sphere is the sphere outside all these $n-1$ spheres 
and yet touching all of them.
Let $\vec{x}_k$ and $r_k$ be the center and radius for the $k^{th}$ hyper-sphere for
$1 \le k < n$.
Let $\vec{x}_{in}$ and $r_{in}$ be the center and radius of the inner Soddy hyper-sphere. If we let 
$$\vec{x}_n = \vec{x}_{in}\quad\text{ and }\quad R_k = \begin{cases}r_k,& 1 \le k < n\\ r_{in},& k = n\end{cases}$$
discussions in Part I tell us
$$\left( \frac{1}{r_{in}} + \sum_{k=1}^{n-1} \frac{1}{r_k} \right)^2 = d 
\left( \frac{1}{r_{in}^2} + \sum_{k=1}^{n-1} \frac{1}{r_k^2} \right)\tag{*6a}$$
We can use this to determine $r_{in}$. If the $n-1$ points $\vec{x}_1, \vec{x}_2, \ldots, \vec{x}_{n-1}$ are in general position, i.e. they are vertices of a non-degenerate $d$-simplex, the $d$ vectors $\vec{x}_2 - \vec{x}_1, \ldots, \vec{x}_{n-1} - \vec{x}_1$ will be linearly independent. This implies there exists $d$ coefficients $\gamma_2, \gamma_3, \ldots, \gamma_{n-1}$ such that
$$\vec{x}_{in} - \vec{x}_1 = \gamma_2 (\vec{x}_2 - \vec{x}_1) + \ldots + \gamma_{n-1} ( \vec{x}_{n-1} - \vec{x}_1 )$$
A consequence of this is $\beta_n \ne 0$. This means we can use $(*5)$ and the relation
$\sum\limits_{k=1}^n \beta_k  \vec{x}_k = \vec{0}$ to compute the center $\vec{x}_{in}$
of the inner Soddy hyper-sphere.
For the outer Soddy hyper-sphere. It is a sphere that contains the original $n-1$ spheres and touching each of them. 
Let $\vec{x}_{out}$ and $r_{out}$ be the center and radius of the outer Soddy hyper-sphere. The touching condition now takes the form:
$$\begin{array}{ccccl}
|\vec{x}_{out} - \vec{x}_j | &=& | r_{out} - r_j |\quad & \text{ for }\quad & 1 \le j < n\\
|\vec{x}_i - \vec{x}_j     | &=& | r_i + r_j | \quad & \text{ for }\quad & 1 \le i < j < n
\end{array}
$$
Once again, if we let 
$$\vec{x}_n = \vec{x}_{out}\quad\text{ and }\quad R_k = \begin{cases}r_k,& 1 \le k < n\\ -r_{out},& k = n\end{cases}$$
we can repeat discussions in Part I to obtain
$$\left( -\frac{1}{r_{out}} + \sum_{k=1}^{n-1} \frac{1}{r_k} \right)^2 = d 
\left( \frac{1}{r_{out}^2} + \sum_{k=1}^{n-1} \frac{1}{r_k^2} \right)\tag{*6b}$$
We can use this to determine $r_{out}$. Once again, if $\vec{x}_1,\ldots,\vec{x}_{n-1}$ are in general position, we will find $\beta_n \ne 0$. As a result, we can use $(*5)$ to compute $\vec{x}_{out}$, the center of the outer Soddy hyper-spheres, from the remaining centers.
If one compare $(*6a)$ and $(*6b)$, they are very similar, $r_{in}$ and $-r_{out}$ are the two roots of the same equation in $R$.
$$\left( \frac{1}{R} + \sum_{k=1}^{n-1} \frac{1}{r_k} \right)^2 = d 
\left( \frac{1}{R^2} + \sum_{k=1}^{n-1} \frac{1}{r_k^2} \right)$$
If the two roots of this equation has different sign, the positive root will be the inner Soddy radius $r_{in}$, the negative root will be $-r_{out}$, the negative of the outer Soddy radius.
