1 Simulacion de sistemas dinamicosPrincipios básicos de la integración numérica

2 Contenido El dominio de estabilidad numéricaCálculo del dominio de estabilidad numérica Iteración de punto fijo Iteración de Newton Conclusiones

3 El dominio de estabilidad numéricaBasado en el libro de Cellier, F.E. and E. Kofman (2006), Continuous System Simulation, Springer-Verlag, New York El dominio de estabilidad numérica

4 Dominio de estabilidad analíticaUn sistema lineal invariante en el tiempo autónomo puede ser representado usando el modelo La solución analítica es La solución es analíticamente estable si Región de estabilidad analítica en el plano λ

5 Dominio de estabilidad numérica para Euler directo

6 Integración numérica usando Euler directoUsando el algoritmo de Euler directo: El sistema original continuo se ha convertido en un sistema de tiempo discreto “equivalente”

7 Valores propios del sistema de tiempo discretoEl sistema continuo autonomo original es El sistema de tiempo discreto autonomo “equivalente” es Demostrar! Si λ es un valor propio de A, entonces λ·h es un valor propio de F

8 Sistema de tiempo discreto equivalenteEn general, el sistema lineal de tiempo discreto “equivalente” puede expresarse usando Región de estabilidad numerica en el plano λ·h La region de estabilidad de un sistema de tiempo discreto es un circulo de radio unitario centrado en el origen

9 Simulación con el algoritmo de Euler directo1 h = 1 Simulacion para valores de a = - 0.1, -1, -2 y -3 con paso de integracion fijo h = 1

10 Simulación con el algoritmo de Euler directoSimulación con el algoritmo de Euler directo del sistema: A fin de obtener un resultado una exactitud del 0.1 %, el tamaño del paso debe ser h = ¡ pasos de integración para simular 10 segundos!

11 El mayor paso de integración posible para Euler directoDado un sistema lineal de segundo orden con dos valores propios complejos λ1 y λ2 Región de estabilidad numerica en el plano λ·h

12 El mayor paso de integración posible para Euler directoDado un sistema lineal con valores propios complejos λi : Ejercicio Demuestre que el tamaño del paso hmarg para el cual el algoritmo directo de Euler dara marginalmente estable es: Ayuda: usar el teorema de Tales

13 El mayor paso de integración posible para Euler directoLa limitacion en el valor de h es particularmente importante cuando los autovalores se encuentran cerca del eje imaginario. d se hace arbitrariamente pequeño Cuando los autovalores estan sobre el eje imaginario no existe ningun paso de integracion que permita obtener una solucion puramente oscilante

14 Dominio de estabilidad numérica para Euler inverso

15 Integración numérica usando Euler inversoUsando el algoritmo de Euler inverso: El sistema original continuo se ha convertido en un sistema de tiempo discreto “equivalente”

16 Valores propios del sistema de tiempo discretoEl sistema continuo autonomo original es El sistema de tiempo discreto autonomo “equivalente” es ¿Si λ es un valor propio de A, cual es el valor propio de F ?

17 Región de estabilidad del Euler inversoDominio de estabilidad numérica del algoritmo de Euler inverso El tamaño del paso de integración es dictado sólo por requerimientos de exactitud, No por estabilidad numérica Región de estabilidad numerica en el plano λ·h

18 Simulación con el algoritmo de Euler inversoSimulación con el algoritmo de Euler inverso del sistema: 1 h = 1 Este tipo de algoritmo es más apropiado que el Euler directo para resolver problemas con valores propios lejanos sobre el real negativo en el plano λ

19 Euler inverso y los sistemas StiffRegión de estabilidad numerica en el plano λ·h Apropiado para sistemas stiff, es decir, sistemas con autovalores cuyas partes reales estan desparramadas a lo largo del eje real negativo

20 Simulación con el algoritmo de Euler inversoSimulación con el algoritmo de Euler inverso del sistema: Con a = +3.0 el sistema original es inestable La simulación sugiere un sistema perfectamente estable Lección a aprender: Puede ser buena idea simular el sistema dos veces, una vez con algoritmo que exhiba un dominio de estabilidad comparable con el algoritmo de Euler directo, y otra con un algoritmo que se comporte como el algoritmo de Euler inverso

