This morning someone asked me if I could “weaponize” his code, i.e. turn his research software into something robust enough for production. I like that term.
Month: March 2012
Approximating galaxies as spheres
A couple days ago I compared different ways of approximating Earth and other ellipsoids by spheres. The earth is so nearly spherical that the difference in the approximations would only matter when you need fairly high accuracy. Elliptical galaxies, however, can be much more eccentric than Earth and so the difference in approximation approaches can matter more.
The Hubble classification of elliptical galaxies uses a scale E0 through E7 where the number following ‘E’ is 10(1 – b/a) where a is the major semi-axis and b is the minor semi-axis. An E0 galaxy is essentially spherical. The most common classification is near E3. The limit is believed to be around E7.
The image above is a photo of Messier 49, an E4 galaxy, taken by the Hubble telescope.
For an E3 galaxy, the minor and major axes are around 7 and 10 in some unit. The average of these is 8.5. A sphere with the same volume would have radius 8.88 and a sphere with the same surface area would have radius 8.98, about 5.7% larger than the average of the axes.
For an E7 galaxy, the minor and major axes would have a ratio of 3 to 10. This gives an average of 6.5. Matching volumes gives a radius of 6.69 and matching surface area gives a ratio of 7.67, about 18% larger than the average of the axes.
Narcissus prime
This morning Futility Closet posted the following.
Repeat the string 1808010808 1560 times, and tack on a 1 the end. The resulting 15601-digit number is prime, and because it’s a palindrome made up of the digits 1, 8, and 0, it remains prime when read backward, upside down, or in a mirror.
I used Mathematica to verify that the number described above is indeed prime.
PrimeQ[ 10*Sum[1808010808*10^(10 i), {i, 0, 1559}] + 1 ]
After a little over two minutes, the function returned True.
Related posts
Approximating Earth as a sphere
Isaac Newton suggested in 1687 that the earth is not a perfectly round sphere but rather an ellipsoid, and he was right. But since our planet is roughly a sphere, it’s often useful to approximate it by a sphere. So if you’re going to do that, what radius do you use? More generally, what radius do you use when approximating any ellipsoid by a sphere?
This post will discuss the more general problem of finding the radius when approximating any ellipsoid by a sphere. We will give the answer for Earth in particular, and we’ll show how to carry out the calculations. Most of the calculations are easy, but some involve elliptic integrals and we show how to compute these in Python.
Ellipsoids and spheroids
First of all, what is an ellipsoid? It is a surface whose (x, y, z) coordinates satisfy
Earth is an oblate spheroid, which means a = b > c. Specifically, a = b = 6,378,137 meters, and c = 6,356,752 meters.
If you wanted to approximate an ellipsoid by a sphere, you could use
r = (a + b + c)/3.
Why? Because the knee-jerk reaction whenever you need to reduce a set of numbers to one number is to average them.
Volume of an ellipsoid
We could do a little better, depending on what property of the ellipsoid we’d like to preserve in our approximation. For example, we might want to create a sphere with the same volume as the ellipsoid. In that case we’d use the geometric mean
r = (abc)1/3.
This is because the volume of an ellipsoid is 4πabc/3 and the volume of a sphere is 4πr3/3.
For the particular case of the earth, we’d use
(a2c)1/3 = 6371000.7
Surface area of an ellipsoid
For some applications we might want a sphere with the same surface area as the ellipsoid rather than the same volume.
The surface area of an ellipsoid is considerably more complicated than the volume. For the special case of an oblate spheroid, like earth, the area is given by
where
The surface area of a sphere is 4 πr2 and so the following code computes r.
from math import sqrt, atanh
e = sqrt(1 - (c/a)**2)
r = a*sqrt(1 + (1-e**2)*atanh(e)/e) / sqrt(2)
This gives r = 6371007.1 for the earth, about 6.4 meters more than the number we got matching volume rather than area.
For a general ellipsoid, the surface area is given by
where
and
Here F is the “incomplete elliptic integral of the first kind” and E is the “incomplete elliptic integral of the second kind.” The names are historical artifacts, but the “elliptic” part of name comes from the fact that these functions were discovered in the context of arc lengths with ellipses, so it shouldn’t be too surprising to see them here.
Computing ellipsoid surface area in Python
In SciPy, F(φ, k) is given by ellipkinc and E(φ, k) is given by ellipeinc. Both function names start with ellip because they are elliptic functions, and end in inc because they are “incomplete.” In the middle, ellipeinc has an “e” because it computes the mathematical function E(φ, k).
But why does ellipkinc have a “k” in the middle? The “complete” elliptic integral of the first kind is K(k) = F(π/2, k). The “k” in the function name is a reminder that we’re computing the incomplete version of the K function.
Here’s the code for computing the surface area of a general ellipsoid:
from math import sin, cos, acos, sqrt, pi
from scipy.special import ellipkinc, ellipeinc
def area(a, b, c):
phi = acos(c/a)
k = a*sqrt(b**2 - c**2)/(b*sqrt(a**2 - c**2))
E = ellipeinc(phi, k)
F = ellipkinc(phi, k)
elliptic = E*sin(phi)**2 + F*cos(phi)**2
return 2.0*pi*c**2 + 2*pi*a*b*elliptic/sin(phi)
The differences between the various approximation radii are small for Earth. See my next post on elliptical galaxies where the differences are much larger.
More geodesy posts
Life off the clock
There was a lot of work to do a few generations ago, but the work wasn’t regulated by a clock.
With the growth of industrial capitalism during the post-Civil War years, more and more Americans were feeling pressure to be “on time.” (The phrase itself was a colloquialism which did not appear until the 1870s.) The corporate drive for efficiency … reinforced the spreading requirement that people regulate their lives by the clock. … And though there was much resistance, especially among workers from a preindustrial background, the triumph of clock time seemed assured by 1890, when the time clock was invented. [1]
Clocks had been around for centuries, but no one punched a time clock until 1890. People had regular schedules, some more so than others, but in general their schedules were not rigid or synchronized.
Increasing numbers of people now enjoy flexible work schedules. This is not something new but a return to something old. Industrialism made synchronization necessary. Post-industrial work is partially returning to pre-industrial norms.
[1] “No Place of Grace: Antimodernism and the Transformation of American Culture, 1880-1920” by T. J. Jackson Lears.