ELAI S.r.l.

Controllo robotico campionato: quando un PD stabile diventa instabile

Derivazione esatta ZOH, criterio di Jury, confronto con Euler e ritardo di un campione: simulazione riproducibile di un asse ideale.

Controllo robotico campionato: quando un PD stabile diventa instabile

Robotica · Analisi e simulazione · 24 settembre 2026

Abstract: la frequenza del controllo fa parte del progetto

Un controllo proporzionale-derivativo può stabilizzare perfettamente un modello continuo e fallire quando viene aggiornato a intervalli discreti. Studiamo una massa ideale comandata in forza, ricaviamo la mappa esatta con comando mantenuto e determiniamo l’intervallo di campionamento che garantisce stabilità asintotica nel modello. Confrontiamo poi questa mappa con Euler in avanti e con un ritardo di un campione. L’obiettivo è capire quale sistema stiamo realmente analizzando, prima di fidarci del grafico di una simulazione.

L’analisi richiede equazioni differenziali lineari, matrici 2×2 e autovalori. Non è un nuovo algoritmo robotico né un benchmark industriale: è una derivazione didattica con codice eseguito, pertinente agli assi di movimento e ai modelli locali di un manipolatore. Non comprende contatti, attrito, elasticità, saturazione o certificazioni. Proprio questa semplicità consente di isolare l’effetto del campionamento senza confonderlo con la cinematica di un robot completo.

1. Massa, forza e guadagni: prima le unità

Sia q lo spostamento rispetto al riferimento, in metri, v la velocità in m/s, m la massa in kg e u la forza in newton. Il riferimento è costante e viene posto a zero. Assumiamo misura esatta di posizione e velocità e un attuatore ideale. Il guadagno k_p ha unità N/m, k_d ha unità N·s/m. L’errore di segno sarebbe decisivo: la forza proporzionale deve opporsi allo spostamento, quella derivativa alla velocità.

dq/dt = v m dv/dt = u u(t) = −k_p q(t) − k_d v(t) m q″ + k_d q′ + k_p q = 0

Nel continuo, m>0, k_p>0 e k_d>0 producono poli con parte reale negativa. Per m=1 kg, k_p=100 N/m e k_d=10 N·s/m, il polinomio è s²+10s+100 e i poli sono −5 ± j√75 s⁻¹. La pulsazione naturale è 10 rad/s e il rapporto di smorzamento è 0,5. Il sistema converge con oscillazioni. Questo risultato riguarda una forza che cambia continuamente con lo stato, non una forza congelata fra due letture.

2. Integrare il plant con il comando mantenuto

Il controllore digitale legge q_k e v_k agli istanti t_k=kh, calcola u_k e mantiene questa forza per h secondi. Si chiama mantenimento di ordine zero, ZOH. Fra i campioni l’accelerazione è costante, quindi si integra senza approssimazioni numeriche: la velocità cresce linearmente e la posizione quadraticamente. Il termine h²u_k/(2m) è essenziale; eliminarlo cambia la mappa. La documentazione MathWorks distingue esplicitamente questa discretizzazione esatta per ingressi a gradini dalle approssimazioni alternative.

u_k = −k_p q_k − k_d v_k q_(k+1) = q_k + h v_k + h²u_k/(2m) v_(k+1) = v_k + h u_k/m x_(k+1) = A_h x_k, x_k = [q_k, v_k]ᵀ A_h = [1−h²k_p/(2m) h−h²k_d/(2m)] [−hk_p/m 1−hk_d/m ]

La matrice non ha tutte le entrate con la stessa unità: l’elemento in alto a destra moltiplica una velocità per produrre una posizione, quindi ha unità di tempo; quello in basso a sinistra ha unità s⁻¹. Gli autovalori sono invece adimensionali e descrivono l’evoluzione da un campione al successivo. Cambiare unità da metri a millimetri modifica la rappresentazione dello stato, non gli autovalori della trasformazione coerentemente riscalata.

3. Una condizione esatta, non una regola empirica

La stabilità asintotica del sistema lineare discreto richiede tutti gli autovalori strettamente dentro il cerchio unitario. Usiamo il criterio quadratico di Jury, riportato nel testo di Kamran Iqbal, e lo applichiamo al nostro esempio. Definiamo τ come traccia e δ come determinante di A_h. Il polinomio caratteristico è z²−τz+δ. Le tre disuguaglianze seguenti sono necessarie e sufficienti per coefficienti reali; l’uguaglianza non dà stabilità asintotica.

