Kinematic support modeling in Sage

Full Text

1. Introduction As one of the types of seismic protection of buildings in seismically active regions of the Earth, among others, auxiliary “seismic suppression” supports are used. A large number of standards of such supports are known with a high level of reliability of operation and a high ability to damp seismic waves of high magnitude. They are high-tech in manufacturing, are sold at high prices and, thanks to efficiency and reliability, are in high demand among large construction companies and in areas with a high standard of living. At the same time, in poor areas of developing countries, citizens and municipalities cannot afford to use them. However, it is precisely in such regions that residents most often suffer from the devastating effects of earthquakes. The Soviet, and then the Russian school of seismic protection architecture was able to offer inexpensive and effective solutions in the form of so-called kinematic supports. In our country and abroad, a large number of active seismic protection systems for buildings have been proposed, developed and applied. Among them are those proposed by A. Kurzanov, S.Yu. Semenov [1], [2], Yu.P. Cherepinsky [3], V.V. Nazina, etc. Some of these systems were practically implemented in real buildings [1], [2], [4], [5], making it possible to assess their workability for building industry. Vibration tests were carried out at many facilities [6], [7], which provided experimental data on the behavior of these systems under dynamic impacts. However, essentially all developed systems need additional analysis under full-scale conditions. Therefore, many aspects of the real behavior of seismic protection systems are difficult to study theoretically or on models due to the very large number of factors affecting the behavior of the structure during an intense earthquake. Kinematic supports are vertically placed cylinders on which the building rests. Neither the place of entry of the support into the ground or concrete slab, nor the contact with the horizontal slab of the building placed on such supports, are fixed rigidly. Supports can be made in the form of short concrete pillars with an outer cage of steel pipe or a reinforcing cage of carbon composite or basalt composite nets. It is promising to use in the construction of concrete racks with dispersed reinforcement basalt fiber, since such concretes have increased resistance to cracking and tensile strength during bending. The essence of the kinematic support is that when the base is displaced by a certain design value, the building slightly rises, receiving some additional kinetic energy. In this case, a returning torque arises, bringing the ”base-building” system to its original state (position before the earthquake). Residential buildings constructed using seismic isolating supports have fullscale confirmation of the reliability of the structure and have proven themselves successful in experimental studies [1], [2]. The aim of our work is to create an adequate mathematical model of the support and its interface with the building, which will help to design kinematic supports taking into account the operational requirements of customers. The solution to this problem can be approached from two sides. First, it is possible to create a model of elastic support in the Ansys, system [8]. With this approach, the main difficulty is the selection of adequate boundary conditions for the place where the support is in contact with the soil and base plate of the building. The second approach proposes to consider the support from the point of view of classical mechanics, assuming it to be a perfectly rigid body, oscillating in a vertical plane above a fixed horizontal absolutely smooth plate, but more adequately describe the interaction of the support with the ground and the building plate. From the point of view of analytical mechanics, the motion of kinematic supports is the motion of a complex system of bodies with non-holding bonds [9]. The mathematical basis of the dynamics of such systems was developed by outstanding mathematicians J. Dalamber, S. Poisson, Yu.A. Arkhangelsky, V.V. Kozlov, A.P. Markeev and others. In our work, we will follow the formalism developed by A.P. Markeev. Models of rigid body mechanics are systems of a large number of ordinary differential equations (ODEs) that are not resolved with respect to derivatives [10], [11]. They, as a rule, do not allow an analytical solution. Therefore, we are going to use the finite difference method. For its successful application, it is important to solve the equations for derivatives by increasing the order. Due to the large number of equations, this procedure turns out to be rather complicated. Therefore, it seems natural to execute it in a computer algebra system (CAS). Computer-aided study of such models is carried out in two stages: 1) symbolic transformations reducing the system to normal form, used in the standard formulation of the Cauchy theorem; 2) numerical solution of this system using the finite difference method (rk4). We chose the Sage system [12], because it can execute both the first and second stages. To test the described approach, we took the simplest model that describes the complete system of equations of motion of a massive rigid body along a fixed horizontal perfectly smooth plane [10], [11]. 2. Description of the mathematical model Let the seismic isolation support be a rigid body and have the shape of a cylinder, in which one of the bases has a spherical shape. The support is installed between the foundation of the building (a horizontal rigid plate) and the building itself. The support touches the foundation of the building always with a spherical end [13]. We assume that when oscillating or when horizontal forces act on the foundation of the structure, the point of contact of the support and the base plate always lies in the plane ????????????. To describe the vibration of the support, we will use the model of motion of a rigid body on the surface [10], [11] and adapt this model to our task. The motion of the body will be considered relative to the fixed laboratory coordinate frame ???????????? with the origin ???? at a certain point of the plane. The axis ???????? is directed vertically, ???? is the unit vector of inner normal to the body surface at point of the axis ????. Let ???????????????? denote the moving coordinate frame rigidly coupled to the body with the origin at its center of gravity ???? and the axes directed along the principal axes of inertia 1. The orientation of the body with respect to the fixed laboratory frame is specified by the Euler angles ????,????,???? or by the matrix of direction cosines ????????????. The unit vector of the ????-axis in the frame ???????????????? is specified by the components ????31, ????32, ????33: ????31 = sin????sin????, ????32 = sin????cos????, ????33 = cos????. Assume ????, ????, ???? to be the principal axes of inertia with respect to the gravity center. Let ???? be the point of contact between the horizontal plane ???????????? and the support (see Figure 1). Its coordinates ????, ????, ???? in the frame ???????????????? will be functions of angles ????,????, determined from the form of equation ????(????,????,????) = 0 that specifies the shape of the body surface [10], [11]. The sign of function ???? is chosen such that ???? = -. Then the quantities ???????????? for the axis ???? are expressed in terms of the Euler angles as 1 ???? ???? = - , |∇????| ∇???? { 1 ???? ⎨????32 = sin????cos???? = -|∇????| ∇????, (1) { { 1 ???? {⎩????33 = cos???? = -|∇????| ∇????. Figure 1. Motion of a perfectly rigid body above the perfectly smooth horizontal plane Following the studies presented in Refs. [10], [11], [14], [15], let us consider the full system of equations of motion for a massive rigid body on a fixed horizontal perfectly smooth plane and introduce the following unknown functions of time ????: center of gravity coordinates ????, ????, ???? of the body in the laboratory frame; Euler angles ????,????,????; components ????, ????, ???? of radius vector ???? of point ???? of contact of the support and the plane (base plate) relative to the gravity center, and magnitude ???? of the normal reaction of the plane. To determine the unknowns listed above the following equations and relations will be used: a) The equations that represent the theorem of momentum variation. Theexternal forces are the reaction of the plane ???? = ???? ⋅ ???? directed vertically (???? ⩾ 0) and the gravity force. The equations are: ⎧{????????̈ = 0, ????????̈ = 0, (2) ⎨{⎩????????̈= -???????? + ????, where ???? is the body mass, ???? is the free fall acceleration. b) The equations that represent the theorem of angular momentum variation: ⎧{???????????????????? + (???? - ????)???????? = ????(????????33 - ????????32), { { ???????? ⎨{???????????? + (???? - ????)???????? = ????(????????31 - ????????33), (3) { ???????? {⎩???????????? + (???? - ????)???????? = ????(????????32 - ????????31), where ????, ????, ???? are the projections of the angular velocity vector ???? on the axes of the coordinate system ????????, ????????, ???????? rigidly bound to the body; ????, ????, ???? are the principal moments of inertia with respect to these axes. c) The relations represented by the Euler kinematic equations: ⎧???? = ????̇ sin ????sin ???? + ????̇ cos ????, { ⎨???? = ???? ????????????????̇ cos ???? - ????̇ sin ????, (4) {⎩???? = ????̇ cos ???? + ????.̇ d) The Poisson equations: ⎧{????̇31 = ????32???? - ????33????, ⎨{⎩????̇????̇3233 == ????????3331???? - ???????? - ????3231????????, (5) indicating the fact that vector ???? defines an invariable direction in the fixed coordinate frame ????????????????; e) The equation of the body surface in the coordinate frame rigidly bound to the body and having the origin at its center of gravity: ????(????,????,????) = 0. (6) In our case, the body surface equation is the equation of a sphere (the support has spherical shape at the point of contact with the plate): ????2 + ????2 + (???? + ????)2 = ????2. (7) f) The constraint equation ???? = -???? ⋅ ???? = -(????31???? + ????32???? + ????33????) (8) means that the support moves contacting with the plate all the time. Equations and relations (1)-(8) determine a complete system of equations of motion of a massive rigid body on a fixed horizontal perfectly smooth plane that can be written in the following form: ⎧????????̈ = 0, {{????????̈ = 0, {{????????̈= -???????? + ????, { ???????? {{???????????? + (???? - ????)???????? = ????(????????33 - ????????32), {{{{???????????????????? + (???? - ????)???????? = ????(????????31 - ????????33), { ???????? {???????????? + (???? - ????)???????? = ????(????????32 - ????????31), { {{????̇31 = ????32???? - ????33????, {{????̇32 = ????33???? - ????31????, {????̇33 = ????31???? - ????32????, (9) ???? + ????32???? + ????33????), { 1 ???????? {{????31 = sin ???? sin ???? = -|∇????| ???????? , {{{????32 = sin ???? cos ???? = -|∇????|1 ????????????????, { {{{????33 = cos ???? = -|∇????|1 ???????????????? , { {{???? = ????̇ sin ????sin ???? + ????̇ cos ????, {{???? = ????̇ sin ????cos ???? - ????̇ sin ????, {⎩???? = ????̇ cos ???? + ????.̇ From relation (1) and equation (7), we get explicit expressions of ???????????? for the spherical base of a kinematic support: ⎧{{????31 = -|∇????|1 ???????????????? = -√4????2 + 4????21+ 4(???? + ????)2 ⋅ 2???? = -????1 ????, { { 1 ???????? 1 1 ⎨{????32 = -|∇????| ???????? = -√4????2 + 4????2 + 4(???? + ????)2 ⋅ 2???? = -????????, (10) 1 ???????? 1 1 ⎩????33 = -|∇????| ???????? = -√4????2 + 4????2 + 4(???? + ????)2 ⋅ 2???? = -????(???? + ????). Using the formulae (10), is it easy to express the components ????, ????, ????, therefore, we can get rid of the Poisson equations in the system (9), since the values of ???????????? and components ????, ????, ???? are already known. According to the two first equations of (9), the center of gravity moves so that its projection on the base horizontal plane moves rectilinearly. Hence, it is obvious that the appropriate equations can be disregarded, too. After a number of executed transformations and assumptions, the system (9) takes the form ⎧???? = ????????̈+ ????????, {{{{???????????????????? + (???? - ????)???????? = ????(????????33 - ????????32), { ???????? {{???????????? + (???? - ????)???????? = ????(????????31 - ????????33), {{{???????????????????? + (???? - ????)???????? = ????(????????32 - ????????31), + (???? + ????)2 = ????2, {{???? = -???? ⋅ ???? = -(????31???? + ????32???? + ????33????), 1 ???????? 1 (11) ⎨????31 = sin ???? sin ???? = -|∇????| ???????? = -????????, { { { 1 ???????? 1 {{????32 = sin ???? cos ???? = -|∇????| ???????? = -????????, { { 1 ???????? 1 {{????33 = cos ???? = -|∇????| ???????? = -????(???? + ????), { {???? = ????̇ sin ????sin ???? + ????̇ cos ????, { {{???? = ????̇ sin ????cos ???? - ????̇ sin ????, {⎩???? = ????̇ cos ???? + ????.̇ Equations (11) are reduced to the form convenient for using the computer algebra system Sage, so that further simplifications and the solution of these equations will be carried out in this CAS. In the next section we investigate the solubility of this system of equations in Sage. 3. Resolving the system with respect to derivatives The system (11) includes both differential equations and algebraic ones (relations). Let us use Sage to reduce it to a simper form. It is obvious that from the system (11) via the values of the direction cosines ????31, ????32, ????33, we can explicitly express variables ????, ????, ???? as follows: ⎧{???? = -???? ⋅ ????31 = -???? ⋅ sin ???? sin ????, ???? = -???? ⋅ ????32 = -???? ⋅ sin ???? cos ????, (12) ???? - ????. Provided that ????, ????, ???? (12) and ????31, ????32, ????33, are known, it is possible to express the values of the function ????. Let us write in Sage the equation ???? = -???? ⋅ ???? = -(????31???? + ????32???? + ????33????) from system (11) and substitute into it the expressions for ????, ????, ???? and ????31, ????32, ????33, then function ???? will be represented by the following symbolic expression: ???? = ???? ⋅ cos????2 ⋅ sin????2 + ???? ⋅ sin????2 ⋅ sin????2 + (???? ⋅ cos???? + ????) ⋅ cos????. We substitute the obtained expression of ???? into equation ???? = ????????̈+ ???????? and calculate the normal reaction of the plane: ???? = (2⋅????⋅cos????2 ⋅cos????2 ⋅????2̇ +2⋅????⋅cos????2 ⋅sin????2 ⋅????2̇ -2⋅????⋅cos????2 ⋅sin????2⋅ ⋅ ????2 - 2 ⋅ ???? ⋅ sin????2 ⋅ sin????2 ⋅ ????2̇ + 2 ⋅ ???? ⋅ cos????2 ⋅ cos???? ⋅ sin???? ⋅ ???? + 2 ⋅ ???? ⋅̈ cos????⋅ ⋅ sin????2 ⋅ sin???? ⋅ ???? - ???? ⋅̈ cos????2 ⋅ ????2̇ + 2 ⋅ ???? ⋅ sin????2 ⋅ ????2̇ - (???? ⋅ cos???? + ????) ⋅ cos????⋅ ⋅ ????2 - ???? ⋅ cos???? ⋅ sin???? ⋅ ???? - (???? ⋅̈ cos???? + ????) ⋅ sin???? ⋅ ????) ⋅ ???? + ???? ⋅ ????.̈ As soon as all quantities ????, ????, ????, ????31, ????32, ????33, ????, and ????, ????, ???? are explicitly expresses, we substitute them into equations 2-4 of the system (11) and of the form ⎧Φ(????,????,????,????,…̇ ????) = 0,̈ { ⎨Ψ(????,????,????,????,…̇ ????) = 0,̈ (13) {⎩Θ(????,????,????,????,…̇ ????) = 0.̈ Explicit expressions for Φ,Ψ,Θ were be found in Sage [12]. The system (13) is linear with respect to ????,̈ ????,̈ ????.̈ Now we resolve the obtained system of differential equations with respect to higher derivatives using the function solve (). Ultimately, we arrive at the system of differential equations of the 6-th order resolved with respect to higher derivatives. The system incorporates three equations of the second order with respect to Euler angles. For further solution of the problem we have to decrease the order of the differential equations to the first one. To reduce the order of the system of differential equations, we perform the changes ???? = ????̇ , ???? = ????̇ , ???? = ????̇ . As a result, we get a system of six first-order differential equations. 4. Numerical experiments in SAGE Let us implement explicit Euler method, in order to confirm the absence of errors related to transformation of data types. The calculations by the Euler method will be performed using the grids with ???? = 400,800,1600 under the following initial conditions: ????0 = 0.1, ????0 = 0.1, ????0 = ????+0.5, ????0 = 1.5, ????0 = 1.5, ????0 = 0, ℎ = 2????/????, ???? = 0, ???? = 304, ???? = 304, ???? = 400, ???? = 1.5, ???? = 0.2, ???? = 500, ???? = 9.8 where ???? is the number of steps; ℎ is the step; ????0, ????0, ????0 are the Euler angles at the initial moment of time; ????0, ????0, ????0 are the initial velocities; ????, ????, ???? are the moments of inertia. The results obtained by the Euler method are presented in Figures 2-7. For ???? = 400 the results are plotted by dash lines, for ???? = 800 by dash-dot lines, and for ???? = 1600 by dot lines. The accuracy of the obtained results is estimated using the Richardson method (see Figures 9, 10). Figure 8 shows graphs. The distance between the center of mass and the plane ground oscillates as it shown on Figure 8. Figure 2. The angle ???? as a function of time Figure 3. The angle ???? as a function of time for the grids with ???? = 400,800,1600 for the grids with ???? = 400,800,1600 Figure 4. The angle ???? as a function of time Figure 5. Rate of angle ???? change with time for the grids with ???? = 400,800,1600 for the grids with ???? = 400,800,1600 Figure 6. Rate of angle ???? change with time Figure 7. Rate of angle ???? change with time for the grids with ???? = 400,800,1600 for the grids with ???? = 400,800,1600 Figure 8. The center of mass height above the plane for the grids with ???? = 400,800,1600 Figure 9. Accuracy of Euler solution Figure 10. Accuracy of Euler solution estimated by Richardson method estimated by Richardson method for the grids with ???? = 400,800,1600 for the grids with ???? = 400,800,1600 The results obtained by the Euler method for the assumed boundary conditions show that the angles ????, ???? linearly increase with time (see Figures 2 and 3), while the angle ???? experiences harmonic-like oscillations with growing amplitude (see Figure 4). Figures 5-7 show that the rates of change of the Euler angles vary according to a harmonic-like law, the amplitude of oscillations increasing with time. 5. Conclusions A crude mathematical model of a kinematic support is constructed, in which the support is considered as a perfectly rigid body oscillating in a vertical plane above a fixed horizontal perfectly smooth plate. The approximate model is rigid and does not take friction into account. A procedure for reducing the system of differential equations of the model to a form suitable for applying the finite difference method is described and implemented in the Sage computer algebra system. An explicit Euler method is implemented for grids with the number of partitions 400, 800, 1600 and the accuracy of the solution is estimated by the Richardson method. For the initial and boundary conditions specified by us, the time dependence of Euler angles and their rates is determined. The results of numerical experiments are consistent with the general idea that small deviations lead to small oscillations of the support. The problems and lines of study that will be addressed at further stages of the research are identified as follows. First, to solve the problem, we have to to determine the correct additional conditions under which the construction is stable or loses stability, i.e., when at strong ground vibrations the center of gravity horizontally shifts beyond the limits of return movement of the support and the support begins to tip over [16]. After finding the correct initial and boundary conditions, to create an adequate mathematical model of the support, we will take friction into account, specifying the reaction of the plane (base plate), and the effect of earthquakes on the movement of the supports.

About the authors

Oleg K. Kroytor

Peoples’ Friendship University of Russia (RUDN University)

Email: kroytor_ok@pfur.ru
Postgraduate of Department of Applied Probability and Informatics 6, Miklukho-Maklaya St., Moscow, 117198, Russian Federation

Mikhail D. Malykh

Peoples’ Friendship University of Russia (RUDN University)

Email: malykh_md@pfur.ru
Doctor of Physical and Mathematical Sciences, assistant professor of Department of Applied Probability and Informatics 6, Miklukho-Maklaya St., Moscow, 117198, Russian Federation

Sergei P. Karnilovich

Peoples’ Friendship University of Russia (RUDN University)

Email: karnilovich_sp@pfur.ru
assistant professor, Ph.d., assistant professor of Institute of Physical Research and Technology 6, Miklukho-Maklaya St., Moscow, 117198, Russian Federation

References