Mathematical modelling of the classic dynamical system of beam and sphere using pneumatic actuator 

María A. Terán, Onexi L. Mendoza, José G. Díaz y Francisco J. Arteaga 

Unidad de Investigación en Automatización Industrial, Escuela de Ingeniería Eléctrica, Centro de Investigación y Tecnología en Automatización, Electrónica y Control CITAEC. Facultad de Ingeniería, Universidad de Carabobo, Valencia, Venezuela E-mail: mteran25@hotmail.com, olmendozam@hotmail.com, farteaga@uc.edu.ve 

Abstract 

This research presents a mathematical model of the classic system of beam and sphere using pneumatic actuator. In the work, models are developed for the plant, a pneumatic double effect cylinder and its proportional valve. These models are obtained from the equations of Lagrange, kinematics of mechanisms, the law for the ideal gas, energy equation for an open system, equation of continuity for the mass flow rate, and equation of the nozzle. The non linear effects produced by the compressed air flow dynamics and the friction between the piston and the cylinder walls are considered to develop the model equations. The model is validated by means of simulation in Matlab Simulink. 

Key words: Beam and sphere system, pneumatic double effect cylinder, proportional valve, dynamic mathematical model. 

Modelación dinámica y simulación del sistema clásico de viga y esfera mediante accionamiento neumático 

Resumen 

Esta investigación presenta un modelo matemático dinámico del sistema clásico de viga y esfera utilizando un actuador neumático. En este trabajo se desarrollan modelos tanto para la planta como para el cilindro neumático de doble efecto y su válvula proporcional. Estos modelos son obtenidos a partir de las ecuaciones de Lagrange, cinemática de mecanismos, ley del gas ideal, ecuación de energía para un sistema abierto, ecuación de continuidad para el flujo másico y ecuación de la tobera. Para desarrollar las ecuaciones se toman en cuenta los efectos no lineales producidos por la dinámica del flujo de aire comprimido y la fricción existente entre el pistón y las paredes del cilindro. El modelo es validado mediante simulación en Matlab Simulink. 

Palabras clave: Sistema viga y esfera, cilindro neumático de doble efecto, válvula proporcional, modelo matemático dinámico. 

Recibido el 27 de Febrero de 2004 

En forma revisada el 20 de Noviembre de 2005 

1. Introducción 

Durante muchos años se ha utilizado el sistema de viga y esfera como ejemplo de un sistema inestable [1]. Por ejemplo, en la Universidad de Michigan se ha estudiado su dinámica, cuando el accionamiento empleado es un motor eléctrico DC con acoplamiento directo [2]. Por otra parte, la modelación de los sistemas con accionamiento neumático ha sido tema de numerosos trabajos durante los últimos quince años. Con la simplificación en mente, los diversos modelos propuestos han considerado un comportamiento lineal [3-6]; sin embargo, en esta investigación se ha tomado en cuenta la alta no linealidad existente en los elementos neumáticos, ocasionada por la dinámica propia del aire comprimido y la fricción en el interior del cilindro. 

2. Sistema de viga y esfera usando accionamiento neumático 

El sistema básico consta de una esfera colocada sobre una viga, a uno de cuyos extremos se fija un actuador, mientras el otro actúa como pivote, tal como se observa en la Figura.1. Cuando la viga cambia su posición en un ángulo q, debido a la acción del actuador, la esfera rueda sin deslizar, estabilizándose en una posición r, diferente de la que tenía inicialmente. 

3. Modelo matemático de la planta 

El modelo matemático de la planta está constituido por dos elementos: el sistema de viga y esfera y el mecanismo de deslizador y manivela. 

3.1. Sistema de viga y esfera 

A partir de la ecuación de la energía cinética de una partícula y de la ecuación de Lagrange para dicha energía, se logra relacionar la posición de la esfera (r) con la posición angular de la viga (q). En coordenadas polares, estas ecuaciones son: 

Sustituyendo la ecuación (1) en la ecuación (2) y desarrollando las derivadas respectivas, se obtiene: 

