Archive for the ‘math’ Category

Unintended Consequences of Math

February 17th, 2010 at 11:30 pm

I thought last week’s episode of In Our Time was particularly good, even though I already knew most of the facts presented and disagreed with much of the discussion of those facts. Actually, maybe I liked it because I disagreed with it. It gave me a lot to think about, which thinking I will now inflict on you.

In Our Time is a weekly radio show broadcast on BBC Radio 4, the hifalutin BBC station. People outside the UK can get the podcast. Each episode covers a fairly restricted topic, with panel discussion moderated and prodded by the host, Melvyn Bragg. Unlike most American shows of this type, In Our Time puts real experts on its panels, not journalists hawking their latest books. Each show’s panel is usually made up of professors of the subject from various British universities. The scope of the show is immense, as you can see from looking at its archive. The most recent three shows, for example, covered

It’s this last one that got me thinking.

Most of the show was given over to stories, told in the style of James Burkes’ Connections series, of concepts developed purely for their own mathematical beauty (that’s the pure math of the episode’s title) that suddenly became useful for practical purposes decades or centuries later.

The first of these stories, told by the delightfully-named Colva Roney-Dougal of the University of St. Andrews, purportedly explains how the solution of cubic equations led to the electric chair. In the 16th century, while investigating the solution of cubic equations, Cardano found square roots of negative numbers in his solutions. This led to the use of complex numbers. It turned out, hundreds of years later, that complex numbers were very convenient in the analysis of alternating current.

(Phasor image blatantly stolen from here.)

If you know anything about the competition between Thomas Edison and George Westinghouse over use of DC and AC, you can see where this is going. As part of his campaign against AC, Edison touted its use in electrocution.1 Therefore, we’re supposed to think, the electric chair was an unintended consequence of complex numbers.

But, and this is pointed out in the show, while complex numbers make the analysis of AC easier, they’re not required. In fact, the electric chair was first designed before anyone used complex numbers to analyze AC. So yes, it’s interesting that complex numbers can be used in AC analysis (as they can be used similarly to analyze mechanical vibrations and problems in plane elasticity), but it really isn’t right to call the connection a consequence.

The use of the term consequence is even less apt in the second example: the application of conic sections to celestial orbits and the trajectories of projectiles. There’s no question this is an interesting connection. The Greeks studied conic sections for their intrinsic beauty, and it was definitely an exercise in pure math. There’s no particular reason to believe they’d be applicable to gravity-controlled motions. But again, this isn’t a consequence, it’s a discovery. Planets didn’t start orbiting in ellipses because Kepler said so.

(Image from Duk at Wikipedia.)

The third story eventually works its way (via Cardano and de Moivre) to Gauss’s use of the normal distribution to predict the position of an asteroid. The asteroid had been observed for only a short time before the view of it was blocked by the sun. Gauss recognized that the observations prior to its disappearance had errors that followed the normal distribution and was able to give the best prediction of where it would be after it reappeared.

The normal distribution had been described by previous mathematicians—de Moivre and Laplace, in particular—studying games of chance and other random phenomena, and Gauss’s use of it to analyze observational errors was definitely a consequence of that earlier work. But here I would disagree with the characterization of that earlier work as pure math. Early studies of probability were intensely practical because of their connection to gambling. There was nothing pure about it.

I would classify the fourth connection, between non-Euclidean geometry and relativity, as I did the conic section/orbit connection: interesting but not a consequence because there’s no cause and effect relation.

The last story, finally, does what all the stories were supposed to do: it describes how the study of prime numbers, once thought of as the purest of pure math topics, now has a very practical and thoroughly unexpected application in cryptography, keeping our credit card numbers safe as we shop online. The historical connection between the pure math and the practical application wasn’t made very strongly—it was basically some hand-waving about Alan Turing—but the story may have been rushed because the show was nearly over.

I was amazed that the show finished without a single mention of “The Unreasonable Effectiveness of Mathematics in the Natural Sciences,” a famous paper written in 1960 by a physicist, Eugene Wigner. (How famous is this paper? Type “unreasonable” into your browser’s Google search field and see what comes up. It was the seventh suggestion when I tried it.)

The theme of the paper, as you can guess from the title, is how amazing it is that mathematics—often developed for completely different purposes—works so well to describe the physical world. While the In Our Time episode’s theme is not exactly the same as this, it’s damned close. Even the episode’s title is an echo of the paper’s title.

Fairly early in the show, one of the guests,2 John Barrow of the University of Cambridge, dives right into the main topic of “Unreasonable Effectiveness,” giving a very reasonable explanation of why many topics in math have applicability beyond their original conception. In his view, the whole of math is a collection of patterns, with the number of patterns increasing as math has progressed. It is, therefore, not surprising that some of these patterns match the patterns found in nature. How Barrow manages to go through this explanation without mentioning Wigner’s paper, I’ll never understand, but this was my favorite part of the show.

