PID Digital (paso a paso)

[krakc]

Hola a todos.

Amigos tengo un proyecto para una materia y me gustaria que me ayudaran, pero la verdad es que no quiero nada hecho, lo que deseo es que me ayuden a hacerlo paso a paso, uds me guian y yo lo hago, quiza algunas cosas ya esten hechas, pero si no hago las cosas no aprendo y esa es mi mas grande ilusion.

Objetivo: Crear un control PID digital usando pic16f877a para controlar la temperatura de un horno.

Datos: La planta sera un cajon de madera y en su interior estara una bombilla de 60W la cual hara las veces de calefactor, esta bombilla se controlara usando pwm, el cual con la ayuda de un circuito adicional controlara la intensidad luminica de la bombilla y a su vez el calor generado. Para sensar la temperatura se usara un LM35 el cual da 10mV/ºC

Nota: Ire editando este post a medida como vaya realizando nuevas cosas y recibiendo sugerencias y correcciones de su parte.


PASO 1:

* Armamos el circuito sensor usando el LM35, probamos de que sense bien la temperatura y lo introducimos en nuestra caja de madera, previamente elaborada; hay que tener en cuenta de que la caja quede bien sellada una vez se cierre, o habra problemas al momento de controlar y sensar la temperatura.


PASO 2:

* Graficamos los datos y hallamos la ecuacion caracteristica de la planta, esta ecuacion es distinta para cada planta, podemos usar un cicruito/programa sensor de temperatura que nos ayude a tomar los datos de la temperatura cada cierto intervalo de tiempo.... o ..... usamos lapiz, papel y cronometro para hacer todo manual... yo lo hare manual.

 -- Con estos datos pasamos a elaborar la grafica, manualmente, usando excel o usando matlab, como aqui lo que queremos es aprender, pues usaremos matlab, para eso ingresamos los datos en matlab en forma de dos matrices (tiempo/temperatura) y que nos la grafique.

Codigo para graficar los datos (esto es un ejemplo, los datos reales son demasiados):

clc;
tiempo[0 1 2 3 4 5 6 7 8 9 10];
temperatura[25 25 26 27 28 29 30 30 31 31 31];
plot(tiempo,temperatura)

Listo ahora tenemos graficada la respuesta de la planta ante una entrada escalon unitario.

Lo siguiente lo pueden hacer como lo explico a continuacion o como lo explica blackcat en un post mas adelante, los dos metodos son validos, ya cada quien escogera.

Para hallar la ecuacion, nos vamos a matlab (donde esta el workspace), alli veremos un boton que dice START (con el logo de matlab), pulsamos ese boton, luego click en toolboxes, luego en curve fitting y por ultimo en curve fitting tool: (figura 1)


Figura 1

Se nos abrira una ventana nueva y alli daremos click en el boton Data... en la ventana que se abrio en X data ponemos la variable tiempo y en Y data la variable temperatura (o sus los nombres que uds les tengan asignados) .. si no les aparece nada es por que tienen que haber corrido su codigo de graficacion de sus datos reales para que las variables del tiempo y temperatura aparezcan en el workspace...

luego daremos click en Create date set. y click en close (figura 2)


Figura 2

La grafica mostrada sera muy similar a la de su figure 1, pero en forma de puntos, lo cual no es ningun problema, por que ahora vamos a hallar nuestra ecuacion, para eso hacemos lo siguiente:

click en el boton fitting, luego en new fit .... aqui es donde les voy a recomendar algo a titulo personal...

Recomendacion: traten de la ecuacion quede en terminos de eulers y no en forma de polinomio, ya que si es asi todos los polos les quedaran en el mismo lugar (en cero) y no es muy buena idea para los calculos, pero en fin ya sabran uds como hacerlo...

Vamos a modificar unicamente donde dice Tye of fitt. ahi podran ver todo el tipo de ecuaciones con las cuales se puede representar nuestra curva caracteristica de la planta... si siguen mi recomendacion entonces buscar donde dice exponential y escogemos el tipo de ecuacion ( a*exp(bx) + c*exp(dx) ) y hacemos click en apply... luego de eso nos aparecera nuestra ecuacion (la que nos describe nuestros datos) ... si tienen dudas de que valores para las constantes escoger, escojan los promedios, osea los valores que estan fuera de los parentesis ... (figura 3)