En el diagrama de cuerpo libre de la esfera (Figura 2.a), al realizar la sumatoria de fuerzas en dirección de r, resulta: 

F_r = -f_re - m_esf g senq.    (4) 

Al sustituir (4) en (3) se obtiene: 

La sumatoria de torques en la esfera (Figura 2.b), produce: 

f_reR = Ja.        (6) 

Las ecuaciones (7) y (8) proporcionan las expresiones para la aceleración angular (a) y el momento de inercia con respecto al centro C (J): 

Sustituyendo las ecuaciones (7) y (8) en la ecuación (4) y despejando fr_e, se obtiene: 

Sustituyendo la ecuación (9) en la ecuación (5) y simplificando, queda finalmente: 

3.2. Mecanismo de deslizador  y manivela 

El mecanismo de deslizador y manivela, que se observa en la Figura 3.a, es el encargado de transmitir el movimiento o desplazamiento del cilindro neumático a la viga. El modelo matemático que se obtiene del análisis de este mecanismo permite obtener una expresión para la posición angular de la viga (q) en función de la posición del vástago del cilindro (y_p ). 

Ubicando a los vectores del mecanismo de deslizador y manivela (Figura 3.b) en el plano complejo, y con el fin de hallar el desplazamiento y, se escribe la ecuación vectorial correspondiente [8]: 

Le^jq + be^jb + ye^jf = 0.        (11) 

Usando la identidad de Euler en la ecuación anterior, resulta: 

y(cosf + j senf) + L(cosq + j senq) - b(cosb + j senb) = 0.   (12) 

Separando la expresión anterior en parte real e imaginaria, con f = 90, se obtienen: 

L cosq = b cosb

y + L senq = b senb.        (13) 

Elevando al cuadrado las expresiones anteriores, sumándolas, relacionando el desplazamiento del punto B (y) con respecto a la posición del pistón (y_p) a través de la expresión y = 0.500 – y_p , se obtiene finalmente la siguiente ecuación: 

Otra expresión que se obtiene del análisis del mecanismo de deslizador y manivela es la fuerza que ejerce la carga de la planta sobre el vástago del cilindro (F_L ). Para ello se toma como referencia la posición del mecanismo mostrada en la Figura 3.a y se aplica el método de fuerza equivalente [8] en el punto A de la viga OA. Este procedimiento se repite para el punto B de la manivela AB y relacionando q con respecto a b, mediante la expresión resulta: 

4. Modelo matemático dinámico del cilindro neumático de doble efecto 

4.1. Cámaras del cilindro 

El cilindro de doble efecto posee dos cámaras separadas. A continuación, se desarrolla un modelo matemático dinámico para cada una de ellas, vinculando el cambio de presión con la tasa de flujo másico y la velocidad de traslación del pistón. En algunos trabajos previos [9-11], al igual que en éste, las ecuaciones se derivan bajo las siguientes consideraciones: 

1. El aire comprimido es un gas perfecto. 

2. La presión de alimentación se considera constante. 

3. La presión y la temperatura dentro de las cámaras son homogéneas. 

4. El proceso se considera isotérmico. 

5. Se desprecian las energías cinética y potencial. 

6. Los procesos de carga y descarga son adiabáticos. 

7. No se consideran escapes de aire en las cámaras. 

El modelo más general para un volumen de gas consta de tres ecuaciones [12] - [13]: la ecuación de estado (ley del gas ideal), la ecuación de la conservación de la masa (ley de continuidad) y la ecuación de la energía. Las ecuaciones se aplican a cada cámara, considerando los volúmenes de control V₁ y V₂ (Figura 4.a). El análisis se hace para una cámara (la otra tiene una expresión análoga). 

El flujo másico se define como en [14]: 

La ecuación de continuidad para la tasa de flujo másico [14] es: 

La ecuación de la energía para un sistema abierto [14] es: 

Con las siguientes expresiones: 

Introduciendo las ecuaciones (19-22) en la ecuación (18), queda: 