21 Implementación en matalbCálculo del dominio de estabilidad numérica Implementación en matalb

22 Definición de la Matrix del sistemaConsiderando un sistema lineal de segundo orden con dos valores propios complejos λ1 y λ2 con la matriz λ1 y λ2 estan localizados sobre el círculo unitario lambda = cos(α) + j*sen(α) cos(α) – j*sen(α)

23 Definición de la Matrix del sistemaConsiderando un sistema lineal de segundo orden con dos valores propios complejos λ1 y λ2 con la matriz Implementación en MATLAB

24 Cálculo de la matriz F Euler directo Euler inverso

25 Cálculo del máximo valor de hCon h = hmax , los valores propios de F se encuentran en el círculo unitario Esta función trabaja para algoritmos con dominios de estabilidad similares al algoritmo Euler

26 Cálculo del dominio de estabilidad numéricaEjercicio: Hacer un grafico del valor de hmax en función de α, para el algoritmo Euler directo en coordenadas polares Utilice las funciones aa.m, ff.m, hh.m , y stabdom.m

27 Iteración de punto fijo

28 Euler inverso para sistemas no linealesEn el caso de los sistemas lineales, en la simulación del algoritmo de Euler inverso la matriz F puede calcularse explícitamente En el caso no lineal esto no puede hacerse En el caso no lineal es necesario resolver un conjunto de ecuaciones algebraicas no lineales implícitas

29 Iteración de punto fijoUna posible aproximación: * Iniciar con una predicción * Continuar con la iteración de varias correcciones predictor: Iteration i corrector:

30 Iteración de punto fijoUna posible aproximación:

31 Iteración de punto fijo: caso linealCuando se aplica la reiteración de punto fijo para un sistema lineal, se tiene Después de un número infinito de interacciones entonces restando las dos ecuaciones La misma matriz F que en el caso del algoritmo de Euler inverso

32 Iteración de punto fijo: caso linealDominio de estabilidad numérica de la técnica predictor-corrector obtenido con el algoritmo que genera los dominios de estabilidad

33 Iteración de punto fijo: caso linealLa aproximación no funciona porque la serie infinita sólo converge si todos los valores propios de A·h se encuentran dentro del círculo unitario Si no es el caso, la resta de las ecuaciones es inválida Dentro del círculo unitario, el dominio de estabilidad numérica del método predictor-corrector es el mismo que para algoritmo de Euler inverso Fuera del círculo unitario, el método es inestable

34 Conclusiones

35 Conclusiones En el análisis de un solver numérico, la estabilidad numérica del algoritmo debe ser siempre tomada en consideración La estabilidad numérica de la mayoría de los solvers puede ser representada por un dominio de estabilidad numérica, dibujado en el plano complejo λ·h

36 Conclusiones En el análisis de un solver numérico, la estabilidad numérica del algoritmo debe ser siempre tomada en consideración La estabilidad numérica de la mayoría de los solvers puede ser representada por un dominio de estabilidad numérica, dibujado en el plano complejo λ·h La estabilidad numérica de los solvers es usualmente analizada solamente para sistemas lineales autónomos invariantes en el tiempo

37 Conclusiones En el análisis de un solver numérico, la estabilidad numérica del algoritmo debe ser siempre tomada en consideración

38 Conclusiones En el análisis de un solver numérico, es también importante considerar la exactitud de la aproximación del algoritmo La exactitud numérica de un solver p está sujeta a diferentes tipos de error, tales como el error por truncado, el error por redondeo, y el error por acumulación

39 Conclusiones En el análisis de un solver numérico, es también importante considerar la exactitud de la aproximación del algoritmo La exactitud numérica de un solver p está sujeta a diferentes tipos de error, tales como el error por truncado, el error por redondeo, y el error por acumulación El más importante de estos errores es el error por truncado, que está caracterizado por el orden de la exactitud de la aproximación del sol ver

40 Fuentes Cellier, F.E. and E. Kofman (2006), Continuous System Simulation, Springer-Verlag, New York

41 FIN