Support for Dymola in the Modeling and Simulation of Physical Systems with Distributed Parameters

1. Support for Dymola in the Modeling and Simulation of Physical Systems with Distributed Parameters Farid Dshabarow ABB Turbo Systems, Switzerland François E. Cellier ETH Zürich, Switzerland Dirk Zimmer ETH Zürich, Switzerland

2. In the 1970s, a number of general-purpose PDE solvers for parabolic and hyperbolic PDEs in one and two space dimensions were developed. The most versatile among those was FORSIM VI, a FORTRAN library for simulating such systems. FORSIM was based on the Method of Lines, a technique that discretizes the spatial axes of a PDE, while leaving the time axis continuous, thereby converting a PDE into a set of ODEs. History of General Purpose PDE Solvers

3. FORSIM VI offered a set of FORTRAN subroutines that the user could call for computing approximations of spatial derivatives and for integrating the resulting ODE system over time. The tool was versatile but not user-friendly. It was versatile, because it offered the user all of the tools that s/he needed for solving a PDE problem, but it wasn’t user-friendly, because it forced the user to manually decompose the problem down to a fairly low level. Other tools, such as PDEL and LEANS III offered higher-level interfaces, but this came at the price of these tools being more specialized tools that could only solve problems belonging to predefined sets of PDEs. History of General Purpose PDE Solvers II

4. In the long run, none of these tools survived, because engineers wanted to solve more and more complex problems, and neither of these tools were versatile and/or efficient enough to do so. Consequently, most engineers and physicists created their own spaghetti codes written for solving one particular PDE problem only in as efficient a manner as possible. Evidently, this caused a lot of duplication of very similar code segments that were rewritten by many programmers over and over again, and ultimately, it must be recognized that this approach is wasteful in terms of human resource utilization. History of General Purpose PDE Solvers III

5. Only quite recently, a renaissance of general-purpose PDE codes could be observed. The best and most advanced among these second-generation codes is FEMLAB (COMSOL). FEMLAB essentially copied the approach taken earlier by PDEL and LEANS III, i.e., it offers solutions to a pre-defined set of PDE problems that are being parameterized to allow a user access to this tool who knows very little about the properties of numerical PDE solvers. Compared to the earlier tools, FEMLAB offers a much improved user interface with neat graphing capabilities. The code has been rather successful and is being widely accepted as the state of the art in numerical PDE technology. History of General Purpose PDE Solvers IV

6. We chose to take the other route and re-implement facets of FORSIM VI in Dymola/Modelica, making use of the built-in matrix manipulation capabilities of Modelica and the graphical user interface of Dymola. In this way, simple PDE problems can be coded graphically without forcing the user to write any programming code directly. Of course, we pay a heavy price in that users of PDELib will have to be somewhat knowledgeable about the algorithms underlying numerical PDE solvers. PDELib

7. Yet, the object-oriented nature of Modelica allows to wrap PDELib codes, offering a “mathematical” user interface, into higher-level abstractions, promoting a “physical” user interface. At that higher level of abstraction, users won’t need to understand the underlying numerical algorithms any longer. Yet, as these higher-level codes will be mapped onto PDELib first, maintaining these higher-level codes will be much easier than if these higher-level codes were programmed in C++ directly, for example. Due to the symbolic formulae manipulation techniques offered by tools such as Dymola, the resulting codes can be made to run as efficiently as any hand-written spaghetti codes of the past. PDELib II

8. Method of Lines • Let us explain by means of a simple example how the method of lines operates. • Given the 1D diffusion equation: • The method of lines discretizes the space dimension: leading to:

9. Method of Lines II • Programmed in PDELib:

10. Simulation Results • The spatial axis had been discretized into 40 segments. • The numerical solution is compared to the analytical solution. • With 40 segments, the results are not yet highly accurate. • Dymola doesn’t offer nifty 3D graphics capabilities yet. • The results could have been exported to MATLAB for graphing, but we chose not to do so.

11. Finite Volume Method • A second method implemented in PDELib is the Finite Volume method. It works as follows. The space gets subdivided into a finite number of cells. The ith cell can be written as: • The average value of the solution u(t) for this cell can be written as: • Mass conservation can be written as: • where f(t) is the flux function.

12. Finite Volume Method II • We integrate: • from t to t+t, and divide by t and by x:

13. Finite Volume Method III • With the abbreviation: • we can rewrite as follows: • which can be reinterpreted as: • In this way, we have once again converted a PDE to a set of ODEs. The only remaining question is how we approximate the average flux values.

14. Advection Equation advection equation initial condition boundary condition Method of Lines Finite Volume Method

15. Simulation Results Finite Volume Method Method of Lines

16. Burger’s Equation Burger’s equation initial condition boundary conditions Method of Lines Finite Volume Method

17. Simulation Results • The FVM solution is less accurate, but remains numerically stable for an arbitrarily large number of segments. • Hence a numerically accurate solution can be obtained using sufficiently many segments. Finite Volume Method Method of Lines • The MoL solution is more accurate, but less stable. • With 10 segments, the solution remains stable, but is very inaccurate. • With 20 segments, it is more accurate, but turns unstable after 0.6 sec. • With 40 segments, it is always unstable.

18. The method of lines works well for parabolic PDEs. It has a tendency to lead to numerically unstable solutions in the case of hyperbolic PDEs, unless a geometric integration algorithm is being used (not currently available in Dymola). The MoL solution of hyperbolic PDEs can frequently be stabilized also by use of an upwind discretization scheme. The finite volume method is not as straightforward. However, it can be stabilized more easily. The FVM was designed for hyperbolic PDEs primarily. Whereas the method can sometimes be abused to solve parabolic PDEs as well, the user has no reason to do so. MoL vs. FVM

19. The first release of PDELib only considers parabolic and hyperbolic PDEs in a single space dimension. Later versions will extend the tool to 2D and 3D problems also. However, this requires finite element approximations in space, which haven’t been implemented yet (similar to what FEMLAB does). Whereas 1D problems lead to ODE sets in a few dozen ODEs, 2D problems quickly lead to thousands of ODEs, and 3D problems almost immediately lead to problems with several hundreds of thousands of ODEs. Hence efficiency matters. Boundary conditions are much more difficult to satisfy in 2D and 3D problems. Conclusions

20. At the current time, PDELib is only a test bed implementation. It can be used for solving toy problems, and also for teaching purposes, but not yet for solving real engineering problems as they need to be solved in practice. Future versions shall enhance the tool, as well as offer higher-level (physical layer) abstractions for dealing with complex problems in a more user-friendly fashion. Conclusions II

21. Dshabarow, F. (2007), Support for Dymola in the Modeling and Simulation of Physical Systems with Distributed Parameters, MS Thesis, Department of Computer Science, ETH Zurich, Switzerland. References