Added support for random sphere packing in a conical container.

[ZoeRichter]

Description

Currently, model.pack_spheres() only supports rectangular prisms, cylinders, and spheres as container. This PR adds a _Cone subclass to _Container that allows the user to randomly pack truncated conical containers using model.pack_spheres(). This is useful for modeling pebble-bed reactors, which often use a cone in the outlet chute.

It uses the _Cylinder container as a guide because it is the closest to the cone, and I stuck as close as I could to the approach/style there as I could. The biggest deviation from the other _Container subclasses is the way it stores its limits.

The radial limit is a linear function dependent upon the z-coordinate of the pebble. It is the equation of a line that is parallel to the side of the cone, but a perpendicular distance further into the cone equal to the radius of the spheres being packed.

Two parallel lines of the form

$$\begin{aligned} ax + by + c_1 &= 0 \\\ ax + by + c_2 &= 0 \end{aligned}$$

Have a perpendicular distance, d, of

$$\begin{aligned} d &= \frac{|c_2 - c_1|}{\sqrt{a^2 + b^2}} \\\ \end{aligned}$$

With the above naming convention, the line for the edge of the cone is

$$\begin{aligned} z_{peb} &= (r_{cone} - r_{minor})\frac{2z}{r_{major} - r_{minor}} + (z0 - z) \\\ 0 &= \frac{2z}{r_{major} - r_{minor}}r_{cone} - z_{sph} - \frac{2z}{r_{major} - r_{minor}}r_{minor} + (z0 - z) \\\ \end{aligned}$$

With $d = r_{sph}$ and $m = \frac{2z}{r_{major}-r_{minor}}$, the constant for the $r_{lim}$ line can be back-calculated:

$$\begin{aligned} r_{sph} &= \frac{|c_{lim} - c_1|}{\sqrt{a^2 + b^2}} \\\ c_{lim} &= r_{sph}\sqrt{a^2 + b^2} \pm c_1 \\\ c_{lim} &= r_{sph}*\sqrt{m^2 + 1} \pm ((z0-z) - r_{minor}) \end{aligned}$$

Plugging the constant back in, solving for $r_{lim}(z_{sph})$, and simplifying, we get to the equation of the radial limit in the cone:

$$\begin{aligned} r_{lim}(z_{sph}) &= \frac{1}{m}z_{sph} - \frac{r_{sph}}{m}\sqrt{m^2 + 1} + r_{minor} - \frac{z0 - z}{m} \\\ \end{aligned}$$

Because the limit depends on the specific sphere you are checking, but I wanted to use a similar data structure to the other containers, the radial limit is broken into a 2 element list, with the first being the coefficient of $z_{sph}$, and the second being the constant in the $r_{lim}$ equation. As far as I have been able to tell, the only place this caused an issue was in updating the domain limits (line 1596 of triso.py). I got the data structure I used to work with a try/except, but there's probably a more elegant/foolproof way to do that.

I have not (yet) checked off the documentation item because I am unsure if I need to add the above equations to the methods section of the documentation. The _Cone class does have docstrings, made by editing the _Cylinder docstring as needed.

Checklist

• I have performed a self-review of my own code
• I have run clang-format (version 15) on any C++ source files (if applicable)
• I have followed the style guidelines for Python source files (if applicable)
• I have made corresponding changes to the documentation (if applicable)
• I have added tests that prove my fix is effective or that my feature works (if applicable)

[ZoeRichter]

I think it is mathematically possible to use a similar approach with equations for lines and the perpendicular distance between them to work with more complex surfaces. I'm not sure of the general interest in expanding the available containers further, though - and I haven't thought as deeply about the implementation.