Mathematics and Computers in Simulation 46 (1998) 535±542
Simulation of modular dynamic systems A. RuÈkgauer*, W. Schiehlen Institute B of Mechanics, University of Stuttgart, Pfraffenwaldring 9, D-70550 Stuttgart, Germany
Abstract A modular structure for dynamic systems is obtained by functional decomposition. This leads to an interconnected set of dynamic modules. A classification for modular dynamic systems is introduced, including constrained multibody systems with explicit joints, and several state space forms. The resulting mechanical differential-algebraic equations are treated by a poststabilizing projection method. For most efficient simulation, a multi-rate-method is incorporated. The problem of algebraic loops is discussed and a method for its solution is proposed. The simulation program NEWMOS supplies all simulation strategies introduced. Several examples from different areas of dynamics show the power of the approach proposed. # 1998 Published by IMACS/Elsevier Science B.V.
1. Introduction Analyzing large dynamic systems as for instance in mechatronics, one usually has to deal with systems including very different time scales and a tremendous potential of complexity. A modular block simulation approach supported by computer programs as SIMULINK [6] among others is well suited to treat standard problems following from the mechatronic system design. However, when the problem size increases, and in particular, when mechanical systems are considered, general purpose block simulators are no longer adequate. An extended strategy supporting modular structures not only for modelling but also for the numerical treatment is preferable. It is the aim of this paper to introduce a simulation approach allowing to exploit the system structure supplied by modular modelling for the needs of efficient simulation by time integration. Modular structures are inherent to any technical problem. A partitioning of the overall problem by engineering intuition leads to a set of interconnected modules. In Fig. 1 the functional decomposition of a platoon of articulated vehicles, Petersen et al. [8], is shown where two vehicles are linked by a rigid tow bar: Modules for the dynamics of the lateral motion of each of the two vehicles, steering, actuators, and control are derived. The modules are interconnected by both constraints and by signal flow from one module to the other. ÐÐÐÐ * Corresponding author. E-mail:
[email protected] 0378-4754/98/$19.00 # 1998 Published by IMACS/Elsevier Science B.V. PII S 0 3 7 8 - 4 7 5 4 ( 9 8 ) 0 0 0 8 2 - 2
A. RuÈkgauer, W. Schiehlen / Mathematics and Computers in Simulation 46 (1998) 535±542
536
Fig. 1. Functional decomposition of a dynamic system.
2. Types of dynamic systems and encapsulation To describe the system behavior, the transfer function approach widely used in control theory is appropriate: the system is described by an output vector function g as y g
u; p; t
(1)
with the nu-input-vector u, the np-parameter-vector p, the time t and the ny-output-vector y. The modular systems structure may contain q modules. Then, all the inputs and Poutputs of the system P are related to each other by the constant ny nu-incidence-matrix P with nu qi1 nui and ny qi1 nyi , as y Pu
(2)
The dynamics of systems is assumed to be an internal quantity and can be described in various manners. The most general dynamic description used is the explicit, nonsingular and nonlinear state space form which reads as x_ f
x; u; t
(3)
where the nx-state-vector x is used. A variation of Eq. (3) is the linear state space form as commonly used in control engineering: x_ Ax Bu
(4)
with the matrices A and B. Experimental data and explicitly discretized dynamic systems can be described using the discrete state space form of Eq. (3) as x
k1 f
x
k ; u
k ; t
k1 ÿ t
k
(5)
with time index k. A very important class of dynamic systems is described by mechanical structures in form of multibody systems (MBS). In the following, the semi-explicit, possibly nonminimal mechanical
A. RuÈkgauer, W. Schiehlen / Mathematics and Computers in Simulation 46 (1998) 535±542
537
descriptor-form, Hairer and Wanner [5], is used, which reads q_ Kz b : a;
(6a)
M_z k qe ÿ
GKT k
(6b)
g
q; t 0
(6c)
and
with the nq-vector q describing the generalized coordinates and the nz-vector z describing the generalized velocities. Eq. (6a) is referred to as kinematics relation, Eq. (6b) describes the dynamics. Additional constraints by closed chains are introduced by Eq. (6c). The set (6) of mechanical differential-algebraic equations (MDAE) is of index k3, Hairer and Wanner [5]. For the numerical treatment, an index reduction must be performed by differentiating the constraint relation (6c) twice: g_ GKz Gb and
@g 0 @t
(7a)
@a @a @2g @2g _ a a 0 __g Ga GK z_ G @qT @t @qT @t @t2
(7b)
To extend the modular block design to mechanical systems, any mechanical system can contain interaction nodes of the form gnode : [r|S] with the node position vector r(q,t) and the node rotation matrix S(q,t), see Fig. 2, and [10]. Now additional joints can be defined using the node data. The tow bar from the above example of Fig. 1 can be defined using a spherical difference joint (spherical joints with given, fixed distance l) as gg : jr1 ÿ r2 j ÿ l 0
(8)
using the node vectors r1 and r2. Appropriate index reductions can be found as well, RuÈkgauer and Schiehlen [9]. The global equations of motion for such modular mechanical systems with q components
Fig. 2. Mechanical systems with joints.
538
A. RuÈkgauer, W. Schiehlen / Mathematics and Computers in Simulation 46 (1998) 535±542
and v joints show a block-diagonal structure:
(9a)
g g1 ; . . . ; gq jgg1 ; . . . ; ggv T
and
(9b)
(9c)
The structure of the global constraint matrix G is diagonal in the upper part, the lower part is sparsely filled according to the orientation of the mechanical components in the global structure. For the actual derivation of dynamic models, a formal and abstract interface allowing for the encapsulation of the above dynamic descriptions is used. This interface is implemented in the programming language C and is based upon a hierarchical structure, see Fig. 3. Each system is described by the structure system_descr which contains the interface functions for initialization (file opening, etc.), termination, output equation, state equation (eval), state events (roots), the type identifier (type) and a corresponding substructure (type_spec). In the case of MBS as shown in Fig. 3, the substructure mbs_descr includes the nb node points bind of type mbs_contact_descr. The special parameter may be used to define special properties of the system, e.g. a constant mass matrix M. A similar interface definition is introduced for joints and other elements needed for simulation, see below. The data structures shown can easily be generated by codeconverters or by hand. Code-converters are currently available from NEWOPT/AIMS [2] and from DSL [7]. The formal interface definition omitting any information about the actual interconnections between the modules allows for a high degree of code-reusability. All interconnections are introduced at run time.
3. Numerical methods For the numerical treatment of the MDAE, a post-stabilizing projection method is used. An extended right hand side by (6a), (6b) and (7b) is discretized by an ODE integrator. After every integration step executed, a projection onto the manifold described by the conditions (6c) and (7a) is performed. An
A. RuÈkgauer, W. Schiehlen / Mathematics and Computers in Simulation 46 (1998) 535±542
539
Fig. 3. Hierarchical set of data-structures describing dynamic systems.
energy-optimal projection can easily be found by the defect relations Gdq dg
(10a)
GKdzdg_
(10b)
and
Relation (10b) can be solved directly as Eq. (7a) is linear in z. Although (6c) can be strongly nonlinear, only a fixed number of 2 iterations must be performed to solve (10a) due to the fact that the performance of the overall projection resulting from the iteration is directly linked to the index of the system to be solved, Ascher [1]. Employing this fact, a quite efficient integration scheme is obtained [9]. To exploit the modular system structure throughout the integration process, a multirate strategy, Gear [4], is used. Such an approach is essential for the efficient treatment of large mechatronic systems that include systems at very different time scales. The multi-rate-approach is described in detail in [9], too. As the multi-rate integration introduces additional errors, a multi-rate error control is supplied that adjusts the multi-rate-stepwidth used for inter-system-communication. For better multi-rate-performance, interpolation procedures are used for the reconstruction of the input data. As Eq. (2) describes dependencies between different dynamic system modules, an order of evaluation must be determined before the simulation can be performed. This task is solved using graph methods.
540
A. RuÈkgauer, W. Schiehlen / Mathematics and Computers in Simulation 46 (1998) 535±542
A directed graph with nodes representing the dynamic system modules, and edges representing the connections, is designed. Searching for the strictly coherent components, Sedgewick [11], leads to the correct order. However, in case of loops, the search algorithm must be modified in order to mark the loop-closing edges. Then, potential loops and actual loops have to be distinguished. Potential loops are those for which a strictly coherent component is found, while actual loops require in addition that all modules involved show throughputs. Potential loops are determined and solved by one evaluation of the output equations of all systems contained. Then, the actual loops are known, and the corresponding marked edges are disconnected and replaced by additional inputs for which a constraint condition holds. These additional constrains cannot be solved by DAE-methods as introduced above. A simplified Newton-iteration using Broyden's vector update secant formula, Dennis and Schnabel [3], is adequate to treat this problem. So even complex dynamic systems with feedback and extensive throughputs can be solved efficiently. 4. Simulation concepts The simulation program NEWMOS was set up allowing for simulation using the above numerical strategies. Any so called services, referred to in the NEWMOS-context for systems, integration codes, joints, and interpolation schemes for multi-rate integration, are loaded at run time as precompiled object code. Also, the interconnections are introduced at run time. To assign the proper simulation topology, an interpreter is supplied allowing the user to introduce the problem in a mnemonic form. In addition, the interpreter allows for the run-time-definition of systems, so that, in case of rather small problems, no compilation must be performed at all. For the articulated vehicles problem, a code fragment for loading the precompiled vehicle models and the tow bar joint, definition of a steering controller and connection of both vehicles by joint and steering control with second vehicle, and some parameter assignment reads as define ss_system (SteeringPrinciple) with inputs (mu1, mu2) outputs (del, delp) parameters (K1, K2) by { y[0] K1 * mu1 K2 * mu2; y[1] 0}; load system (Vehicle1) of type (spatial_car_model) from file (``car.so''); load system (Vehicle2) of type (spatial_car_model) from file (``car.so''); load joint (TowBar) of type (spher_diff_joint) from file (``spher_diff_joint.so''); link(Vehicle1.rear_towbar_node, Vehicle2.front_towbar_node) by (TowBar); connect(Vehicle2.del, SteeringPrinciple.del); TowBar.1 1.8;
Several dynamic systems are combined to a so called set which is treated by an integrator as assigned by the user. Several sets are allowed concurrently, resulting in a multi-rate-simulation. Using this approach, one can easily switch between different simulation strategies and dynamic system structures, module tests can be performed independently and finally the whole dynamic problem can be assembled quite easily.
A. RuÈkgauer, W. Schiehlen / Mathematics and Computers in Simulation 46 (1998) 535±542
541
Fig. 4. Autonomous (left) and controlled (right) 1-body-pendulum, crossed circles mark initial position, zero initial velocity.
5. Examples A 10-body-pendulum is set up using mass points with an interaction node at the center and spherical difference joints. The simulation results for the autonomous pendulum are shown in Fig. 4 (left). To study the multi-rate-simulation strategy, an actor with very high frequency eigendynamics is added to hieve the first arm from the low equilibrium to the horizontal position, as shown in Fig. 4 (right). In closed form using a stiff integrator, this simulation task for tsim10 s requires 207 s. Applying the multi-rate-strategy, the computational cost can be dropped down to 60 s. The articulated vehicles problem as introduced in Fig. 1 is also treated by a multi-rate-integration. Here, the fast and stiff components as hydraulic servo, electric actuator, and control device are placed in
Fig. 5. Articulated vehicles undergoing ISO-lane-change.
Fig. 6. Belt-gear, implementation containing algebraic loop, and simulation results.
542
A. RuÈkgauer, W. Schiehlen / Mathematics and Computers in Simulation 46 (1998) 535±542
one set, the vehicles are placed in another, as described in detail in [9]. The results from simulation of an ISO-lane-change at a velocity of v 20m=s are shown Fig. 5. It turns out that, for an optimal integration configuration, the results for tsim8 s can be obtained within 80 s. Without multi-ratestrategy, a simulation of the given manoeuvre takes several days. Fig. 6 (left) shows a gearbox with two pulleys connected by a belt. The given model, Fig. 6 (right) includes one potential and one actual loop. As input to the system, step functions acting as torque at the pulley are used. The simulation program detects the loop and supplies the correct solution after one iteration due to the fact that all throughputs contained are linear. 6. Summary A modular modelling and simulation approach has been proposed for the treatment of large dynamic problems as arising in mechatronics. Mechanical systems with additional joints and constraints that allow for a modular description are supported. The resulting mechanical differential-algebraic equations can efficiently be treated by a post-stabilizing projection method. Most efficiently, the modular dynamic system can be simulated by a multi-rate integration. Due to the modular structure, algebraic loops can occur. A strategy for loop-detection and -treatment is shown. Examples from multibody dynamics, vehicle dynamics, and system dynamics show the validity of the proposed approach. Employing the multi-rate-strategy, big savings in computation time can be obtained. References [1] U.M. Ascher, Stabilization of invariants of discretized differential systems, Report, Department of Computer Science, University of British Columbia, Vancouver, 1996. [2] D. Bestle, P. Eberhard, NEWOPT/AIMS 2.2. Ein Programmpaktet zur Analyse und Optimierung von mechanischen Systemen. Anleitunng AN-35, Institut B fuÈr Mechanik, UniversitaÈt Stuttgart, 1994. [3] J.E. Dennis, R.B. Schnabel, Numerical Methods for Unconstrained Optimization and Nonlinear Equations, Prentice Hall Inc., Englewood Cliffs, 1983. [4] C.W. Gear, Multirate methods for ordinary differential equations, Report UIUCDCS-F-74-880, Department of Computer Science, University of Illinois at Urbana-Champaign, Urbana, 1974. [5] E. Hairer, G. Wanner: Solving Ordinary Differential Equations II. Stiff and Differential-Algebraic Problems, Springer, Berlin, 1991. [6] The Mathworks Inc., Nattick: SIMULINK User's Manual, 1992. [7] Mechatronik Laboratorium Paderborn (MLaP), UniversitaÈt-Gesamthochschule Paderborn: DSL und CAMeL Handbuch, 1994. [8] U. Petersen, A. RuÈkgauer, W. Schiehlen: Lateral Control of a Convoy Vehicle System, in: L. Segel (Ed.), Proceedings of the 14th IAVSD Symposium held at Ann Arbor, Michigan, USA, August, 21±25, 1995, volume 25 of THE DYNAMICS OF VEHICLES on roads and tracks, pp. 519±532, Lisse, 1996. Swets & Zeitlinger. [9] A. RuÈkgauer, W.O. Schiehlen, Efficient Simulation of the Behaviour of Articulated Vehicles, in: H. Wallentowitz (Ed.), Proceedings of the International Symposium on Advanced Vehicle Control, pp. 1057±1070, Aachen, 1996, AVEC'96 Organizing Committee. [10] W. Schiehlen, Technische Dynamik, Teubner Verlag, Stuttgart, Germany, 1986. [11] R. Sedgewick, Algorithmen in C, Addison-Wesley, Bonn, 1992.