Resolver Ecuaciones Diferenciales en Matlab

En el siguiente vídeo se explica de manera breve cómo utilizar el comando "dsolve" para resolver ecuaciones diferenciales en Matlab.

Criterio de estabilidad de Routh-Hurwitz en MATLAB

Enseguida adjunto un "script" mediante el cual podemos analizar la estabilidad de un sistema con el método de Routh-Hurwitz, para ello sólo es necesario correr el "script" y posteriormente introducir el vector de coeficientes del polinomio característico que nos interese. El script fue desarrollado por Farzad Sagharchi y está disponible en el File Exchange de MathWorks.

%The Routh-Hurwitz stability criterion is a necessary (and frequently
%sufficient) method to establish the stability of a single-input,
%single-output (SISO), linear time invariant (LTI) control system. 
%More generally, given a polynomial, some calculations using only the 
%coefficients of that polynomial can lead us to the conclusion that it
%is not stable.
%in this program you must give your system coefficents and the
%Routh-Hurwitz table would be shown
%
%   Farzad Sagharchi ,Iran
%   2007/11/12

clc
clear
r=input('input vector of your system coefficents: ');
m=length(r);
n=round(m/2);
q=1;
k=0;
for p = 1:length(r) 
    if rem(p,2)==0
        c_even(k)=r(p); 
    else
        c_odd(q)=r(p); 

        k=k+1;
        q=q+1;
    end
end
a=zeros(m,n); 

if m/2 ~= round(m/2)
    c_even(n)=0;
end
a(1,:)=c_odd;
a(2,:)=c_even;
if a(2,1)==0
    a(2,1)=0.01;
end
for i=3:m
    for j=1:n-1
        x=a(i-1,1);
        if x==0
            x=0.01;
        end

        a(i,j)=((a(i-1,1)*a(i-2,j+1))-(a(i-2,1)*a(i-1,j+1)))/x;

    end
    if a(i,:)==0
        order=(m-i+1);
        c=0;
        d=1;
        for j=1:n-1
            a(i,j)=(order-c)*(a(i-1,d));
            d=d+1;
            c=c+2;
        end
    end
    if a(i,1)==0
        a(i,1)=0.01;
    end
end
Right_poles=0;
for i=1:m-1
    if sign(a(i,1))*sign(a(i+1,1))==-1
        Right_poles=Right_poles+1;
    end
end
fprintf('\n Routh-Hurwitz Table:\n')
a
fprintf('\n Number Of Right Poles =%2.0f\n',Right_poles)

reply = input('Do You Need Roots of System? Y/N ', 's');
if reply=='y'||reply=='Y'
    ROOTS=roots(r);
    fprintf('\n Given Polynomials Coefficents Roots :\n')
    ROOTS
else
end

Comunicación Serial MATLAB - Arduino

El objetivo de esta entrada es proporcionar una guía para establecer comunicación mediante el puerto serial entre Matlab y la placa Arduino, graficando los valores que habrán de leerse y además guardando en un archivo de datos los valores recogidos u obtenidos durante la lectura, los cuales podrían servir para un posterior análisis. Para ello sólo necesitaremos obviamente disponer de una placa Arduino y el software correspondiente, además de Matlab.

Primeramente cargaremos en nuestra placa Arduino el código que se muestra enseguida, el cual como puede apreciarse es muy sencillo, lo único que se hace es leer un valor del pin analógico A0 del arduino e imprimirlo con la instrucción "Serial.println". El valor que estaremos leyendo en A0 puede provenir de cualquier sensor analógico que nosotros dispongamos o bien el que nos interese para un fin determinado, e inclusive podemos hacer uso de un potenciómetro que simplemente cambiará su valor a capricho nuestro.

/*
Comunicación Serial Matlab Arduino
Sketch Arduino
@Jorge De Los Santos
*/

int val; 

void setup(){
Serial.begin(9600);
}

void loop(){
val=analogRead(A0);
Serial.println(val);
delay(100);
}

Nosotros haremos uso de un potenciómetro para variar la lectura que habremos de tomar. Recordemos que el valor que obtendremos del pin analógico será un número entero entre 0 y 1023. La conexión del potenciómetro a la placa Arduino se asume que es algo de trámite para vosotros.

Una vez que hemos cargado el código anterior en nuestra placa, cerramos el entorno de desarrollo de Arduino para evitar conflictos con el puerto serial, puesto que si lo mantenemos abierto, Matlab nos enviará un mensaje de error cuando intentemos ejecutar la comunicación serial.

Considerando lo anterior, ahora abriremos Matlab y crearemos un nuevo "script" al cual le asignaremos un nombre a conveniencia propia, y en este escribiremos las siguientes instrucciones:

clear all;clc;

delete(instrfind({'Port'},{'COM11'}));
pserial=serial('COM11','BaudRate',9600);
fopen(pserial);

figure('Name','Gráfica de valores obtenidos')
title('LECTURA ANALOGICA CON ARDUINO');
xlabel('Muestra');
ylabel('Voltaje de Salida');
val=zeros(1,1000);

for i=1:1000
    ylim([0 5.1]); 
    xlim([i-100 i+10]);
    lectura=fscanf(pserial,'%d');
    val(i)=lectura*(5/1023);
    hold on
    plot(i,val(i),'x');
    drawnow