I doubt I’ll have as much to say about this week’s episode, but I bet I’ll enjoy it.

  1. No, I’m not going to link to the video of the electrocution of Topsy the elephant. I’m sure you can find it on your own if you really want to see it. ↩

  2. Occasional listeners of In Our Time know that the guest panel consists of three experts and may be wondering who the third guest is. Regular listeners will have guessed that it’s Marcus du Sautoy of Oxford, who seems to appear on every BBC show on math. ↩

Tags: , , ,
Post #1180 in entertainment, math, media | No Comments

That LaTeX jQuery plugin

December 22nd, 2009 at 10:54 pm

I’ve been seeing a link to this article, “A jQuery plugin to directly embed LaTeX into your website or blog,” on delicious.com/popular today. I have a strong interest in this topic—most recently expressed here and here—so I followed the link to have a look. It was disappointing.

Basically, the plugin

You can use other online LaTeX generators in place of Code Cogs, but they all generate images.

How is this disappointing? Let me count the ways:

  1. It won’t do inline equations. Reliance on <div>s forces each equation into its own block.
  2. You have to write the LaTeX and call the plugin for each equation.
  3. It relies on the LaTeX generator site to be up and running when someone visits your page.
  4. Although the equation is just a static image, that image has to be generated every time the page is visited.1
  5. Images don’t scale with the font size of the page.
  6. Images look like shit when printed because they’re lo-res bitmaps.

Update 12/24/09
Andreas Grech, the plugin’s author corrects me in a comment below.

Point 2 is wrong. You need only call the plugin once per page. I read the source code of his example page too quickly and came to the wrong conclusion. He calls the plugin several times on that page because he’s demonstrating various ways of using it.

Point 1 is at least half wrong. Because the plugin works off classes, you can put equations in <span>s. I’m not sure there’s a way to get LaTeX’s inline style, which is different from its display style.

His comment is a good defense of his plugin. It doesn’t persuade me to use the plugin, but it’s worth reading.

The jsMath system, which has been around for quite a while, has none of these disadvantages. OK, the equations are “typeset” on the page dynamically, but the work is done locally on the browser’s computer, so it’s well distributed and not being performed by some innocent third-party server.

You could argue that images have more universal support than the CSS tricks jsMath uses. That’s true, but apart from Mobile Safari, is anyone really using a browser that isn’t on jsMath’s supported list?

Maybe I shouldn’t be so harsh; after all, I’ve written plenty of scripts that must seem foolish to others. But if you want to use images for your equations, it’s far more efficient to generate them statically. You could even write a script that does what this plugin does, but saves the generated images to your site and inserts the appropriate <img> tags into the HTML. Such a script would be run just once—or, at most, a few times as you’re editing the page.

And if you want equations that actually look good on your page, use jsMath.

  1. I can see where this would be an advantage as you’re writing and editing the page, but once it’s done you’re just using the generator site’s bandwidth for no good reason. ↩

Tags: , , , , ,
Post #1140 in internet, math | 6 Comments

MathJax: equations on the web

December 15th, 2009 at 10:15 pm

When I started reading the MathJax homepage, I wondered why I would ever want to switch to it from the tried and true jsMath, which I’ve used for years to get good looking equations in my web page. Then I saw that Davide Cervone, the creator of jsMath, is one of the MathJax designers, and decided this is a project worth watching.

MathJax continues the jsMath tradition of accepting TeX/LaTeX notation to generate the equations, and intends to extend that to cover MathML as well. The MathJax preview page has many examples. Here’s the MathJax code for showing the binomial distribution on a line of its own:

\[P(E)   = {n \choose k} p^k (1-p)^{ n-k} \]

The result (apart from a difference in fonts; see below) looks like this:

P(E) = {n \choose k} p^k (1-p)^{ n-k}

That code should be familiar to anyone who’s used LaTeX. It is, in addition, the same way you’d generate the equation if you were writing an HTML page using jsMath and the tex2math plugin.

It’s also how I generated the equation for this post. I don’t use the tex2math plugin, but I modified Michel Fortin’s PHP Markdown Extra to take input like that and create the equivalent in “classic” jsMath:

<div class="math">P(E)   = {n \choose k} p^k (1-p)^{ n-k}</div>

MathJax will have the ability to use fonts likely to already be on the user’s computer and to download fonts as needed via webfonts. This will be a big improvement over jsMath, which really prefers that the user have the TeX font set already installed.

Update 12/17/09
Looking at the MathJax preview page again, I see that it is using the TeX fonts (i.e, Computer Modern), at least on systems which, like mine, have them installed. I wonder if it uses @font-face to bring those fonts to systems that don’t. The fonts are free, so I can’t imagine there being a licensing problem with them.

