The computation of the circumsphere center is rather easy (and IMHO the following is simpler than the current proposals under [this thread](https://math.stackexchange.com/q/2585392)). Since (for good reasons, of course) equations do not look good when produced with the text field, the derivation is coded with LaTeX. The TeX implementation of the `3d/circumsphere center` key, which was added for the purpose of this answer, can be found in the [3d tools source](https://github.com/marmotghost/tikz-3dtools).
![Screen Shot 2020-08-22 at 12.29.59 AM.png](/image?hash=d8ee0b033ebfbfafa9fa754978d9b208bb5dc5362736dc3069136fcbf83894e4)
```
\documentclass[fleqn]{article}
\usepackage[margin=1in]{geometry}
\usepackage{amsmath}
\usepackage{tikz}
\usetikzlibrary{3dtools}% https://github.com/marmotghost/tikz-3dtools
\begin{document}
Given four noncoplanar points $A$, $B$, $C$ and $D$, the center $I$ of a sphere
through these points is at the intersection of three planes in the middle
between two given points. For instance, the plane in the middle between $A$ and
$B$ has the normal
\[ \vec n_{AB}=\vec B-\vec A=:\vec B\,'\]
and has the equation
\[ \vec n_{AB}\cdot\vec
I=\vec n_{AB}\cdot(\vec A+\vec B)/2=:b_{AB}\;.\]
To simplify the discussion a bit, it helps to shift the coordinate system
so that the origin is $\vec A$. In these coordinates
\[\vec A\,'=0\;,\quad
\vec B\,'=\vec B-\vec A\;,\quad
\vec C\,'=\vec C-\vec A\quad\text{and}\quad
\vec D\,'=\vec D-\vec A\;.
\]
So the center can be obtained by solving a simple matrix equation
\[ N\, \vec I\,'=\vec b
\]
with $N=(\vec B\,',\vec C\,',\vec D\,')^T$ and $\vec
b=(|\vec B\,'|^2/2,|\vec C\,'|^2/2,|\vec D\,'|^2/2)$. The matrix is
invertible if $\det N\ne0$, where
\[
\det N=\vec B\,'\cdot\left(\vec C\,'\times\vec D\,'\right)
=\vec D\,'\cdot\left(\vec B\,'\times\vec C\,'\right)
=\vec C\,'\cdot\left(\vec D\,'\times\vec B\,'\right)\;.
\]
The condition that the points be non coplanar is equivalent to $\det
N\ne0$. This matrix equation can be solved by projecting on $\vec
C\,'\times\vec D\,'$, $\vec D\,'\times\vec B\,'$ and $\vec B\,'\times\vec
C\,'$. The result is
\[
\vec I\,'
=\frac{1}{2\det N}
\left(
|\vec B\,'|^2\,\vec C\,'\times\vec D\,'
+|\vec C\,'|^2\,\vec D\,'\times\vec B\,'
+|\vec D\,'|^2\,\vec B\,'\times\vec C\,'
\right)\;.
\]
In the end we have to shift back to the original coordinate system,
\[
\vec I=\vec I\,'+\vec A\;.
\]
An example is shown in figure~\ref{fig:circumsphere}.
\begin{figure}[h]
\centering
\begin{tikzpicture}[declare function={a=3;},
dot/.style={circle,inner sep=1pt,fill}]
\begin{scope}[3d/install view={phi=100,psi=0,theta=70}]
\path (a,0,0) coordinate (A)
(0,a,0) coordinate (B)
(0,0,a) coordinate (C)
(0,0,0) coordinate (D);
\path[3d/circumsphere center] coordinate (I);
\end{scope}
% radius
\pgfmathsetmacro{\csr}{sqrt(TD("(A)-(I)o(A)-(I)"))}
\draw (I) circle[radius=\csr];
\path foreach \X in {A,B,C,D}
{(\X) coordinate[dot,label=above:{$\X$}]}
(I) (I) coordinate[dot,label=above:{$I=\pgfmathparse{TD("(I)")}%
(\pgfmathprintvector\pgfmathresult)^T$}];
\end{tikzpicture}
\caption{Example of a circumsphere computation.}
\label{fig:circumsphere}
\end{figure}
\end{document}
```
**OLD ANSWER**: This is a possible way of definining the coordinates via arrays, and making things work with `fpu`. The resulting complexity is a real mess, in particular compared to what one can do with `3dtools`. However, developing the latter further only makes real sense if there is a committment by the maintainers of Ti*k*Z not to randomly drop some internal routines that allow us to keep track of the 3d coordinates.
On the bright side, your computation seems to get confirmed, and the following code provides one with the computation of the center of sphere of four points which do not lie all in one common plane. It also has a discussion of one possible way to compute this center, which the code is based on.
```
\documentclass[border=2mm,12pt,tikz,fleqn]{standalone}
\usepackage{tikz-3dplot-circleofsphere}
\usetikzlibrary{fpu}
\makeatletter
\pgfmathdeclarefunction{detN}{4}{%
\begingroup
\pgfkeys{/pgf/fpu=false}%
\pgfmathsetmacro\pgfutil@tempxA{{#1}[0]}%
\pgfmathsetmacro\pgfutil@tempyA{{#1}[1]}%
\pgfmathsetmacro\pgfutil@tempzA{{#1}[2]}%
\pgfmathsetmacro\pgfutil@tempxB{{#2}[0]}%
\pgfmathsetmacro\pgfutil@tempyB{{#2}[1]}%
\pgfmathsetmacro\pgfutil@tempzB{{#2}[2]}%
\pgfmathsetmacro\pgfutil@tempxC{{#3}[0]}%
\pgfmathsetmacro\pgfutil@tempyC{{#3}[1]}%
\pgfmathsetmacro\pgfutil@tempzC{{#3}[2]}%
\pgfmathsetmacro\pgfutil@tempxD{{#4}[0]}%
\pgfmathsetmacro\pgfutil@tempyD{{#4}[1]}%
\pgfmathsetmacro\pgfutil@tempzD{{#4}[2]}%
\pgfkeys{/pgf/fpu,/pgf/fpu/output format=fixed}%
\pgfmathparse{-(\pgfutil@tempxB*\pgfutil@tempyC*\pgfutil@tempzA) + \pgfutil@tempxB*\pgfutil@tempyD*\pgfutil@tempzA + \pgfutil@tempxA*\pgfutil@tempyC*\pgfutil@tempzB - \pgfutil@tempxA*\pgfutil@tempyD*\pgfutil@tempzB +
\pgfutil@tempxB*\pgfutil@tempyA*\pgfutil@tempzC - \pgfutil@tempxA*\pgfutil@tempyB*\pgfutil@tempzC + \pgfutil@tempxA*\pgfutil@tempyD*\pgfutil@tempzC - \pgfutil@tempxB*\pgfutil@tempyD*\pgfutil@tempzC +
\pgfutil@tempxD*(-(\pgfutil@tempyB*\pgfutil@tempzA) + \pgfutil@tempyC*\pgfutil@tempzA + \pgfutil@tempyA*\pgfutil@tempzB - \pgfutil@tempyC*\pgfutil@tempzB - \pgfutil@tempyA*\pgfutil@tempzC +
\pgfutil@tempyB*\pgfutil@tempzC) - \pgfutil@tempxB*\pgfutil@tempyA*\pgfutil@tempzD + \pgfutil@tempxA*\pgfutil@tempyB*\pgfutil@tempzD - \pgfutil@tempxA*\pgfutil@tempyC*\pgfutil@tempzD +
\pgfutil@tempxB*\pgfutil@tempyC*\pgfutil@tempzD + \pgfutil@tempxC*(-(\pgfutil@tempyD*\pgfutil@tempzA) - \pgfutil@tempyA*\pgfutil@tempzB + \pgfutil@tempyD*\pgfutil@tempzB +
\pgfutil@tempyB*(\pgfutil@tempzA - \pgfutil@tempzD) + \pgfutil@tempyA*\pgfutil@tempzD)}%
\pgfmathsmuggle\pgfmathresult\endgroup}%
\pgfmathdeclarefunction{fourpointcenterx}{4}{%
\begingroup
\pgfkeys{/pgf/fpu=false}%
\pgfmathsetmacro\pgfutil@tempxA{{#1}[0]}%
\pgfmathsetmacro\pgfutil@tempyA{{#1}[1]}%
\pgfmathsetmacro\pgfutil@tempzA{{#1}[2]}%
\pgfmathsetmacro\pgfutil@tempxB{{#2}[0]}%
\pgfmathsetmacro\pgfutil@tempyB{{#2}[1]}%
\pgfmathsetmacro\pgfutil@tempzB{{#2}[2]}%
\pgfmathsetmacro\pgfutil@tempxC{{#3}[0]}%
\pgfmathsetmacro\pgfutil@tempyC{{#3}[1]}%
\pgfmathsetmacro\pgfutil@tempzC{{#3}[2]}%
\pgfmathsetmacro\pgfutil@tempxD{{#4}[0]}%
\pgfmathsetmacro\pgfutil@tempyD{{#4}[1]}%
\pgfmathsetmacro\pgfutil@tempzD{{#4}[2]}%
\pgfmathsetmacro\pgfutil@tempe{detN("#1","#2","#3","#4")}%
\pgfmathsetmacro\pgfutil@tempf{abs(\pgfutil@tempe)>0}%
\ifdim\pgfutil@tempf pt>0pt\relax
\pgfkeys{/pgf/fpu,/pgf/fpu/output format=fixed}%
\pgfmathparse{(-(\pgfutil@tempxB^2*\pgfutil@tempyC*\pgfutil@tempzA) - \pgfutil@tempyA*\pgfutil@tempyB*\pgfutil@tempyC*\pgfutil@tempzA - \pgfutil@tempyB^2*\pgfutil@tempyC*\pgfutil@tempzA +
\pgfutil@tempyA*\pgfutil@tempyC^2*\pgfutil@tempzA + \pgfutil@tempyB*\pgfutil@tempyC^2*\pgfutil@tempzA + \pgfutil@tempxB^2*\pgfutil@tempyD*\pgfutil@tempzA +
\pgfutil@tempyA*\pgfutil@tempyB*\pgfutil@tempyD*\pgfutil@tempzA + \pgfutil@tempyB^2*\pgfutil@tempyD*\pgfutil@tempzA - \pgfutil@tempyA*\pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzA -
\pgfutil@tempyB*\pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzA + 2*\pgfutil@tempyA^2*\pgfutil@tempyC*\pgfutil@tempzB - 2*\pgfutil@tempyA*\pgfutil@tempyC^2*\pgfutil@tempzB -
2*\pgfutil@tempyA^2*\pgfutil@tempyD*\pgfutil@tempzB + 2*\pgfutil@tempyA*\pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzB + \pgfutil@tempyC*\pgfutil@tempzA^2*\pgfutil@tempzB -
\pgfutil@tempyD*\pgfutil@tempzA^2*\pgfutil@tempzB - \pgfutil@tempyC*\pgfutil@tempzA*\pgfutil@tempzB^2 + \pgfutil@tempyD*\pgfutil@tempzA*\pgfutil@tempzB^2 +
\pgfutil@tempxB^2*\pgfutil@tempyA*\pgfutil@tempzC - \pgfutil@tempyA^2*\pgfutil@tempyB*\pgfutil@tempzC + \pgfutil@tempyA*\pgfutil@tempyB^2*\pgfutil@tempzC -
\pgfutil@tempyA^2*\pgfutil@tempyC*\pgfutil@tempzC + \pgfutil@tempyA*\pgfutil@tempyB*\pgfutil@tempyC*\pgfutil@tempzC - \pgfutil@tempxB^2*\pgfutil@tempyD*\pgfutil@tempzC +
2*\pgfutil@tempyA^2*\pgfutil@tempyD*\pgfutil@tempzC - \pgfutil@tempyA*\pgfutil@tempyB*\pgfutil@tempyD*\pgfutil@tempzC - \pgfutil@tempyB^2*\pgfutil@tempyD*\pgfutil@tempzC -
\pgfutil@tempyA*\pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzC + \pgfutil@tempyB*\pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzC - 2*\pgfutil@tempyB*\pgfutil@tempzA^2*\pgfutil@tempzC +
\pgfutil@tempyC*\pgfutil@tempzA^2*\pgfutil@tempzC + \pgfutil@tempyD*\pgfutil@tempzA^2*\pgfutil@tempzC + \pgfutil@tempyA*\pgfutil@tempzA*\pgfutil@tempzB*\pgfutil@tempzC -
\pgfutil@tempyC*\pgfutil@tempzA*\pgfutil@tempzB*\pgfutil@tempzC + \pgfutil@tempyA*\pgfutil@tempzB^2*\pgfutil@tempzC - \pgfutil@tempyD*\pgfutil@tempzB^2*\pgfutil@tempzC -
\pgfutil@tempyA*\pgfutil@tempzA*\pgfutil@tempzC^2 + 2*\pgfutil@tempyB*\pgfutil@tempzA*\pgfutil@tempzC^2 - \pgfutil@tempyD*\pgfutil@tempzA*\pgfutil@tempzC^2 -
\pgfutil@tempyA*\pgfutil@tempzB*\pgfutil@tempzC^2 + \pgfutil@tempyD*\pgfutil@tempzB*\pgfutil@tempzC^2 +
\pgfutil@tempxC*\pgfutil@tempxD*(-(\pgfutil@tempyB*\pgfutil@tempzA) + \pgfutil@tempyC*\pgfutil@tempzA + \pgfutil@tempyA*\pgfutil@tempzB - \pgfutil@tempyC*\pgfutil@tempzB - \pgfutil@tempyA*\pgfutil@tempzC +
\pgfutil@tempyB*\pgfutil@tempzC) + \pgfutil@tempxA*(\pgfutil@tempxC - \pgfutil@tempxD)*(\pgfutil@tempyC*\pgfutil@tempzA + \pgfutil@tempyA*\pgfutil@tempzB - \pgfutil@tempyC*\pgfutil@tempzB -
\pgfutil@tempyA*\pgfutil@tempzC + \pgfutil@tempyB*(-\pgfutil@tempzA + \pgfutil@tempzC)) - \pgfutil@tempxB^2*\pgfutil@tempyA*\pgfutil@tempzD +
\pgfutil@tempyA^2*\pgfutil@tempyB*\pgfutil@tempzD - \pgfutil@tempyA*\pgfutil@tempyB^2*\pgfutil@tempzD + \pgfutil@tempxB^2*\pgfutil@tempyC*\pgfutil@tempzD -
\pgfutil@tempyA^2*\pgfutil@tempyC*\pgfutil@tempzD + \pgfutil@tempyB^2*\pgfutil@tempyC*\pgfutil@tempzD + \pgfutil@tempyA*\pgfutil@tempyC^2*\pgfutil@tempzD -
\pgfutil@tempyB*\pgfutil@tempyC^2*\pgfutil@tempzD + 2*\pgfutil@tempyB*\pgfutil@tempzA^2*\pgfutil@tempzD - 2*\pgfutil@tempyC*\pgfutil@tempzA^2*\pgfutil@tempzD -
\pgfutil@tempyA*\pgfutil@tempzA*\pgfutil@tempzB*\pgfutil@tempzD + \pgfutil@tempyC*\pgfutil@tempzA*\pgfutil@tempzB*\pgfutil@tempzD - \pgfutil@tempyA*\pgfutil@tempzB^2*\pgfutil@tempzD +
\pgfutil@tempyC*\pgfutil@tempzB^2*\pgfutil@tempzD + \pgfutil@tempyA*\pgfutil@tempzA*\pgfutil@tempzC*\pgfutil@tempzD - 2*\pgfutil@tempyB*\pgfutil@tempzA*\pgfutil@tempzC*\pgfutil@tempzD +
\pgfutil@tempyC*\pgfutil@tempzA*\pgfutil@tempzC*\pgfutil@tempzD + \pgfutil@tempyA*\pgfutil@tempzB*\pgfutil@tempzC*\pgfutil@tempzD - \pgfutil@tempyC*\pgfutil@tempzB*\pgfutil@tempzC*\pgfutil@tempzD +
\pgfutil@tempxC^2*(\pgfutil@tempyB*\pgfutil@tempzA - \pgfutil@tempyD*\pgfutil@tempzA - \pgfutil@tempyA*\pgfutil@tempzB + \pgfutil@tempyD*\pgfutil@tempzB + \pgfutil@tempyA*\pgfutil@tempzD -
\pgfutil@tempyB*\pgfutil@tempzD) + \pgfutil@tempxA^2*(-(\pgfutil@tempyD*\pgfutil@tempzB) - \pgfutil@tempyB*\pgfutil@tempzC + \pgfutil@tempyD*\pgfutil@tempzC +
\pgfutil@tempyC*(\pgfutil@tempzB - \pgfutil@tempzD) + \pgfutil@tempyB*\pgfutil@tempzD))/(2*\pgfutil@tempe)}%
\else
\edef\pgfmathresult{0}%
\fi
\pgfmathsmuggle\pgfmathresult\endgroup}%
\pgfmathdeclarefunction{fourpointcentery}{4}{%
\begingroup
\pgfkeys{/pgf/fpu=false}%
\pgfmathsetmacro\pgfutil@tempxA{{#1}[0]}%
\pgfmathsetmacro\pgfutil@tempyA{{#1}[1]}%
\pgfmathsetmacro\pgfutil@tempzA{{#1}[2]}%
\pgfmathsetmacro\pgfutil@tempxB{{#2}[0]}%
\pgfmathsetmacro\pgfutil@tempyB{{#2}[1]}%
\pgfmathsetmacro\pgfutil@tempzB{{#2}[2]}%
\pgfmathsetmacro\pgfutil@tempxC{{#3}[0]}%
\pgfmathsetmacro\pgfutil@tempyC{{#3}[1]}%
\pgfmathsetmacro\pgfutil@tempzC{{#3}[2]}%
\pgfmathsetmacro\pgfutil@tempxD{{#4}[0]}%
\pgfmathsetmacro\pgfutil@tempyD{{#4}[1]}%
\pgfmathsetmacro\pgfutil@tempzD{{#4}[2]}%
\pgfmathsetmacro\pgfutil@tempe{detN("#1","#2","#3","#4")}%
\pgfmathsetmacro\pgfutil@tempf{abs(\pgfutil@tempe)>0}%
\ifdim\pgfutil@tempf pt>0pt\relax
\pgfkeys{/pgf/fpu,/pgf/fpu/output format=fixed}%
\pgfmathparse{(\pgfutil@tempxC*\pgfutil@tempyB^2*\pgfutil@tempzA - \pgfutil@tempxD*\pgfutil@tempyB^2*\pgfutil@tempzA - \pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempyC*\pgfutil@tempzA +
\pgfutil@tempxD*\pgfutil@tempyC^2*\pgfutil@tempzA + \pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempyD*\pgfutil@tempzA - \pgfutil@tempxC*\pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzA -
\pgfutil@tempxC*\pgfutil@tempyA^2*\pgfutil@tempzB + \pgfutil@tempxD*\pgfutil@tempyA^2*\pgfutil@tempzB + \pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempyC*\pgfutil@tempzB -
\pgfutil@tempxD*\pgfutil@tempyC^2*\pgfutil@tempzB - \pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempyD*\pgfutil@tempzB + \pgfutil@tempxC*\pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzB -
\pgfutil@tempxC*\pgfutil@tempzA^2*\pgfutil@tempzB + \pgfutil@tempxD*\pgfutil@tempzA^2*\pgfutil@tempzB + \pgfutil@tempxC*\pgfutil@tempzA*\pgfutil@tempzB^2 -
\pgfutil@tempxD*\pgfutil@tempzA*\pgfutil@tempzB^2 - \pgfutil@tempxD*\pgfutil@tempyA^2*\pgfutil@tempzC + \pgfutil@tempxD*\pgfutil@tempyB^2*\pgfutil@tempzC -
\pgfutil@tempxC*\pgfutil@tempzA^2*\pgfutil@tempzC - \pgfutil@tempxD*\pgfutil@tempzA^2*\pgfutil@tempzC + \pgfutil@tempxC*\pgfutil@tempzA*\pgfutil@tempzB*\pgfutil@tempzC +
\pgfutil@tempxD*\pgfutil@tempzB^2*\pgfutil@tempzC + \pgfutil@tempxD*\pgfutil@tempzA*\pgfutil@tempzC^2 - \pgfutil@tempxD*\pgfutil@tempzB*\pgfutil@tempzC^2 +
\pgfutil@tempxB^2*(\pgfutil@tempxD*(-\pgfutil@tempzA + \pgfutil@tempzC) + \pgfutil@tempxC*(\pgfutil@tempzA - \pgfutil@tempzD)) + \pgfutil@tempxC*\pgfutil@tempyA^2*\pgfutil@tempzD -
\pgfutil@tempxC*\pgfutil@tempyB^2*\pgfutil@tempzD + 2*\pgfutil@tempxC*\pgfutil@tempzA^2*\pgfutil@tempzD - \pgfutil@tempxC*\pgfutil@tempzA*\pgfutil@tempzB*\pgfutil@tempzD -
\pgfutil@tempxC*\pgfutil@tempzB^2*\pgfutil@tempzD - \pgfutil@tempxC*\pgfutil@tempzA*\pgfutil@tempzC*\pgfutil@tempzD + \pgfutil@tempxC*\pgfutil@tempzB*\pgfutil@tempzC*\pgfutil@tempzD +
\pgfutil@tempxB*(-(\pgfutil@tempyC^2*\pgfutil@tempzA) + \pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzA + \pgfutil@tempxC*\pgfutil@tempxD*(\pgfutil@tempzA - \pgfutil@tempzC) +
\pgfutil@tempyA*(\pgfutil@tempyC - \pgfutil@tempyD)*(\pgfutil@tempzA - \pgfutil@tempzC) - \pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzC + 2*\pgfutil@tempzA^2*\pgfutil@tempzC -
2*\pgfutil@tempzA*\pgfutil@tempzC^2 + \pgfutil@tempyA^2*(\pgfutil@tempzC - \pgfutil@tempzD) + \pgfutil@tempyC^2*\pgfutil@tempzD -
2*\pgfutil@tempzA^2*\pgfutil@tempzD + 2*\pgfutil@tempzA*\pgfutil@tempzC*\pgfutil@tempzD + \pgfutil@tempxC^2*(-\pgfutil@tempzA + \pgfutil@tempzD)) +
\pgfutil@tempxA*(-(\pgfutil@tempyA*\pgfutil@tempyC*\pgfutil@tempzB) + \pgfutil@tempyC^2*\pgfutil@tempzB + \pgfutil@tempyA*\pgfutil@tempyD*\pgfutil@tempzB - \pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzB +
\pgfutil@tempxB*(\pgfutil@tempxC - \pgfutil@tempxD)*(\pgfutil@tempzA - \pgfutil@tempzC) - \pgfutil@tempyB^2*\pgfutil@tempzC + \pgfutil@tempyA*\pgfutil@tempyC*\pgfutil@tempzC -
\pgfutil@tempyA*\pgfutil@tempyD*\pgfutil@tempzC + \pgfutil@tempyC*\pgfutil@tempyD*\pgfutil@tempzC - \pgfutil@tempzA*\pgfutil@tempzB*\pgfutil@tempzC - \pgfutil@tempzB^2*\pgfutil@tempzC +
\pgfutil@tempzA*\pgfutil@tempzC^2 + \pgfutil@tempzB*\pgfutil@tempzC^2 + \pgfutil@tempxC*\pgfutil@tempxD*(\pgfutil@tempzA - 2*\pgfutil@tempzB + \pgfutil@tempzC) +
\pgfutil@tempyB^2*\pgfutil@tempzD - \pgfutil@tempyC^2*\pgfutil@tempzD + \pgfutil@tempzA*\pgfutil@tempzB*\pgfutil@tempzD + \pgfutil@tempzB^2*\pgfutil@tempzD -
\pgfutil@tempzA*\pgfutil@tempzC*\pgfutil@tempzD - \pgfutil@tempzB*\pgfutil@tempzC*\pgfutil@tempzD - \pgfutil@tempxC^2*(\pgfutil@tempzA - 2*\pgfutil@tempzB + \pgfutil@tempzD) +
\pgfutil@tempxB^2*(-\pgfutil@tempzC + \pgfutil@tempzD)) + \pgfutil@tempxA^2*(2*\pgfutil@tempxD*(\pgfutil@tempzB - \pgfutil@tempzC) +
\pgfutil@tempxB*(\pgfutil@tempzC - \pgfutil@tempzD) + \pgfutil@tempxC*(-2*\pgfutil@tempzB + \pgfutil@tempzC + \pgfutil@tempzD)))/(2*\pgfutil@tempe)}%
\else
\edef\pgfmathresult{0}%
\fi
\pgfmathsmuggle\pgfmathresult\endgroup}%
\pgfmathdeclarefunction{fourpointcenterz}{4}{%
\begingroup
\pgfkeys{/pgf/fpu=false}%
\pgfmathsetmacro\pgfutil@tempxA{{#1}[0]}%
\pgfmathsetmacro\pgfutil@tempyA{{#1}[1]}%
\pgfmathsetmacro\pgfutil@tempzA{{#1}[2]}%
\pgfmathsetmacro\pgfutil@tempxB{{#2}[0]}%
\pgfmathsetmacro\pgfutil@tempyB{{#2}[1]}%
\pgfmathsetmacro\pgfutil@tempzB{{#2}[2]}%
\pgfmathsetmacro\pgfutil@tempxC{{#3}[0]}%
\pgfmathsetmacro\pgfutil@tempyC{{#3}[1]}%
\pgfmathsetmacro\pgfutil@tempzC{{#3}[2]}%
\pgfmathsetmacro\pgfutil@tempxD{{#4}[0]}%
\pgfmathsetmacro\pgfutil@tempyD{{#4}[1]}%
\pgfmathsetmacro\pgfutil@tempzD{{#4}[2]}%
\pgfmathsetmacro\pgfutil@tempe{detN("#1","#2","#3","#4")}%
\pgfmathsetmacro\pgfutil@tempf{abs(\pgfutil@tempe)>0}%
\ifdim\pgfutil@tempf pt>0pt\relax
\pgfkeys{/pgf/fpu,/pgf/fpu/output format=fixed}%
\pgfmathparse{(\pgfutil@tempxC*\pgfutil@tempyA^2*\pgfutil@tempyB - \pgfutil@tempxD*\pgfutil@tempyA^2*\pgfutil@tempyB - \pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempyB^2 + \pgfutil@tempxD*\pgfutil@tempyA*\pgfutil@tempyB^2 +
\pgfutil@tempxC*\pgfutil@tempyA^2*\pgfutil@tempyC + \pgfutil@tempxD*\pgfutil@tempyA^2*\pgfutil@tempyC - \pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempyB*\pgfutil@tempyC -
\pgfutil@tempxD*\pgfutil@tempyB^2*\pgfutil@tempyC - \pgfutil@tempxD*\pgfutil@tempyA*\pgfutil@tempyC^2 + \pgfutil@tempxD*\pgfutil@tempyB*\pgfutil@tempyC^2 -
2*\pgfutil@tempxC*\pgfutil@tempyA^2*\pgfutil@tempyD + \pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempyB*\pgfutil@tempyD + \pgfutil@tempxC*\pgfutil@tempyB^2*\pgfutil@tempyD +
\pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempyC*\pgfutil@tempyD - \pgfutil@tempxC*\pgfutil@tempyB*\pgfutil@tempyC*\pgfutil@tempyD +
\pgfutil@tempxB^2*(\pgfutil@tempxD*(\pgfutil@tempyA - \pgfutil@tempyC) + \pgfutil@tempxC*(-\pgfutil@tempyA + \pgfutil@tempyD)) -
\pgfutil@tempxA^2*(2*\pgfutil@tempxD*(\pgfutil@tempyB - \pgfutil@tempyC) + \pgfutil@tempxB*(\pgfutil@tempyC - \pgfutil@tempyD) +
\pgfutil@tempxC*(-2*\pgfutil@tempyB + \pgfutil@tempyC + \pgfutil@tempyD)) + \pgfutil@tempxC*\pgfutil@tempyB*\pgfutil@tempzA^2 -
\pgfutil@tempxD*\pgfutil@tempyB*\pgfutil@tempzA^2 + \pgfutil@tempxD*\pgfutil@tempyC*\pgfutil@tempzA^2 - \pgfutil@tempxC*\pgfutil@tempyD*\pgfutil@tempzA^2 -
\pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempzB^2 + \pgfutil@tempxD*\pgfutil@tempyA*\pgfutil@tempzB^2 - \pgfutil@tempxD*\pgfutil@tempyC*\pgfutil@tempzB^2 +
\pgfutil@tempxC*\pgfutil@tempyD*\pgfutil@tempzB^2 + \pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempzA*\pgfutil@tempzC - \pgfutil@tempxC*\pgfutil@tempyB*\pgfutil@tempzA*\pgfutil@tempzC -
\pgfutil@tempxD*\pgfutil@tempyA*\pgfutil@tempzC^2 + \pgfutil@tempxD*\pgfutil@tempyB*\pgfutil@tempzC^2 - \pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempzA*\pgfutil@tempzD +
\pgfutil@tempxC*\pgfutil@tempyB*\pgfutil@tempzA*\pgfutil@tempzD + \pgfutil@tempxC*\pgfutil@tempyA*\pgfutil@tempzC*\pgfutil@tempzD - \pgfutil@tempxC*\pgfutil@tempyB*\pgfutil@tempzC*\pgfutil@tempzD +
\pgfutil@tempxA*(-(\pgfutil@tempxB*(\pgfutil@tempxC - \pgfutil@tempxD)*(\pgfutil@tempyA - \pgfutil@tempyC)) + \pgfutil@tempyA*\pgfutil@tempyB*\pgfutil@tempyC +
\pgfutil@tempyB^2*\pgfutil@tempyC - \pgfutil@tempyA*\pgfutil@tempyC^2 - \pgfutil@tempyB*\pgfutil@tempyC^2 -
\pgfutil@tempxC*\pgfutil@tempxD*(\pgfutil@tempyA - 2*\pgfutil@tempyB + \pgfutil@tempyC) + \pgfutil@tempxB^2*(\pgfutil@tempyC - \pgfutil@tempyD) -
\pgfutil@tempyA*\pgfutil@tempyB*\pgfutil@tempyD - \pgfutil@tempyB^2*\pgfutil@tempyD + \pgfutil@tempyA*\pgfutil@tempyC*\pgfutil@tempyD + \pgfutil@tempyB*\pgfutil@tempyC*\pgfutil@tempyD +
\pgfutil@tempxC^2*(\pgfutil@tempyA - 2*\pgfutil@tempyB + \pgfutil@tempyD) + \pgfutil@tempyC*\pgfutil@tempzB^2 - \pgfutil@tempyD*\pgfutil@tempzB^2 +
\pgfutil@tempyB*\pgfutil@tempzA*\pgfutil@tempzC - \pgfutil@tempyC*\pgfutil@tempzA*\pgfutil@tempzC - \pgfutil@tempyB*\pgfutil@tempzC^2 + \pgfutil@tempyD*\pgfutil@tempzC^2 -
\pgfutil@tempyB*\pgfutil@tempzA*\pgfutil@tempzD + \pgfutil@tempyC*\pgfutil@tempzA*\pgfutil@tempzD + \pgfutil@tempyB*\pgfutil@tempzC*\pgfutil@tempzD - \pgfutil@tempyC*\pgfutil@tempzC*\pgfutil@tempzD) +
\pgfutil@tempxB*(\pgfutil@tempxC*\pgfutil@tempxD*(-\pgfutil@tempyA + \pgfutil@tempyC) + \pgfutil@tempxC^2*(\pgfutil@tempyA - \pgfutil@tempyD) -
2*\pgfutil@tempyA^2*(\pgfutil@tempyC - \pgfutil@tempyD) + \pgfutil@tempyA*(2*\pgfutil@tempyC^2 - 2*\pgfutil@tempyC*\pgfutil@tempyD -
(\pgfutil@tempzA - \pgfutil@tempzC)*(\pgfutil@tempzC - \pgfutil@tempzD)) - (\pgfutil@tempzA - \pgfutil@tempzC)*
(-(\pgfutil@tempyD*(\pgfutil@tempzA + \pgfutil@tempzC)) + \pgfutil@tempyC*(\pgfutil@tempzA + \pgfutil@tempzD))))/(2*\pgfutil@tempe)}%
\else
\edef\pgfmathresult{0}%
\fi
\pgfmathsmuggle\pgfmathresult\endgroup}%
\tikzset{coorf/.style={insert path={({{#1}[0],{#1}[1],{#1}[2]})}}}
\begin{document}
\tdplotsetmaincoords{65}{70}
\begin{tikzpicture}[tdplot_main_coords,scale=1,line join = round,
line cap = round,
declare function={a=5;
A={a,0,0};B={0,a,0};C={0,0,a};D={0,0,0};}]
\path
[coorf=D] coordinate (O)
[coorf=A] coordinate (A)
[coorf=B] coordinate (B)
[coorf=C] coordinate (C);
\pgfmathsetmacro{\myx}{fourpointcenterx("A","B","C","D")}
\pgfmathsetmacro{\myy}{fourpointcentery("A","B","C","D")}
\pgfmathsetmacro{\myz}{fourpointcenterz("A","B","C","D")}
\path (\myx,\myy,\myz) coordinate (T);
\typeout{The center is at T=(\myx,\myy,\myz).}
\foreach \p/\g in {O/90,A/60,B/-90,C/-90,T/90}
{\path (\p) node[inner sep=1pt,fill,circle]{}
++(\g:3mm) node{$\p$}; }
\pgfmathsetmacro{\R}{sqrt(\myx*\myx+\myy*\myy+\myz*\myz)}
\begin{scope}[tdplot_screen_coords]
\draw (T) circle circle[radius=\R] ;
\path (T) node[above={\R*1cm+2mm},text width=pi*\R*1cm] {%
Given four noncoplanar points $A$, $B$, $C$ and $D$, the center of a
sphere through these points is at the intersection of three planes in the
middle between two given points. For instance, the plane in the middle
between $A$ and $B$ has the normal
\[ \vec n_{AB}=\vec B-\vec A\]
and has the equation
\[ \vec n_{AB}\cdot\vec
r=\vec n_{AB}\cdot(\vec A+\vec B)/2=:b_{AB}\;.\]
So the center can be obtained by solving a simple matrix equation
\[ N\, \vec T=\vec b\]
with $N=(\vec n_{AB},\vec n_{AC},\vec n_{AD})^T$ and $\vec
b=(b_{AB},b_{AC},b_{AD})$, i.e.\
\[\vec T=N^{-1} \vec b\;.\]
This equation can be solved with \texttt{Mathematica}, say. The condition
that the points be non coplanar is equivalent to $\det N\ne0$, which is a
necessary condition for the last step.
};
\end{scope}
\end{tikzpicture}
\end{document}
```
![Screen Shot 2020-08-12 at 2.30.12 AM.png](/image?hash=61a7c9269ae20a5ce487c6b023ea5aaaaf8809bcca54848f066de50044cdb747)