end

dlmwrite('Valores_Obtenidos.dat', val, 'delimiter', '\n', 'precision', '%.2f')

fclose(pserial); 
delete(pserial);
clear all;

Ahora veamos una breve descripción de cada una de las instrucciones dadas en el "script" de Matlab que hemos creado.

Comenzamos borrando variables y limpiando la pantalla del "command window" para evitar cualquier tipo de incoveniente, para ello usamos "clear all" y "clc".

Luego proseguimos a buscar si existe algún objeto creado en Matlab asociado al puerto serial que usaremos y en caso de haberlo lo borramos para evitar que Matlab aborte la operación, para ello usaremos la instrucción de la primera línea del código mostrado enseguida. Con la segunda línea inicializamos el puerto serial al cual le hemos llamado "pserial", indicamos el puerto serial que usaremos en este caso "COM11" (este valor variará dependiendo de cada usuario) y definimos los bps a los que se ejecutará la conexión o comunicación serial, en este caso 9600 bps, tal como se definió en el sketch del Arduino. Finalmente con "fopen" abriremos el puerto serial que hemos inicializado.

delete(instrfind({'Port'},{'COM11'}));
pserial=serial('COM11','BaudRate',9600);
fopen(pserial);

Posteriormente definimos una serie de paramétros estáticos que contendrá el objeto "figure" en el cual graficaremos los datos leídos, tales como título y etíquetas en los ejes, las siguientes cuatro lineas muestran la parte de código implicada.

figure('Name','Gráfica de valores obtenidos')
title('LECTURA ANALOGICA CON ARDUINO');
xlabel('Muestra');
ylabel('Voltaje de Salida');

Ahora creamos un vector de ceros llamado "val" de dimensión 1x1000 y en el cual guardaremos todos los datos que obtengamos de la lectura que habremos de realizar.

val=zeros(1,1000);

Enseguida creamos un ciclo "for" que se ejecutará desde 1 hasta 1000. Definimos los límites superior e inferior de la gráfica mostrada con "xlim" y "ylim", nótese que el valor para "y" permanece estático de 0 a 5.1 (valores en volts), mientras que el valor de los límites de "x" dependen de la ejecución del ciclo. Ahora creamos una variable llamada "lectura" la cual tomará el valor entero (eso nos indica el '%d' incluido en esa linea) que leamos del puerto serial con el comando "fscanf". Acto seguido, en el vector "val(i)" iremos guardando cada uno de los valores obtenidos en la lectura en cada ejecución del ciclo, además reconvertimos los valores obtenidos multiplicando por (5/1023) para con ello cambiar el rango de 0 a 1023 (Valor analógico) a uno proporcional de 0 a 5 (Voltaje de Salida). Con "hold on" nos aseguramos que los valores graficados en cada ciclo se mantegan, y evidente con "plot" graficamos los vectores "i" y "val(i)" usando como símbolo una "x". Finalmente el comando "drawnow" (literalmente) nos permite ir graficando los valores tan pronto como se vayan leyendo del puerto serial.

for i=1:1000
    ylim([0 5.1]); 
    xlim([i-100 i+10]);
    lectura=fscanf(pserial,'%d');
    val(i)=lectura*(5/1023);
    hold on
    plot(i,val(i),'x');
    drawnow
end

Siguiendo con el código, con la línea que se muestra enseguida escribimos en un fichero con extensión *.dat los valores que hemos guardado en el vector "val" durante la ejecución de todo el ciclo. Usamos saltos de línea como separadores (\n) y además guardamos los datos en modo de coma flotante con precisión de dos lugares decimales (%.2f)

dlmwrite('Valores_Obtenidos.dat', val, 'delimiter', '\n', 'precision', '%.2f')

Finalmente cerramos (fclose) y borramos (delete) el puerto serial que hemos inicializado. Además limpiamos las variables que hemos creado a lo largo del proceso (clear all).

fclose(pserial); 
delete(pserial);
clear all;

Enseguida os muestro las imágenes de la gráfica generada y de los datos que se han recolectado durante el proceso.

Finalmente sólo mencionar que todo lo anterior es adaptable para cualquier sensor analógico, para ello sólo bastará con adecuar los valores de salida leídos del puerto serial. Y claro, hay muchos mejoras que se pueden implementar, o digamoslo de otro modo, decorarlo un poco, pero obviamente eso depende en gran manera del uso para el cual se destine. 

Raíces de polinomios en MATLAB

Para calcular raíces de polinomios de cualquier grado, sean las raíces reales o complejas, se utiliza la función roots. La estructura que ha de escribirse en el "command window" es la siguiente:

>> roots([Vector de coeficientes del polinomio])

Donde el vector de coeficientes es un arreglo de los coeficientes ordenados de acuerdo al grado de mayor a menor.

Por ejemplo, supongamos que queremos calcular las raíces de $3x^3-4x^2+7x-10$, entonces el vector de coeficientes para este polinomio sería [3 -4 7 -10], y la instrucción correspondiente quedaría tal como se muestra:

>> roots([3 -4 7 -10])

ans =

  -0.0261 + 1.5508i
  -0.0261 - 1.5508i
   1.3856         

MATLAB nos devolverá un vector que contiene las raíces del polinomio.