Since there isn’t a MathJax download link, it appears that the authors don’t consider it ready for prime time just yet. But with talented programmers and some significant sponsors behind it, MathJax has a good chance to become the standard way people write equations for the web.

Tags: , , , ,
Post #1133 in internet, math, programming | No Comments

An application of Castigliano’s Second Theorem with Octave

October 24th, 2009 at 5:53 pm

This week I did an analysis of a set of leaf springs and found myself using a favorite technique that I seldom get to apply anymore: Castigliano’s Second Theorem. There is some real, though not especially difficult, math ahead, so if you come here for the little Mac-based scripts I typically write you may want to skip this one. If you do want to follow along, I suggest you visit this page and download one of the TeX font packages there. This page, and the other pages in which I use Davide Cervone’s wonderful jsMath system, will load faster and print better if you have those fonts on your computer.

First, let’s talk about leaf springs. Leaf springs are long, flat steel strips bundled together to provide part of the suspension system of heavy vehicles like trucks. They sit between the axle and the frame and smooth out the ride by flexing. You can think of them as allowing the body of the vehicle to bounce above the axle (roll your mouse over the images to see the springs flex),

or you can think of them as allowing the axle to bounce up toward the frame as the vehicle goes over bumps in the road.

In fact, the usual behavior is a combination of the two, but since we will be investigating the spring stiffness—a static property—the analysis is the same no matter which frame of reference you choose.

In some leaf spring packages, the individual leaves are curved to nestle against one another; in others, like the kind I dealt with this week, the leaves are curved to meet only at the center and at the tips, as in the simplified schematic drawings above.

To determine the stiffness of a spring package like that in the drawing, each leaf can be treated as a simply-supported beam with a central point load, and the stiffnesses of the package is the sum of the stiffnesses of the individual leaves. Because the initial curvatures of the leaves are relatively small (the ratio of leaf thickness to initial radius of curvature is much less than one), we can treat them as straight beams with no significant error.

The deflection under the load of a simply-supported beam of uniform thickness and width loaded at the center is

\Delta = \frac{F L^3}{48EI}

where E is the modulus of elasticity of the steel (29,000,000 psi in US Customary units), and I is the area moment of inertia of the cross section,

I = \frac{b t^3}{12}

for a rectangle of width b and thickness t. The equivalent spring stiffness, k, is the ratio of load to deflection, so

k = \frac{F}{\Delta} = \frac{48EI}{L^3} = \frac{4Ebt^3}{L^3}

Repeat this calculation for each leaf, add them up, and you get the equivalent spring stiffness for the package. Very simple if you have leaves of uniform thickness.

Unfortunately, the leaves I was analyzing were tapered, about twice as thick in the center as out at the tips. It’s because there is no simple formula for a tapered beam that I went to Castigliano’s Second Theorem.

Carlo Castigliano was a 19th-century Italian mathematician/physicist/engineer who studied at the Polytechnic of Turin. His dissertation presented the two theorems that bear his name. The First Theorem is interesting but doesn’t have a lot of practical application. Castigliano’s Second, though, is very helpful in solving many interesting real-world problems in structural and mechanical engineering. Including, as I said at the top, the problem of the tapered beam.

Castigliano’s Second Theorem states that the derivative of the complementary strain energy of a body with respect to a point load acting on that body is equal to the deflection, in the direction of the load, of the point on the body to which the load is applied. In mathematical terms, this is

\frac{\partial U^*}{\partial F} = \Delta

A trivial demonstration of this can be found in the behavior of a simple linear spring. The complementary strain energy in the spring, U^* is the area above the usual force-deflection diagram, written in terms of the applied force and the spring stiffness. Recall that the regular strain energy, U, is the area below the force-deflection diagram, and is written in terms of the deflection and the spring stiffness.

According to Castigliano’s Second, the deflection of the spring is

\Delta = \frac{\partial U^*}{\partial F} = \frac{\partial}{\partial F}\left(\frac{F^2}{2k}\right) = \frac{F}{k}

which is exactly what we expect. The genius of Castigliano was his generalization of this trivial case to any type of structure.

Update 10/25/09
I meant to add here that the strain energy and the complementary energy are equal for a linear elastic spring, but would not be equal for a nonlinear spring.

To apply Castigliano’s Second to our tapered beam problem, we’re going to first take advantage of the symmetry of the beam and analyze just the left half. The right half will behave as a mirror image of the left, so once we’ve analyzed the left half, we have the solution for the right half.

The left half acts like a cantilevered beam, as shown in the diagram below. The complementary strain energy of a beam in bending is governed by this equation:

U^* = \int_0^L \frac{M^2}{2EI} dx

where M is the bending moment in the beam due to the applied force. Both M and I are functions of x,