Figura 3

Ahora damos click en close y confirmamos que la ecuacion reproduce lo mas fielmente posible la grafica de nuestros datos si es asi ya tienemos nuestro y(t) = a*exp(bt) + c*exp(dt) y tambien los valores de las constantes... y si no s asi buscamos la que mejor lo haga.. pero estoy un 90% seguro de que la funcion exponential lo hara de maravilla.

Por lo tanto obtenemos una ecuacion y(t) = a*exp(b*t) + c*exp(d*t) aplicando laplace obtenemos:

                32.98 s + 0.06097
Y(S) =   ----------------------------
         s^2 + 0.00141 s - 9.773e-009

Para mi caso ... a uds les dara una ecuacion similar .. matlab tambien lo hace usando el codigo

syms t
a = 42.98;
b = -6.9x10^-6;
c = -10.15;
d = 0.001417;
yt = a*exp(b*t) + c*exp(d*t)
YS = laplace (yt)

pero la verdad es que se ve muy raro y desordenado, asi que recomiendo lapiz y papel que hacerlo no es dificil y ademas sale directo, luego se aplica un poco de algebra, se reemplazan valores y listo.

Ahora sabemos que dentro de Y(S) tenemos implicito nuestra entrada escalon unitario, po lo tanto nuestra G(S) sera:

           32.93 s^2 + 0.06097 s
G(S) = ----------------------------
           s^2 + 0.00141 s - 9.773e-009

Ahora tenemos que determinar el tiempo de muestreo (T) para eso podemos usar el metodo que describe blackcat mas adelante o podemos asumirlo, ya cada quien decidira, por mi parte lo sumire como T=0.1

Pasamos a la parte en que tenemos que aplicar la transformada Z para poder discretizar nuestro G(S) y obtener nuestro G(Z) .. pero recuerden que G(Z) no es igual a la transformada Z de G(S) si no que lleva incorporado un retenedor de un orden cualquiera, para este caso usaremos uno de orden cero... y lo llamaremos Bo(S) asi que:

G(Z) = Z/ Bo(S) G(S) donde Z/ sera la tranformada Z de todo eso... bueno para los amantes del lapiz y papel, podeis hacerlo, por que por mi parte no lo voy a hacer, para eso voy a usar matlab y el siguiente codigo.

n1=32.98;
n2=0.06097;
d1=1;
d2=0.0014101;
d3=-(9.773/1000000000);
GS=tf( [n1 n2 0], [d1 d2 d3])

T=0.1;
GZ=C2D(GS,T,'zoh')

donde zoh indica a matlab que se va a usar un retenedor de orden cero. osea nuestro Bo(S)

y el resultado que nos entrega matlab sera el siguiente:

             3.293 z - 3.292
G(Z) =   ------------------
           z^2 - 2 z + 0.9999
 
Sampling time: 0.1


Bueno hasta aqui todos deberiamos llevar algo similar, ahora es donde yo me he decidido por hacer un controlar por el metodo de ragazzini .. asi que por ahi es donde lo haremos....

lo primero es hallar los polos deseados teniendo en cuenta las siguientes condiciones de lazo cerrado.
essp = 0 ; Test 5% = 10 seg ; T = 0.1 ; Zcita = 0.5

usando formulas ya conocidas en nuestros estudios realizados obtenemos los dos polos deseados Z(1,2):

Z(1,2) = 0.9691 +/- j 0.0504    lo cual llamaremos (alfa +/- j beta) para comodidad

El metodo ragazzini nos dice que:

F(Z) = K / P(Z)

ahora P(Z) = (Z - Z1) (Z - Z2) como ya sabemos los valores de Z1 y Z2 entonces reemplazamos... pero podemos hacer algo mejor y es reemplazar por las letras y despejar para que nos quede algo asi:

P(Z) = Z^2 - 2 alfa Z + (alfa^2 + beta^2) ahora se reemplazan los valores y listo.