También se aplican las siguientes ecuaciones [15]: 

Introduciendo las ecuaciones (24-26) en la ecuación (23) y despejando : 

Debido a que la transferencia de calor en el proceso de entrada y salida es diferente, se sustituye la relación de calor específico (k) por un coeficiente de característica térmica (a) [16] en el segundo término de la ecuación (27): 

Si se escoge el origen del desplazamiento del pistón en la mitad de la carrera del vástago, el volumen de cada cámara viene dado por la siguiente ecuación: 

Despejando en la ecuación (22) y sustituyendo por PA_iy_p, queda: 

Sustituyendo las ecuaciones (29) y (30) en la ecuación (28) y considerando además que en el proceso de carga de la cámara del cilindro los flujos másicos están dados por: y que en el proceso de descarga de la cámara, los flujos vienen dados por: , se obtienen las ecuaciones que describen el comportamiento de las cámaras del cilindro: 

4.2. Émbolo del cilindro 

El comportamiento del émbolo del cilindro puede ser modelado usando la segunda ley de Newton. En la Figura 4.b se muestra su diagrama de cuerpo libre, del cual resulta (F_f es la fuerza de fricción): 

Sustituyendo la ecuación (34) en la (33) se obtiene: 

La fuerza de fricción presenta un comportamiento impredecible, originado por diversos factores, como son: la condición de lubricación del aire, las condiciones de operación y, principalmente, de la velocidad. Esta fuerza de fricción ha sido objeto de estudio por algunos investigadores con el fin de obtener su modelo matemático dinámico [17] y [18]. Con estos modelos se ha desarrollado una descripción y representación más detallada y precisa del comportamiento de la fricción. Éstos van desde los más simples, pasando por otros más complejos como el modelo de la fuerza de fricción de Coulomb y el modelo de fricción viscosa, hasta los más elaborados, como lo es el modelo de LugGre, propuesto en [17], donde la fuerza de fricción puede expresarse como la sumatoria de una fuerza de fricción de Coulomb y otra viscosa: 

F_f = F_C + F_V.        (36) 

El modelo de fuerza de fricción de Coulomb ha sido validado en [19-21], con la validación experimental de Kang [22]; por consiguiente, éste es el modelo utilizado en este trabajo. 

La fuerza de fricción de Coulomb (F_C ) y la fuerza de fricción viscosa (Fn) se expresan como: 

Sustituyendo las ecuaciones (37) y (38) en la ecuación (36), se obtiene la siguiente expresión: 

5. Modelo matemático dinámico de la válvula reguladora de caudal 

En este sistema se utiliza una servoválvula de cinco vías y tres posiciones (ver símbolo en la Figura 5), la cual se caracteriza por ser de posicionamiento continuo entre sus dos ubicaciones extremas. 

A continuación se desarrolla un modelo del flujo de aire que se desplaza a través de un orificio de área variable de una válvula proporcional, en función de las presiones de las cámaras (P_u, P_d ) y de la presión de alimentación (P_s ) y sin tomar en cuenta los escapes de aire. Para este propósito se utiliza como ecuación básica la fórmula de la tobera [23], la cual expresa el flujo másico ideal isentrópico (para un solo orificio) de la siguiente manera: 

Con respecto al área efectiva del orificio de la válvula, su modelación se hace con la ayuda de curvas obtenidas experimentalmente en el trabajo de Bashir [24], cuya expresión analítica es: 

Para la función y, se identifican dos formas, cualitativamente iguales, pero cuantitativamente diferentes, las cuales corresponden al proceso de carga (y_C) y descarga (y_d), respectivamente. Sus expresiones analíticas se derivan de curvas obtenidas experimentalmente [24] y son: 

La ecuación 40 tiene dos expresiones, dependiendo del proceso que se esté realizando, ya sea de carga de aire del cilindro o de descarga de aire al cilindro , y las expresiones son: 