M = Fx
I = \frac{b (t_0 + \alpha x)^3}{12}

where \alpha = (t_1 - t_0)/L is the rate at which the thickness increases as we move from the tip to the wall. Thus,

U^* = \int_0^L \frac{F^2 x^2}{2 E \frac{b (t_0 + \alpha x)^3}{12}} = \frac{6 F^2}{E b} \int_0^L \frac{x^2}{(t_0 + \alpha x)^3} dx

Putting this into the Castigliano equation, we get

\Delta = \frac{\partial U^*}{\partial F} = \frac{12F}{Eb} \int_0^L \frac{x^2}{(t_0 + \alpha x)^3} dx

Back when I was a student, we didn’t have all these fancy personal computers, so we’d have to do that integral analytically, probably through integration by parts. Now I can’t be bothered. I just fire up Octave, enter the values, and let it do the integration numerically. Here’s my Octave session for a particular set of tapered beam dimensions (I used Octave’s diary command to save the session as a text file):

octave-3.2.3:2> global E=29e6 b=3 t0=.5 t1=1 L=25;
octave-3.2.3:3> function retval = A(x)
> global t0 t1 L;
> retval = x^2/(t0 + (t1-t0)/L*x)^3;
> endfunction
octave-3.2.3:4> [int,icode,nfun,err] = quad("A",0,L)
int =  8518.4
icode = 0
nfun =  21
err =  9.4573e-11
octave-3.2.3:5> k = E*b/12/int
k =  851.10

So the stiffness of one-half of a leaf is 851 lb/in, which makes the stiffness of the full leaf 1702 lb/in. Compare this to

\frac{4\cdot29,000,000\cdot3\cdot1^3}{50^3} = 2785\: \textrm{lb/in}

for a 1″ flat leaf, and

\frac{4\cdot29,000,000\cdot3\cdot0.5^3}{50^3} = 348\: \textrm{lb/in}

for a 0.5″ flat leaf. How close would we have come if we’d used the flat-leaf formula with the average thickness? We’d get

\frac{4\cdot29,000,000\cdot3\cdot0.75^3}{50^3} = 1174.5\: \textrm{lb/in}

which underestimates the stiffness by more than 30%. Another way to think of this is that by tapering the beam from 1″ to ½″ we get 30% more stiffness for the same amount of steel as a flat beam of ¾″.

For the purposes of this post, I made the taper linear, which allowed me to use Octave’s quad function for integration. In fact, the leaf thickness tends to vary in a more complicated way, so the more general way to solve the problem is to

  1. Express the thickness as a vector of values uniformly-spaced along the length of the spring.
  2. Use Octave’s many matrix/vector operators to generate a vector of the integrand.
  3. Use the trapz function to perform the integration.

It takes a bit longer to do it this way because there are more thickness values to enter, but conceptually it’s the same problem.

The tapered leaf spring problem is tailor-made for Castigliano’s Second Theorem.

When you run into a problem with these characteristics, it’s good to have Castigliano in your toolbox.

Tags:

Post #1073 in math, mechanics, programming | No Comments

AIDS vaccine and statistical significance

October 20th, 2009 at 1:01 pm

This morning came word that the results of the AIDS vaccine trial, reported as a modest success last month, were, on further analysis, not statistically significant. I’ve been waiting a couple of weeks for this shoe to drop, as even the originally-reported results were pretty weak. I’m not claiming any expertise in vaccine research, but I do know how to do elementary statistics.

According to the LA Times article, the switch from statistically significant to not statistically significant came when seven of the test subjects were reclassified. That’s seven subjects out of over 16,000. How can such a small change flip the results from success to failure? Because the original results were on the edge of statistical significance to begin with, as anyone with a little statistics background could have calculated.

Here are the results as they were first given:

Observed Vaccine No vaccine
HIV positive 51 74
HIV negative 8,146 8,124

This is what’s known as a 2×2 contingency table. There are two treatments, vaccine and no vaccine, and two outcomes, HIV positive and HIV negative. The basic test of statistical significance poses the following question:

Assume that the vaccine has no effect, that there would have been 125 HIV-positive test subjects even if all 16,395 subjects had been untreated. This “pretend there’s no effect” assumption is called the null hypothesis. If one had then arbitrarily split the test population into two groups, 8,197 in one group and 8,198 in the other, what is the probability that the distribution of HIV-positive subjects would be at least as skewed as they turned out?

By “at least as skewed” I mean a distribution of the HIV-positive subjects between the vaccinated and unvaccinated that’s at least as far from an equal distribution as the test results were. Since the test results were 51-74, distributions like 50-75, 49-76, 48-77, etc., all the way to 0-125, would all fall into the “at least as skewed” category. (So, by the way, would results skewed the opposite way: 74-51, 75-50, and so on. Results like this would suggest the vaccine has a negative effect, but that’s also a violation of the null hypothesis.)

