Open access peer-reviewed chapter

# Appell-Gibbs Approach in Dynamics of Non-Holonomic Systems

By Jiří Náprstek and Cyril Fischer

Submitted: November 15th 2017Reviewed: March 6th 2018Published: July 18th 2018

DOI: 10.5772/intechopen.76258

## Abstract

Hamiltonian functional and relevant Lagrange’s equations are popular tools in the investigation of dynamic systems. Various generalizations enable to extend the class of problems concerned slightly beyond conventional limits of Hamiltonian system. This strategy is very effective, particularly concerning two-dimensional (2D) and simpler three-dimensional (3D) systems. However, the governing differential systems of most non-holonomic 3D systems suffer from inadequate complexity, when deduced using this way. Any analytical investigation of such a governing system is rather impossible and its physical interpretation can be multivalent. For easier analysis, particularly of systems with non-holonomic constraints, the Appell-Gibbs approach seems to be more effective providing more transparent governing systems. In general, the Appell-Gibbs approach follows from the Gaussian fifth form of the basic principle of dynamics. In this chapter, both Lagrangian and Appell-Gibbs procedures are shortly characterized and later their effectiveness compared on a particular dynamic system of a ball moving inside a spherical cavity under external excitation. Strengths and shortcomings of both procedures are evaluated with respect to applications.

### Keywords

• Appell-Gibbs function
• Lagrangian approach
• non-holonomic systems
• engineering applications

## 1. Introduction

The energy contained in a dynamic system is given by a scalar potential Et. It is a function of time and system response components (displacement, velocity, and acceleration vectors). Moreover, Etis a function of system parameters, position in a field of forces (potential or not), internal sources of energy and of the system evolution including a residual energy. The total energy of the system increases or decreases accordingly with external excitation and dissipation of energy. The form of energy contained within the system can have a deterministic or stochastic character and similarly also excitation and dissipation.

Considering the mechanical energy only, the total energy increase/decrease of the system with respect to time should be in equilibrium with the energy supplies and energy losses due to dissipation. This relation can be outlined by the following equilibrium:

ddtEt=Pt+St,E1

where Ptis power supply (excitation energy per unity time) and Stthe specific dissipation of energy also per unity time (supposed to be independent on accelerations x¨). Functions Pt,Stcan dispose in special cases with a superior potential, which, however, cannot be incorporated into the potential part of total energy. Eq. (1) has a scalar character.

The energy is a primary value characterizing the system state and its evolution in time. The function Etand external influences are a background for the derivation of a governing differential system characterizing the system response with respect to initial and boundary conditions. The governing differential system is then deduced from the equivalence of Eq. (1) type using an adequate variational principle. It claims that the form of the system response corresponds with the minimum of energy spent among all admissible shapes of the system reaction. Take a note that many important settings of external forces and dissipation mechanisms do not admit the formulation by means of potentials. In such cases, they should be incorporated separately into the governing differential system using complementary principles and theorems, for example, virtual works, and so on.

We can find in monographs, for example, [1, 2, 3, 4, 5] and many others, various formulations of potentials Etand functions Pt,Stcombining the system parameters (physical and geometric) and the system response vectors x-displacements, ẋ-velocities, and x¨-accelerations. They can be selected in individual cases with respect to physical or geometric complexity of the system, components of the response, which are to be found, deterministic or stochastic character of the system and its excitation, and so on.

## 2. Basic considerations

Approaches commonly applied to construct mathematical models of dynamic systems with multiple degrees of freedoms (MDOF) follow mostly from principles symbolically outlined by Eq. (1). The equation of this type can be deduced using, for instance, a procedure of virtual displacements. They balance the energy flow in every step and subsequently applied minimization steps try to select such response trajectories, which represent a minimum of energy consumption among all admissible shapes. Let us get briefly through Lagrangian and Appell-Gibbs procedures in order to compare their basic properties. Later, we recognize that most of these properties can be regarded as positive or negative in dependence on a particular problem. Therefore, the solution method should be selected in every particular case very sensitively.

Let us remember that the aim of this chapter is a comparison of Lagrangian and Appell-Gibbs approaches effectiveness to process dynamic systems in holonomic and non-holonomic settings and to help estimate which one is more suitable to be employed in a particular case. Despite that the most important features of non-holonomic systems themselves are briefly treated as well, but for thorough evaluation of their properties, special literature should be addressed. Except five monographs cited in introductory section containing a large number of additional relevant references, a vast number of papers have been published concerning the investigation of various properties of non-holonomic systems.

Motion of an MDOF system with ndegrees of freedom can be described by a system of ndifferential equations and lconstraints:

msx¨s=Xs+r=1lλrArs,s=1,..,n,x=xs,X=Xs,x,XRn,as=1nArsẋs+Br=0,r=1,..,l,λ=λr,B=Br,λ,BRl,bA=Ars,ARl×n.E2

Vector Xrepresents external forces, while λare unknown multipliers. The summation in Eq. (2a) characterizes influence of constraints (holonomic and non-holonomic) related with constraints (Eq. (2b)). These constraints reduce the number of the original degrees of freedom from nto k=nl. The system (Eq. (2)) includes n+ldifferential equations for xand λunknown functions t, which can be determined, provided x,ẋare given in an initial point t0. If the system (Eq. (2b)) is fully integrable, it provides lfunctions fr=frxt,r=1,..,land constraints can be formulated as fr=frxt=cr. They are exclusively of a geometric character and the system is holonomic. Corresponding constraints are formulated in displacements only. In principle, lcomponents of xcan be eliminated and then remains to analyze the system with nlunknowns. Then, it can be considered λ0, and the second part on the right side of Eq. (2a) vanishes. The system with holonomic constraints takes the form:

msx¨s=Xs,fr=frxt=cr,s=1,..,k,r=1,..,l,k=nl.E3

However, frankly speaking, such an operation is possible rather exceptionally. In general, the full form of Eq. (2) should be treated, despite the system is holonomic. If some (or all) of constraints (Eq. (2b)) are not integrable, then the system is non-holonomic. In practice, we encounter these cases when the formulation of constraints includes velocities (more often velocities only).

