[VA] SRC #011 - April 1st, 2022 Bizarro Special

[Valentin Albillo]

 
Hi, all,

Thanks a lot to those of you who posted some comments, namely Albert Chan, EdS2, ijabbott, Massimo Gnerucci, rawi, Ren, rprosperi, vaklaff, Werner, and last but certainly not least, J-F Garnier, much appreciated. Now, for my own final Comments:

My Comments

1) On my solution's algorithm and poor-man's random number generator

As I've said, there's no way to implement a deterministic cubature algorithm for a non-trivial definite 6-D integral like the one featured here in just 79 bytes of RAM (including both program and required data storage) and with no subroutines, so I had to use the poor-man's RNG which I advocated (and sent to Richard Nelson for publication in PPC CJ, where it was indeed eventually published 40+ years ago) as the fastest & smallest one which was still capable of producing decent results.

As you can see in the above linked vintage letter (page 7), it can be used in the HP-41, HP-67/97 and any calculators featuring a built radians-to-degrees conversion, and is particularly useful for games and simulations, due to its speed and simplicity.

This RNG generates uniformly distributed pseudo-random numbers in the interval 0-1 and the seed can be any integer or real number except 0, Pi or its multiples. It essentially uses a multiplier equal to 180/Pi, and as you can see in the vintage letter, a trial test generating and analyzing 1,000 values produced a decent uniform distribution with mean = 4.4 and standard deviation = 2.9 where the theoretical values are 4.5 and 3.0, respectively, close enough. I also couldn't find its period back then after generating 3,000 values.

J-F Garnier had the correct idea and almost succeeded in duplicating my program, as seen in his post, but he couldn't fit his program in the HP-10C and more important, he didn't discover the possibility of using the ->DEG instruction available in the HP-10C, which essentially uses 180/Pi = 57.2957795+ as the multiplier, and so he used Pi = 3.14159265+ instead. The problem with such a low-valued multiplier is that it's very prone to ascendent runs. For instance, if the seed ever becomes a low value such as 0.01, you'll get a long ascendent run:

0.0100 -> 0.0314 -> 0.0987 -> 0.3101 -> 0.9741

that is, 5 consecutive values in increasing order, which means a linear dependency from the previous value and damages the overall randomness, which is probably why J-F's program couldn't achieve a sufficiently accurate answer. If the seed eventually gets smaller than 0.01, you'll get even longer ascending runs (7 values, 10 values, ...).

My method has this problem too, but to a far lesser extent because the ->DEG multiplier (57.30+) is much bigger than J-F's 3.14+ and thus any ascending runs are far shorter. For instance, with the same seed:

0.0100 -> 0.5730

and the next number generated may or may not be in ascending order. Matter of fact, the worst that can happen to a method which uses just a multiplier (unlike linear congruential generators, which use multiplication, addition and modulus operations) is that if ever the seed becomes exactly 0, then the method will be stuck on an indefinite loop, always producing 0, and as I've explained above, whenever the seed becomes very small, then you'll get a long ascendent run if the multiplier is also small, like J-F's Pi.


2) On non-deterministic cubature methods

When numerically computing definite integrals in a single variable, it's quite common and efficient to use deterministic methods that evaluate the function being integrated at a number of well-chosen arguments (such as 16-point Gaussian quadrature, say), which for reasonably well-behaved integrands are both fast and accurate.

However, computing a double integral in two variables to the same level of accuracy would roughly need 16² integrand evaluations and in general computing a multiple integral in D variables will require grosso modo about 16^N evaluations, which for the 6-D integral featured here would be 16 ⁶ ~ 17 million evaluations.

When the number of evaluations needed to compute the integral grows exponentially with the dimension D, then the integration method suffers from the so-called curse of dimensionality, which means that deterministic methods are utterly inefficient for high-dimensional integrals (such as the ones appearing in mathematical/computational finance, where integrals having hundreds (D > 100) and even thousands (D > 1000) of variables aren't uncommon) and in practice it is mandatory to resort to non-deterministic Monte Carlo (MC) methods, which do not suffer the curse of dimensionality but converge very slowly.

Matter of fact, simple MC methods converge as slowly as 1/√D, which means that to get one additional correct digit (10x accuracy) we must use 100x the number of evaluations, but we can resort instead to Quasi-Monte Carlo (QMC) methods, which attempt to speed up the convergence to 1/D (i.e. to increase 10x the accuracy you have to increase 10x the number of evaluations, not 100x) by using low-discrepancy sequences (aka quasi-random sequences) instead of sequences of (pseudo-)random numbers, as MC uses.

The gains in speed and accuracy that QMC methods afford over simple MC (let alone deterministic methods) for multidimensional integration is extremely noticeable, e.g. requiring as little as 500 integrand evaluations to compute a 25-D test integral within 0.01, as compared to 220,000 evaluations using MC.

Mind you, none of this would fit in 79 bytes at all, so I did the best I could given the circumstances !!


3) On the gravitational force F between two cubical planets