para que ragazzini funcione se dice que F(1) = 1 .. lo cual es de gran ayuda.. sabemos por formula que:

F(Z) = K / P(Z) --- reemplazamos P(Z) y nos quedara asi:

F(Z) = K / Z^2 - 1.9382 Z + 0.94169 --- hacemos F(1) = 1 --- quedando: 1 = K / 1^2 - 1.9382 + 0.94169 -- despejando K

K = 0.00349 ahora pasamos a la parte interesante...

D(Z) = [ 1/G(Z) ] [ F(Z) / 1 - F(Z) ]

pero como todo se puede simplificar entonces lo hacemos.

D(Z) = [ 1 / G(Z) ] [ K / P(Z) - K ]

como ya conocemos todas las variables y ecuaciones, simplemente hallamos D(Z)

D(Z) = [ 0.001058Z^2 - 0.002116Z + 0.001058 ] / [ Z^3 - 2.93769Z^2 - 0.9379 ]

Bueno la parte final es hallar la ley de control, para eso cambiamos D(Z) por u(Z)/ e(Z) y multiplicamos al otro lado por Z^-3 / Z^-3

quedandomos asi:

u(Z)     0.001058Z^-1 - 0.002116Z^-2 + 0.001058Z^-3
----- = ---------------------------------------------------------
e(z)     1 - 2.93769Z^-1 + 0.9382Z^-2 - 0.9379Z^-3

ahora el denominador pasa a multimplicar a u(Z) y e(Z) a multiplicar al numerador. luego de lo cual despejamos u(Z) y de paso usamos una un teorema que dice en resumen que u(Z)Z^-1 = u(Z-1) .. lo mismo con e(Z) .... ahora cambiamos Z por k quedando finalmente la ley de control asi:

u(k) = 0.001058 e(k-1) - 0.002116 e(k-2) + 0.001058 e(k-3) + 2.93769 u(k-1) - 0.9382 u(k-2) + 0.9379 u(k-3)

esta ley de control u(k) es la que programamos para poder controlar nuestro sistema



SE EDITARA A MEDIDA QUE SE AGREGUEN NUEVAS COSAS....

[blackcat]

Hola ...

Para poder estimar un sistema con matlab se utiliza la herramienta IDENT .. esta es muy facil de usar, lo que se necesita es lo siguiente:

-> un vector de entradas (x) ... la señal de estimulo que aplicaste al sistema
-> un vector de salidas (y) ... la señal de respuesta que obtuviste del sistema al aplicar la señal de entrada
-> un vector de tiempo (t) ... los tiempos de cada muestra, deben ser constantes, es decir, 1, 2, 3, 4 segundos etc ... y no 1, 1.5, 2, 3, 5, 8 .. etc..

Con eso IDENT identificará el sistema que se encuentra adentro.



Lo que necesitas es una señal de prende o apaga el bombillo .. o bien una señal de volaje RMS que le aplicas al bombillo, luego debes medir en tiempo real el comportamiento de la temperatura ante la señal. Un escalon como lo estas haciendo, es posible, sin embargo solo tenes el comportamiento de calentar y no de enfriar.

Para poner a prueba el IDENT me inventé un sistema tipico y muy facil de entender .. un filtro pasobajos RC de primer orden, representado en Laplace como:

\begin{equation}
G(s) = \frac{1}{1 + 0.001s}
\end{equation}

Es decir R = 1k y C = 1uF ... siendo tau = 1ms;


>> R = 1000
>> C = 1e-6
>> G = tf( 1, [R*C 1] )

Lo que voy hacer es inventar un vector de entradas, estimulare el sistema con esas entradas y utilizaré IDENT para identificar el sistema; lo ideal es que obtenga un resultado como en de la ecuacion anterior. En el caso practico, no conozco la ecuación del sistema, pero lo que se acostumbra es estimular el sistema con diferentes vectores de entrada y estimar varios modelos y comparar el resultado con los otros vectores.