τ = 2−hk_d/m−h²k_p/(2m) δ = 1−hk_d/m+h²k_p/(2m) 1−τ+δ = h²k_p/m > 0 1+τ+δ = 4−2hk_d/m > 0 1−δ = hk_d/m−h²k_p/(2m) > 0 0 < h < min(2m/k_d, 2k_d/k_p)

Il primo vincolo è automaticamente soddisfatto per h positivo e k_p positivo. Il secondo limita il periodo rispetto al guadagno derivativo; il terzo rispetto al rapporto fra guadagni. Entrambi i limiti hanno unità di secondi. Con i parametri scelti valgono entrambi 0,2 s. Nel modello senza ritardo, quindi, la stabilità asintotica si ha per 0<h<0,2 s. Non stiamo raccomandando di comandare un robot a 5 Hz: abbiamo calcolato una frontiera per questa massa ideale e questi guadagni.

Una conseguenza meno intuitiva riguarda k_d: aumentarlo non estende indefinitamente il periodo ammissibile. Il limite 2k_d/k_p cresce, ma 2m/k_d diminuisce. A m e k_p fissati, il massimo della loro funzione minimo si trova quando coincidono, cioè k_d=√(mk_p). Qui vale 10 N·s/m. È un’ottimizzazione del margine di campionamento del modello, non del transitorio, del rumore o della precisione: criteri diversi possono preferire guadagni diversi.

4. La frontiera e il falso conforto di un autovalore

Esattamente a h=0,2 s, la matrice è [[−1,0],[−20,−1]]. Ha un autovalore doppio −1 ma non è diagonalizzabile. Scrivendo A=−I+N con N²=0, risulta A^k=(−1)^k(I−kN). Compare un termine che cresce linearmente con k. Un raggio spettrale pari a uno non basta dunque a dire che tutte le traiettorie restano limitate. Per q iniziale non nullo, la velocità alterna segno e cresce in modulo.

Vicino a questa frontiera, aritmetica floating point e autovalori quasi coincidenti rendono fuorviante un confronto automatico senza tolleranza. Nel run, il valore calcolato al confine è circa 1,000000024. Non lo interpretiamo come una nuova soglia fisica: il risultato analitico determina il caso limite. Lo script esclude proprio il confine dalla verifica binaria sulla griglia e confronta le altre 299 posizioni con la disuguaglianza derivata.

5. Perché Euler dà una risposta diversa

Euler in avanti sul sistema già chiuso nel continuo usa E_h=I+hA_c. Per il nostro esempio, i poli discreti sono 1−5h ± j√75h; il loro modulo al quadrato è 1−10h+100h². La condizione diventa 0<h<0,1 s, metà del limite ZOH. A h=0,15 s, Euler è instabile mentre la massa con forza mantenuta è stabile. Non è una contraddizione: l’integratore approssimato e il sistema fisico campionato hanno matrici diverse.

A_c = [0 1 ] [−k_p/m −k_d/m] E_h = I+hA_c exp(A_c h) ≠ A_h

Anche discretizzare esattamente l’intero anello continuo con exp(A_c h) risponde a un’altra domanda: restituisce campioni di una retroazione che continua ad aggiornarsi fra le letture. Poiché i poli continui sono stabili, i loro esponenziali restano nel cerchio unitario per qualunque h positivo. Concluderne che il controllore digitale sia stabile per ogni periodo sarebbe sbagliato. Occorre discretizzare il plant con il mantenimento dell’ingresso, poi chiudere la retroazione campionata.

6. Un campione di ritardo cambia l’ordine del modello

Supponiamo che il comando appena calcolato sia applicato soltanto al campione successivo. La forza nell’intervallo corrente è quindi u_(k−1). Lo stato deve ricordare questa forza: y_k=[q_k,v_k,u_(k−1)]ᵀ. Il modello diventa di ordine tre. Non basta sottrarre una costante dalla frequenza o riutilizzare la frontiera precedente: si deve analizzare la nuova matrice. Qui il ritardo è esattamente h, senza jitter né ritardi frazionari.

y_(k+1) = D_h y_k D_h = [1 h h²/(2m)] [0 1 h/m ] [−k_p −k_d 0 ]

