The basic MOC algorithm for a system

The Method of Characteristics - MOC

To get the full understanding, MOC should be explained in more detail. That's for a later post. For the time being I just refer to Wylie and Streeter and Wikipedia:

In mathematics, the method of characteristics is a technique for solving particular partial differential equations. Typically, it applies to first-order equations, though in general characteristic curves can also be found for hyperbolic and parabolic partial differential equation. The method is to reduce a partial differential equation (PDE) to a family of ordinary differential equations (ODEs) along which the solution can be integrated from some initial data given on a suitable hypersurface.

Explained below is the basic system algorithm for MOC in a system, not MOC itself. It's how it is programmed in LVTrans, and how it will be programmed in any eventual future versions.

System of hydraulic pipes

Each hydraulic pipe is modelled with MOC. A Pipe has a length and two ends. Each of these ends must be connected to a "non pipe" element, called NP. This NP element is typically a valve, a T-joint, a turbine, a reservoir but can also simply be an open or closed end of the Pipe. Thus for each Pipe, the NP represent the spatial boundary elements at each end.

Each NP can have one or more boundaries. These boundaries are typically Pipes, but can also be other NP elements. An example is a PID for a turbine. Many of these NP elements also have natural boundaries. For instance, a reservoir can be connected to one or several pipes but also have a surface level that is relatively constant seen in relation to waterhammer transients. This surface level is a natural boundary condition. The surface level can also vary, and can be measured and given as input for a digital twin.

The figure above shows a simple system of Pipes and NPs. Each Pipe has an arrow towards the NP elements. Why becomes apparent soon.

Basic element algorithm

Each pipe is divided into smaller sections. For each time step, the H and Q (head and flow) for each section is calculated based on the values for the previous time step. The boundaries are yet unknown, but the C+ and C- characteristics can be calculated. This is how far each Pipe can be calculated alone. All the internal values and the C+ and C+ characteristics at the boundaries. What is missing is the H and Q at each end.

Each NP element do almost the same. They start with the values from previous time step and the C+ and C- characteristics from the connecting pipes that just were calculated. All the inner variables and the H and Q at the boundaries can be fully resolved. These H and Q are then sent back to each pipe and make their calculations complete. Now the time step is fully resolved everywhere. The procedure is repeated for the next time step and so on.

Note that NP elements normally use Runge-Kutta or Newton-Raphson methods.

NP elements also can have manual inputs. For instance the opening of a valve. This is also updated for each time step.

MOC System algorithm

From the figure the full arrows are the C+ and C- characteristics at each end of the pipe. The back propagation of H and Q at the boundaries are the dotted lines. The PID does not send H, Q or C, but rather sends A and receives B, where A and B typically but not necessarily, are guide vane opening and turbine speed.

The MOC system algorithm therefore becomes:

1. For all Pipes and (PID).
1. Calculate the internal variables; H Q and predict the A in (PID) for T1 based on T0
2. Send C+, C- and A to each connected NP element
2. For all NP elements
1.  Calculate all internal variables and the boundary values for T1 based on T0 and values received from Pipe and (PID)
2. Iterate each NP-PID pair once more to correct the predicted A if needed*.
3. Send H, Q and B to each Pipe and (PID)
3. Update the boundary H and Q for each pipe and B for the (PID)
4. Start at 1 again for the next time step.

Programmatically there are many ways to to this depending on choice and the programming language. In "native" LVTrans these connections happens almost by itself. It's the way LabVIEW works by default. In a text based programming language this can be done by lists or arrays of Pipes and NPs. A major point is that the inner workings and transfer of values of each element is the same no matter how this is done. An independent library of MOC can thus be made. Here independent means independent from the actual implementation of the MOC system algorithm. 

*Note that the (PID) is not a Pipe with MOC inside. It does not have analytically correct C+ and C- characteristics for T1. The initial A it sends is a forward in time predicted value at T1 based on the time step T0 for all variables in (PID) including B. The B which is sent back from NP, is based on this predicted value of A. This B is used to update all variables inside (PID) and can be used to correct the predicted value of A at T1. This corrected value can be sent once more for a new iteration in the connected NP. Default is a predictor method. One extra iteration becomes a predictor-corrector method. An error analysis of this is found in the LVTrans manual, and is a part of choosing the correct time step.

It must be noted that the initial "error" comes from the fact that the B used when calculating the PID at T1 comes from the previous time step T0. This is also very much what happens in real life. The update frequency of a real PID is seldom better than 0.01 Hz for practical purposes. There's no reason it should be faster than this, because the physical frequencies for which it shall control is already several decades lower. What all this means is that with a time step in the simulation of roughly 0.05 s or lower, none of this matters - at all. This is also explained in detail in the LVTrans manual (I may also write a blog post at some time).  

In Wylie and Streeter the procedure is to gather all connected NP elements together and iterate them as one single block with regard to connected Pipe elements. In practice this is equivalent to an NP-NP pair iteration. It's a matter of choice, but the method of Wylie and Streeter may be more practical for blocks of NPs with more than two NPs. Such blocks never happens in LVTrans. With the exception of PID elements, it's probably more practical to create separate elements where larger collections of NP blocks are are pre-made.  

Conclusion

