Saturday, March 23, 2013

Good set of Bounds for Colebrook-White Numerical Solution

The value for fmax from my last post, it turns out, isn't as good as what it could be.  On further investigation, it turns out that for some combinations of roughness values and Reynolds values it still fails.  I found a solution, which is actually more simple, so I'm making this post to describe it.  To begin with, I will keep the expression for the minimum bound I gave in the last post.  The following is rearranged a little bit.

Now, to find a better expression for fmax, I realized there was a better approach attainable by rearranging the Colebrook-White equation.

The plot:

  • LHS limits to 1
  • RHS limits to Rr/3.7
This looks a little more complicated, but the critical difference is that the left hand side asymptotes to a constant value, just like the right side.  That's helpful, because the fundamental error in my last post was trying to fit a square peg in a round hole by replacing a non-converging function with a converging function.  Here, I don't have to do that.

In order to get an upper-bounded expression for f, the LHS only needs to be replaced with a function that is absolutely lower than itself.  Then not only will that produce a value that is always necessarily larger than the correct solution for f, but it will also be very well bounded, because it's replacing one asymptotic function with another similarly asymptotic function.  The obvious replacement is to replace the LHS with a function of the form of 1-1/sqrt(f), to match the RHS.  The ultimate form I went with was:

Plugging into the last form of the Colebrook-White equation gives:

It follows that the solution to this equation is an absolute upper limit for the equation.  The actual algebra is fairly trivial, and gives:

The similarity to the equation for the minimum is obvious.  In fact, this is a fairly elegant solution to the problem of setting bounds.  It turns out to work fairly well.  The edited code is:
Function Colebrook(Rr As Double, Re As Double) As Double
  Dim f_min As Double, f_max As Double, f As Double
  Dim g_lower As Double, g_middle As Double
  Dim n As Integer

  f_min = (2.51 / Re) ^ 2 * (1 - Rr / 3.7) ^ (-2)
  f_max = ((2.51 / Re + Log(10) / 2) / (1 - Rr / 3.7)) ^ 2
  n = 0

    n = n + 1
    f = (f_min + f_max) / 2
    g_middle = -2 / Log(10) * Log(Rr / 3.7 + 2.51 / Re * f ^ (-1 / 2)) - f ^ (-1 / 2)
    g_lower = -2 / Log(10) * Log(Rr / 3.7 + 2.51 / Re * f_min ^ (-1 / 2)) - f_min ^ (-1 / 2)
    If (g_middle = 0 Or (f_max - f_min) / 2 < 0.000001) Then
      Exit Do
      If (g_middle * g_lower > 0) Then
        f_min = f
        f_max = f
      End If
    End If
    Colebrook = f
    If (n > 100) Then
      Colebrook = CVErr(xlErrNA)
      Exit Do
    End If
End Function
This begs one more question.  Could you use this general form to estimate the friction coefficient itself?  Yes, you could.  Notice the only difference between fmin and fmax is only the one constant there.  I wanted to figure out a value of this costant that would approximate the answer, as opposed to setting a minimum or maximum.  To do this, I required that the limit as Reynolds number goes to infinity to approximate the real answer.  This could be used to make a more efficient Newton's method for solution.  I didn't do this, but I will share my expression for the expression.

This isn't a great approximation by the way.  In the transition areas where it really matters, it can easily be off by 60%.  In the laminar and turbulent limits it approaches the correct solution.  Still, that gives it some utility for numerical approaches.

No comments:

Post a Comment