If the calculated probability is sufficiently low—in other words, if it’s quite unlikely that the test results were due to chance alone—we say that the results are statistically significant. It’s fairly common to use 5%, one chance in twenty, as the upper limit for statistical significance, but lower values are sometimes used. In fact, it’s often thought to be good practice to provide that probability instead of just reporting whether the results are significant or not.

Which leaves us with the problem of calculating that probability. Fortunately, the 2×2 contingency table is a very well-studied problem, and the procedure is straightforward. First, we rewrite the table of test results (often called the observed values) including the row and column sums.

Observed Vaccine No vaccine Sum
HIV positive 51 74 125
HIV negative 8,146 8,124 16,270
Sum 8,197 8,198 16,395

If the null hypothesis is true, then the probability of any subject in the test group becoming HIV positive is

p = \frac{125}{16,395} = 0.007624

or about three-quarters of one percent. We then make a similar table filled, not with the test results, but with the expected values based on the above probability

Expected Vaccine No vaccine Sum
HIV positive 62.496 62.504 125
HIV negative 8134.504 8135.496 16,270
Sum 8,197 8,198 16,395

These are the expected results under the null hypothesis. Note that the row and column sums remain the same.

Now that we have the expected values, we calculate the deviation of the observed values from the expected values, which is just a subtraction. We’ll leave out the row and column sums this time, because they don’t deviate.

Deviation Vaccine No vaccine
HIV positive -11.496 11.496
HIV negative 11.496 -11.496

These deviations are usually called errors, even though no mistakes have been made. Statisticians spend a a lot of time studying the properties of errors. For a 2×2 contingency table with a sufficiently high count in each table cell, statisticians have shown that if we square the errors,

Square error Vaccine No vaccine
HIV positive 132.16 132.16
HIV negative 132.16 132.16

standardize the square errors by dividing by the expected values,

Standardized Vaccine No vaccine
HIV positive 2.11473 2.11447
HIV negative 0.01625 0.01625

and add all the values in the table,

\chi^2 = 2.11473 + 2.11447 + 0.01625 + 0.01625 = 4.2617

we get a number known as the chi-squared statistic, so called because, under the null hypothesis, the sum of the standardized square errors is a random variable with a chi-squared distribution.

The chi-squared distribution is actually a family of distributions, parameterized by the number of degrees of freedom. For a 2×2 contingency table, the number of degrees of freedom is 1. Why is it called the degrees of freedom and why is it 1 in this case? Consider our 2×2 contingency table with none of the observed values filled in, but with the row and column sums fixed.

Observed Vaccine No vaccine Sum
HIV positive 125
HIV negative 16,270
Sum 8,197 8,198 16,395

If I were to give you just one of the missing values, you’d be able to fill in the rest of the table through subtraction from the row and column sums. In a sense, then, I am free to fill in only one of the values; all the others are contingent on that one. That’s why there’s only one degree of freedom in a 2×2 contingency table.

(If you’re wondering why the row and column sums are fixed, it’s because those sums are bound up in the null hypothesis, and we are doing our calculations under that hypothesis.)

We’re now just one step away from our answer. The probability, under the null hypothesis, that the distribution of HIV-positive test subjects would be at least as skewed as the observed results is equal to the probability that a chi-squared random variable with one degree of freedom would be larger than 4.2617. In the old days, when I was a student, we’d look this number up in a chi-squared table, but today we have more convenient options.

I used the chi2cdf function in Octave to get an answer of 0.03898. You could also use this Wolfram Alpha page to get the same answer.1 So there’s a 4% chance of getting results at least as skewed as those observed in the study. This is close to the usual 5% boundary, so it’s not surprising that reclassifying a few subjects would make the results not statistically significant.

Addendum
Here’s a transcript of the Octave session in which I made all the calculations. It takes a lot less time to do it than to explain it. Note that in some cases I’m using regular matrix operators (-, *), and in other cases I’m using element-by-element operators (.^, ./).

octave-3.2.3:1> obs = [51, 74; 8146, 8124]
obs =

     51     74
   8146   8124

octave-3.2.3:2> sum(obs,1)
ans =

   8197   8198

octave-3.2.3:3> sum(obs,2)
ans =

     125
   16270

octave-3.2.3:4> p = sum(obs,2)(1)/sum(sum(obs))
p =  0.0076243
octave-3.2.3:5> exp = [p; 1-p]*sum(obs,1)
exp =

     62.496     62.504
   8134.504   8135.496

octave-3.2.3:6> err = obs - exp
err =

  -11.496   11.496
   11.496  -11.496

octave-3.2.3:7> sqerr = err.^2
sqerr =

   132.16   132.16
   132.16   132.16