We should remember that the non-holonomic constraints introduced in Eq. (2b) represent the most simple version of such constraints, as they are linear and in velocity. Many applications, for example, robotics, wind engineering, automotive systems, plasma physics, and so on, present more complicated types of non-holonomic constraints. Notifications to nonlinear non-holonomic constraints in velocity are given in elderly monographs [2, 3]. Later, many papers have appeared presenting results of systematic research at this field originating from particular physical or engineering problems, for example, [6, 7, 8], where higher derivatives of velocities in non-holonomic constraints are discussed. These attributes have been reflecting also in pure mathematical studies with respect to control theory and systems with a delayed feedback, see, for example, series [9, 10, 11] dealing with generalized Lagrange-d’Alembert-Poincaré equations or other studies devoted to non-holonomic reduction and related problems, see, for example, [12] and many others.

Let us realize now that the virtual work of every constraint force should vanish in the meaning as follows:

λrs=1nArsδxs=0,r=1,..,l.E4

Therefore, we have with respect to Eq. (2a):

s=1nmsx¨sXrδxs=0.E5

This equation holds for any arbitrary virtual displacements and represents a generalization of the principle of virtual works in statics and of the d’Alembert principle. The important issue is that it does not include any reactions of constraints. It has been well investigated in the study. For many details, see monographs, for example, [1, 3] and many others.

Let us consider that velocities in constraints Eq. (2b) are increased by virtual increments δẋ, so that they read

s=1nArsẋs+δẋs+Br=0,r=1,2,..,l,E6

Deducting from Eq. (6), the initial state (Eq. (2b)) holds

s=1nArsδẋs=0,r=1,2,..,l.E7

Virtual increments of velocities δẋfit into constraints requested for constraints (Eq. (4)) and, consequently, in Eq. (5) the δxscan be replaced by δẋs:

s=1nmsx¨sXsδẋs=0.E8

We revisit Eq. (2b) and perform differentiation with respect to t:

s=1nArsx¨s+dArsdtẋs+dBrdt=0,r=1,..,l,E9

where d/dtrepresents the operator /t+i=1nẋi/xi. Considering two possible movements of the system in identical initial state and velocities in time t, but with different accelerations x¨and x¨+δx¨, then we obtain with respect to Eq. (9):

s=1nArsδx¨r=0,r=1,..,l.E10

It means that virtual accelerations δx¨satisfy constraints requested similarly like virtual displacements or velocities following Eqs. (4) or (7). Therefore, we can write:

s=1nmrx¨sXsδx¨s=0.E11

Some authors call relations (Eqs. (5), (8), and (11)) as the first, second, and third form of the system equation, see, for example, [3] and others.

Let us multiply each equation in the system (Eq. (2a)) by velocities ẋs. Summing them together, one obtains

s=1nmsx¨sẋs=s=1nXsẋs+r=1ls=1nλrArsẋs,E12

which can be rewritten in the form

dTdt+dVdt=s=1nX˜sẋs,E13

where T,Vare kinetic and potential energy, respectively, and X˜sare forces, which cannot be included into the potential energy V. In other words, relation Eq. (13) indicates that the change of the full energy (kinetic and potential) is equivalent to power (work on velocities) of all forces X˜s, which do not contribute to the potential energy V. Relation Eq. (13) corresponds to equilibrium condition Eq. (1), where functions of excitation and dissipation Pt,Stcorrespond with the influence of non-potential forces X˜s.

## 3. Lagrange’s equations

The original coordinates xshould be replaced with respect to Lagrangian coordinates q. The reason is that they represent the most inherent coordinates respecting the real movement of the system and configuration of external forces. Let us write basic coordinates as functions of Lagrangian ones:

xr=xrq1..qnt,r=1,..,n.E14

It can be easily shown

ẋrq̇s=xrqs,ẋrqs=ddtxrqs.E15

We reconsider Eq. (5) where the virtual displacement δxis replaced by

δxr=s=1nxrqsδqs,r=1,..,ns=1nr=1nmrx¨rXrxrqsδqs=0.E16

The last Eq. (16) can be modified using Eqs. (15), which implies

s=1nr=1nmrddtẋrẋrq̇sẋrẋrqsr=1nXrxrqsδqs=0,E17

This equation can be rewritten now in the form:

s=1nddtLq̇sLqsQsδqs=0,E18

where it has been denoted:

L=r=1nẋrẋrqs,Qs=r=1nXrxrqs.E19

Inspecting the polynomial L, we recognize that it consists of the polynomial of the second and first degrees of components q̇(coefficients are still functions of displacements qand time t) and the absolute part without any velocity components q̇. We can now assign the first part to the kinetic energy T, while the part without velocities to the potential energy V. So that Lcan be understood as the Lagrange function as usually defined

L=TV,E20

provided the dynamic system studied is holonomic and no constraints are applied. In such a case, all variations δqsare independent and Eq. (18) can be fulfilled only if every coefficient in curly brackets vanishes individually. Consequently, we obtain Lagrange’s equations in the form:

ddtTq̇sTqs+Vqs=Qs,s1,..,n,E21

where Qsare generalized external forces as functions of qand t. These forces are basically linear transforms of original forces Xr, see Eq. (19). If holonomic constraints are inserted, then the number of remaining degrees of freedom is lower k=nl. Nevertheless, if there is possibility to define the system after elimination of inactive DOFs, then we can consider formally k=nagain and Eq. (21) remains in force.

Let us suppose now that the system includes lnon-holonomic constraints and those holonomic, which cannot be eliminated. Whatever is the reason for that, it still holds k+l=n. These constraints are described by constraints in Lagrange’s coordinates (analogous with Eq. (2b)) as follows:

s=1nCrsq̇s+Dr=0,r=1,..,l,D=Dr,DRl,C=Crs,CRl×n,E22

This time, the variations δqsare not fully independent and only those components, which satisfy conditions:

s=1nCrsδqs=0r=1,2,..,l,E23

can be regarded as independent.

In such a case, the right side of Eq. (21) should be completed:

ddtTq̇sTqs+Vqs=Qs+r=1lλrCrs,s=1,..,n.E24

