Efraín Domínguez1, Felipe Ardila2 y Santiago Bustamante3
1 Ingeniero Hidrólogo. M.Sc., en ecología hidrometeorologia. Ph.D., en ciencias técnicas. Universidad Javeriana, Bogotá, Colombia. e.dominguez@javeriana.edu.co
2 Ingeniero de Ssitemas, Universidad Distrital. Estudiante de M.Sc., en Hidrosistemas, Universidad Javeriana, Bogotá, Colombia. Instituto Geofísico. felipeardilac@gmail.com.
3 Pregrado, Universidad Javeriana, Bogotá, Colombia. santiagobus@gmail.com
RESUMEN
]]> En este artículo se presenta una herramienta de código abierto para la simulación de sistemas dinámicos, representados a través de ecuaciones diferenciales de distinto orden o de sistemas de ecuaciones diferenciales ordinarias (EDO). La idea principal que incentivó el desarrollo de esta herramienta consiste en que muchos procesos físicos, biológicos, ecológicos, económicos, químicos y sociales, y un sinnúmero de problemas ingenieriles, pueden ser expresados de esta forma. Adicionalmente, la solución de ese tipo de problemas exige experticia en métodos numéricos y programación. Este tipo de experticia no es muy común para los expertos de las áreas mencionadas; por ello, se requiere de una herramienta que permita sortear esta falta de experticia e incrementar la productividad en trabajos relacionados con tópicos de modelación. System-Solver es una aplicación que facilita la formulación de problemas de valor inicial en EDO y su solución numérica, y mediante un programa en Visual Basic la automatiza. El programa obtenido se puede compartir fácilmente con otros investigadores, facilitando la reproducción del ejercicio de modelación incluso en diferentes sistemas operativos. La presentación del software está acompañada con ejercicios de diversas áreas temáticas, incluyendo un ejemplo de modelación estocástica.Palabras claves: ODE Solver, diferenciación por computadoras, modelación, sistema dinámicos.
ABSTRACT
The following paper presents a freeware modelling tool simulating dynamic systems that can be represented by either an ordinary differential equation (ODE) or a set of differential equations of different orders. The main idea leading to this software development is related to the fact that many physical, biological, ecological, economical, chemical, social and engineering problems can be expressed in this way. Furthermore, the solution to these problems requires some expertise in numerical methods and programming. Such knowledge is uncommon in some of the experts in such scientific domains. A tool to fill in this knowledge gap, increase productivity within modelling-related research and support the teaching of mathematical modelling topics is thus needed. This paper introduces System Solver, a computer application that facilitates the formulation of initial value problems for ODE systems, numerically solves these problems and provides a user with not only the solution but also debugged Visual Basic source code for the application. The obtained code can be easily shared among researchers, which facilitates the replication of numerical experiments even across different operating systems. This software's introduction is accompanied by examples from different domains, including one example from stochastic modelling.
Keywords: ODE solver, computational differentiation, modelling, dynamical system.
Recibido: junio 23 de 2009
Aceptado: noviembre 15 de 2010
Introducción
Muchos problemas del mundo real pueden ser estudiados utilizando modelación matemática. El operador diferencial es una herramienta matemática ampliamente utilizada para describir procesos físicos, biológicos, ecológicos, etcétera. Usualmente el problema en estudio se describe a través de una EDO de orden r o a través de un sistema de EDO. Y con frecuencia la solución de esa ecuación, o sistema de ecuaciones, se establece con la aplicación de métodos numéricos. Los interesados en solucionar numéricamente un modelo matemático deben estar familiarizados con lenguajes de programación y algoritmos numéricos. Sin esta experiencia el proceso de investigación del modelo matemático se retrasa hasta que haya disponibilidad de un especialista en las áreas mencionadas.
La integración automática de sistemas de EDO siempre ha sido un tema de interés para la comunidad científica. El desarrollo de esta temática comenzó en 1970 (Moore, 1979). Tolsma y Barton (1998) sugieren cinco aproximaciones para la derivación computacional: 1) codificación manual, 2) aproximación por diferencias finitas, 3) derivación simbólica, 4) aplicación de la notación polaca inversa (RPN, por su sigla en inglés) (Bardsley y Prasad, 1997), 5) diferenciación automática. Iri (1991) postuló los prerrequisitos de un buen método de integración: debe ser rápido, libre de error de truncamiento y preferiblemente aplicado en forma automática. Adicionalmente, se requiere añadir otro prerrequisito: el método a utilizar debe permitir su generalización a sistemas de cualquier complejidad. De los métodos enunciados, el quinto y el sexto cumplen con casi todos los requerimientos planteados. El error de truncamiento se mantiene presente en el segundo método. El quinto método requiere alta intervención del usuario, quien debe programar un procedimiento inicial para iniciar la integración de un sistema de EDO concreto. Adicionalmente, la diferenciación automática puede requerir algún esfuerzo de adaptación de código para resolver problemas planteados en forma discreta (tabulada). Los métodos mencionados se aplican en la solución de problemas de modelación para variadas áreas científicas, y usualmente todas estas técnicas son estudiadas en cursos estandarizados de modelamiento matemático. En el mercado hay una cantidad considerable de programas especializados en la solución de sistemas de EDO. La mayoría de éstos son de licencia privativa y se pueden acceder sólo comprándolos y generalmente no son transparentes para el usuario en el sentido de los algoritmos utilizados, ya que no proveen libre acceso al código fuente de los algoritmos preprogramados en ellos.
]]> Este artículo ofrece los fundamentos teóricos, métodos, algoritmos, y la implementación de un programador automático de soluciones numéricas para sistemas de EDO. Dicho sistema de modelación implementa los métodos explícitos de Euler y Runge–Kutta y exporta la solución numérica para visualización en Scilab o Matlab. Tal aplicación representa la primera versión de un sistema de modelamiento que implementa la generación automática de programas que resuelven ecuaciones deferenciales para el estudio de sistemas complejos.System–Solver es una herramienta de distribución gratuita, recomendada para personas que quieren resolver EDO en forma numérica pero que no son programadores. Su aplicación es más fácil de usar que programas como Matlab porque no requiere conocimientos de programación, usa pocos comandos y está diseñada específicamente para resolver sistemas de EDO. System-Solver es un software pedagógico de alta calidad para la introducción al mundo de la modelación. La ventaja de usarlo en lugar de los bien conocidos Matlab, Scilab u Octave consiste en su diseño, orientado a la enseñanza. Este diseño permite desarrollar habilidades programadoras y de modelación matemática a la vez. Otra ventaja potencial reside en la posibilidad de escoger el lenguaje de programación que el usuario desea aprender y utilizar para compartir los programas de los sistemas modelados. Dicha habilidad será habilitada en la siguiente versión del paquete.
System-Solver se creó en Java utilizando una librería especializada llamada Qt Jambi, una tecnología con soporte oficial para programadores que desean crear interfaces gráficas muy ricas aplicando código en lenguaje Java.
Comparado con herramientas de modelamiento conceptual como Stella, Powersim, Simile o Vensim, System-Solver explica en forma más clara los conceptos básicos en modelación. De hecho, un usuario puede ver las EDO que explican el sistema que se modela, y también tiene el control total del código numérico desarrollado para el sistema en estudio. El diseño de problemas n-dimensionales en modelamiento conceptual puede ser problemático a través de la representación visual del diagrama para el problema propuesto. La estructura de System–Solver fue diseñada para permitir la configuración amigable de sistemas n-dimensionales, de modo que se facilite la modelación de sistemas distribuidos.
Estructura de System-Solver
Se propone la siguiente notación:
El sistema (1) sugiere un algoritmo cíclico que aplique un método numérico (Euler o Runge-Kutta en nuestro caso) n-veces para encontrar el vector y (Боглаев, 1990). System-Solver implementa un esquema de diferencias finitas explícitas de Euler como el siguiente:
]]>
También es posible seleccionar un esquema explícito de Runge-Kutta de cuarto orden de la forma:
En ambas ecuaciones [(2) y (3)], nij representa un error inducido por el truncamiento de las derivadas de alto orden en la serie de Taylor durante la deducción de las diferencias finitas para la EDO j. De hecho, ni→0 si Δx→0. Los esquemas (1) y (2) son explícitos, y su desempeño depende fuertemente del tamaño del paso Δx. El índice i fija el nodo x que está siendo evaluado.
Para un programador experimentado es fácil reconocer un comportamiento genérico en las ecuaciones (1), (2) y (3) para la solución de un sistema de EDO. También se debe tener una rutina para evaluar yji+1 (procedimiento 3 en la Figura 1), y la parte derecha de la ecuación (1) debe resolverse con un procedimiento externo (2 en la Figura 1). La configuración de condiciones iniciales, del intervalo de solución, los coeficientes de la EDO y el paso de integración (Δx deben ser especificados. Esto puede realizarse al comienzo de los cálculos (procedimiento 1 en la Figura 1). Para resolver la función de la parte derecha se implementó un procedimiento parcelador de ecuaciones que permite identificar las variables requeridas y los parámetros de la EDO. Estos parámetros también pueden ser función del tiempo o de las variables de estado del sistema (1).
En la Figura 2 se muestra la estructura básica de System-Solver. Esta estructura propone una interface para configurar la cantidad de ecuaciones en el sistema de EDO y campos para las expresiones que definen las funciones de la parte derecha de las EDO. Un método llamado "Solve" activa la traslación de la notación algebraica hacia la notación de registro polaco inverso (RPN). De la cadena de caracteres RPN se deducen y clasifican los parámetros, las variables independientes y variables de estado, que son insertados en una plantilla de código preparado previamente. La primera plantilla (RPEvaluator) evalúa las funciones de la parte derecha.
]]>
Una vez se rellena la plantilla RPEvaluator con el número correcto de variables y los parámetros requeridos por la EDO, se completa otra plantilla con el intervalo de discretización, vectores de condiciones iniciales y el número de pasos en el tiempo. EDOSystem es la implementación del diagrama de flujo presentado en la Figura 1. Este procedimiento es el subprograma principal que invoca una instancia de RPEvaluator, el cual retorna un vector solución para las funciones de la parte derecha requerido por los subprocedimientos Runge-Kutta o Euler y se crea una nueva evaluación de yji+1 para cada nodo de tiempo x.
Como ejemplo, se examinan algunos sistemas dinámicos y se compara el código generado por System-Solver contra el código escrito a mano por un programador. Considérese un sistema predador-presa del siguiente tipo:
Donde r y q son los coeficientes de crecimiento, mientras que ay b son las tasas de beneficio por ataque del predador y la tasa de perjuicio de la presa, respectivamente.
Dadas las siguientes condiciones iniciales,
Un modelador experimentado encontraría la solución numérica de (4) utilizando el método de Euler.
]]> Siguiendo la estructura propuesta en la Figura 1, el procedimiento Sistema_ODE_Euler es equivalente al procedimiento 1 en la Figura 1; el procedimiento de la subderecha (…) equivale al procedimiento 2, y el subprograma Euler (…) evalúa xt+1 y yt+1 (procedimiento 3 en la Figura 1).Para resolver el mismo problema en System-Solver se debe configurar el sistema EDO añadiendo cada ecuación usando la notación algebraica comúnmente aplicada en hojas de cálculo para definir fórmulas. El orden del sistema EDO está limitado sólo por las capacidades computacionales del ordenador. La Figura 2 describe la ventana principal de System-Solver. Esta aplicación ofrece una interface amigable incluso para investigadores con mínima experiencia en modelación, métodos numéricos y programación. La configuración de un modelo requiere no más de un par de minutos. El investigador debe invertir su tiempo encontrando las ecuaciones necesarias para su sistema, definiendo los valores de los parámetros en estas ecuaciones e interpretando resultados, en lugar de gastarlo definiendo, codificando y depurando el código de un algoritmo numérico.
Para el caso del modelo predador-presa (ecuación 4), System-Solver genera el código en Visual Basic, que se puede consultar en http://www.mathmodelling.org/home/system-solver-1.
El código Visual Basic generado por System-Solver puede ser copiado directamente en un módulo Visual Basic de Microsoft Excel para ser utilizado como un macro guardado directamente en un libro de Excel. El código generado también puede ser preparado para el compilador de Microsoft Visual Basic (versión Stand Alone). Una versión gratuita de este compilador se puede obtener desde el sitio web http://www.microsoft.com. Adicionalmente, el proyecto de código abierto "Mono" financiado por Novell incluye un compilador de Visual Basic, con entorno de desarrollo integrado, que puede ser ejecutado en diferentes sistemas operativos, incluyendo: Linux, Mac OS X, Sun Solaris, BSD - OpenBSD, FreeBSD, NetBSD y Microsoft Windows. (Más información sobre el proyecto "Mono", en: http://www.mono-project.com).
Si el código generado por System-Solver es preparado para la versión Stand Alone del compilador de Microsoft Visual Basic, se deberá editar el código para implementar procedimientos de pre y posprocesamiento. Independientemente de la plataforma para la cual fue generado el código, las herramientas Scilab y Matlab ofrecen una opción de posprocesamiento. System-Solver puede guardar los resultados del experimento numérico en formato ASCII y generar automáticamente el script que realizará la salida gráfica en los entornos de Scilab o Matlab.
Muchos problemas de modelación incluyen sistemas de tres o más EDO para describir los sistemas reales que se estudian. Estos sistemas multidimensionales requieren gráficas 3D para el posprocesamiento. En este caso System-Solver ofrece una opción para exportar resultados, la cual los prepara para ser procesados con Scilab o Matlab. La solución discreta se guarda en un archivo de texto y luego se importa en Scilab o Matlab con un script generado por System-Solver. Para examinar esta opción se resuelve el siguiente sistema predador-presa 3-dimensional (Sáez et al., 2007):
Con las siguientes condiciones iniciales:
]]>
y utilizando los siguientes parámetros: a=1.5, b=5.0, c=8.0,q=1.0, f=0.16, g=0.1, w=1.1 , l=2.0, m=0.16, Δt=0.1, tmin=0,tmax=500. Desarrollando esta configuración en System-Solver se obtiene el código Visual Basic presentado en la Figura 2.
La función "export" de System-Solver guarda los resultados y genera un script de Scilab en una ubicación seleccionada por el usuario. El script generado por System-Solver se presenta en la Figura 3.
Cargando y ejecutando este script en el entorno de Scilab, se produce la gráfica presentada en la Figura 4.
Conclusiones
Los ejemplos presentados aquí muestran que el primer prototipo de System-Solver es capaz de manejar problemas EDO de diferente complejidad. Esta versión puede generar programas de soluciones numéricas eficientes y tiene herramientas para exportar los resultados de los experimentos numéricos a Scilab o Matlab, habilitando posprocesamiento 3D de los resultados. Esa aplicación ofrece una herramienta útil para el modelamiento de sistemas representados por conjuntos de EDO, acelerando la investigación relacionada a sistemas de EDO. System-Solver también puede acompañar los esfuerzos de enseñanza en curso de modelación matemática, programación y métodos numéricos. Usar System-Solver cierra el ciclo de definición del modelo matemático, desarrollo del algoritmo numérico y codificación del programa para la simulación de un sistema dinámico (Samarsky y Mikhailov, 1997).
]]> La facilidad de System-Solver para aplicar el método de Euler, o el de Runge-Kutta de cuarto orden, permite la exploración de sus propiedades para el sistema de EDO investigado.SystemSolver se desarrolla como herramienta multiplataforma. Se proyecta que las futuras versiones generen código para diversos compiladores como: C#, C++, Java, Object Oriented Pascal, Component Pascal y Fortran. En este momento el código generado se escribe para compiladores de Visual Basic (como el incluido en las aplicaciones de MS Office). Un compilador de Visual Basic está disponible en el proyecto "Mono", por ello el código Visual Basic generado por System-Solver también puede ser ejecutado en Linux y otros sistemas operativos. Una versión beta de System-Solver puede descargarse gratuitamente desde: http://www.mathmodelling.org/home/system-solver-1.
Se proyecta que distintos algoritmos numéricos de solución de sistemas EDO, incluidos algoritmos para ecuaciones rígidas, serán implementados en System-Solver en el futuro cercano. Algoritmos para la solución de ecuaciones diferenciales estocásticas están siendo asimilados. Las sugerencias con requerimientos de los usuarios de System-Solver también son bienvenidas.
Bibliografía
Bardsley, W.G., Prasad, N., Using ASCII text files in post-fix notation (reverse Polish) to define mathematical models and systems of differential equations for simulation and non-linear regression., Computers & Chemistry, 21(2), 1997, pp. 71-82. [ Links ]
Iri, M., History of automatic differentiation and rounding error estimation., Automatic differentiation of algorithms: Theory, implementation and application. SIAM, Philadelphia, PA., 1991. [ Links ]
Moore, R.E., Methods and applications of interval analysis., SIAM Publications, Philadelphia, 1979. [ Links ]
Saez, E., Stange, E., Szanto, I., Chaotic Dynamics and Coexistence in an Interaction Model Between Three Species., Universidad Técnica Federico Santa María - Departamento de Matemáticas, 2007, pp. 5. [ Links ]
Samarsky, A., Mikhailov, A., P., Matematicheskoie modelirovanie: Idei, metodi, primieri., Nauka, Moscow, 1997, 316 pp. [ Links ]
Tolsma, J.E., Barton, P.I., On computational differentiation., Computers & Chemical Engineering, 22(4-5), 1998, pp. 475-490. [ Links ] ]]>