Implementando una función de transferencia de tiempo continua en código de PC

4

Tengo una función de transferencia como esta:

\ $ H (s) = \ frac {1} {Ts + 1} \ $

Lo convertí en una ecuación de diferencias para resolverlo de manera iterativa:

\ $ H (s) = \ frac {Y (s)} {X (s)} = \ frac {1} {Ts + 1} \\ TsY (s) + Y (s) = X (s) \\ \ mbox {Suponiendo que las condiciones iniciales son cero,} \\ T \ frac {dy (t)} {dt} + y (t) = x (t) \\ y (t) = x (t) - Ty '(t) \ $

Escribí un código como este;

const double dt = 0.001;
double x;       // Input
double y;       // Output
double yp=0;    // Previous output value
double yd;      // Derivative of output
for(double t=0; t<TIMEMAX; t+=dt)
{
    // ...
    x = ReadNewInput();
    yd = (y - yp) / dt;
    yp = y;
    y = x - T*yd;       // y(t) = x(t) - Ty'(t)
    StoreNewOutput(y);
    //...
}

T simboliza el retraso del sensor (el sensor que lee x (t)) en mi código. Si T = 0 (sin retraso del sensor) mi código funciona muy bien. Pero si configuro un retardo de sensor muy pequeño (por ejemplo, T = 0.1s), la salida se vuelve inestable (y se acerca al infinito después de unas pocas iteraciones).

¿Estoy haciendo algo mal? Esta es la primera vez que implemento una función de transferencia en un algoritmo de computadora como este. ¿La realización de las funciones de transferencia continua de tiempo se realiza de esta manera o las personas utilizan un método diferente? Por favor, ¿puede confirmar la validez de mi método o hacer correcciones en él?

    
pregunta hkBattousai

2 respuestas

8

Está preguntando sobre el mundo de soluciones numéricas para ecuaciones diferenciales ordinarias . Y creo que lo que intentas obtener en tu código es esencialmente Método de Euler para resolver este tipo de ecuación. Sin embargo, creo que de alguna manera te has volcado en analizar el problema.

En lugar de su ecuación, \ $ y (t) = x (t) - Ty '(t) \ $, es más común analizar este problema reorganizando los términos:

\ $ y '(t) = \ frac {1} {T} (x (t) -y (t)) \ $

Para aplicar el método de Euler, estimamos la evolución de \ $ y (t) \ $ by

\ $ y (t + h) = y (t) + hy '(t) \ $,

donde h es un paso de tiempo para la solución numérica (el dt en su código). Podemos poner esto en notación de tiempo discreto como:

\ $ y_ {n + 1} = y_n + hy'_n \ $

Esta es una formulación general para cualquier ecuación diferencial ordinaria. Para tu problema, tendrías

\ $ y_ {n + 1} = y_n + h \ frac {1} {T} (x_n - y_n) \ $.

Que puedes traducir fácilmente en código C.

Dicho esto, este método tiende a no ser particularmente preciso. El error general en la solución crece en proporción al tamaño de paso de la simulación, y el método numérico en sí mismo puede ser inestable (los errores pueden ir al infinito incluso si la solución real es finita) en algunas circunstancias.

El método habitual (un método que funciona para muchos problemas, algo que se debe probar primero antes de buscar métodos más avanzados en caso de que su problema resulte especialmente difícil) es integración Runge-Kutta de cuarto orden . Las ecuaciones para esto son algo más complicadas, pero 1) hay muchas bibliotecas escritas para hacerlo, y 2) si trabajas las ecuaciones con cuidado, verás que todas encajan muy bien en su lugar. También hay numerosas elaboraciones disponibles en las bibliotecas, por ejemplo, para automatizar la búsqueda del tamaño de paso más grande (y por lo tanto, el menor tiempo de procesamiento) para obtener la precisión requerida en la solución.

    
respondido por el The Photon
1

Su filtro sería más estable calculado como un filtro IIR con dos coeficientes.

La función de transferencia

H (s) = 1 / (T * s + 1)