Durante el proceso de carga, la presión de alimentación se considera la presión aguas arriba, mientras que la presión en la cámara del cilindro es la presión aguas abajo; mientras que en el proceso de descarga, la presión en la cámara es la presión aguas arriba y la presión de la atmósfera es la presión aguas abajo. 

6. Simulaciones del modelo matemático dinámico del sistema 

Con la finalidad de validar el modelo matemático dinámico obtenido, se simuló el mismo en Matlab Simulink, utilizando tensiones de entrada variables para la servoválvula. 

En la Figura 6 se resume un primer experimento de simulación, en donde se presentan cambios en la tensión de alimentación de la servoválvula. Como se puede observar, tanto la posición del pistón como la inclinación de la viga siguen los cambios en la tensión de entrada; es decir, cuando en t = 1 seg, la tensión de entrada aumenta de 5V a cerca de 10V, el vástago comienza a avanzar hacia su posición de máxima extensión (0.10 m), alcanzándola justo antes de t = 2 seg. A esta extensión corresponde a la máxima inclinación de la viga (13°). En t = 2 seg ocurre otro cambio brusco en la tensión de entrada, la cual disminuye a 0V, con lo que el vástago retrocede y la viga se inclina hasta –10.56°, para un tiempo algo menor a 4 seg. 

Una vez en movimiento, la esfera necesita tiempo para detenerse y cambiar de dirección, debido a su inercia propia. Esto se observa en la gráfica, ya que, cuando la viga pasa por la posición de equilibrio (0°), en t = 3 seg, la esfera continúa rodando en dirección a r = 0.5 m, deteniéndose y cambiando la dirección de movimiento en t = 3.2 seg, aproximadamente. En cuanto a las presiones y los flujos másicos, se puede decir que, cuando el pistón avanza hacia su extremo, la presión 1 aumenta y la presión 2 disminuye, mientras que el flujo másico 1 es positivo y el flujo másico 2 es negativo. 

Cuando se aplica una tensión de entrada que va desde la tensión de equilibrio (5 voltios) hasta un valor cercano al máximo (8 voltios) y luego se disminuye a un valor inferior al anterior, pero mayor de 5 voltios (Figura 7), el pistón avanza al extremo máximo y luego retrocede, sucediendo lo mismo con la inclinación de la viga. La esfera inicia su movimiento desde el reposo hasta el extremo fijo de la viga y continúa su recorrido a pesar de los cambios. Esto se debe a que la inclinación de la viga permanece positiva durante un largo tiempo, con lo que la esfera se saldría de la viga. 

7. Conclusiones 

El modelo matemático dinámico del sistema de viga y esfera con accionamiento neumático está constituido por las ecuaciones (10), (14) y (15) (para la planta), las ecuaciones (31-32), (35) y (39) (para el cilindro) y las ecuaciones (42-46) (para la válvula). Se pueden hacer modificaciones al modelo presentado, mediante cambios en algunas de las consideraciones asumidas, convirtiéndolo así en un modelo aún más complejo, pero más exacto. El modelo desarrollado resulta ser de una gran utilidad para el análisis y simulación de otros sistemas que incluyan accionamientos neumáticos, porque el mismo ha sido obtenido tomando en cuenta la no linealidad existente en los elementos neumáticos, proveniente de la dinámica del aire comprimido y la fricción interna del cilindro. De esta manera, el modelo puede utilizarse con sistemas automáticos de control que incluyan actuadores neumáticos, permitiendo el ensayo de técnicas de control para sistemas no lineales. 

Como validación del modelo obtenido se simuló el sistema completo de viga y esfera empleando diferentes valores de tensión de entrada para la servoválvula. En los dos experimentos realizados, los resultados fueron totalmente satisfactorios con respecto al comportamiento esperado. 

8. Lista de términos 

J, R, m_esf : Momento de inercia para una esfera maciza, radio y masa de la esfera. 

r, g: Posición de la esfera y aceleración de gravedad, respectivamente. 

q, L: Ángulo de inclinación de la viga y longitud de la viga, respectivamente. 

d: Ángulo que forma un extremo del brazo de palanca con la horizontal. 