In the case of spherical planets in contact (m₁=1, m₂=1, d=1, G=1), the gravitational force F is 1, but if the planets are cubical and in contact, we have F<1 because they have part of their mass in the corners, which are farther away.

If instead of being in contact the cubes were at a distance d > 1 or even d >> 1, then F would quickly approach 1/d² and the cubes would act more and more like spheres in that their mass could be considered as a point mass at their centers, like in the spherical case. Indeed, by the time the centers of the cubic planets are 4 or more units apart (d >= 4), they can be practically considered spherical as far as gravity is concerned.


4) On real-life applications of computing the gravitational field of a cubical object

In the past few years, a number of spacecraft have been sent to visit diverse astronomical objects, from 1 Ceres (dwarf planet, 939 km mean diameter, visited by Dawn) to 101955 Bennu (asteroid, 490 m mean diameter, visited by OSIRIS-REx, which first orbited, then successfully touched down on its surface and later departed for Earth). Ceres is big enough to have a reasonably spherical shape, but Bennu is markedly "squarish":

and other irregular objects also visited by spacecraft include two-lobed comet 67P/Churyumov–Gerasimenko, which is much further away from sphericity:

As another such instance, the asteroid 433 Eros (~17 km mean diameter) also has a highly irregular shape and was visited by the NEAR Shoemaker spacecraft, which was initially put on a relatively distant ~320-360 km elliptical orbit. At that distance, Eros' gravitational field could be considered as if the mass of the asteroid were concentrated in its center but later, when NEAR was moved to a much closer orbit and eventually landed on the asteroid, it was necessary to compute a more accurate gravitational field, least the spacecraft would impact the asteroid at a potentially dangerous speed.

The problem is compounded if, as it's usually the case, the object not only has an irregular shape but it's also rotating. In the case of a cubic planet with a side equal to Earth's diameter, can an artificial satellite orbit it ? For starters, it will feel a stronger gravitational attraction near the cube's corners and there will be additional resonances as the planet rotates. Moreover, again due to the corners, the satellite won't follow a closed elliptical orbit but will instead be subject to rapid precession and in general the orbit won't close at all.

If both the planet and the satellite rotate in the same direction, with the satellite orbiting a few planet's radii high, the differential rotation will perturb the orbit so much when the satellite is near the cube's corners that eventually it will collide with the planet, like this:

Thus, launching satellites in low orbits around non-spherical bodies requires careful calculation to overcome the perturbations, and there's a number of academic publications on the gravity field of a cube, with the resulting formulae being used in real life to compute the gravitational field of a body of irregular shape by superimposing on the object a 3D grid of cubic blocks, iteratively reducing their size until the desired accuracy is attained, like this (approximately, you get the idea) :

Such methods can be tested as I did here, by applying them to a case whose solution is known, i.e. a spherical planet, whose shape is approximated by iteratively filling up its volume with cubic blocks of diminishing size as per the algorithm, then integrating the gravitational force over all the cubic blocks. This can also be applied to non-homogeneous bodies by using cubic blocks small enough for them to be considered individually homogeneous, then having their individual densities vary as needed.

And just in case you'd think that cubical planets wouldn't be taken seriously by anyone ...  

---

Well, that's more than enough !

If you want to comment something about my OP, my Original Solution and/or my ersatz RNG you're welcome to post it to this very thread. But for comments on general Monte Carlo or quasi-Monte Carlo methods or general space exploration please create another thread so that this one may remain on-topic and focused on my OP. Thanks !

This will be my last SRC for a long while, hope you enjoyed it. Thanks for your interest and

Best regards.
V.

P.S.: A final question:  knowing that the Bizarro given name "Nitnelav" is unisex (like the English given names Morgan, Cameron or Hayden, say), what do you think ?  Is Rd. Albizarro a Bizarro-man or a Bizarro-woman ?