To the system, Eq. (24) should be attached lconstraints Eq. (22). So that, finally we have the system of n+lequations with unknowns qand λ. Multipliers λare linearly related with forces in constraints. In particular cases, multipliers λcan be physically interpreted, for instance, they can have a meaning of reactions of a body moving along a given trajectory. Very knowledgeable explanation about manipulation and interpretation of Lagrange’s multipliers from the viewpoint of a general theory as well as of employment in particular cases can be found in the monograph concerning non-holonomic systems, see [2]. For additional information and a large overview of additional literature resources, see [5].

The real dynamic system is always influenced by energy dissipation. Some simple models can be introduced using Rayleigh function R, see, for example, [3]. This way is typically applicable, if linear viscous damping is considered and the Rayleigh function has a quadratic form in velocities q̇s. We can include this factor symbolically into Eq. (24), which reads now:

ddtTq̇sTqs+Vqs=Qs+∂ℛq̇s+r=1lλrCrs,s1,..,n.E25

Hence, the completed system Eqs. (22) and (25) with n+lunknowns can be considered. However, we should be aware that this supplement is rather intuitive and does not follow from any rigorous derivation, although in practice it is widely and successfully used. Nevertheless, comparison of this system with the general relation Eq. (1) introduced in Section 1 is obvious. Take a note that more sophisticated versions of Lagrange’s equations have been developed inspired by physical problems; see, for instance, generalized Lagrange-d’Alembert-Poincaré equations discussed in [11].

Let us add that many details (internal mechanisms and inclusion into governing system) concerning more sophisticated models of damping can be found in monographs of the rational dynamics, for example, [1, 4]. See also papers oriented to practical aspects of the damping either of natural, for example, rheological, aeroelastic origin, or intentionally included in order to achieve the highest damping effectiveness, for instance [13].

## 4. Appell-Gibbs function and equation system

Although the Appell-Gibbs approach is not referred so often in the study as the Lagrangian procedure, there are some monographs treating the analytical dynamics, for example, [1, 3], where detailed features of this method are explained. Moreover, journal papers can be found where special aspects of the Appell-Gibbs approach are discussed. A close relation of the fifth Gaussian form and the Gibbs equations from the viewpoint of Dynamics is studied, for example, [14, 15], important remarks for application are concerned in [16, 17], as well as possibilities of extension for systems with time-dependent masses [18] are indicated.

Let us briefly outline principal steps leading to the Appell-Gibbs differential system with respect to essentials ascertained and introduced in Section 2. We should be aware that generalized external forces Qs, introduced in Eq. (19), follow in principle only kdegrees of freedom, which remained free; thereafter, lconstraints have been applied and the original number nof DOFs has been reduced to k=nl,0<ln. However, due to complicated relations inside the dynamic system, this fact is rather impossible to be employed in basic coordinates xs,s=1,..,nand Lagrange’s coordinates qs,s=1,..,nshould be addressed, as we have also seen in previous Section 3. Nevertheless, it is worthy to involve only such coordinates qs, which correspond to kremained DOFs. It can be easily expressed in Lagrange’s coordinates, unlike basic coordinates xs. So that, as the first step, we reformulate some expressions of Section 2 concerning the transform from basic to Lagrange’s coordinates.

Velocities ẋr,r=1,..,nshould be evaluated with respect to the fact that coordinates xrare functions of all Lagrange’s coordinates qs,s=1,..,kand time t, see Eq. (14):

ẋr=s=1kαrsq̇s+αr,r=1,..,n,where:αrs=xrqs,αr=xrt,E26

which also implies

δẋr=s=1kαrsδq̇s,r=1,..,n,E27

Differentiation of Eq. (26) with respect to time gives

x¨r=s=1kαrsq¨s+s=1kdαrsdtq̇s+dαrdt,ddt=t+m=1kq̇mdt,r=1,..,n.E28

The incremented acceleration vector, when keeping velocities and displacements, can be formulated as follows:

x¨r+δx¨r=s=1kαrsq¨s+δq¨s+s=1kdαrsdtq̇s+dαrdt,r=1,..,n.E29

Deducting Eq. (28) from Eq. (29), one obtains

δx¨r=s=1kαrsδq¨s,r=1,..,n.E30

Hence, it can be written

r=1nXrδx¨r=s=1kr=1nXrαrsδq¨s=s=1kQsδq¨s,E31

where Qsare identical generalized forces, as they have been defined in Eq. (19). With reference to Eq. (11), we can reformulate this equation as follows:

r=1nmrx¨rδx¨rs=1kQsδq¨s=0.E32

This relation will be used later, see Eq. (36).

As a principal step of this section, we define now the Gibbs function Gconcentrating “acceleration energy” included in all nDOFs as follows:

G=12r=1nmrx¨r2,E33

When we pass from basic to Lagrange’s coordinates, only kactive coordinates remain in force and so the expression Eq. (33) can be rewritten:

G=12r=1kmrq¨r2,E34

Expressions Eqs. (33) and (34) differ only in terms independent from accelerations.

Let us introduce the function H:

H=Gs=1kQsq¨s,E35

and evaluate its virtual increment:

δH=δGs=1kQsq¨s=12r=1nmrx¨r+δx¨r212r=1nmrx¨r2s=1kQsδq¨s=12r=1nmrδx¨r2+r=1nmrx¨rδx¨rs=1kQsδq¨s.E36

The last parenthesis vanishes due to the relation Eq. (32). Therefore, if δx¨0, then the function δHis always positive:

δH=δGs=1kQsq¨s>0,E37

which implies that accelerations q¨s,s=1,..,kshould lead to a minimum of the function H, which means:

Gq¨r=Qr,r=1,..,k.E38

The energy dissipation terms Rx,Ry,Rzshould be added to the right side of Eq. (38). At this moment, the conformity of Eq. (38) with the equivalence Eq. (1) is well pronounced, similar like in the previous section. The system Eq. (38) should be completed by geometric constraints:

q̇r=s=1kβrsq̇s+βrr=k+1,..,n.E39

Equations (38) and (39) are the Gibbs-Appell differential system including nequations, which can be written in the normal form and hence it is suitable to be immediately investigated using common methods.

