More Perspective on Software Unit Testing for Scientific Programming

Each unit test is a fork in the road - and your job isn't done until the right turn is taken every time.

Each unit test is a fork in the road – and your job isn’t done until the right turn is taken every time.

A few days ago I published a post on software unit testing in MATLAB. I’d figure I’d share more perspective on what I have learned so far, in terms of strategy and tactics for software testing, especially as it pertains to scientific and engineering applications. For reference, this presentation, this presentation, and this paper have been informative. I also direct the reader’s attention to this paper in PLOS Biology from January 2014, titled “Best Practices for Scientific Computing“.

The Tougher Path

You see friends, we scientists and engineers have it tougher than the average programmer. With typical programming, the inputs and outputs to a function are well defined. Expected outputs are known well in advance. Heck, there is even a software engineering paradigm called “test-driven development“, where the actual unit tests for the code are written before the actual code is written! We also are not exposed to much software engineering methods in school, and most of us just “picked up” programming and tried to roll with the punches. I would complain we are also babied too much with MATLAB, and that engineering schools would better serve their students by switching to teaching Python and its scientific computing add-ons.

In complex scientific and engineering applications, it can be very difficult to know what your program is supposed to be spitting out, without becoming bogged down by huge numbers of over-specific test cases. We also encounter small (but unavoidable) errors due to the numerical precision of our machine. Solution algorithms also have a given amount of error attached to their solution (for example, Runge-Kutta integration possesses error O(h^4)). Sometimes, we even use a stochastic algorithm to solve a problem, making the results, well… random! How do you debug an algorithm with random inputs? These hurdles make it difficult to know if a program is operating correctly but is simply overcome with numerical imprecision, or is genuinely bugged.

More Ideas on Unit Testing Scientific and Engineering Simulation Codes: The example of a mixing point

As an example, consider a mixing point that is mixing two streams of the same material at two different temperatures and two different mass flowrates, and the assumption of constant specific heats cannot be made. What is the outlet temperature? To solve this problem, a nonlinear algebraic equation needs to be solved for T_{out}:

m_1\int_{T_{ref}}^{T_1} C_pdT + m_2\int_{T_{ref}}^{T_2} C_pdT = (m_1+m_2)\int_{T_{ref}}^{T_{out}} C_pdT

Where C_p is a heat capacity polynomial.

From what we know about how heat transfer works, there is an upper and lower bound on what T_{out} can be. The lower bound is min(T_1,T_2), and the upper bound is max(T_1,T_2). Bracketing the root between these bounds is guaranteed to find the root, e.g. using the regula falsi method. The procedure to solve this problem is fairly straightforward: use fzero() in MATLAB and bracket the root.

But let’s think a little further. What should happen if we alter the two mass flowrates, m_1 and m_2? Logically, we should observe the following asymptotic behavior:

  1. As m_1 tends to infinity, T_{out} should asymptotically approach T_1.
  2. As m_2 tends to infinity, T_{out} should asymptotically approach T_2.
  3. If m_1=0, T_{out} = T_2 exactly.
  4. If m_2=0, T_{out} = T_1 exactly.
  5. If m_1=m_2, over small temperature differences, T_{out} \approx (T_1+T_2)/2 (if we had constant specific heats, this equation would be exact).
  6. The choice of T_{ref} should not change the results of cases 1-5 whatsoever.

These 6 items represent 6 test cases that a code would need to pass. Passing all these tests wouldn’t prove the code right, but failing any of these tests would prove the code wrong. Furthermore, the code would need to be dummy-proofed to ensure m_1 \geq 0 and m_2 \geq 0. It would also need to be dummy-proofed to accept the correct range of temperatures over which the C_p expression has good accuracy.

Some Questions for Yourself

  1. How might you generalize the solution of the problem for the case of N inlet streams, with N specific temperatures and mass flowrates?
  2. What should the bracketing bounds be? What should the asymptotic behavior be with regards to the inputs?
  3. How should we handle the case of all flowrates being equal to zero?
  4. Now let’s make the problem even more complicated: What if we have N different substances being mixed, each with a different heat capacity? What are the new bounds on T_{out}? Do they change?

More Ideas on Testing Scientific Codes

Unfortunately, the testing of scientific and engineering computer codes is something of an art, and as such, I can only brainstorm general advice. The details of your specific code and specific problem will probably generate the bulk of your unit test cases. Again, passing your test cases doesn’t really prove your code is working correctly – but failing them proves it is working incorrectly.

  1. Should something ALWAYS be equal to something else?
  2. Should a particular input be INVARIANT with certain output? In the mixing point example above, T_{out} should be totally invariant with T_{ref}.
  3. Think up a simple approximation to your problem. Does your code (which solves the exact problem) come “close” to the approximation, or are you way off? See point 5 in the mixing point example above.
  4. Do you have access to experimental data? Can you reproduce it to within experimental uncertainty?
  5. If you are testing an optimization algorithm, can your algorithm identify the global optima for the numerous toy problems encountered in the literature?
  6. What should happen as one of your inputs grows extremely large or extremely small? What is the asymptotic behavior of the inputs to your code? Does it line up with physical reality, reason, and logic?
  7. Are there upper and lower bounds?
  8. Does your code produce physically unreasonable behavior, e.g. does it violate the 1st or 2nd laws of thermodynamics, or does not balance momentum?
  9. Is mass and energy always conserved in the system you are simulating, to within numerical accuracy?
  10. What if your user is a bonehead, and enters zero in for ALL the inputs?
  11. Can the user enter inputs that violate the assumptions of your model? For example, if you are solving a laminar flow problem, can the user insert inputs that make the Reynolds number unreasonably high?
  12. Does a closed-form solution exist for a particular case? Do your numbers match the closed-form solution?

Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s