octave-3.2.3:8> stdsqerr = sqerr./exp
stdsqerr =

   2.114726   2.114468
   0.016247   0.016245

octave-3.2.3:9> chi2stat = sum(sum(stdsqerr))
chi2stat =  4.2617
octave-3.2.3:10> 1 - chi2cdf(chi2stat,1)
ans =  0.038981

Tags:

  1. Note that I used 10,000 as the “right endpoint” on the Wolfram Alpha page. What I really wanted as the right endpoint was infinity, but I couldn’t figure out how to tell Alpha that (neither “infinity” nor “Infinity” worked). Because the chi-squared density function drops off rapidly with increasing values, 10,000 was a good proxy for infinity. Even 100 would have given four digits of accuracy in the calculated probability. ↩

Post #1067 in math, programming | No Comments

Octave’s str2num function

October 10th, 2009 at 8:12 am

In the Octave script for computing the binomial distribution that I posted a few days ago, I used the str2num function to covert command-line arguments, which the argv function returns as strings, into numbers. It turns out that str2num can do more than that.

Here’s a brief Octave session demonstrating str2num’s power:

octave-3.2.3:1> str2num('5')
ans =  5
octave-3.2.3:2> str2num('5.33')
ans =  5.3300
octave-3.2.3:3> str2num('5+7')
ans =  12
octave-3.2.3:4> str2num('(5+7)/5')
ans =  2.4000
octave-3.2.3:5> str2num('sin(pi/4)')
ans =  0.70711
octave-3.2.3:6> str2num('fred')
ans = [](0x0)
octave-3.2.3:7> str2num('[1, 2].^2')
ans =

   1   4

Rather than just converting a string to a number, str2num evaluates the string as an Octave expression and then returns the result as a number. More generally, if the result is an array of numbers, it will return that array. Only if the result cannot be interpreted as a number will str2num fail, and even in that case it won’t issue an error message; Perl-like, it will just return an empty array.

This makes my binomial script even easier to use. Recall that to demonstrate it, I ran

binomial 10 .166667

to get the distribution of ones for 10 rolls of a standard six-sided die. I could have used the simpler

binomial 10 1/6

and taken advantage of the power of str2num to get the same result

 0  0.161506
 1  0.323011
 2  0.290710
 3  0.155045
 4  0.054266
 5  0.013024
 6  0.002171
 7  0.000248
 8  0.000019
 9  0.000001
10  0.000000

Strangely, the Octave documentation for str2num doesn’t mention that it does an expression evaluation, but the documentation for a similar function, str2double says:

str2double can replace str2num, but avoids the use of eval on unknown data.

which is a pretty good hint that str2num does an evaluation, but it would be nice if the docs did more than hint.

Tags:

Post #1050 in math, software | No Comments

Octave 3.2.3 for Mac

October 9th, 2009 at 9:38 am

Just a quick sidebar related to yesterday’s post: the latest version of Octave, version 3.2.3, has just been compiled for the Mac by the Octave-Forge people and made available as an easy-to-install standard Macintosh application that you just download and drag into your Applications folder. A a similar precompiled Gnuplot application—version 4.2.6, also the latest—is included in the download.

Last year, I wrote a post about configuring Octave and Gnuplot and using them from the command line. Everything I wrote in that post can still be applied to the new versions, but I have changed one aspect of my configuration. The actual executables are buried deeply within the app packages, and I used to have shell aliases,

alias octave=/Applications/Octave.app/Contents/Resources/bin/octave
alias gnuplot=/Applications/Gnuplot.app/Contents/Resources/bin/gnuplot

to make them easy to launch from the Terminal. I’ve since gotten rid of the aliases and made symbolic links to the executables,

ln -s /Applications/Octave.app/Contents/Resources/bin/octave ~/bin/octave
ln -s /Applications/Gnuplot.app/Contents/Resources/bin/gnuplot ~/bin/gnuplot

in my ~/bin/ directory, which is where I put all my self-written scripts and is in my $PATH. This means that tab completion on the octave and gnuplot names works in the bash shell and I can get scripts to run with a shebang line of

#!/usr/bin/env octave -q

or

#!/usr/bin/env gnuplot

at the top of the script.

Tags:

Post #1047 in mac, math, software | No Comments

Octave and command-line arguments

October 8th, 2009 at 12:21 pm

I needed a script that would print the values of the probability mass function for the binomial distribution for any number of trials and any probability of success in an individual trial. I wanted to pass the number of trials and the probability in on the command line and get the PMF values back in a form compatible with Gnuplot; that is, I wanted to be able type this

binomial 10 .166667

and get back this

 0  0.161505
 1  0.323011
 2  0.290710
 3  0.155046
 4  0.054266
 5  0.013024
 6  0.002171
 7  0.000248
 8  0.000019
 9  0.000001