The differential system (Eqs. (38) and (39)) represents the simplest and in the same time the most general form of equations of the dynamic system movement. The form of this system is very simple, and it can be used with the same effectiveness to the investigation of holonomic as well as non-holonomic systems, as the constraints can represent non-holonomic but also holonomic type of constraints. Unlike the Lagrangian approach, the non-holonomic or non-eliminable constraints do not augment the number of differential equations.

Procedure of the Appell-Gibbs equations employment in particular cases is obvious, looking back at this section. In the first step, the “so called kinetic energy of accelerations” 12r=1Nmrx¨r2is composed using nacceleration components of the vector x¨. It represents the Appell-Gibbs function G. In a general case, this function includes also all coordinates xand velocities ẋ. Nevertheless, it is important that Gin Lagrange’s coordinates contains only kselected components of accelerations q¨. Anyway, all ncomponents of q̇and qare still included as a result of a transformation from basic to Lagrange’s coordinates.

It is worthy to remind that the differentiation outlined in Eq. (38) is very easy in a particular case. Indeed, let us realize that Gcan be symbolically expressed as a sum of quadratic function of accelerations q¨s,s=1,..,kG2, linear function of these components G1and function without accelerations G0. Differentiating G2, one obtains the relevant acceleration component in a linear form, which will be moved onto the left side together with a coefficient, which can be a function of all velocities and displacements q̇s,qs,s=1,..,n. Differentiation of G1leads to acceleration-free coefficients and G0can be omitted leading to zeroes. Sometimes, the so-called reduced Appell-Gibbs function Gis defined where G0is a priori omitted.

In the second step, the work of kgiven forces Qon kvirtual displacements qis carried out. It has the form s=1kQsδqs. We substitute now back into Eqs. (38) and add l=nkgeometric constraints following Eqs. (39). So we obtain k+l=ndifferential equations for ncomponents of the vector q¨t. Take a note that no unknown multipliers λemerge here, which on the other hand increases the number of unknowns in a Lagrangian approach.

The procedure working with accelerations instead with velocities provides much simpler governing differential system. Unlike velocities, the acceleration components in the Appell-Gibbs function are included only in a few parts of energy expression. Therefore, all parts including only velocity and displacement components disappear during the differentiation of the Appell-Gibbs function with respect to q¨r,r=1,..,k, and therefore they can be considered beforehand as unimportant.

Investigating problems with rotations, we work with Lagrange’s coordinates ω, which represent in fact velocities. So that by solving the abovementioned differential system, the displacements and velocities ωemerge as results. Rotations themselves remain unattended. May be, it is a forfeit for a relative simplicity of the governing system in comparison with the Lagrangian approach. However, this shortcoming is mostly apparent only. The main part of the result represents usually displacement components, which are obtained without restrictions. Together with velocities ω, they represent a full set of information needed to get through the shape of trajectories of the system response including rotation (illustrative example will be presented later in Section 6). If detailed rotations (not only velocities) are still needed, a subsequent integration can be performed independently using differential relations between rotation velocity vector ωand (for instance) Euler angles, see monographs [1, 3, 4] and others. They provide a detailed description of time history of a body orientation as a function of time t. This step can be useful, for instance, when a detailed animation is needed for presentation purposes.

## 5. Planar movement of a ball in a spherical cavity

### 5.1. Engineering motivation

Passive vibration absorbers of various types are very widely used in civil engineering. TV towers, masts, and other slender structures exposed to wind excitation are usually equipped by such devices. Conventional passive absorbers are of the pendulum type. Although they are very effective and reliable, they have several disadvantages limiting their application.

These shortcomings can be avoided using the absorber of the ball type. The basic principle comes out of a rolling movement of a metallic ball of a radius rinside of a rubber-coated cavity of a radius R>r. This system is closed in an airtight case, see, for instance, Figure 1. First papers dealing with the theory and practical aspects of ball absorbers have been published during the last decade, see [13, 19].

### 5.2. Planar layout of the system, Lagrangian procedure

The version, when the ball is forced to move solely in a vertical plane, has been thoroughly studied using Lagrangian approach in [20, 21] and other detailed papers dealing not only with theoretical aspects but also with experimental verification in the laboratory and in situ examining absorbers installed on real structures.

The cavity is fixed to a vibrating structure. Their dynamic character is represented by a linear single degree of freedom (SDOF) system represented by a mass M. Inside of the cavity, the ball min a vertical plane is moving, that is two degrees of freedom (TDOF) system should be investigated, as it is outlined in Figure 2. It follows from geometric relations:

Rφ=rψ+φ=ϱrφ,E40

where ρr=Rr. It holds for vertical or horizontal components of a displacement and velocity of the internal ball center:

horiz.:u+ϱrsinφu̇+ϱrφ̇cosφ,vert.:ϱrcosφϱrφ̇sinφ.E41

Kinetic energy of a moving system of the ball mand the cavity Mcan be written in a form:

T=12mu̇+ϱrφ̇cosφ2+ϱr2φ̇2sin2φ+12Jψ̇2+12Mu̇2=12m+Mu̇2+mϱru̇φ̇cosφ+m2κϱr2φ̇2,E42

where m/κ=m+J/r2κ=5/7, while the potential energy is given by an expression:

V=mgϱr1cosφ+12Cu2.E43

The damping should be introduced in a form of a simple Rayleigh function:

R=12Mbuu̇2+mbφϱr2φ̇2.E44

m,M– mass of the ball m, mass of the cavity, Mrepresenting the protected structure;

J– inertia moment of the ball m;

bu,bφ– damping coefficients (logarithmic decrements, linear viscous damping);

Expressions Eqs. (42), (43), and (44) should be put into Lagrange’s equations of the second type, see Eqs. (24) or (25) and monographs, for example, [1, 3, 4] and others:

r=1nddtTq̇rTqr+Vqr+∂ℛq̇rδqr=Prt,q1=u=ζϱr,q2=φ,Put=ptMϱr,Pφt=0,E45

which give the governing equations of the system:

φ¨+κbφφ̇+κωm2sinφ+κζ¨cosφ=0,aμφ¨cosφμφ̇2sinφ+1+μζ¨+buζ̇+ωM2ζ=pt,bμ=m/M,ωM2=C/M,ωm2=g/ϱr.cE46