Ec, fr_e : Energía cinética de la partícula y fuerza de roce entre la esfera y la viga. 

: Velocidad y aceleración de la esfera 

: Velocidad y aceleración angular de la viga (OA) 

b: Ángulo formado entre la manivela y la horizontal 

y, b: Desplazamiento del punto B y longitud de la manivela, respectivamente. 

F_L, m_vi, m_ma: Fuerza de carga en el punto B, masa de la viga (OA) y de la manivela (AB). 

: Flujo másico, densidad y volumen. 

: Flujo másico de entrada y de salida de la cámara del cilindro. 

: Derivada con respecto al tiempo del volumen, de la densidad y del calor. 

h, h_ent, h_sal: Entalpía, entalpía de entrada y de salida de la cámara. 

: Energía total en la cámara y tasa de cambio del trabajo realizado. 

k, C_v: Relación de calor específico y calor específico a volumen constante. 

T, P, q_ent,q_sal: Temperatura, presión, calor que entra y calor que sale de la cámara. 

: Cambio de la energía interna y constante del gas ideal. 

: Derivada de la presión y coeficiente de característica térmica. 

V_i, V_oi : Volumen de la cámara i y volumen muerto de la cámara i. 

A_i, L_c: Área efectiva de la cámara i y carrera máxima del vástago. 

åF_y: Sumatoria de todas las fuerzas aplicada al émbolo en el eje y. 

: Masa, desplazamiento, velocidad, aceleración y diámetro del émbolo. 

F_f , F_Lp, i: Fuerza de fricción, fuerza de carga sobre el pistón y número de la cámara del cilindro (1,2). 

m_va,d_va: Masa y diámetro del vástago. 

F_C, F_sC, F_dC : Fuerza de fricción de Coulomb, estática de Coulomb y dinámica de Coulomb. 

F_v, j: Fuerza de fricción viscosa y coeficiente de fricción viscosa. 

P_u, P_d: Presión absoluta aguas arriba del orificio y aguas abajo del orificio, respectivamente. 

A_ef, P_crit: Área efectiva del flujo del orificio y relación de presión crítica (P_crit = 0.528 = b). 

k, Y_T: Constante isentrópica, función factor de flujo de aire. 

T_u: Temperatura absoluta aguas arriba. 

A_m, C, D, e, g, A_o: Constantes. 

u: Voltaje de mando de la válvula proporcional. 

u_s: Voltajes de la válvula proporcional correspondiente a A_ef = 0. 

u_o: Voltajes de la zona muerta de la válvula proporcional. 

Referencias Bibliográficas 

1. Dorf R., B.R. (2002). Ball and Beam System, Modern Control Systems, Addison Wesley.          [ Links ]

2. Carnegie M. (2003). Control tutorials for Matlab modeling the ball and beam experiment. Universidad de Michigan, http:/engin. umich.edu/group/ctm 425.         [ Links ]

3. Bobrow J., M.W. (1998). Modeling, Identification, and Control of a Pneumatically Actuated, Force Controllable Robot. IEEE Transactions on Robotics and Automation, Vol.15, No. 5, pp. 732-742.         [ Links ]

4. Gal, T. (2001). Measuring the air use of the servo-pneumatic drive. Periodica Polytechnica ser. mech. eng. Vol. 45, No. 1, pp. 11-16 (2001). Budapest University of Technology. Department of Precision Mechanics and Optics.H–1111 Budapest, Hungary.         [ Links ]

5. Moreno H., D.S. (2000). Control of a Pneumatic Servosystem Using Fuzzy Logic. Proc. of 1er FPNI-PhD Simposio. Hamburgo 2000, pp. 189-201. Departamento de Mecánica de Fluidos - UPC. Campus de Terrassa, Colón 7-11 08225 Terrassa. Barcelona, España.         [ Links ]

6. Perondi E., G.R. (2000). Control of a servopneumatic drive with friction compensation. robotics. Laboratorio de Ingeniería Mecánica. Universidad Federal de Santa Catarina - UFSC 88040-900-Florianópolis, SC Brasil.         [ Links ]