Il confronto numerico mostra un effetto sostanziale: a h=0,1 s il raggio spettrale del modello senza ritardo è circa 0,707, mentre con il ritardo sale a circa 1,441. Una scelta stabile nel primo modello diventa instabile nel secondo. Il numero non è una misura del ritardo di un vero bus o di un runtime AI; è una sensibilità controllata ottenuta cambiando una sola ipotesi strutturale.

h (s)ZOH esattoEulerRitardo h
0.020.9055390.9165150.914089
0.10.7071071.0000001.440965
0.150.7905691.3228761.848672
0.180.9055391.5620502.082966
0.21.0000001.7320512.236068
0.222.0143441.9078782.387195

7. Esperimento, transitori e riproducibilità

Il codice usa Python 3.14.0 e NumPy 2.5.3, senza casualità. La griglia contiene 300 periodi da 0,001 a 0,300 s. Per le traiettorie abbiamo q_0=0,01 m e v_0=0, simulazione fino all’ultimo campione entro tre secondi, e h in {0,02; 0,15; 0,22} s. Non viene integrata un’ODE con un passo nascosto: ogni aggiornamento usa direttamente la mappa esatta. Le figure sono calcolate con Matplotlib 3.11.2.

Sopra: raggio spettrale rispetto al periodo; la linea a uno delimita la stabilità asintotica. Sotto: modulo della posizione ai campioni, scala logaritmica; valori sotto 10⁻¹² m sono limitati graficamente. Solo simulazioni del modello ideale.
Sopra: raggio spettrale rispetto al periodo; la linea a uno delimita la stabilità asintotica. Sotto: modulo della posizione ai campioni, scala logaritmica; valori sotto 10⁻¹² m sono limitati graficamente. Solo simulazioni del modello ideale.

Il caso instabile può raggiungere spostamenti enormi perché il modello permette forze e corse illimitate. Quei valori non sono una previsione di movimento reale: indicano che la ricorrenza diverge e che le ipotesi locali cesserebbero presto di rappresentare una macchina. Nei casi stabili possono comunque comparire transitori importanti. Il raggio spettrale decide l’esito asintotico, non il massimo errore né il massimo sforzo. Anche confrontare raggi per periodi diversi richiede attenzione: un campione non rappresenta la stessa durata fisica.

import numpy as np
m, kp, kd = 1., 100., 10.
for h in [.02, .1, .15, .18, .2, .22]:
    A = np.array([[1-h*h*kp/(2*m), h-h*h*kd/(2*m)],
                  [-h*kp/m, 1-h*kd/m]])
    print(h, max(abs(np.linalg.eigvals(A))))

8. Che cosa cambia quando entra l’AI

Un modello visivo che aggiorna il riferimento a frequenza bassa non coincide con un modello che chiude direttamente l’anello di forza. Nel primo caso si possono separare una pianificazione lenta e una regolazione veloce; nel secondo la latenza d’inferenza diventa parte della dinamica del controllo. Per interpretare un esperimento bisogna indicare dove si trova il modello AI, quando il suo output diventa disponibile e quale comando resta attivo durante l’attesa. Il solo tempo medio di inferenza non descrive questa sequenza.

Restano da analizzare rumore di velocità, filtri, jitter, perdita di campioni, variazioni di massa, saturazioni e contatto. Per esempio, filtrare la velocità riduce certe componenti del rumore ma introduce stati e ritardo: non si può applicare automaticamente la formula a due stati. Analogamente, un sistema con h variabile è descritto da prodotti di matrici diverse; controllare ogni matrice isolatamente non dimostra in generale stabilità del prodotto. Questi sono problemi distinti da affrontare con modelli ed esperimenti dedicati.

La robotica industriale e collaborativa è una direzione che EL-AI intende esplorare. Questo articolo contribuisce al metodo di analisi, senza attribuire all’azienda robot installati o risultati hardware. La conclusione è operativa: scrivere la sequenza misura-calcolo-applicazione, ricavare la mappa che la rappresenta e verificare stabilità e transitori di quella mappa. Una simulazione convincente del sistema sbagliato resta una risposta alla domanda sbagliata.

Fonti, codice e limiti editoriali

Kamran Iqbal, Stability of Sampled-Data Systems (2023). MathWorks, Continuous-Discrete Conversion Methods.

Fonti consultate il 24 settembre 2026. Derivazione ed esempio elaborati per questo articolo con assistenza AI; nessuna peer review o certificazione dichiarata. Codice, risultati e istruzioni. Dati JSON. Copertina ImageGen illustrativa: non è una macchina o sede EL-AI reale.