Desplazamiento de una estructura en oscilación amortiguada El desplazamiento de una estructura en una oscilación amortiguada viene dada por la ecuación: y
=
9e − kt cos ω t
Si k = 0.7;
=
4,
determinar por métodos numéricos numéricos y programas en MATLAB MATLAB el tiempo para que el desplazamiento disminuya a 3.5 con una una precisión de cuatro dígitos decimales significativos.
SOLUCIÓN: Se trata de hallar la raíz de la ecuación: 9e − kt cos ω t − 3.5 = 0
Vamos a utilizar dos métodos distintos: el de Newton-Raphson y el de la secante. a) Newton-Raphson f (t ) = 9e − kt cos ω t − 3.5
cuya derivada es: f ′(t ) = −9ke − kt cos ω t − 9ω e − kt senω t = −9e − kt ( k cos ω t + ω senω t ) La fórmula recurrente de Newton-Raphson es: t i +1
=
t i
−
f (t i ) f ′(t i )
Escribimos la siguiente función (genérica que nos va a servir para aplicar el método de Newton-Raphson a otros problemas diferentes) ( Programa 1 ):
function[x]=NR(x,p) ant=-inf; while abs(x-ant)>p ant=x; x=x-f(x)/fprima(x); end Observamos que simplemente le pasamos un valor inicial, x, y le ordenamos que aplique iterativamente la fórmula de Newton-Raphson hasta alcanzar la precisión que p, y nos devuelve el valor buscado . Hay que definir dos funciones, f le hemos fijado, p y fprima que son, respectivamente, la función de la que buscamos los ceros y su derivada , correspondientes a cada problema. La variable ant , en donde vamos a guardar el valor calculado en la última iteración para compararlo con el nuevo, la
inicializamos a –inf (representación en MATLAB del menor número negativo que puede almacenar el ordenador) para asegurarnos que siempre entramos la primera vez en el bucle while-end Para este problema en concreto, las funciones f y fprima son, en este caso:
function[y]=f(x) y=9*exp(-0.7*x)*cos(4*x)-3.5; y
function[y]=fprima(x) y=-6.3*exp(-0.7*x)*cos(4*x)-36*exp(-0.7*x)*sin(4*x); Tomamos como valor inicial t 0
=
0.5 y como p =0.00005 (nos piden precisión de
cuatro dígitos decimales significativos) En el Command Windows hacemos la llamada:
>> t=NR(0.5,0.00005) y obtenemos el resultado:
t= 0.2704 que es el resultado buscado. Podemos comprobar también en el Command Windows:
>> f(0.2704) ans = -6.4316e-005 >> que nos confirma que 0.2704 es un valor aproximado, con la precisión exigida, al cero exacto buscado de la función. Recordamos que el método de Newton-Raphson es un método abierto en el que, en principio , no está asegurada la convergencia. Sin embargo en este caso sí podíamos haber asegurado que iba a converger por el valor inicial escogido (0.5) y cumplirse la relación de convergencia: f (t ) f ′′(t )
<
( f ′(t )) 2 en un punto próximo al valor buscado.
Podríamos estar también interesados en el número de iteraciones que han sido necesarias para conseguir la precisión deseada. Utilizaríamos, entonces, una variante de nuestra función de Newton-Raphson:
function[x,nit]=NR(x,p) ant=-inf; nit=0; while abs(x-ant)>p ant=x; x=x-f(x)/fprima(x); nit=nit+1; end Y así , al llamar de nuevo a la función en el Command Windows:
>> [t,nit]=NR(0.5,0.00005) obtenemos:
t= 0.2704 nit = 4 Hemos necesitado únicamente 4 iteraciones (¡El método de Newton-Raphson converge realmente rápido si escogemos adecuadamente el valor inicial!). Podríamos también estar interesados en limitar el número de iteraciones para evitar entrar en un proceso no convergente o excesivamente lento en converger. Una variante de nuestra función NR sería:
function[x,nit]=NR(x,p,nmaxit) ant=-inf; nit=0; while abs(x-ant)>p ant=x; x=x-f(x)/fprima(x); if nit>nmaxit disp(‘Se han superado el numero maximo de iteraciones permitidas’); return end nit=nit+1; end
b) Método de la Secante La fórmula recurrente e el método de la secante es: t n +1
=
t n
−
f (t n )(t n
− t n −1 )
f (t n ) − f (t n −1 )
que viene a ser la fórmula de Newton-Raphson sustituyendo la primera derivada por la diferencia finita dividida de primer orden. Escribimos una función MATLAB genérica para el método de la secante ( Programa 2):
function[x3]=secante(x1,x2,p) while abs(x2-x1)>p x3=x2-f(x2)*(x2-x1)/(f(x2)-f(x1)); x1=x2; x2=x3; end En caso de que nos interesara contar el número de iteraciones o limitar el número máximo de ellas a realizar emplearíamos variantes del tipo de las funciones de Newton-Raphson que desarrollamos en el apartado anterior. La función f es la misma, obviamente, que utilizamos para el método de NewtonRaphson. Una de las ventajas (quizás la única) del método de la secante frente al de Newton-Raphson es que no utiliza la derivada, por tanto en este caso no necesitamos la función fprima. Recordamos que ambos métodos son del tipo abierto por lo que no está asegurada la convergencia y es necesario efectuar convenientemente la elección de los valores iniciales, lo que se puede hacer por tanteos o aplicando criterios físicos del problema. Tomamos como valores iniciales: t 0
=
0.3; t 1
=
0.7 y la misma función f que
utilizamos en el método anterior:
>> t=secante(0.3,0.7,0.00005) t= 0.2704 obteniendo, el mismo resultado que con el método de Newton-Raphson.