Como mencioné, es necesario que tenga un vector de datos de entrada y ver la respuesta del sistema (salida) ante esa entrada, el mejor vector de entrada es aquel que tenga un buen contenido de frecuencia, ruido por ejemplo; sin embargo, dependiendo del sistema esto se complica, por ejemplo, en un motor; lo recomendable es estimar el modelo matemáticamente y ver que frecuencias son las aptas para crear una señal de estímulo; por ejemplo, el modelo de un motor es parecido a un filtro pasobajos, entonces utilizamos una señal que contenga componentes de frecuencia dentro del ancho de banda del motor.

Como en este caso, el ejemplo es simulado, utilizo la función RAND para generar un vector de entrada aleatorio de 100 muestras:

>> x = rand([100, 1])



Esta es mi señal de entrada, ahora, necesito el vector de salida, estimularé el sistema con la función LSIM, entre los datos de esta funcion está el vector de tiempo. Aqui voy a utilizar tiempo discreto y elejiré que el tiempo de muestreo sea un décimo de la constante de tiempo del sistema (tau = RC); en este caso, ts = 0.1ms. Genero un vector de tiempo de 100 datos desde 0 hasta 99*0.1e-3.

>> n = (0:99)'
>> ts = 0.1e-3
>> t = n*ts

Ahora utilizo LSIM ....

>> y = lsim(G, x, t)

La grafica que obtengo en el tiempo es:



donde la linea azul es la entrada y la linea verde es la salida, veran que esta filtrada.

Ya tengo listo mi vector de entrada, salida y tiempo; esto es lo que tendria en la practica, una señal de entrada, el comportamiento de la salida y el tiempo en que fue tomado cada dato.  En el caso de la temperatura, el vector de entrada seria pulsos rectangulares, que me dicen cuando se prendio el bombillo y cuando permanecio apagado; la salida seria el aumento de temperatura y la disminucion de temperatura, el tiempo seria cada cuanto tome la muestra, como es temperatura puede ser de 1 o 3seg, constantes y no interrumpidos ... en este caso, tomar los datos a reloj en mano no es conveniente.

Yo conozco el sistema; sin embargo, la idea es estimar el sistema en tiempo discreto utilizando el vector de entrada, salida y tiempo; luego transformar el sistema a tiempo continuo y comparar.

Ejecuto IDENT

>>ident

Y aparece una pantalla como la siguiente:



Lo que hacemos es en IMPORT DATA seleccionamos TIME DOMAIN DATA ... esto porque nuestros datos de entrada, salida y tiempo estan en el tiempo. Aparece un cuadro como asi:



Y lo lleno con la informacion que tengo, entrada (x), salida (y), tiempo de inicio (0), tiempo entre muestras (ts) como lo hice en la figura. Damos IMPORT ..



En la columna de cuadros izquierda me aparece el bloque de datos que cargue en IDENT; ahi puedo cargar la cantidad de datos que desee y puedo combinar datos como de estimulo y verificacion, esto es que con un conjunto de datos puedo estimar el sistema y con otro conjunto puedo verificar si la estimacion está correcta.

Ahora, sobre ESTIMATE damos LINER PARAMETRIC MODELS ... aqui escogemos el método de estimación que mas nos guste o que mejor de resultados, sobre el tipo de estimador se puede encontrar información en internet:

http://www.ie.itcr.ac.cr/einteriano/control2/Laboratorio/3.Models.pdf

Yo elejí ARX tanto de orden 1 como de orden 2 ... entre mas orden es mejor la estimación; sin embargo, hay un compromiso, pues estimar un compensador para un sistema que tiene un orden alto es mas dificil que para uno bajo, entonces escogemos el modelo de menor orden. Seleccionando la casilla MODEL OUTPUT aparecerá una gráfica que nos dice cual modelo es el que mejor se ajusta a los datos de salida. Vemos que con el modelo de 1º orden tenemos un ajuste de 100.



Ahora arrastramos el cuadro que dice ARX111 de la columna derecha sobre el cuadro TO WORKSPACE ... eso hara que los datos del modelo lo tengamos en la linea de comando. Para ver el modelo ejecutamos:

>> H = tf(arx111)
 
Transfer function from input "u1" to output "y1":
 0.09516
----------
z - 0.9048
 
Transfer function from input "v@y1" to output "y1":
1.665e-016 z
------------
 z - 0.9048
 