10  0.000000

(For 10 rolls of a standard six-sided die, the probability of getting 0 ones is 16.15%, the probability of getting 1 one is 32.30%, the probability of getting 2 ones is 29.07%, and so on.)

I decided to use Octave as my scripting language because it has a binomial PMF function built in (binopdf1). I’ve written scripts in Octave in the past, but I’ve never tried to pass arguments to those scripts and wasn’t sure how to do it.

It won’t surprise anyone with C/Shell/Perl/Python/Ruby experience to learn that the Octave function that grabs the command-line arguments is called argv. It takes all the space-separated arguments and puts them in a cell array, a heterogeneous data structure comparable to a list in Perl/Python/Ruby. Cell arrays are distinguished from regular (homogeneous) arrays by the use of curly braces ({}) for indexing instead of straight brackets ([]). The arguments come in as strings, so if you want to pass in numbers (as I did), you’ll have to convert them.

Here’s the script I came up with

1:  #!/usr/bin/env octave -q
2:  
3:  n = str2num(argv(){1});
4:  p = str2num(argv(){2});
5:  x = 0:n;
6:  pdf = [x; binopdf(x, n, p)];
7:  
8:  printf("%3d  %f\n", pdf)

The shebang line includes the -q option to prevent Octave from printing its version and copyright information when it starts up. Lines 3 and 4 get the arguments and turn them into numbers. Line 5 creates the 1D array of success counts, and Line 6 uses that to create the 2D array of success counts and associated PMF values. Here you see one of Octave’s strengths: its ability to apply functions and operators to entire arrays at once—no need for map or list comprehensions. Line 8 formats and prints the results.

With the script written—and made executable via chmod +x—I could call it from within Gnuplot to generate data for plots like this

an example of how the binomial distribution tends toward the normal distribution for large numbers of trials.

Tags:

  1. Yes, Octave calls it binopdf even though it might be more appropriate to call it binopmf. I think the Octave people wanted naming consistency across all distributions, whether continuous or discrete. So all the cumulative distribution functions end with “cdf,” all the inverse CDFs end with “inv,” and all the density/mass functions end with “pdf.” ↩

Post #1045 in math, programming | No Comments

A simple PCalc tip

October 6th, 2009 at 2:51 pm

If you’re the type of person who likes to figure out how to use your calculator to make upside-down naughty words, you probably heard about PCalc’s latest revision. I didn’t write about PCalc 1.8 when it was released because—apart from its little dig at the App Store’s age-rating system1—there weren’t many new features worth writing about. Plus, upside-down words only work when if you use PCalc’s LCD font, and I always use Easier To Read.

The new release did remind me, though, that I have for some time combined two of PCalc’s features to make its extensive function set more readily available. I haven’t seen anyone else mention this simple trick, so I figured it was worth a quick post.

PCalc has several different keyboard layouts with different sets of functions available in each layout. Normally, I work with the Engineering layout, which has the power, trig, and exponential functions usually associated with a scientific calculator. Occasionally, though, I need the hex conversion and bit functions found in the Programming layout. In the past, I’d go into the PCalc preferences to change the layout when I needed a Programming feature, then change back to Engineering when I was done. Eventually, I realized this was stupid and set up PCalc to save me some time.

Like many iPhone apps, PCalc can be used in either the vertical or horizontal position, and the keyboard layouts for these two positions are independent of one another. So now I have the vertical layout set to Engineering

and the horizontal layout set to Programming

and whenever I need to do something in hex, I just turn the phone sideways.

Am I embarrassed that I didn’t see this opportunity right away? Yes. I was probably distracted by the 5318008.

Tags:

  1. To James Thomson’s horror, this guy didn’t get the joke and made it look like Thomson himself was the bluenose. ↩

Post #1041 in iphone, math, software | No Comments

PCalc 1.7

July 19th, 2009 at 12:50 am

PCalc 1.7 just hit the App Store, and it is—despite the meager 0.1 increment in version number—a major update to the iPhone’s premier calculator application. The new and/or improved features (in order of importance1) are:

  1. Easier to Read digits that are truly easy to read at all sizes.
  2. Multiple memory locations.
  3. Buttons for Stack Clear (SC) and Memory Clear (MC).
  4. A display that can show up to four lines of data.
  5. A faster startup.
  6. More logical Enter behavior.
  7. More unit conversions.
  8. A setting that prevents the calculator from going to horizontal mode when the phone is rotated.

We’ll look at each of these features in turn. Let me first point out that I’ve written several earlier posts on PCalc. Most of what I wrote in my series of Definitive Posts on PCalc (part 1, part 2, part 3, and part 4) and my descriptions of the version 1.5 and version 1.6 updates still holds and won’t be repeated here. This will be just the new stuff. I should also mention that I use PCalc’s RPN mode exclusively—improvements to its algebraic mode fly under my radar, and I can’t comment on them.

