How to create nice-looking nuclei in TikZ?

Question

In the responses to Draw Bohr atomic model with electron shells in TeX?, there are nice drawings of atoms. However, the nuclei don't look very appealing or realistic.

I was wondering: Can anyone think of an algorithm to (semi-)automatically (for example in a randomized fashion) create large nuclei that look more realistic as for example the in the image attached?

It seems to be key that the balls are sufficiently spaced and that the spherical look requires more centered balls to be on top. Both requirements are not met with my code:

\documentclass{standalone}
\usepackage{tikz}

\begin{document}
\begin{tikzpicture}
\path (-2,-2) rectangle (2,2);
\pgfmathdeclarerandomlist{color}{{red}{white}}
\foreach \a in {1,...,200} {
    \pgfmathsetmacro{\r}{rnd}
    \pgfmathsetmacro{\a}{random(0,360)}
    \pgfmathrandomitem{\c}{color}
    \shade[ball color=\c] (\a:-\r) circle (5pt);
    }
\end{tikzpicture}
\end{document}

The result is:

EDIT:

In case anyone is interested, here is what I am quite happy with: Based on the answer, I have defined slightly modified versions of the suggested nucleus in three different sizes and with the option to feed a random seed to get different species.

\documentclass{standalone}
\usepackage{tikz}
\usetikzlibrary{calc}

\begin{document}