Input groups:           
      Name      Channels
    Measured       1   
     Noise         2   
                       
Sampling time: 0.0001

Entonces nuestro modelo estimado en tiempo discreto es:

\begin{equation}
G(z) = \frac{0.09516}{z - 0.9048}
\end{equation}

Usamos:

>> G2 = d2c(H)

y obtenemos:

\begin{equation}
G(s) = \frac{1000}{s + 1000}
\end{equation}

Vemos que IDENT identificó el sistema correctamente ...  Utilizando este modelo podemos usar SISOTOOL y calcular de manera sencilla un compensador 

Saludos y espero que les sirva de algo!

PD: Espero sus comentarios y estamos listos para la implementación PID en un PIC!

[krakc]

Bueno amigo me di a la tarea de buscar lo que no entendia y segun como voy armando las cosas paso a paso y encontre la primera pregunta y la segunda.

Lo que has escrito esta completisimo y muy interesante, que cosas las que sabes, hombre yo ni lo habia pensado asi, pues lo que queria era algo sencillo, pero eso seguro que servira, pero se necesita algo de tiempo para poder aprender a usarlo no es facil...

pero eso no seria usando un circuito para leer los datos y que matlab los interprete?????
si es asi a mi no me sirve jejeje pues yo tome los datos manualmente....


quisiera que veas lo que he hecho y veas que era lo que yo queria, haber si me entiendes:

****************
usando datos experimentales (no son los reales) se intrducen a matlab, luego se grafican, una vez graficados lo que se hace es usar BASIC FITTING hallamos la ecuacion de esa grafica y este es el resultado... (se puede ver la ecuacion que mas se asemeje o modele mas fielmente la grafica)

--- aclaro que la grafica de los datos reales de la planta seran muy similares a esta grafica, pero con mucho mas datos ---





una vez teniendo esta ecuacion (la ecuacion de la grafica) se puede usar tf() para hallar la funcion de transferencia y empezar el pid matematicamente ??????''

[blackcat]

El detalle está en que estas estimando un modelo en el tiempo ante un escalon no podes usar TF asi por asi ya que TF crea el sistema deacuerdo a los coeficientes de la ecuacion diferencial ...

Podes hacer un tiempo de 10 segundos calentando (bombillo encendido), 5 seg enfriando (bombillo apagado), 3 segundos calentando, 5 enfriando. Notaras un comportamiento en la temperatura como de carga o descarga. Con esta secuencia e IDENT podes estimar un buen modelo.

Sin embargo, recuedo que unos compañeros estimaron un modelo termico similar a lo que estas haciendo, el siguiente PDF quiza te ayude:

http://www.ie.itcr.ac.cr/einteriano/control2/Clase2008/EjemploControldeunSistemaTermico.pdf

Saludos!

[krakc]

Gracias amigo, creo que voy entendiendo, pues lo que yo he hecho es como me lo han enseñado, no me han dicho nada de lo que tu me dices, pero lo hare y te dire como me va, lo mas seguro es que lo haga de mi forma y de tu forma...

por cierto el pdf esta completisimo, lo guarde y lo revisare mañana con mi compañero del pryecto y discutiremos los dos metodos.

revise rapidamente el pdf y me fije que el siguiente paso a seguir es analizar la grafica....

bueno mañana revisare mas cosas, por que esto si va en serio XD

[blackcat]

Bueno ... metodos hay muchos .. sin embargo, segun mi experiencia, el metodo que te explique es el que mejor resultados me ha dado, mi interes no se enfoca en solo los sistemas termicos sino en muchos como mecanicos (motores, etc ), magneticos, fluidos y hasta quimicos ...

El detalle que tienes que tener cuidado es que el PDF considera que el calentamiento y el enfriamiento es homogeneo, es decir, como la carga y descarga de un capacitor ... sin embargo, hay factores que pueden provocar que este comportamiento no sea homogeneo, como la disipacion termica, que tan aislante es la caja del exterior .. etc etc ...

