Foundations for Continuous Time Hierarchical Co-simulation Cláudio Gomes University of Antwerp claudio.gomes@uantwerpen.be Abstract having to be made early in the development of each Complex systems have to decomposed into sub-systems component. which are developed by specialized teams. Modeling Modeling and simulation techniques are used to mit- and simulation techniques help each team achieve a igate these issues: models of components are created locally optimal solution but they fail to help all teams and simulated before any physical prototype is built. achieve a globally optimal solution. This is due to each Simulation can also be used to analyze the behavior of team having its own models made in its own tools, interacting models of components, created with differ- and external suppliers having intellectual property. The ent languages [10]. The simulation of interactions be- result is that it is difficult to simulate the coupled tween models specified in different languages is an open system. challenge [16], mostly done with a small number of sim- Co-simulation is proposed as a way to enable such ulators with similar features, making it difficult to gen- simulation. Simulators communicate and collaborate as eralize to other scenarios. To aggravate, specialized sup- black boxes. The technique has been used in a number pliers of components are interested in protecting their of domains improving the overall cost of engineered sys- Intellectual Property (IP) leading to the situation were tems. In these use cases, a minimum common denomina- the simulation needs to be made without access to the tor is assumed for the capabilities of simulators. Lever- models [6]. aging the optional capabilities of simulators improves Co-simulation is a technique to couple multiple sim- the performance/accuracy of a co-simulation, but the ulators, each simulating a single component, often seen number of decisions the orchestration mechanism has as a black box, in order to simulate the whole system. to make grows exponentially. Different, conflicting con- An orchestration algorithm is responsible for ensuring cerns have to be addressed in an optimal way. the communication between the simulators. The lack of We propose a way to deal with each concern inde- information about the simulators makes co-simulation a pendently and possibly reuse existing orchestration al- hard challenge which must be addressed in a systematic gorithms that perform better with respect to the con- way. Despite this, co-simulation shortens development flicting concerns. Our approach leverages Model Driven time and improves quality, as reported by the industrial Development techniques to process a co-simulation sce- partners of the DESTECS project [7, 11, 26]. For other nario into a canonical/trivial version, where fewer deci- application domains where co-simulation has been ap- sions remain to be made. Along the process, conflicting plied, see [19, 23, 25, 31], for example. concerns will be addressed as an optimization problem, Recognizing the potential of co-simulation, the Func- and an appropriate cost function identified. If success- tional Mockup Interface (FMI) standard [6] was created ful, the result of this research allows co-simulation sce- to provide a common interface for different simulation narios that may offer real-time guarantees, bounds in tools to communicate. Obviously, these tools have differ- the maximum error made, or packet size in the commu- ent capabilities and the FMI standard acts as a common nication between simulators. denominator, requiring the basic features that enable co-simulation. Many other features are available in the 1. Introduction state of the art. See Fig. 1 (a) for some of these features. The optional capabilities of simulators make for a Integration – the interconnection of the components combinatorial explosion when these are coupled in a co- that comprise a system – is identified as a major source simulation, making the development of an orchestration of problems [28, 31] in the concurrent development of mechanism a hard challenge, with many different, often complex engineered systems. It is usually caused by conflicting, concerns. In this project, we address this assumptions about other components of the system Simulator Causality Requirements Rollback Support Information Causal A-causal Exposed None Single Multiple Availability Model Solver Time Detailed Deadreckoning Micro-step Constraints Local Remote Model Model Outputs None Scaled RT Nominal Values Discontinuity Indicator State Order of Static Dynamic Outputs State Step-size Accuracy Values Serialization Preferred Next WCET Input Frequency I/O Signal Kind Extrapolation I/O Abstract Feature Causality Outputs State Feature FMI CS 2.0 Derivative Propagation Feedthrough Mandatory Exclusive Or Delay Time Jacobian Optional Or Dependency Kind Outputs State Outputs State Non-Linear Linear (a) (b) Figure 1. (a) Feature model mapping the capabilities of the simulators encountered in the state of the art. (b) A 2-DOF Oscillator as a coupling of two sub-systems. challenge and we propose a way to take advantage of the Experimental Frame The experimental frame de- optional features provided by simulators, instead of just scribes a set of assumptions in which the behavior trace using the common denominator, without resorting to a of the dynamical system can be compared with the one complex orchestration algorithm. This possibly allows of the original system [3, 29, 32, 33, 35]. for the reuse of existing orchestration algorithms while Validity In order to be used successfully as models at the same time providing a better control over the of the systems, dynamical systems have to be valid accuracy/performance tradeoff. within the experimental frame in which the they are defined. The validity of the dynamical system is then the 2. Background difference between the behavior trace of the dynamical In order to properly define co-simulation scenarios, system and the behavior trace of the original system, some background terms have to be introduced. measured under the assumptions of the experimental Dynamical System A dynamical system is a model frame. For example, the Hooke’s law in the mass-spring- characterized by a state and a notion of evolution rules. damper system can only be used to predict the reaction An example is a mass-spring-damper system, depicted force of the spring for small deformations. in Fig. 1 (b): Solver A solver is an algorithm capable of obtaining the approximate behavior trace of a dynamical system. m1 ·x¨1 = −c1 · x1 − d1 · x˙1 + Fe It is typically an iterative procedure that advances (1) x1 (0) = p1 ; ẋ1 (0) = s1 the simulated time and approximates the values of the variables at that point in time. For a dynamical system where c1 is the spring stiffness constant and d1 the in the form of Eq. (2), the Forward Euler solver is given damping coefficient; m1 is the mass; p1 and s1 the initial as: position and velocity; and Fe denotes the external force x̃(t + h) := x̃(t) + f (x̃(t), u(t)) · h over time acting on the mass. We consider dynamical (3) x̃(0) := x(0) systems that can be written in the state space form: where x̃ is the approximated state vector, u(t) the ẋ = f (x, u) ; x(0) = x0 (2) input, and h > 0 is the micro-step size. Accuracy Accuracy is the difference between an ap- where x(t) is the state vector, u(t) the input vector, proximate behavior trace and an exact one. In most and x0 is the initial state. practical cases, the correct behavior trace is difficult to Behavior Trace The trajectory followed by the state obtain. However, it is possible to get a worst case es- over time is called the behavior trace of the dynamical timate in the order of accuracy of a solver, provided system. Behavior traces can be exact (also called an- that the system obeys certain, physically meaningful, alytical) or approximations. In the above example, the assumptions (e.g., state continuity, Lipschitz conditions function x1 (t) that satisfies Eq. (1) is the behavior trace. [2], and conservation laws [24]). Simulator The composition of a solver with a specific have: model is called a simulator. For example, to get a sim- CS = hR, {S1 , S2 } , Li   ulator of the mass-spring-damper system, write Eq. (1) xc − v1 (7) in state space form, combine with Eq. (3) to get: L =  ẋc − x1  Fe − Fc x̃1 (t + h1 ) := x̃1 (t) + v1 (t) · h1 where: (−c1 x̃1 (t) − d1 ṽ1 (t) + Fe (t)) S1 is the simulator defined in Eq. (4) and the definition ṽ1 (t + h1 ) := ṽ1 (t) + · h1 m1 of S2 is omitted; xc , ẋc are the inputs of S2 , and Fe is (4) the input of S1 ; x1 , v1 are outputs of S1 and Fc is the where h1 is the micro-step size, x̃1 (0) := p1 , and output of S2 ; ṽ1 (0) := s1 In general a simulator can be represented as: Trivial Co-simulation Scenario A co-simulation scenario is trivial when the coupling conditions can be Si = hXi , Ui , Yi , δi , λi , xi (0), φUi i transformed into a set of assignments from outputs to δi : R × Xi × Ui → Xi inputs. To achieve this: (1) no input/output is a func- tion of itself; (2) for each input, there is an output that λi : R × Xi × Ui → Yi or R × Xi → Yi (5) provides its value. The co-simulation scenario described xi (0) ∈ Xi by Eq. (7) is trivial. φUi : R × Ui × . . . × Ui → Ui Orchestration Algorithm Given a co-simulation scenario, an orchestration algorithm coordinates the where: simulators ensuring that each progresses in time and Xi is the state set, typically Rn ; Ui is the input set, receives inputs. A trivial scenario can be simulated with typically Rm ; Yi is the output set, typically Rp ; xi (0) Algorithm 1. is the initial state; δi (t, xi (t), ui (t)) = xi (t + H) is the function that instructs the simulator to compute a be- ALGORITHM 1: Generic orchestration mechanism for havior trace from t to t + H, making use of the input trivial co-simulation scenarios. extrapolation function φUi ; H ∈ R is the communica- Data: Stop time Tf , a co-simulation scenario hS, Li, and a tion step size; and λi (t, xi (t), ui (t)) = yi (t) is the output communication step size H. Result: A co-simulation trace. function. The input extrapolation function φUi plays an t := 0; important role in guaranteeing that the simulator does while t < Tf do not read values from the environment while computing Solve the following system for the unknowns the behavior trace in the interval t → t + H. Often- y1 (t), . . . , yn (t), u1 (t), . . . , un (t): yi (t) = λi (t, xi (t), ui (t)), for i = 1, . . . , n times, constant extrapolation from the last known input ; L(y1 (t), . . . , yn , u1 (t), . . . , un (t)) = 0̄ is used, that is, φUi (τ, ui (t)) = ui (t), for τ ∈ [t, t + H]. The values [y1 (t), . . . , yn (t), u1 (t), . . . , un (t)]T denote a point at time t of the co-simulation trace; Co-simulation Scenario Simulators can have in- Instruct each simulator to advance to the next puts and outputs, which capture the environment in communication step: which the original system operates. They can be com- xi (t + H) := δi (t, H, xi (t), ui (t)), for i = 1, . . . , n; bined by specifying how their inputs/outputs are con- Advance time: t := t + H; nected. A co-simulation scenario is a specific arrange- end ment of simulators and their I/O coupling conditions. An autonomous scenario requires at least the following information: Hierarchical Co-simulation A co-simulator is ob- tained when an orchestration mechanism is coupled CS = hS, Li with a co-simulation scenario. According to our nomen- S = {S1 , . . . Sn } (6) clature, a co-simulator is a simulator and can be speci- m fied as in Eq. (5). This means that a co-simulation sce- L : Y1 × . . . × Yn × U1 × . . . × Un → R nario can be comprised of simulators, which can them- S is the set of causal simulators, each defined as in selves be co-simulation scenarios with suitable orches- Eq. (5); and L induces the following coupling constraint: trators. This is important because hierarchical systems are best described by hierarchical co-simulation scenar- L(y1 , . . . , yn , u1 , . . . , un ) = 0̄ ios. In the following sections, we describe non-trivial co- As an example, for the co-simulation scenario cor- simulation scenarios and, instead of describing how responding to the multi-body system of Fig. 1 (b), we these can be solved with more complex orchestration algorithms, we describe how they can be translated as is shown in [2, 18], and empirically in [4], breaking into trivial co-simulation scenarios. This approach, sup- an algebraic loop instead of solving it can lead to a ported by modeling the co-simulation scenarios, allows high error in the co-simulation. A better way is to use for clear separation of concerns and provides flexibility a fixed point iteration technique, where in the general in choosing how to deal with each of them. case, simulators are asked to compute the interval t → t + H many times, with improved inputs, until some 3. Concerns in Co-simulation convergence criteria is met. 3.1 Accuracy Concern Given a co-simulation scenario with algebraic loops, extra information is necessary to be able to identify the The accuracy of a co-simulation trace is the degree loops, as pointed out in [2, 8, 30]. Assuming that this to which it conforms to the real trace. Error – the information exists, the simulators that are involved in difference between the co-simulation trace and the real an algebraic loop can be “lifted out” and replaced by trace – is then a measure of accuracy. Obtaining the a single simulator SSC whose δSC implements the it- real trace, for most dynamical systems, is currently eration techniques that solves the loop. The result is impossible. However, there is an important result in a trivial scenario that can be simulated by the orches- simulation – convergence – which allows the order of tration mechanism of Algorithm 1. This adaptation is the worst case deviation from the real trace made by summarized in Fig. 2 (b) a numerical method to be controlled by adjusting the micro-step size hi of the solver. The same result can be Inaccurate Scenario: Scenario with Loop: applied to certain co-simulation scenarios [2], allowing ... ... the communication step size H to control the global error. Accurate Scenario: Scenario without Loop: Adjusting the communication step size H is then an accuracy concern. Given a co-simulation scenario, H ... Fixed point iterate: ... can be controlled by an extra simulator SH , introduced artificially, whose outputs are the time variable t and H, and inputs are relevant outputs of other simulators. (a) (b) A new independent variable s with a communication step size of 1 is introduced. Variables t and H become Figure 2. (a) Transformation that deals with accuracy functions of s. A similar translation has been proposed concern. (b) Dealing with algebraic loop concern. in [22] and the simulator SH can implement a well known PI-Controller. See [2, 9, 13, 24, 34] for error control alternatives in co-simulation. Note that these 3.3 Communication Concern can also be applied in our approach because SH can If simulators in a co-simulation scenario execute in dif- be a co-simulation scenario (e.g., a copy of the original ferent computers, a small H incurs a too high commu- scenario running at a communication step size of H 2 , for nication cost. On the other hand, using a large H places Richardson extrapolation). Fig. 2 (a) summarizes this the burden in the functions φ to accurately extrapolate approach. the inputs of each simulator across a large interval. In many cases – in particular, for the FMI Standard –, 3.2 Algebraic Loops Concern each simulator Si is the one responsible for implement- Algebraic loops occur whenever there is a variable that ing φUi . Therefore a problem exists where H should be is a function of itself. The state and output of each high to reduce the communication cost, but the accu- simulator Si in a co-simulation can be written as: racy of functions φ cannot be improved to compensate. xi (t + H) = δi (t, H, xi (t), ui (t)) To show how this problem can be addressed, con- (8) sider the scenario shown in Fig. 3 (a), where the simu- yi (t + H) = λi (t, xi (t + H), ui (t + H)) lators communicate every H units of time. The purpose Taking into account the coupling conditions, it is easy is to increase H and mitigate the accuracy loss. For to see that an output of a simulator may depend on that, replace each group of simulators in the same com- itself. These kinds of algebraic loops in the output puter by a single simulator. Then, the new simulator equations can be avoided by replacing ui (t + H) in encapsulates a co-simulation scenario where the inter- Eq. (8) by the corresponding extrapolation φui (H, ui (n· nal communication step size is Hsmall < H. An artificial H), ui ((n−1)·H), . . .) which does not depend on ui ((n+ simulator is introduced to provide approximated values 1) · H), thus breaking the algebraic loop 1 . However, of the outputs of the simulators in the other computers. These can be extrapolations from the values collected 1 See [18] for the other kind of algebraic loops. at the other computers. In the example, S20 collects the outputs of S1 and sends them over the wire to S10 at ev- that ensure equal outputs. For the details of how those ery H time units. The smaller Hsmall , the finer grained calculations can be done, see [1, 15, 27]. the extrapolation of S1 will be. This approach can be applied whenever an input ex- 4. Related Work trapolation function needs to be provided, regardless The aim of this work is to generalize the work done by of the computer in which the simulators execute. For Van Acker et al. [30], where a language is proposed to instance, when the scenario is comprised of simulators configure the co-simulation scenario with extra informa- whose outputs evolve at very different rates, as happens tion identifying the optimal rates for each simulator and in circuit simulation [21], better extrapolation functions algebraic loops. In our work, we recognize that, for each can be provided to save computation on the “slow” com- concern, there are multiple solutions, with differing or- ponents. Furthermore, if simulators provide rollback ca- ders of “cost”, that depend on the sensitivities between pabilities, an iterative predictor correct method can be simulators. The work in Kajtazovic et al. [17] is similar made, yielding a generalized waveform relaxation itera- to ours in the sense that a generative approach is fol- tion [20]. lowed, but there is no focus into identifying and solving 3.4 Modularity Concern the multiple concerns involved in devising an orches- tration algorithm. The work in Benedikt and Holzinger It is possible that, even without algebraic loops, the [5] presents some initial steps toward an orchestration coupling conditions do not yield a set of assignments. To mechanism that adapts at run-time. Similarly to our show how this can happen, consider the co-simulation work, it recognizes that the orchestration mechanism scenario that represents the coupled system on top of is highly dependent on the co-simulation scenario and Fig. 3 (b). The input to the first simulator is the external that it should be tuned automatically. However, we dif- T force Fe and the outputs are [x̃1 , ṽ1 ] (see Eq. (4)). fer in the approach: we do it statically, as opposed to at The input to the second simulator is the external force run-time, like they do. T Fc and the outputs are [x̃3 , ṽ3 ] . Clearly, there is a T mismatch: the outputs [x̃1 , ṽ1 ] of the first simulator 5. Conclusion cannot be coupled directly to the input Fc of the second This project aims at dealing separately with the many simulator, and vice versa. However, the massless link concerns that originate in continuous time co-simulation. restricts the outputs of the two sub-systems to be the Our approach is to stick to a simple orchestration al- same and Fe = Fc , whatever that force may be. gorithm, and transform the scenarios, by introducing Coupled System: artificial simulators. Fig. 4 summarizes our overall ap- proach. Starting Scenario Causality Conflicts Out-to-In Assignments Non-Trivial Scenario: Scenario Algebraic Loops Slow Distributed Scenario: Trivial Scenario Computer A Computer B Communication Trivial Scenario: Optimization Optimized Step: Trivial Scenario Fast Distributed Scenario: Figure 4. Overview of the main transformation stages. Step: Step: Step: The advantage is that there is a clear set of pre- (a) (b) conditions and post-conditions for each transformation, showing the separation of concerns. Based on anecdo- Figure 3. (a) Dealing with distribution concern. (b) tal evidence, we propose the order of the stages to be Transformation that solves causality conflicts. the one in the figure. This order ensures that no concern resurfaces in later stages of the transformation. Because Our solution to this concern is similar to that of of the performance/accuracy tradeoff, the communica- [14] and is summarized in Fig. 3 (b). The essence is to tion concern can only be addressed as an optimization add an artificial simulator to the co-simulation, which problem and we allow for the application of optimiza- calculates the appropriate inputs to the simulators, tion techniques to find an optimal co-simulation sce- nario. One disadvantage of our approach is that co- national MODELICA Conference, pages 173–184, Mu- simulation scenario can quickly become unreadable, due nich, Germany, nov 2012. Linköping University Elec- to the injected artificial simulators. Further evaluation tronic Press; Linköpings universitet. doi: 10.3384/ is necessary to measure the how complex non-trivial ecp12076173. co-simulation scenarios can become after being trans- [7] J. F. Broenink and Y. Ni. Model-driven robot-software formed. design using integrated models and co-simulation. In The co-simulation scenarios used in the current work Embedded Computer Systems (SAMOS), 2012 Interna- were created artificially. In the future, we aim at testing tional Conference on, pages 339–344, 2012. ISBN VO -. these approaches with real co-simulation scenarios, such doi: 10.1109/SAMOS.2012.6404197. as the ones developed in [12] and in the INTO-CPS [8] D. Broman, C. Brooks, L. Greenberg, E. A. Lee, project 2 . Furthermore, the order of convergence has to M. Masin, S. Tripakis, and M. Wetter. Determinate be studied for the approach described in Section 3.3. composition of FMUs for co-simulation. In Eleventh ACM International Conference on Embedded Software, Montreal, Quebec, Canada, 2013. IEEE Press Piscat- Acknowledgment away, NJ, USA. ISBN 978-1-4799-1443-2. This work was partially funded with PhD fellowship [9] M. Busch. Continuous approximation techniques for from the Agency for Innovation by Science and Tech- co-simulation methods: Analysis of numerical stability nology in Flanders (IWT). The author is grateful for all and local error. ZAMM - Journal of Applied Mathemat- the pertinent comments provided by the reviewers and ics and Mechanics / Zeitschrift für Angewandte Math- other contestants of the Student Research Competition. ematik und Mechanik, 96(9):1061–1081, sep 2016. ISSN 00442267. doi: 10.1002/zamm.201500196. References [10] C. W. De Silva. Mechatronics: an integrated approach. [1] M. Arnold and M. Günther. Preconditioned Dynamic CRC press, 2004. Iteration for Coupled Differential-Algebraic Systems. [11] U. Eliasson, R. Heldal, J. Lantz, and C. Berger. BIT Numerical Mathematics, 41(1):1–25, 2001. doi: Agile Model-Driven Engineering in Mechatronic Sys- 10.1023/A:1021909032551. tems - An Industrial Case Study. In J. Dingel, [2] M. Arnold, C. Clauß, and T. Schierz. Error Analy- W. Schulte, I. Ramos, S. Abrahão, and E. Insfran, edi- sis and Error Estimates for Co-simulation in FMI for tors, Model-Driven Engineering Languages and Systems Model Exchange and Co-Simulation v2.0. In S. Schöps, SE - 27, volume 8767 of Lecture Notes in Computer A. Bartel, M. Günther, W. E. J. ter Maten, and Science, pages 433–449. Springer International Pub- C. P. Müller, editors, Progress in Differential-Algebraic lishing, 2014. ISBN 978-3-319-11652-5. doi: 10.1007/ Equations, pages 107–125. Springer Berlin Heidelberg, 978-3-319-11653-2 27. Berlin, Heidelberg, 2014. ISBN 978-3-662-44926-4. doi: [12] J. Fitzgerald and K. Pierce. Collaborative Design for 10.1007/978-3-662-44926-4 6. Embedded Systems. Springer Berlin Heidelberg, Berlin, [3] P. I. Barton and C. C. Pantelides. Modeling of combined Heidelberg, 2014. ISBN 978-3-642-54117-9. doi: 10. discrete/continuous processes. AIChE Journal, 40(6): 1007/978-3-642-54118-6. 966–979, jun 1994. ISSN 0001-1541. doi: 10.1002/aic. [13] V. Galtier, G. Plessis, and L. Renardi. FMI-Based Dis- 690400608. tributed Multi-Simulation with DACCOSIM. In Spring [4] J. Bastian, C. Clauß, S. Wolf, and P. Schneider. Master Simulation Multi-Conference, pages 804–811. Society for Co-Simulation Using FMI. In 8th International Mod- for Computer Simulation International, 2015. elica Conference, pages 115–120, Dresden, Germany, jun [14] B. Gu and H. H. Asada. Co-simulation of algebraically 2011. Fraunhofer Institute for Integrated Circuits IIS. coupled dynamic subsystems. In American Control doi: 10.3384/ecp11063115. Conference, 2001. Proceedings of the 2001, volume 3, [5] M. Benedikt and F. R. Holzinger. Automated configu- pages 2273–2278 vol.3, 2001. ISBN 0743-1619 VO - 3. ration for non-iterative co-simulation. In 17th Interna- doi: 10.1109/ACC.2001.946089. tional Conference on Thermal, Mechanical and Multi- [15] B. Gu and H. H. Asada. Co-Simulation of Algebraically Physics Simulation and Experiments in Microelectron- Coupled Dynamic Subsystems Without Disclosure of ics and Microsystems (EuroSimE), pages 1–7, Mont- Proprietary Subsystem Models. Journal of Dynamic pellier, apr 2016. IEEE. ISBN 978-1-5090-2106-2. doi: Systems, Measurement, and Control, 126(1):1, apr 2004. 10.1109/EuroSimE.2016.7463355. ISSN 00220434. doi: 10.1115/1.1648307. [6] T. Blockwitz, M. Otter, J. Akesson, M. Arnold, [16] C. Hardebolle and F. Boulanger. Exploring Multi- C. Clauss, H. Elmqvist, M. Friedrich, A. Junghanns, Paradigm Modeling Techniques. SIMULATION, 85(11- J. Mauss, D. Neumerkel, H. Olsson, and A. Viel. Func- 12):688–708, nov 2009. doi: 10.1177/0037549709105240. tional Mockup Interface 2.0: The Standard for Tool in- dependent Exchange of Simulation Models. In 9th Inter- [17] S. Kajtazovic, C. Steger, A. Schuhai, and M. Pistauer. Automatic Generation of a Coverification Platform. In 2 http://into-cps.au.dk/ A. Vachoux, editor, Applications of Specification and Design Languages for SoCs, pages 187–203. Springer 188, 2007. ISSN 00078506. doi: 10.1016/j.cirp.2007.05. Netherlands, Dordrecht, 2006. ISBN 978-1-4020-4997-2. 044. doi: 10.1007/978-1-4020-4998-9 11. [29] M. K. Traoré and A. Muzy. Capturing the dual rela- [18] R. Kübler and W. Schiehlen. Two Methods of Simulator tionship between simulation models and their context. Coupling. Mathematical and Computer Modelling of Simulation Modelling Practice and Theory, 14(2):126– Dynamical Systems, 6(2):93–113, jun 2000. ISSN 1387- 142, feb 2006. ISSN 1569190X. doi: 10.1016/j.simpat. 3954. doi: 10.1076/1387-3954(200006)6:2;1-M;FT093. 2005.03.002. [19] M. Kudelski, L. M. Gambardella, and G. A. Di Caro. [30] B. Van Acker, J. Denil, P. D. Meulenaere, RoboNetSim: An integrated framework for multi-robot H. Vangheluwe, B. Vanacker, and P. Demeulenaere. and network simulation. Robotics and Autonomous Generation of an Optimised Master Algorithm for FMI Systems, 61(5):483–496, may 2013. ISSN 09218890. doi: Co-simulation. In Proceedings of the Symposium on 10.1016/j.robot.2013.01.003. Theory of Modeling & Simulation-DEVS Integrative, pages 946–953. Society for Computer Simulation [20] E. Lelarasmee, A. E. Ruehli, and A. L. Sangiovanni- International, 2015. Vincentelli. The Waveform Relaxation Method for [31] H. Van der Auweraer, J. Anthonis, S. De Bruyne, and Time-Domain Analysis of Large Scale Integrated Cir- J. Leuridan. Virtual engineering at work: the challenges cuits. In Computer-Aided Design of Integrated Cir- for designing mechatronic products. Engineering with cuits and Systems, IEEE Transactions on, volume 1, Computers, 29(3):389–408, 2013. ISSN 0177-0667. doi: pages 131–145, 1982. ISBN 0278-0070 VO - 1. doi: 10.1007/s00366-012-0286-6. 10.1109/TCAD.1982.1270004. [32] H. Vangheluwe. Foundations of Modelling and Simu- [21] W. J. McCalla. Fundamentals of computer-aided circuit lation of Complex Systems. EASST, 10, 2008. doi: simulation, volume 37. Springer Science & Business 10.14279/tuj.eceasst.10.162.148. Media, 1987. ISBN 1461320119. [33] H. Vangheluwe, J. De Lara, and P. J. Mosterman. An [22] P. J. Mosterman, J. Zander, G. Hamon, and B. Denckla. introduction to multi-paradigm modelling and simula- A computational model of time for stiff hybrid sys- tion. In Proceedings of AIS2002 (AI, Simulation & tems applied to control synthesis. Control Engineer- Planning), pages 9–20. SCS, 2002. ing Practice, 20(1):2–13, 2012. ISSN 09670661. doi: [34] A. Verhoeven, B. Tasic, T. G. J. Beelen, E. J. W. 10.1016/j.conengprac.2011.04.013. ter Maten, and R. M. M. Mattheij. BDF compound- [23] N. Pedersen, T. Bojsen, J. Madsen, and M. Vejlgaard- fast multirate transient analysis with adaptive stepsize Laursen. FMI for Co-Simulation of Embedded Control control. J. Numer. Anal. Ind. Appl. Math, 3(3-4):275– Software. In The First Japanese Modelica Conferences, 297, 2008. May 23-24, Tokyo, Japan, number 124, pages 70–77. [35] B. P. Zeigler, H. Praehofer, and T. G. Kim. Theory of Linköping University Electronic Press, may 2016. ISBN modeling and simulation: integrating discrete event and 1650-3740. doi: 10.3384/ecp1612470. continuous complex dynamic systems. Academic press, [24] S. Sadjina, L. T. Kyllingstad, E. Pedersen, and 2 edition, 2000. S. Skjong. Energy Conservation and Power Bonds in Co-Simulations: Non-Iterative Adaptive Step Size Control and Error Estimation. arXiv preprint arXiv:1602.06434, 2016. [25] R. A. Saleh, S.-J. Jou, and A. R. Newton. Mixed- mode simulation and analog multilevel simulation, vol- ume 279. Springer Science & Business Media, 2013. ISBN 1475758545. [26] S.-A. Schneider, J. Frimberger, and M. Folie. Signifi- cant Reduction of Validation Efforts for Dynamic Light Functions with FMI for Multi-Domain Integration and Test Platforms. In 10th International Modelica Confer- ence, 2014. [27] B. Schweizer, D. Lu, and P. Li. Co-simulation method for solver coupling with algebraic constraints incorpo- rating relaxation techniques. Multibody System Dy- namics, 36(1):1–36, jan 2016. ISSN 1384-5640. doi: 10.1007/s11044-015-9464-9. [28] T. Tomiyama, V. D’Amelio, J. Urbanic, and W. El- Maraghy. Complexity of Multi-Disciplinary Design. CIRP Annals - Manufacturing Technology, 56(1):185–