Es un filtro de paso bajo de primer orden con una frecuencia de corte 1 / T.

Para calcular el resultado y [n] utilizando el paso de tiempo discreto dt, en cada momento n * dt, como filtro de tiempo discreto IIR, una forma general es ...

y [n] = Beta * y [n - 1] + Alpha * x [n]

y [n] es la salida del filtro en el momento n * dt
y [n-1] es el valor del filtro en el momento (n - 1) * dt, el paso de tiempo anterior.
x [n] es la entrada del filtro en el momento n * dt
Beta = exp (-dt / T)
Alfa = dt

Tenga en cuenta que esta forma es incondicionalmente estable siempre que Beta sea menor que 1, lo que será para todos los pasos de tiempo positivo y las frecuencias de corte. Cuando se calcula en un dispositivo digital con precisión finita, existe la restricción adicional de que es posible que los números de punto flotante deban redondearse hacia 0 para la estabilidad.

La derivación es la siguiente ...

H (s) = 1 / (T * s + 1) = 1 / T * 1 / (s + 1 / T)

Tomando la transformada inversa de laplace para encontrar la respuesta de impulso de dominio de tiempo que obtenemos ...
h (t) = 1 / T * e ^ (- t / T)

En general ...
Para calcular la respuesta de la función de transferencia h (t) a una señal de entrada x (t) en tiempo continuo, tome la convolución ...

y (t) = integral (h (t - tau) * x (tau) * dtau), de tau = -infinito a t

en el caso de tiempo discreto, usando el paso de tiempo dt, la integral se aproxima como una suma, y la convolución es ...

y [n] = suma (h [n-i] * x [i] * dt), sumado de i = -infinito a n, el tiempo de simulación actual es n * dt

y [n] es la salida del filtro en el momento n * dt.
h [n-i] es la respuesta de impulso del filtro en el momento (n-i) dt
x [i] es el valor de la señal de entrada en el momento i
dt

Dada su función de transferencia específica, la convolución de tiempo discreto es ...

y [n] = 1 / T * sum (e ^ - ((n-i) * dt / T) * x (i * dt) * dt), sumado de i = -infinito a n
Tenga en cuenta que existe la siguiente relación de recursión ...
y [n] = y [n-1] * e ^ (- dt / T) + 1 / T * x [n] * dt


Para generalizar aún más, casi cualquier función de transferencia polinomial de la forma ...
Nu [n] * s ^ n + Nu [n-1] * s ^ (n-1) + ... + Nu [1] * s + Nu [0] / (De [n] * s ^ n + De [n-1] * s ^ (n-1) + ... + De [1] * s + De [0])

Donde Nu [] son coeficientes de numerador constantes.
De [] son coeficientes constantes del denominador.

Se puede tener en cuenta en la forma ...

(Nu [n] * s ^ n + ... + Nu [0])) * (C [n] / (sR [n]) + ... + C [0] / (sR [0]))

Donde C [n] son constantes y R [n] son las raíces del denominador en la ecuación.

De la forma de la ecuación es claro que en el dominio del tiempo esto es solo una serie de exponenciales y sus derivadas. Para una ecuación que tenga X términos exponenciales, obtienes ...

y [0] [n] = Beta [0] * y [0] [n-1] + Alpha [0] * x [n]
y [1] [n] = Beta [1] * y [1] [n-1] + Alpha [1] * x [n]
...
y [X] [n] = Beta [X] * y [1] [n-1] + Alfa [X] * x [n]

Donde y [i] [n] es la respuesta de filtro parcial y [i] en el tiempo n * dt debido a uno de los términos exponenciales.
Beta [i] es el factor de disminución constante para la respuesta de filtro parcial y [i].
Alfa [i] es el factor de ponderación de entrada constante para la respuesta de filtro parcial y [i]
x [n] es la entrada del filtro en el momento n * dt.

El resultado total es la suma del resultado de cada respuesta de filtro parcial.

y [n] = y [0] + y [1] + ... + y [X]

    
respondido por el user4574

Lea otras preguntas en las etiquetas