Veo que eres nuevo por el foro y propusiste un tema que a mi parecer es bastante interesante y con mucha ingenieria ... (Aunque en la practica industrial el problema se soluciona con un regulador de temperatura autosintonizado  .. ) que ayuda a explotar la imaginacion y ciencia de todos los del foro. Te felicito pues me gustan este estilo de proyectos con mucho contenido matemático; sin embargo, aprender hacer un PID y todo el modelado involucra mucho conocimiento, en el que el fin es que todos aportemos ideas y aprendamos.

Saludos!

[Renatox_]

La explicación que te ha dado blackat es buena a mi me ayudó para hallar la TF de mi motor, te permite hallar la función de transferencia con la transformada Z, con esta diseñas tu controlador discreto que lo llevas a tu pic en forma de ecuacion en diferencias.

Si usas ziegler-nichols tu respuesta al escalón tiene que tener la forma de S bien definida, para el caso de un motor dc este método no me sirvió porque la S no estaba bien definida al comienzo y el PID no era nada real porque trataba de establecer al motor en pocos milisegundos cosa nada real.

El Lugar geométrico de las raices me resultó mucho más seguro y preciso, manejando el estado permanente y transitorio con mayor aproximación a como salia en la simulación.
También podrías analizar tu sistema en el plano s y tu controlador analógico lo conviertes a uno digital.

saludos,

[Mario]

Información muy completa Blackcat.


Quizás les sirva estoesto


Tendré que realizar un modelado como lo hizo Blackcat para comprender y después hacer una identificación de parámetros para un motor de DC (tomando RPMs vs tiempo).

[krakc]

Bueno segun he leido y leyendo sus aportes he concluido que el metodo de sintonizacion de ziegler-nichols no es muy exacto, pues se basa en algo un poco ideal....

entonces como ya tengo la ecuacion que define mi sistema, hallada con el maravilloso matlab, tengo que hallar ahora el G(s) de esa funcion.

para eso usamos matlab y hallamos la funcion de transferencia, para esto recurrimos a la ecuacion encontrada con basic fitting, lo que tenemos que hacer es crear una matriz y llenarla con los datos linealizados, si no saben el codigo aqui se los dejo:

X=0:1:10;
Y=[];
for jj=1:1
    for ii=1:11
            ecu=(0.002*(ii-1)^4) - 0.053*(ii-1)^3 + 0.4*(ii-1)^2 - 0.1*(ii-1) + 25;
            Y(jj,ii)=(ecu);
    end
end
Y
plot(X,Y)

ahora vemos que matlab nos muestra la matriz "Y" con los datos de la ecuacion "ecu" y si los comparamos con los datos experimentales, nos fijaremos que son muy, pero muy similares, ahora al sacar la funcion de transferencia de "X" vs "Y" y la de "tiempo" vs "temperatura" y compararlas vemos que es aun mas similares.

se hace asi:

tf(X,Y)
tf(tiempo,temperatura)

AQUI es donde quiero que me ayuden a aclarar o a corregir lo que este mal (o todo, si esta mal todo) ... teniendo la funcion de transferencia

Transfer function:
 
     s^9 + 2 s^8 + 3 s^7 + 4 s^6 + 5 s^5 + 6 s^4 + 7 s^3 + 8 s^2 + 9 s + 10
--------------------------------------------------------------------------------------------------------------------------------------------------
25 s^10 + 25.25 s^9 + 26.01 s^8 + 27.03 s^7 + 28.12 s^6 + 29.13 s^5 + 29.94 s^4 + 30.52 s^3 + 30.86 s^2 + 30.98 s + 31

le aplico transformada Z para discretizarla y poder hacer el metodo de lugar geometrico de las raices y diseñar el controlador ?????????????????

o me salte algun paso, seria mejor usndo algun otro metodo o como hago para seguir adelante....

gracias.

PD: por el metodo de ziegler nichols tambien voy a probar paralelamente (para quienes solo les han enseñado ese metodo) asi que teniendo una grafica en forma de S basta, donde se vea claramente un tiempo "muerto". Como en la del ejemplo.

[blackcat]

MMMmmmm ... hay algo que no cuadra ... parece que estas confundiendo las peras con los tomates ... voy a retornar a los cuadernos y repasar un poco de matematica avanzada para ingenieria.