Easier to Read digits

PCalc can display digits in a manner reminiscent of older physical calculators, imitating either LCD elements or chunky pixels. Because I am in my late 40s and have the faltering eyesight of middle-age, I have always prefered the Easier to Read digits, even though they weren’t as easy to read as i would have liked.

Actually, the digits themselves were pretty easy to read, it was the tiny comma and decimal point that I had trouble with. James Thomson, PCalc’s developer, bulked up the punctuation several point releases ago, which worked pretty well; but when he later introduced the two-line display, the punctuation shrunk back down to a size I found hard to see.

But with 1.7, Thomson has introduced fonts with big, robust commas and decimal points that aren’t squeezed tightly in between the numerals.

The new display can pass the most stringent of visibility tests: that of a tired middle-aged man trying to do a few calculations on a late night flight.

Multiple memory locations

PCalc has always had a memory slot, and it’s always had a stack that can be loaded up with lots of values. But most serious scientific calculators have serveral memory slots that are quickly addressable. Now PCalc joins them with ten memory slots, labeled 0–9. To store a value in memory you

  1. Type in the number.
  2. Press the M in button. Note how most of the buttons are dimmed; pressing them has no effect.
  3. Press one of the digit buttons.

The value is now stored in that memory slot. To retrieve a value from memory, you

  1. Press the M re button.
  2. Press one of the digit buttons.

The value will be copied from memory and pushed onto the stack.

NB: if you retrieve from a memory slot that hasn’t had anything stored in it, zero will be pushed onto the stack. As I consider zero to be distinctly different from empty, I see this as a flaw. Fortunately, it’s a flaw that does no damage. No items on the stack are destroyed by errant zero, and you can always return to the state before the inadvertent retrieval by Dropping the zero.

New Clear buttons

After several calculations, my stack is usually cluttered with old results, and it would be nice to blow them away with a single command. PCalc has the AC (All Clear) button, but pressing it has the unfortunate effect of clearing the memory in addition to the stack. This is an even more destructive operation now that there are 10 memory slots. Fortunately, PCalc 1.7 comes with a new SC (Stack Clear) button that will clean your stack while leaving your memories intact. It’s accessed through the 2nd key.

For you Phil Dick fans, 1.7 also has an MC (Memory Clear) button that does just the opposite: erases all your memories while leaving your stack in place.

Multi-line display

PCalc has had an optional two-line display since version 1.2. Now you can chose to have up to four lines in the display, and the lines don’t have to be values on the stack. The settings let you choose from many options.

Right now I’m using a three-line display of the stack.

Faster startup

If you follow James Thomson on Twitter (@jamesthomson), you know that he spent a lot of time working on PCalc’s code to get it to launch more quickly. He succeeded. While PCalc certainly doesn’t start as quickly as Apple’s standard iPhone calculator, it’s become quite peppy. And it presents a cute aphorism to amuse (i.e., distract) you as it gets going.

There are many of these, all with a math or calculator theme.

New Enter behavior

Turning on the HP48 Style RPN option (it’s in the Advanced settings) eliminates a behavior that’s been bothering me for ages: the duplication of the x-register into the y-register when the Enter key is pushed. Here’s the scenario:

This is batty. Pressing Enter should just enter it (hence the name) in the x-register—there shouldn’t be a second copy of the value in the next stack position. I know that some calculators act this way, but that doesn’t make it right.

If you turn on the HP48 Style RPN setting, PCalc 1.7 does the right thing: it enters the number in the x-register and waits for you to type another number before pushing the entered value into the y-register.

New unit conversions

There are now US and UK versions of the teaspoon, tablespoon, and fluid ounce so US cooks aren’t confused by a cup with 8.11537 fluid ounces.

Locked orientation

The new Do Not Rotate setting will lock PCalc in portrait mode, which can be handy if you find yourself holding the calculator in a variety of orientations as you move around, write, adjust an instrument, etc. Locking isn’t a big deal for me, but I can see where it would be very useful for some people.

It should be noted that you cannot lock PCalc in the horizontal mode. Although you can turn on the Do Not Rotate setting while the calculator is horizontal, it will shift to portrait mode the next time the calculator is turned vertical and will then be locked in that position.

Finally

I bought PCalc at version 1.0.x, and it’s been my calculator ever since. I use it more than my physical HP 35s calculator, even when I’m at my desk and the HP is within reach. It’s annoying that the 35s refuses to mutate when I press the ↰ and ↱ keys.

Tags:

  1. In order of importance to me, that is, which is all that matters. ↩

Post #864 in iphone, math, software | No Comments

« Older Entries