January 19, 2011

Interpolation

To check for maximum bending stress in flat plates — this for a hanger clamp, Roark’s formulas are handy. Here’s a table that offers a recipe.

Manually, this is straightforward: First, I calculate the aspect ratio (a / b) of the flat plate, and then pick a value of β1 that corresponds to the ratio. If I do not find a ready match, then I perform a linear interpolation between the two values that form the intermediate lower and upper boundaries of (a / b).

For example, if my plate size is 200×100, then (a / b) is 2.0 (in the sixth column), whose corresponding β1 would then be equal to 1.226.

But if my plate size is 220×100, then (a / b) becomes 2.2. Now, the nearest lower and upper boundary values in the table are 2.0 and 3.0 respectively. Since Roark’s table does not readily offer a corresponding β1 value, this further requires interpolation — let’s stick to linear, for simplicity.

  p / q = P / Q
      p = (P . q) / Q
        = (2.106 - 1.226) . (2.2 - 2.0) / (3.0 - 2.0)

     β1 = 1.226 + p

How to automate this in python? Here’s one way:

from bisect import bisect
# Roark's Table 11.4, Case 10, with three edges fixed
# a / b, and b1 are defined as lists below:
raob = [0.25, 0.50, 0.75, 1.00, 1.50, 2.00, 3.00]
b1 = [0.020, 0.081, 0.173, 0.321, 0.727, 1.226, 2.105]

a, b = input("Enter a, b: ")
a, b = float(a), float(b)
aob = a / b

try:
    i = raob.index(aob)
    beta1 = b1[i]
except ValueError:
    # aob did not match any value in raob
    i = bisect(raob, aob)
    j = i - 1
    P = b1[i] - b1[j] 
    q = (aob - raob[j])
    Q = (raob[i] - raob[j])
    beta1 =  b1[j] + (P * q / Q)
finally:
    print "a/b               : ", round(aob, 3)
    print "Index of a/b ratio: ", i
    print "beta1             : ", round(beta1, 3)
pass

And here’s how the above code works:

  1. First, I import bisect standard library for referencing across two lists (raob, and b1).
  2. Get user input for plate size in terms of a (length) and b (width).
  3. Then, in the try statement, I attempt to find the index, i (column number) that matches the aspect ratio (a / b).1 Finding the index is key, which lets me refer to the corresponding value in the second list.
  4. If an exact numerical match of (a / b) ratio is not found in the raob list, then the control jumps to process the except ValueError part.
  5. Using bisect, I try and find the nearest next numerical value and its corresponding index number. So to take the aforementioned example, if my (a / b) ratio turns out to be 2.2, then the next numerical value that bisect finds for me is 3.0, and its corresponding i value, which is 6 in this case.2
  6. Rest of the code thereafter deals with linear interpolation between two limits (2.0 and 3.0 in this case, and their corresponding β1 values) to get the β1 value (corresponding to 2.2 in this example); and then print those values.

The only shortcut I’ve taken in this crude — doesn’t-catch-all-errors — code yet I think, is in the exception part. I’ve used that to move the control from simply throwing up ValueError to the next part of the code where the interpolation occurs. (Note: I haven’t addressed the IndexError that would occur — when the ratio is either lower than 0.25 or greater than 3.0 — in this code.) Until I work out a more elegant way to do the above, this I think will do for now.


  1. For instance, the index of 0.75 in the raob list is 2. (Index numbers in python start with 0, not 1.) 

  2. For tracing index, and the subsequent interpolation, I preferred using bisect from the standard library (instead of using a custom package like numpy that comes with better interpolation tools), so I could keep the number of dependencies (and package installation requirements) to a minimum in order to run this code. Without having to resort to using a sledge hammer to kill a mosquito, I’m happy with it so far.