\begin{tikzpicture}
\tikzset{
    pics/proton/.style={code={\shade[ball color=red] circle (3pt);}},
    pics/neutron/.style={code={\shade[ball color=white] circle (3pt);}},
    pics/nucleussmall/.style={code={%
        \pgfmathdeclarerandomlist{nucleon}{{proton}{proton}{neutron}{neutron}{neutron}}
        \pgfmathsetseed{#1+1}
        \foreach \A/\R in {8/0.2, 5/0.13, 1/0}{
        \pgfmathsetmacro{\S}{360/\A}
            \foreach \B in {0,\S,...,360}{
                \pgfmathrandomitem{\C}{nucleon}
                \pic at ($(\B+2*\A+5*rnd:\R)$) {\C}; } }} },
    pics/nucleusbig/.style={code={%
        \pgfmathdeclarerandomlist{nucleon}{{proton}{proton}{neutron}{neutron}{neutron}}
        \pgfmathsetseed{#1+1}
        \foreach \A/\R in {24/0.4, 24/0.3, 24/0.2, 13/0.35, 11/0.27, 6/0.15, 1/0}{
        \pgfmathsetmacro{\S}{360/\A}
            \foreach \B in {0,\S,...,360}{
                \pgfmathrandomitem{\C}{nucleon}
                \pic at ($(\B+2*\A+5*rnd:\R)$) {\C}; } }} },
    pics/nucleusbiggest/.style={code={%
        \pgfmathdeclarerandomlist{nucleon}{{proton}{proton}{neutron}{neutron}{neutron}}
        \pgfmathsetseed{#1+1}
        \foreach \A/\R in {24/0.5, 24/0.4, 24/0.3, 24/0.2, 13/0.47, 15/0.44, 13/0.37, 11/0.27, 6/0.15, 1/0}{
        \pgfmathsetmacro{\S}{360/\A}
            \foreach \B in {0,\S,...,360}{
                \pgfmathrandomitem{\C}{nucleon}
                \pic at ($(\B+2*\A+5*rnd:\R)$) {\C}; } }} },
    }
\pic at (0,0) {nucleussmall};
\pic at (2,0) {nucleusbig=1};
\pic at (4,0) {nucleusbiggest=1};    
\end{tikzpicture}
\end{document}

Answer 1

Here is a proposal that makes the nucleus look more like a compact ball. It works by building up circular rings starting from the outside in. By adjusting the number of protons/neutrons in each ring and its distance from the center, you can create a ball effect.

\documentclass{standalone}
\usepackage{tikz}
\usepackage[version=4]{mhchem}
\begin{document}
\begin{tikzpicture}
\path (-2,-2) rectangle (2,2);
\pgfmathdeclarerandomlist{color}{{red}{white}}
\pgfmathsetseed{1}
\foreach \A/\R in {25/1,12/0.9,15/0.8,20/0.7,12/0.5,7/0.3,1/0}{
      \pgfmathsetmacro{\S}{360/\A}
           \foreach \B in {0,\S,...,360}{
               \pgfmathrandomitem{\C}{color}
               \shade[ball color=\C] (\B+\A:\R) circle (5pt);
           }
}
\node at (-1,1.3) {\ce{^{226}_{88}Ra}};
\end{tikzpicture}
\end{document}

Answer 2

Based on your code, I first draw protons/neutrons following a circular pattern three times, at radius 1, 0.5 and 0.2. I also draw random protons/neutrons in between.

\documentclass{standalone}
\usepackage{tikz}

\begin{document}
\begin{tikzpicture}
\path (-2,-2) rectangle (2,2);
\pgfmathdeclarerandomlist{color}{{red}{white}}

\foreach \a in {0,10,...,360}{
    \pgfmathrandomitem{\c}{color}
    \shade[ball color=\c] (\a:1) circle (5pt);
}

\foreach \a in {0,20,...,360}{
    \pgfmathrandomitem{\c}{color}
    \shade[ball color=\c] (\a:0.5) circle (5pt);
}

\foreach \a in {1,...,350} {
    \pgfmathsetmacro{\r}{rnd}
    \pgfmathsetmacro{\a}{random(0,360)}
    \pgfmathrandomitem{\c}{color}
    \shade[ball color=\c] (\a:\r) circle (5pt);
    }

\foreach \a in {0,60,...,360} {
    \pgfmathrandomitem{\c}{color}
    \shade[ball color=\c] (\a:0.2) circle (5pt);
}
\end{tikzpicture}
\end{document}

The results is:

Answer 3

Here is another version in which the spheres are put on the root lattice of A_3 and allowed to wiggle a bit. More explanations can be found here.

\documentclass[tikz,border=3.14mm]{standalone}
\usepackage{tikz-3dplot}
\usetikzlibrary{3d}
\tikzset{declare function={posx(\x,\y,\z)=\x-\y/2;
posy(\x,\y,\z)=\y/sqrt(2);
posz(\x,\y,\z)=-\y/2+\z;
}}
\newsavebox\Proton
\newsavebox\Neutron
\sbox\Proton{\tikz{\shade[ball color=red] circle({0.85/sqrt(2)});}}
\sbox\Neutron{\tikz{\shade[ball color=gray!20] circle({0.85/sqrt(2)});}}
\begin{document}
\xdef\Lst{{-1, 0, 2}, {-2, -1, 1}, 
 {0, 0, 2}, {-1, -1, 1}, 
 {-2, -2, 0}, {-1, 1, 2}, 
 {-2, 0, 1}, {1, 0, 2}, {0, -1, 1}, 
 {-1, -2, 0}, {-2, -3, -1}, 
 {0, 1, 2}, {-1, 0, 1}, 
 {-2, -1, 0}, {1, -1, 1}, 
 {0, -2, 0}, {-1, -3, -1}, 
 {1, 1, 2}, {0, 0, 1}, {-1, -1, 0}, 
 {-2, -2, -1}, {0, 2, 2}, 
 {-1, 1, 1}, {2, 1, 2}, {-2, 0, 0}, 
 {1, 0, 1}, {0, -1, 0}, 
 {-1, -2, -1}, {-2, -3, -2}, 
 {1, 2, 2}, {0, 1, 1}, {-1, 0, 0}, 
 {2, 0, 1}, {-2, -1, -1}, 
 {1, -1, 0}, {0, -2, -1}, 
 {-1, -3, -2}, {2, 2, 2}, 
 {1, 1, 1}, {0, 0, 0}, 
 {-1, -1, -1}, {-2, -2, -2}, 
 {1, 3, 2}, {0, 2, 1}, {-1, 1, 0}, 
 {2, 1, 1}, {-2, 0, -1}, {1, 0, 0}, 
 {0, -1, -1}, {-1, -2, -2}, 
 {2, 3, 2}, {1, 2, 1}, {0, 1, 0}, 
 {-1, 0, -1}, {2, 0, 0}, 
 {-2, -1, -2}, {1, -1, -1}, 
 {0, -2, -2}, {2, 2, 1}, {1, 1, 0}, 
 {0, 0, -1}, {-1, -1, -2}, 
 {1, 3, 1}, {0, 2, 0}, {-1, 1, -1}, 
 {2, 1, 0}, {1, 0, -1}, 
 {0, -1, -2}, {2, 3, 1}, {1, 2, 0}, 
 {0, 1, -1}, {-1, 0, -2}, 
 {2, 0, -1}, {1, -1, -2}, 
 {2, 2, 0}, {1, 1, -1}, {0, 0, -2}, 
 {2, 1, -1}, {1, 0, -2}}
\tdplotsetmaincoords{-90+109.471}{-90+70}
\foreach \X in {1,...,10}
{\begin{tikzpicture}
\path[use as bounding box] (-3.5,-3.5) rectangle (3.5,3.5);
\draw (0,0) circle ({1}); % /sqrt(2)
\begin{scope}[tdplot_main_coords]
 \draw[-latex] (0,0,0) coordinate (O) -- (1,0,0) node[right]{$\alpha_1$};
 \draw[-latex] (O) -- (-1/2,{1/sqrt(2)},-1/2) node[right]{$\alpha_2$};
 \draw[-latex] (O) -- (0,0,1) node[right]{$\alpha_3$};
 \draw[red,-latex] (O) -- (1/2,{1/sqrt(2)},1/2) node[right]{$-\theta$};
 \foreach \Z in \Lst
  {\pgfmathsetmacro{\myx}{{\Z}[0]}
  \pgfmathsetmacro{\myy}{{\Z}[1]}
  \pgfmathsetmacro{\myz}{{\Z}[2]}
  \pgfmathsetmacro{\mydeltax}{0.1*(rnd-0.5)}
  \pgfmathsetmacro{\mydeltay}{0.1*(rnd-0.5)}
  \pgfmathsetmacro{\mydeltaz}{0.1*(rnd-0.5)}
  \pgfmathtruncatemacro{\mycol}{int(2*rnd)}
  \ifnum\mycol=1
  \node at ({posx(\myx+\mydeltax,\myy+\mydeltay,\myz+\mydeltaz)},
  {posy(\myx+\mydeltax,\myy+\mydeltay,\myz+\mydeltaz)},
  {posz(\myx+\mydeltax,\myy+\mydeltay,\myz+\mydeltaz)}) {\usebox\Neutron};
  \else
  \node at ({posx(\myx+\mydeltax,\myy+\mydeltay,\myz+\mydeltaz)},
  {posy(\myx+\mydeltax,\myy+\mydeltay,\myz+\mydeltaz)},
  {posz(\myx+\mydeltax,\myy+\mydeltay,\myz+\mydeltaz)}) {\usebox\Proton};
  \fi}
\end{scope}
\end{tikzpicture}}
\end{document}

Just for fun: more nuclei. And no, it does not look like sphere, but like a set of sphere which are packed with maximum density. This is of course not the same as demanding that the nuclei should fill out a sphere. The latter might translate in the requirement that the sum of distances gets minimized or something like that, which obviously is not the same requirement as maximal packing. I do not know if there is a simple algorithm that minimizes the sum of distances while making sure the spheres do not overlap.