Lo que tenes es una ecuacion que representa algo como funcion del tiempo ante una entrada escalon ... sin embargo, no estas representando a un sistema por completo ... esto es algo mas complicado (y la peor pesadilla de muchos!  )

Aqui tenemos que ser claros .. la entrada del sistema no es el tiempo . sino la funcion escalon! es decir [0 0 0 1 1 1 1 1 ... ] donde 1 es el momento en que se encendio el bombillo. Para describir un sistema en el tiempo y en palabras matematicas tenemos algo como asi:

\begin{equation}
y(t) = h(t) * x(t)
\end{equation}

Donde y(t) es la salida o mas bien la temperatura, x(t) es la entrada o la tension del bombillo, si es fija, puede usarse 1 para encendido o 0 para apagado. h(t) es la funcion de transferencia del sistema descrita en el tiempo (muchas veces como ecuaciones diferenciales y exponenciales)

Y el operador * no es una multiplicacion!!!; si no, la temible, terrorifica y tragica funcion de convolucion que ha destrozado miles de sueños de mucha gente!. Quiza la han oido mencionar ...

Ahora ... ¿¿¿porque todo el cuento de la convolucion??? ... si averiguas que es la funcion de convolucion veras que es una formula un poco complicada (tanto de entender como de resolver ) ... Dios trabaja de manera misteriosa pero aveces deja portillos abiertos que simplifican mucho las cosas, pero sin antes haberse roto el coco .. en este caso el que se rompio el coco fue un señor de apellido Laplace (un genio y casi padre de la ingenieria, tambien como Heaviside y Fourier), la ventaja es que Laplace "transforma" la convolucion en una simple multiplicacion. Entonces, se veria algo como asi:

\begin{equation}
Y(s) = H(s)X(s)
\end{equation}

Donde Y(s) es la transformada Laplace de y(t) que es la señal de salida, X(s) es la transformada de Laplace de x(t) o la señal de entrada y H(s) es nuestro querido y desconocido sistema a encontrar! ... que lo averiguamos simplemente haciendo:

\begin{equation}
H(s) = \frac{Y(s)}{X(s)}
\end{equation}

Las cosas van mejorando, solo tenemos que encontrar Y(s) y X(s) .... la transformada de Laplace de una funcion escalon o x(t) = 1 (para t > 0) es:

\begin{equation}
X(s) = \frac{1}{s}
\end{equation}

Ahora solo tenes que averiguar la transformada de Laplace de la funcion que obtuviste osea Y(s) y hacer:

\begin{equation}
H(s) = sY(s)
\end{equation}

Pero te advierto que: Te deseo muy buena suerte transformando a Laplace (a lapiz y papel) esa funcion (mas si nunca lo has hecho antes pero me imagino que si, entonces no hay problema  ) ... sin embargo, hay muchas herramientas que te facilitan el proceso ... creo que en MATLAB tiene funciones para eso ... nunca las he usado ... la HP49G me lo ha hecho todo! 

Para entender mejor como esta todo este cuento utilizaré el ejemplo del filtro RC ... tenemos que:

\begin{equation}
H(s) = \frac{1}{RCs + 1}
\end{equation}

Si aplicamos un escalon de voltaje a un circuito RC ... el capacitor se empieza a cargar .... ya todos conocemos la funcion de carga de un capacitor ... si multiplico (en el dominio s) H(s) por X(s) donde X(s) es la funcion escalon, entonces tendria:

\begin{equation}
Y(s) = H(s)X(s) = \frac{1}{RCs^2 + s}
\end{equation}

Ahora si aplico Laplace inverso .. obtengo la funcion de carga de un capacitor en el tiempo:

\begin{equation}
y(t) = 1 - e^{\frac{-t}{RC}}
\end{equation}

Si la matematica no miente, yo no miento y todo esta comprobable!

El otro detalle importantisisimo es que la funcion TF de MATLAB no hace eso que estas pensando ... TF crea el sistema conociendo ya los coeficientes, ejemplo, teniendo el sistema:

\begin{equation}
G(s) = \frac{s^2 + 2s + 3}{4s^2 + 5s + 6}
\end{equation}