7. Wells, Dare (1988). Shaum´s outline of theory and problems of Lagragian Dynamics, Editorial McGraw-Hill, U.S.A.         [ Links ]

8. German A., G.S. (1998) Diseño de Mecanismos, Análisis y Síntesis. Editorial Prentice-Hall, 3^ra Edición, México.         [ Links ]

9. Ben-Dov, D., S. (1995). A Force-Controlled Pneumatic Actuator. IEEE Transactions on Robotics and Automation, Vol.11, No. 6, pp. 906-911.         [ Links ]

10. Liu S., B.J. (1988). An analysis of a pneumatic servo system and its application to a computer-controlled robot. Journal of Dynamic Systems, Measurement, and Control, Vol. 110, pp. 228-235.         [ Links ]

11. Richard E., S.S. (1996). Comparison between linear and nonlinear control of an electropneumatic servodrive. Journal of Dynamic Systems, Measurement and Control, Vol.118, pp. 245-118.         [ Links ]

12. Hullender D., Woods. R. (1985). Modeling of fluid control components. 1^ra Conferencia Control y Medida de Fluido, FLUCOME 1985. Tokio, Londres: Pergamon Press.         [ Links ]

13. Richer E., H.Y. (2000). A high performance pneumatic force actuator system. (Part1: Nonlinear Mathematical Model). Articulo de Septiembre 2000. Publicación del Periódico ASME. Journal of Dynamic Systems Measurement and Control, Vol.122, No. 3, pp. 416-425. Universidad Metodista del Sur, Escuela de Ingeniería y Ciencias Aplicadas a la Ingeniería Mecánica, Dallas, TX75275.         [ Links ]

14. Van Wylen G., S.R. (1982). Fundamentos de Termodinámica. Editorial Limusa. México.         [ Links ]

15. Rotty, Ralph (1968). Introducción a la Dinámica de los gases. Editorial R.A. Hermanos Sucesores S.A. Editores. 1ra Edición. México.         [ Links ]

16. Al-Ibrahim A. y Otis, D. (1992). Transient air temperature and pressure measurements during the charging and discharging processes of an actuating pneumatic cylinder. 45ta Conferencia Nacional Fluid Power.         [ Links ]

17. Armstrong, B. y Helouvry, H. (1994). Control of Machines with Friction. Publicaciones Académicas Kluwer, Londres.         [ Links ]

18. Canudas De Wit, C., O.H., A. K., L.P. (1995). A new model for control of systems with friction. IEEE Transactions on Automatic Control, Vol. 40, No. 5, pp. 419-425.         [ Links ]

19. Hatipoglu C., O.U. (1998). Robust control of systems involving non-smooth nonlinearities using modified sliding manifolds. Proc. Amer. Contr. Conf., Vol. June, pp. 2133-2137.         [ Links ]

20. Dupont P., D. P. (1995). Friction modeling and proportional-derivative compensation at very low velocities. Trans. ASME Vol. 117, No. 3, pp. 8-14.         [ Links ]

21. HAESSING D., F. B. (1990). On the modeling and simulation of friction. Proc. Amer. Contr. Conf. Vol. 2, pp. 1256-1261.         [ Links ]

22. Kang Y., K.K. (1997). Friction identification in a sight stabilisation system at low velocities. Mechanical Systems and Signal Processing, Vol. 11. No. 3, pp. 491-505.         [ Links ]

23. Fox R., M.A. (1994) Introduction to fluid mechanics. 4^ta Edición, Editorial John Willey & Sones, INC.         [ Links ]

24. Bashir M., Y.N., F. A., S.J., V.P., V.H. (2000). Modeling a pneumatic servo positioning system with friction. (Procedente de la Conferencia Americana de Control) Junio 2000. Universidad Católica Leuven. Departamento de Ingeniería Mecánica. Div. de Ing. de Producción, Diseño de Máquina y Automatización Chicago, Illinois, EE.UU.         [ Links ]