Equation (46) describes 2D movement of a ball absorber under excitation by the force Ptat any arbitrary deviation amplitudes including incidental transition through a limit cycle toward an open regime.

### 5.3. Illustration of some planar system features

Analysis of the governing system (Eqs. (46)) has been done in a couple of papers, for example, [20, 21]. Investigation has been carried out using the harmonic or multi-harmonic balance method, see, for example, [22, 23], respectively.

The system is auto-parametric, see, for example, [24] and other resources. Very rich overview of a theoretical basis of auto-parametric systems can be found in [25]. Expecting a single mode response, the Harmonic balance-based methods are applicable. Following approximate expressions for excitation and response can be written (cf., e.g., [22]):

pt=p0sinωt,φt=αsinωt+βcosωt,ζt=γsinωt+δcosωt.}E47

Having four new variables α=αt,β=βt,γ=γt,δ=δtinstead of two original unknowns φt,ζt, two additional conditions can be freely chosen:

α̇sinωt+β̇cosωt=0,γ̇sinωt+δ̇cosωt=0.E48

After substituting Eqs. (47) and (48) into Eqs. (46) and substituting the sinϖand cosφfunctions by two terms of Taylor expansion, the harmonic balance procedure gives the differential system for unknown amplitudes Z=αβγδT, see, for example, [21, 23]:

MZŻ=FZ.E49

System (49) for amplitudes Ztis meaningful if they are functions of a “slow time,” in other words, if their changes within one period 2π/ωare small or vanishing and individual steps of the harmonic balance operation are acceptable. The matrix Mand the right-hand side vector Fhave the following form:

M=0ω14αβκω18κωAαω018κωAβ14αβκω18μωAβ14αβμωμ+1ω014αβμω18μωAα0μ+1ω,E50

F=1486A0κ3γω2αωm2+12ω2καβδ+8β2γ4α48βκωbφ6A0κδω2βωm2+12ω2αγκ+βδκ4β+48ακωbφω2A0A0+22βμ163δμ+14βμ+48γωbu+δωM2ω2A0A0+22αμ163γμ+14αμ48δωbuγωM2+p0,E51

where A0=α2+β28,Aα=3α2+β28,Aβ=α2+3β28.

Let us consider stationary response of the system. In this case, the derivatives dZ/dtvanish and the right-hand side has to vanish too. Eq. (49) degenerates to the form of

FZ=0E52

Thus, to identify the stationary solutions, the zero solution points of F, depending on the excitation frequency and amplitude, should be traced. In the same time, the signum and the zero points of the Jacobian detJFhave to be checked. The negative value of the Jacobian for a particular point indicates that the corresponding solution is stable, whereas when Jacobian vanishes a bifurcation could occur.

The curve Fαβγδω=0, projected into the planes ωRor ωS(for S2=γ2+δ2), forms the resonance curves known from the analysis of linear oscillators. However, correspondence of this curve to the original Eq. (46) is limited to the case of stationary response. It is necessary to remind that limits of stationarity of the response cannot be determined from properties of Eq. (52) itself. The complete Eq. (49) has to be taken into account for this purpose.

With respect to actual experiences regarding passive vibration absorbers and some interesting properties of system (46), the following reference input data have been introduced:

M=10.0;m=2.0;ϱr=0.71;bφ=0.1;bu=0.2;C=140;po=0.5÷2.5.E53

Utilizing Eqs. (52) and (51), the nonlinear resonance curves describing the stationary response of system (46) can be obtained. A set of such curves for excitation amplitudes p0=0.25,0.5,1,1.5,2.5is shown in Figure 3. It is obvious for the first view the nonlinear character manifesting oneself by a dependence of a position of extreme points on an amplitude of excitation force. This effect is visible predominantly in a neighborhood of a conventional “linear” natural frequency of the absorber, although also the second natural frequency corresponding to the original natural frequency of the structure is affected. The resonance curves are typical for a system with “softening” nonlinearities.

## 6. Spatial version of the system, Appell-Gibbs procedure

### 6.1. Gibbs function

The spatial version of the ball absorber on the basis of rational dynamics has been widely investigated by authors of this chapter, see, for example, [26, 27]. Lagrangian approach and Appell-Gibbs procedures have been discussed in these papers combining analytical and numerical methods. Some important issues will be roughly outlined and for details see cited papers.

Unlike the planar version discussed in the previous section, the Appell-Gibbs approach is used to formulate the governing nonlinear differential system. The authors tried to formulate the spatial version using the Lagrangian procedure as well, see [28]. Although the governing system of the respective holonomic system has been successfully assembled, the further analysis appeared very cumbersome, and therefore, it has been given up to follow this way. Thus, the Appell-Gibbs approach is used to formulate the governing system. Its structure is much more transparent and represents a wider option of analytical-numerical investigation of detailed properties of the ball trajectories within the cavity.

With respect to Sections 2 and 4, the first step represents to construct the Appell-Gibbs function (often referred to as an energy acceleration function) defined as follows:

G=12Mu¨Gx2+u¨Gy2+u¨Gz2+12Jω̇x2+ω̇y2+ω̇z2,E54

where Mis the mass of the ball, Jis central inertia moment of the ball with respect to point G, ωthe angular velocity vector of the ball with respect to its center G, uGthe displacement of the ball center with respect to absolute origin O, Ccontact point of the ball and cavity, Amoving origin related with the cavity in its bottom point, see Figure 4. Coordinates x=xyzare Cartesian coordinates with origin in the point O. Hence, it holds:

uG=uA+uC+un,un=rnu̇G=u̇A+u̇C+ru̇n,ρu̇C,ρ=1r/R,cf.Eq.40:ϱr=Rr,E55

where uAis the displacement of the moving origin Awith respect to absolute origin O, uCthe displacement of the contact point Cwith respect to moving origin A, unthe displacement of the ball center Gwith respect to contact point C, nthe cavity normal unit vector in point C.

Geometry of the cavity (radius R) with respect to moving origin Ais given by equation:

xA2+yA2+zAR2=R2,E56

where xA=xAyAzAare Cartesian coordinates with origin in the moving origin A.

Using Pfaff theorem and adopting a conjecture of non-sliding contact between the ball and the cavity, the respective non-holonomic constraints of “perfect” rolling can be deduced after a longer manipulation:

u̇Gx=u̇Ax+ρωyuCzRωzuCy,u̇Gy=u̇Ay+ρωzuCxωxuCzR,u̇Gz=+ρωxuCyωyuCx,E57

where ρ=1r/R.

In order to substitute for accelerations uointo the Appell function (Eq. (54)), let us differentiate constraints Eqs. (57).

Several manipulations provide expressions for components of the ball center acceleration u¨G, which consist of acceleration in the moving origin A:u¨Arepresenting the given external kinematic excitation and acceleration related to the point Abeing given by an expression: ρu¨C:

u¨Gx=u¨Ax+ρω̇yuCzRω̇zuCy+ρωyu̇Czωzu̇Cy,u¨Gy=u¨Ay+ρω̇zuCxω̇xuCzR+ρωyu̇Czωzu̇Cy,u¨Gz=u¨Az+ρω̇xuCyω̇yuCx+ρωxu̇Cyωyu̇Cx,E58

Because the kinematic excitation is supposed to be horizontal, u¨Az=0into Eqs. (58) should be substituted.

Expressions Eqs. (58) are to be substituted into Eq. (54). Thereby, we obtain the Appell-Gibbs function Gfor the system investigated. The function Gcan be significantly simplified keeping only terms including second-time derivatives u¨Gand ω̇, which represent second-time derivatives of respective rotations. This step provides the reduced Appell-Gibbs function Gr. Using Gr, one can write the Appell-Gibbs differential system:

Gr/ω̇x=FGx,Gr/ω̇y=FGy,Gr/ω̇z=FGz,E59

where FGis the external force vector acting in ball center G. Vector FGis determined subsequently using the virtual displacements principle. Let us introduce the quasi-coordinates φx,φy,φzwhere ωx=φ̇x,ωy=φ̇y,ωz=φ̇z. The only external force acting in the ball center is the gravity. Therefore, the elementary work performed can be expressed as

δFG=mgδuGz.E60

Virtual displacement δuGzcan be determined using the third non-holonomic constraint in Eqs. (57). It holds

δuGz=ρuCyδφxuCxδφy,E61

and therefore

δFG=mgρuCyδφxuCxδφy.E62

At the same time, the elementary work can be expressed in terms of quasi-coordinates:

δFG=FGxδφx+FGyδφy+FGzδφz.E63

Comparing coefficients at respective virtual components δφx,δφy,δφz, in Eqs. (61) and (63), one obtains

FGx=ρmguCy,FGy=ρmguCx,FGz=0.E64

The damping will be introduced later in Section 6.3 in order to separate energy conservative approach and enable to discuss various stationary regimes with respect to parameter and excitation settings.

### 6.2. Governing system

Carrying out the differentiation outlined in Eqs. (59), and respecting Eqs. (64), it can be written after some adaptations:

Jsω̇xuCxΩ̇s=u¨Ay+ρωzu̇Cxωxu̇CzuCzRuCyg+ρωxu̇Cyωyu̇Cx,
Jsω̇yuCyΩ̇s=u¨Ax+ρωyu̇Czωzu̇CyuCzR+uCxg+ρωxu̇Cyωyu̇Cx,E65
Jsω̇zuCzRΩ̇s=u¨Ax+ρωyu̇Czωzu̇CyuCyu¨Ay+ρωzu̇Cxωxu̇CzuCx,

where

Ω̇s=uCxω̇x+uCyω̇y+uCzRω̇z,Js=J+mρ2R2/mρ2.mass interia moment of the ballwithrespecttocenterofthecavityE66

It can be shown that Ω̇s=0and therefore the second column on the left side of the system Eqs. (65) should be omitted. External excitations are specified by movement or acceleration in the point A. Hence, kinematic excitation in Ais given as follows:

EAxt=u¨Ax/ρ,EAyt=u¨Ay/ρ,EAzt=0,E67

as it can be seen in Eqs. (65). Provided we need to investigate the response processes in a vertical plane, only one component remains non-zero and the second vanishes as well.

In order to obtain the system Eqs. (65) in the form with first-time derivatives concentrated on the left side, the first derivatives u̇Cin its right sides should be expressed in displacements uCusing non-holonomic constraints Eqs. (57):

u̇Cx=ωyuCzRωzuCy,u̇Cy=ωzuCxωxuCzR,u̇Cz=ωxuCyωyuCx.E68

Therefore, we obtained the system of six non-linear ODEs (Eqs. (65) and (68)) in a normal form with six unknown functions of time: uCx,uCy,uCz,ωx,ωy,ωz. Vector uCdepicts displacements of the contact point and can be used to study the movement of the ball from a global point of view. Detailed behavior of the ball as a rotating body is given by angular velocities ω. If the time history of rotation should be traced, then a subsequent run is necessary to obtain rotations by means of Euler angles as solution of the system of three ODEs with an input of angular velocities ω.

### 6.3. Influence of the damping

Influence of the damping will be taken into account. Basically, two sources of the energy dissipation are ruling in the system: (1) dissipation due to air dynamic resistance and (2) energy loss in contact of the cavity and rolling ball. The former one can be neglected with respect to obvious geometric configuration of the device and relative velocity ball/cavity. Concerning the latter one, complicated energy dissipating processes are ruling in contact of the ball with cavity. Nevertheless, supposing that no slipping arises in the contact, the dissipation process can be approximated as proportional to relevant components of the angular velocity vector ωand the quality of the cavity/ball contact. Considering the obvious setting, the respective material coefficients characterizing the rolling movement of the ball can be considered constant regardless of the direction in the tangential plane to the cavity in the point C, see Figure 4. The coefficient determining the rotation resistance around the normal vector nin the contact point Cis different as a rule. Therefore, the resistance moment vector Dcan be expressed in moving coordinates p,q,n, see Figure 4, as follows:

D=DpDqDnT.E69

Components of the above vector can be written in a form as follows:

Dp=κrωp,Dq=κrωq,Dn=κsωn,E70

where κr,κsare coefficients of “viscous resistance” of rolling and spinning. Their meaning is: the moment for a unity rotation per second, that is (Nms/rad).

