Aspects of derivative causality in bond-graph models

The bond-graph method used in the analysis of system dynamics problems leads to a system containing a number of differential equations equal to the number of energy storing elements in integral causality and a number of algebraic equations equal to the number of energy storing elements in derivative causality. When the system has many energy storing elements in derivative causality, the system of equations becomes a very difficult differential-algebraic one. The paper proposes the simplest procedure for finding the appropriate form for numerical simulation of the equations system. The efficiency of the proposed solution is exemplified by using the bond-graph model of a mechanical system which contains nine energy storing elements in derivative causality.


Introduction
The modelling of system dynamics by using the bond-graph method presumes a profound analysis of the manner of transmitting and transforming the power from the driving element to the working element [1][2][3].It represents a valuable alternative to other methods for solving dynamics problems [4,5].More difficult problems in obtaining the motion equations and in performing the numerical simulation occur when the bond-graph model contains energy storing elements in derivative causality [6,7].The system of equations which results in such circumstances is a differential-algebraic one, difficult to be resolved.In case of a numerical simulation with block diagrams, derivative blocks occur, which, because of the anticipative character of the derivative, can lead to computing errors and simulation stopping [8].In order to render evident these difficulties and to propose and demonstrate the simplest solution, the paper presents the bond-graph model of a mechanism which contains two elements in integral causality and nine elements in derivative causality.

System description
The mechanical system, whose dynamics will be further modelled, is presented in Fiure.1.It consists of o bar of mass M 1 and length 2l, whose ends can move away along two fixed mutual perpendicular directions.At the end of the bar which moves along the vertical direction, a slide of mass m is pinned, which is attached to a spring of stiffness k.When the bar is in the vertical position, the spring is undeflected.At the end of the bar which moves along the horizontal direction, a disk of mass M 2 and radius R is pinned.The disk rolls without slip along the horizontal direction.The two weights of mass M 4 are supported by ideal wires, which are attached at the centre of the disk of mass M 2 .Each wire passes over a pulley of mass M 3 and radius R. The friction is ignored.The system oscillates if, for example, the bar is inclined with angle with respect to the vertical direction at the initial moment and then it is let free.

The bond-graph model
The bond-graph model is presented in Figure 2. The element I that is assigned with integral causality, corresponds to the bar rotation (bond 10), after the causal assignment of sources.This fact leads to the causal assignment of all bonds.The element C, corresponding to the spring, results in integral causality and all the other nine elements I, corresponding to element translation and rotation inertia (bonds 9, 1, 13, 33, 16, 19, 22, 26 and 29) result in derivative causality.This means that the system of equations which describes the mechanical system dynamics will contain two differential equations and nine algebraic equations.These equations are: 2 sin 2 sin 2 sin cos  4 The use of equations system resulted from the bond-graph model In order to perform a numerical simulation, it is necessary to obtain a system of differential equations in explicit form.They are found by differentiating the previous nine algebraic equations ( the equations ( 3) to (11)) with respect to time.The derivatives have the form: Parameter is determined from the constitutive equation of element I, which is connected by the bond 10 and is in integral causality: (21) The system of 12 differential equations containing the equations (1), (2) (12), (13), ( 14), (15), ( 16), (17), (18), ( 19), (20) and (21), is in implicit form and has the unknowns , , , system of equations should be resolved like an algebraic system with respect to derivatives of the unknowns, in order to obtain a system of differential equations in explicit form, that can be numerically integrated.This operation involves a difficult calculus.A simpler solution is to substitute the derivatives given by equations ( 12), ( 13), ( 14), (15), ( 16), (17), (18), ( 19), (20) in equation (1).In the final shape of the resulted equation, the variable and its derivative are replaced by: where denotes the angular velocity of the bar.In this manner, it is much easier to obtain a system of three differential equations in explicit form:  The following numerical values are considered: m=4kg, M 1 =8kg, g=9.81m/s, k=1000N/m, L=0.5m, M 2 =6kg, M 3 =6kg and M 4 =10kg.Equations ( 24), ( 25) and ( 26) are numerically integrated in MATLAB, by using different initial conditions.
The variation of the angular velocity in terms of vertical bar position, given by angle  is shown in Figure 3.
The curve represented by the continuous line is obtained when the initial conditions are 0 0.5rad and 0 0 .The initial potential energy is not high enough to produce the passing of the bar to the left side of the vertical position, that is to have negative values of the angle .
The curve represented by the dashed line is obtained when the initial conditions are 0 0.6rad and 0 0 .The initial potential energy is now enough big to produce the passing of the bar to the left side of the vertical position, that is to have negative values of the angle .The angular velocity varies in a very interesting manner.It increases at the beginning and then decreases, when the bar approaches the vertical position.The evolution is then symmetric with respect to axis Oy.
Here there are some comments referring to the two curves represented in Figure 3.The system does not move if the initial conditions are zero.When the initial conditions are slightly different from zero, for example they are 0 0.01 rad and 0 0 , the maximum value for angle is 0.56779 rad.If the initial value of the angle is 0.56779 rad, then the initial potential energy will determine the bar to come very close to the vertical position, but it will not pass to the left side of the vertical position.The initial conditions for the angle chosen for the two diagrams are less and higher than the value 0.56779 rad.The bar remains at the right side of the vertical position in the first case because the initial potential energy is not high enough, so that the bar to pass to the left side of the vertical position.The bar passes to the left side of the vertical position in the second case, because the initial potential energy is high enough to make the bar to pass at the left side of the vertical position.

Conclusions
The bond-graph method used in the analysis of system dynamics problems leads to a system of equations, containing a number of differential equations equal to the number of energy storing elements in integral causality and a number of algebraic equations equal to the number of energy storing elements in derivative causality.After differentiating the algebraic equations, there are two possibilities for finding a system of differential equations in explicit form, suitable for numerical simulation.To get this explicit form represents a very tedious task in case of systems which contain a great number of elements in derivative causality.The easier way, revealed by the present paper, is to substitute the derivatives of the algebraic equations in the differential equations corresponding to the energy storing elements in integral causality.In this manner, the final system has a low number of differential equations, equal to the number of energy storing elements in integral causality.
This procedure is exemplified by using the bond-graph model of a mechanical system which contains two energy storing elements in integral causality and nine energy storing elements in derivative causality.There are also presented the results of the numerical simulation and some comments about the motion of the mechanism in different initial conditions, in order to demonstrate the efficiency of the proposed procedure.

Fig. 3 .
Fig. 3.The angular velocity of the bar with respect to the angle .