La introduccion en MATLAB de ese sistema es:

>>G = tf( [1 2 3], [4 5 6] )

Es decir, TF solo sirve para agregar un sistema ... no te facilita muchas cosas!

Si ya pudiste introducir el sistema, podes verificarlo haciendo la funcion STEP o escalon,

>>step(G)

Esto te genera un grafico con la respuesta del sistema ante un escalon, deberias ver algo similar a lo que obtuviste experimentalmente.

Bueno eso es un resumen ... Espero que te ayude y no te moleste, pues esa es la matematica del asunto, muchos la evitan diariamente ¿¿Como?? usando reguladores PID ya fabricados ... pero lamentablemente los que estudian ingenieria no tienen como salvarse de eso ... Saludos!

[Renatox_]

hola, un sistema de orden 10 es demasiado para  tu planta, si tu respuesta tiene la forma de S o es como la carga de un condensador ante un escalón, entonces aproximalo a un sistema de primer orden, incluso eso lo puedes hacer manualmente, hallas el tiempo de establecimiento T,,, entonces tau=T/5, y tu ganancia seria la temperatura final entre el voltaje de entrada al foco, osea:

Gplanta=A/(tau*s+1);

donde:  A=Tss/U      ; Tss=temp.final, U=voltaje de entrada.
tau=T/5;

Con Gp puedes usar él método LGR (Lugar geométrico de las raices para hallar tu controlador),

Agregando algo más al concepto de convolución que dio Blackat, h(t) es la respuesta que tiene un sistema ante una entrada impulso. Entonces conociendo esta respuesta se puede calcular la respuesta que tendría este sistema ante cualquier otra entrada, convolucionando h(t) con la entrada x(t), que es un juego de sumas y productos entre ellos.

saludos,

[Renatox_]

Blackat tu última funcion G(s) no existe fisicamente, debería estar invertida.

saludos,

[blackcat]

Blackat tu última funcion G(s) no existe fisicamente, debería estar invertida.

saludos,

A si si ... fue algo que me inventé en el momento ... voy a corregirlo


El método de Renatox es correcto ... es el mismo que está propuesto en el pdf

[krakc]

Excelente, como me habeis abierto los ojos, definitivamente nunca debi pensar que esas aburridas clases con ese aburrido profesor de matematicas especiales eran aburridas, pero en fin, lo de convolucion, la place, fourier y Z lo tengo algo claro, tendre que repasar y no me es dificil volver a entenderlo.

Leyendo lo que ha escrito blackcat me perdi en una cosa:

haber
H(s) = Es quien me describe la planta y lo que tengo que hallar
X(s) = Entrada al sistema, en este caso seria el escalon unitario, por que solo inicio el sistema. (solo bombilla y sensor) ¿verdad?
Y(s) = Aqui esta mi gran duda, este es mi Y(t) (hallado con basic fitting) pasandolo al plano S con ayuda de LaPlace ?????

Si las cosas son asi, no resta si no pasar Y(t) -> Y(s) y luego hallar mi anhelado H(s) ????

---------------------

Renatox, lo que dices es cierto, pero si lo hacemos por ese lado, simplemente uso el metodo de sintonizacion de ziegler nichols, el cual con ayuda de la grafica ya me da todo, ya tengo mi h(t) etc... pero como se sabra, no es muy exacto y toca empezar a variar Kp, Ti y Td, luego de variarlos, se escoje el que mejr controle y sale, aunque no se aplica mucha matetica ni mucha intuicion pues todo viene prefabricado XD

-------------------------

renatox tienes razon, blackcat se ha equivocado, el denominador y numerador deberian ir invertidos...

[blackcat]

Si ... y(t) es una ecuacion que varia con el tiempo y esa es la respuesta ante un escalon .. asi como al obtuviste con basic fitting ... deberias convertir y(t) -> Y(s) ... te doy una pista:

La transformada Laplace de una señal descrita como:

\begin{equation}
\frac{t^{n-1}}{(n-1)!}
\end{equation}

es

\begin{equation}
\frac{1}{s^n}
\end{equation}