The MOC system algorithm is very simple. MOC is also superfast. Orders of magnitude faster than anything else while being analytically correct and with no numerical artifacts. With an object oriented approach for each elements, the same core library can be used in several different applications. This is very useful when different use cases are needed.

Ponderings about the future of open source transient simulation software

Intro

It's more than 25 years ago since LVTrans started. Lots have happened during that time. Neither AI or Python existed in some generally usable form 25 years ago, but today AI and Python are everywhere along with many other new and modern programming languages.

My first programming was done with interpreted BASIC on a Commodore 64. At computer courses at the university I learned Pascal and then FORTRAN. I liked them both. Then came C and C++ and ruled the world. I wrote a FEM software in C++ for my PhD, linking in FORTRAN matrix library. Compilers were expensive stuff. I used Watcom, which was the only one that really could combine C/C++ and FORTRAN. Today compilers are free and open source and can link in whatever other language you want.

C/C++ and FORTRAN

Even though FORTRAN has received more attention lately due to vector programming and AI/machine learning, it's still very specialized. Perfect for "hard numerical computing" but not particularly suited for anything else. C/C++ is still popular and can be used for everything. Compared with newer languages, C/C++ is complicated and not particularly suited for the casual programmer or engineers cheating in the field. I really like the cheer power, the brute force, of C/C++ but for the programming stuff I actually do, I always go in another direction, a simpler direction that gets the job done faster.

High level specialization

Another kind of programming has developed during the years, partly from academia, partly from the industry and goes way back. This is Matlab/Simulink, LabVIEW, Modelica, Scilab/Xcos and so on. Of these it's Labview first and foremost I have used but I have also used Scilab/Xcos a lot. Then there's Excel with Visual Basic/C or whatever "Visual" it has these days. It's a very powerful tool for a wide variety of tasks. All of them are perfect tools for the casual programmer and engineer who don't want to make an app, but rather want to get solutions to engineering problems.

Python

In comes Python. It has existed since long, but remained in the shadows for a long time. The largest effect in engineering and academia is that it has made Matlab etc. much less relevant. Matlab has still some considerable momentum due to its sheer size (and perhaps it's similarities with Python), but LabVIEW has gone backwards. It's no longer in the top 50 at the TIOBE index. Ladder Logic (PLC) and COBOL are much more popular. I never envisioned LabVIEW to become anything like C/C++ in popularity, but didn't really envision the ill fate it has received either. An open source program on a programming language that hardly exists anymore is a real tough hill to climb.

What to do?

The question is what to do. Does anything need to be done? LabVIEW is not likely to simply disappear from the face of the earth, but it will perhaps shrink in generality and focus more on the stuff it's actually designed to do. The choice of LabVIEW for LVTrans was more of a coincidence anyway. I was doing some LabVIEW work in the lab, and was asked to do a transient analysis of acoustics in gas risers. I said yes, and got myself a license of the commercial transient simulation software Flowmaster. It turned out that Flowmaster couldn't handle the particular problem, so I made a quick MOC solver for one pipe in LabVIEW that could. The Method of Characteristics (MOC) is very easy to program, for one pipe at least. The rest is history.

There have been several good points about this:

• LabVIEW is compiled to machine code, and is relatively fast. About 1/3-1/2 the speed of C in my experience. Sometimes even directly comparable to C.
• All graphics, graphs and charts are already there, professional looking and ready to use.
• The LabVIEW programming environment is used directly in setting up systems. It's very functional and easy to use, and additionally sets and organize the topology of the system.

This has enabled me to work as an engineer, doing lots of other stuff, while still maintaining and developing LVTrans, and using it in my work. There is no way this would have been possible in C/C++ for instance. Using Matlab (or Python for that matter) would be too slow execution vise. Last, but not least, LVTrans is also used as the numerical core in digital twins at Aker Solutions This is due to its speed and the native embeddedness of LabVIEW. It was just a matter of compiling it down to one of National Instrument's boxes, very simplistically speaking. 

It's safe to say that LabVIEW has served me well. It has enabled me to get the job done just perfectly. I have never been super happy with the solution though. The open source part has never taken off to put it mildly. Closed proprietary binary source code is no boost for open source (have to put in a "duh" here :-) ) 12-13 years ago I started making a C++ library. I never finished that, but I still have it for the core MOC solver. In the beginning, LVTrans actually had the core MOC solver programmed in (pure) C for 3-4 years. The long term plan was always to finish the C++ library and move on from there. Due to being busy with work, it never happened. Linking in the C++ library as DLL, LVTrans will be the same. LabVIEW would be used as the HMI only, while the numerical core could move on to other projects. It would be free and truly open source.

C++ or? 

Is C++ the right way forward today? In some ways yes, in other ways no. It's perfect for the job, but C++ is bloated and overly complicated for this simple job. There are other languages today that could be much better suited. Even professional programmers find C++ to be overkill for most of the work they do, or simply not suited and not efficient enough. I don't need efficiency as such, but there is a strong correlation between efficiency and having the right tool for the job.

This is where I'm right now. Is C++ the right tool for the job today, and if not, which other options are there. One of my sons is a lead programmer in IT now. He has opened my eyes to other stuff. Lots have happened in the last 25 years, lots more than just Python and AI. That's for another post.