Turning of the vector DG=DGxDGyDGzTexpressed in xyzcoordinates into the vector Dcan be written as

D=TCDG,E71

TC=xCzC+R,yCzC+R,νRyCν,xCν,0xCR,yCR,zC+RRE72

where ν2=xC2+yC2. The matrix TCis orthogonal and, therefore, the inverse transformation goes using matrix TC1=TCT, in particular:

DG=TCTD,E73

Components of the vector DGshould be incorporated onto the right side of Eqs. (59), where right sides should be completed. It means that the elementary work δFGfollowing Eq. (60) must be completed by a negative dissipating work due to DG.

δFG=mgδuGzDGδφ,E74

Repeating the further derivation like in Section 6.2, one can revisit the system Eqs. (65) and (68), where the right sides are completed and instead (Eqs. (65)) they read:

Jsω̇x=(u¨Ay+ρωzu̇Cxωxu̇CzuCzRuCyg+ρωxu̇Cyωyu̇Cx)DGx/m,Jsω̇y=(u¨Ax+ρωyu̇Czωzu̇CyuCzR+uCxg+ρωxu̇Cyωyu̇Cx)DGy/m,Jsω̇z=(u¨Ax+ρωyu̇Czωzu̇CyuCyu¨Ay+ρωzu̇Cxωxu̇CzuCxDGz/m.E75

Terms DGx/m,DGy/m,DGz/mwhich are linear functions of ωx,ωy,ωzdetermine the viscous type of the damping, although intensity in individual coordinates is variable depending on the position of the ball within the cavity.

### 6.4. Ball trajectories within the fixed cavity due to initial conditions

A large program of a ball trajectory investigation within a spherical cavity has been performed using the differential system (Eqs. (68) and (75)). Basically, it consists of two groups which are briefly illustrated in this and the next subsections. The first group concerns the fixed cavity (no excitation is applied). The only source of energy introduced is given by the initial deflection of the ball from equilibrium position in the point A(“southern pole”), or in other words by non-homogeneous initial conditions.

Differential system (Eqs. (68) and (75)) admits a number of singular solutions which can serve as separating limits of zones within which regular solutions exhibit certain character of trajectory shape. Some of them can be found analytically from the differential system taking into account their special properties concerning individual response component along the trajectory as a whole or in certain points of these curves. For details, special papers should be referred. Take a note that most of them emerge when no damping is considered. The reason is that the trajectory should be quasi-periodic (or cyclic-stationary), which is impossible when damping is respected and no external energy supply is considered. Trajectories start in a certain point on a meridian into which the ball is elevated. Then, it is thrown horizontally along the cavity parallel circle. Let us mention a few of the most important:

1. circular trajectory in horizontal plane. No initial spin is considered ωn0=0. The impulse applied corresponds with the initial velocity ω=ωps,0,0, where it holds for ωps:

ωps=gMρuCz02RuCz0J+Mρ2R2RuCz0.E76

This case is the most important and can be called separating circle (SC).

1. 2. circular trajectory in inclined plane, see Figure 5, ω0=100.0,0.0,0.0. This state is exactly valid for ωp0. The space spiral type trajectory changes from SC upwards successively into the upper hemisphere. Before the limit state for ωp0is reached, the osculation plane of the trajectory can be recognized. It rotates around the vertical axis with a descending angular velocity as far as it vanishes and osculation and operating planes coincide.

2. 3. trajectory of “kings crown form,” see Figure 5, ω0=5.817,0.0,5.0537. Cases, when the initial spin is considered. For a special value of ωn0=5.0537takes a shape visible in the picture. The apexes of this curve correspond to ω=0.0,0.0ωn0and u=0.0,0.0,0.0, which is a clue to find forms and parameters of this special case. This trajectory is reached from SC, increasing the initial spin velocity until the limit value. If it is lower, the trajectory has the spiral form. For a higher value, it became a curly form, see Figure 5, ω0=5.817,0.0,10.0. The limit state for infinite initial spin represents the ball apparently fixed in the initial point and not moving neither horizontally nor vertically.

Let us have a look at the bottom two pictures in Figure 5. They respect the influence of the damping. Coefficients κr,κsare different as it corresponds to conditions in the real system. The left demonstrates trajectory for positive initial spin and the right for negative initial spin. The transition through limit cases mentioned earlier is visible. The trajectory obviously finishes in the bottom “southern pole” of the cavity.

### 6.5. Ball trajectories within kinematically excited cavity

The second group of tests deals with the cavity which is undergone to kinematic excitation in a horizontal plane (only one-direction excitation is reported here).

Two extensive series of tests demonstrate the auto-parametric character of the system. In the first series, the response has been evaluated separately for every excitation frequency ωstarting from homogeneous initial conditions. Figure 6 shows some selected results of numerical simulations which follow from the differential system (Eqs. (68) and (75)). We briefly point out a couple of features visible in Figure 6. In the picture (a), we can see the maximal horizontal amplitude of the ball trajectory, when the cavity is kinematically excited in the horizontal plane in xdirection. The solid curve represents max uCxand the dashed curve is max uCyas functions of the exciting frequency ω. We can see that in the interval ω02.84, the semi-trivial solution is stable and so uy=0. The point ω=2.84is a beginning of the resonance zone, which spans in ω2.842.99, where auto-parametric resonance occurs and amplitudes of both response components are commeasurable. For ω>2.99, the semi-trivial solution is regained. Samples of the trajectory shape are plotted in picture (b) for four frequencies ω=2.84,2.88,2.92,2.96. Their vertical views demonstrate the character of the semi-trivial and the auto-parametric resonance states. Take a note that the trajectory since ω=2.94is a simple ellipse-like curve, which does not exhibit any symptom of a chaotic process. Compare this finding with analysis concerning the sweeping up and down excitation frequency for ωaround and above B2bifurcation point (BP) (see Figure 7 and explanation later in this subsection).

The second series has been controlled by sweeping the excitation frequency up and down in a large interval and in several detailed regimes in the area of the auto-parametric resonance zone. A few of the results are visible in Figure 7. Picture (a) demonstrates amplitudes max uCx(solid curves) and max uCy(dashed curve) and the total amplitude uCrin the interval ω1.08.5. Picture (b) is the magnified detail of picture (a) within the interval ω2.803.05in order to make visible the resonance zone.

Let us pay attention to bifurcation points (BPs). There are obviously concentrated in the resonance zone. In principle, they can be classified into two categories. The most important reveal B1and B2. In the latter one, two branches start. The lower one bl2approaches zero for ωwhich indicates the non-moving ball in the vertical view. This branch takes place in the vertical plane and basically has a form of semi-trivial solution. Its stability increases with rising ω>ωB2as it follows from decreasing negative values of the Lyapunov exponent and of inspection of the relevant stability basins. The upper branch bu2is spatial. It follows from the resonance zone where the spatial response type has a chaotic character. The relevant attractor reveals as an annular concentric area with diminishing width with increasing ω. The trajectory very quickly approaches a circular form in the horizontal plane. Its level with respect to the vertical axis rises and approaches “equatorial” position. However, the stability of this trajectory decreases, and we can see in Figure 7 that around ω=8.0even numerical perturbations of the integration process can overcome the stability limit (despite very small integration step) and the response trajectory falls down to the lower branch in the point D2. Its position is not fixed. If hypothetically zero perturbation occurs, it could shift to infinity and approach together with the branch bu2the asymptote at the level R=1. Observing black max uCxand red max uCyparts of bu2, we can see that they are getting coincide with increasing ω. It means that trajectory approaches the circle with radius R=1.

Let us briefly discuss the shape of the response amplitudes for ωbelow BP B1and B4. The BP B1is reached sweeping up along the branch bl1, when it loses planar character passing through B4. In such a case, the spatial response type emerges, exhibiting a chaotic response since B1. This fact is obvious also looking at the dashed red curve representing uCy, which is trivial as far as B4and can bring the system from the semi-trivial solution into the auto-parametric resonance starting B4. Take a note that passing BP B4, planar response can remain in force, if any perturbation is avoided. It meets in B3the branch bl2following also a planar path for swept up ω. Branch bu1starts in B1. Its stability rapidly decreases with descending ω. The point D1illustrates its limited extent in sub-resonance zone. This feature is visible observing the curve bs, which represents a limit separating the area of attraction to bu1and to bl1. Take a note that bsstarts in B4and approaches D1, despite hypothetically it goes together with D1as far as the vertical axis in the point R=1. The bscan be earned from stability basins for ωin the adequate interval for initial value ωp=0. It corresponds to amplitude of uCxas a testing value for decision about affiliation to bu1or bl1attractiveness.

The interval between B1and B2includes the spatial response, see non-trivial amplitude max uCy. The spatial response has a chaotic character, as it has been already outlined in the previous paragraph, when commenting the branch bu2.

## 7. Conclusion

The common physical origin of Lagrangian and Appell-Gibbs approaches has been shown. It originates from the equilibrium of energy-level evolution in time on one side and power supply together with energy dissipation on the other side. Various formulations of this principle lead finally to different variational principles, although they follow from the same minimization of the energy spent to system response portrait. Comparing individual sections of the chapter, we can see that each one of commonly used procedures based on particular energy formulations is preferable for a certain type of problems. It can be concluded that there does not exist a single universal approach which should be recommended.

Some detailed properties of both approaches have been demonstrated in Sections 5 and 6. Both of them discuss non-holonomic problem of the ball movement within the spherical cavity under external excitation. The former one deals with a simple planar problem and shows that the Lagrangian approach is easily applicable to obtain reasonable results as far as a wide parametric discussion, which enable to earn a detailed insight into the system dynamic properties. The latter alternative represents the full space problem with six DOFs and three non-holonomic constraints. Some earlier studies tried to formulate this problem also in Lagrangian style using Lagrange’s multipliers. Finally, it proved that the relevant governing differential system is too complex and does not enable appropriate detailed analysis of dynamic properties of the system. Therefore, the space problem outlined in Section 6 has been formulated using Appell-Gibbs approach. Transparent results have been obtained as needed for practical purposes in a device design and in further study of multi-body system dynamics. Take a note regarding the classification of singular solutions and their applicability for detailed analysis, stability of various regimes of the system under kinematic excitation, transitions among semi-trivial, auto-parametric, chaotic, and other states typical for nonlinear system. Let us add that both 2D and 3D problems have been investigated respecting the full nonlinearity without any simplifications of transcendent functions and thus enabling to study all effects without any limitations in amplitudes.

A certain shortcomings which apparently follow from the knowledge of rotation velocities only (no rotations themselves are calculated) can be disregarded, when displacements have been obtained. The rotation velocities represent mostly satisfactory information. Nevertheless, if rotations are still needed, there exist several variants of a simple differential system (following rotation vector definition) relating velocity and rotation vector components. This system can be subsequently easily solved, when necessary. A hidden complexity of the Lagrangian approach follows from an implicit connection of both parts, which are independent when Appell-Gibbs procedure is applied.

## Acknowledgments

The kind support of the Czech Science Foundation project No. 17-26353 J and of the RVO 68378297 institutional support are gratefully acknowledged.

## More

© 2018 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution 3.0 License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

## How to cite and reference

### Cite this chapter Copy to clipboard

Jiří Náprstek and Cyril Fischer (July 18th 2018). Appell-Gibbs Approach in Dynamics of Non-Holonomic Systems, Nonlinear Systems - Modeling, Estimation, and Stability, Mahmut Reyhanoglu, IntechOpen, DOI: 10.5772/intechopen.76258. Available from:

### chapter statistics

2Crossref citations

### Related Content

Next chapter

#### Canonical Generalized Inversion Form of Kane’s Equations of Motion for Constrained Mechanical Systems

By Abdulrahman H. Bajodah and Ye-Hwa Chen

First chapter

#### On Nonoscillatory Solutions of Two-Dimensional Nonlinear Dynamical Systems

By Elvan Akın and Özkan Öztürk

We are IntechOpen, the world's leading publisher of Open Access books. Built by scientists, for scientists. Our readership spans scientists, professors, researchers, librarians, and students, as well as business professionals. We share our knowledge and peer-reveiwed research papers with libraries, scientific and engineering societies, and also work with corporate R&D departments and government entities.

View all Books