come valutare la derivata della funzione in MATLAB?

Question

Questo dovrebbe essere molto semplice. Ho una funzione f(x) e voglio valutare f'(x) per un dato x in MATLAB.

Tutte le mie ricerche sono arrivate con matematica simbolica, che non è quello che mi serve, ho bisogno di differenziazione numerica.

E.g. se io definisco: fx = inline('x.^2')

voglio trovare dire f'(3), che sarebbe 6, non voglio trovare 2x

Answer 1

Per ottenere una differenza numerica (differenza simmetrica), si calcola (f(x+dx)-f(x-dx))/(2*dx)

fx = @(x)x.^2; 
fPrimeAt3 = (fx(3.1)-fx(2.9))/0.2; 

In alternativa, è possibile creare un vettore di valori della funzione e applicare DIFF, cioè

xValues = 2:0.1:4; 
fValues = fx(xValues); 
df = diff(fValues)./0.1; 

noti che diff prende la differenza in avanti e che presuppone che dx sia uguale a 1.

Tuttavia, nel caso, potrebbe essere meglio definire fx come polynomial e valutare la derivata della funzione, piuttosto che i valori della funzione.

Answer 2

hai provato diff (calcola le differenze e si avvicina un derivato), gradient, o polyder (calcola la derivata di un polinomio) funzioni?

È possibile leggere ulteriori informazioni su queste funzioni utilizzando help <commandname> nella console MATLAB oppure utilizzare il browser delle funzioni nel menu?.

Answer 3

Per una data funzione in forma analitica, è possibile valutare la derivata in un punto desiderato con il seguente codice:

syms x 
df = diff(x^2); 
df3 = subs(df, 'x', 3); 
fprintf('f''(3)=%f\n', df3); 

Per derivate numeriche puri utilizzare le soluzioni già fornite da Jonas e posdef.

Answer 4

Manca la casella degli strumenti simbolici, nulla ti impedisce di utilizzare Derivest, uno strumento per la differenziazione numerica adattiva automatica.

derivest(@sin,pi) 
ans = 
      -1 

Per il tuo esempio lo fa molto bene. In effetti, fornisce persino una stima dell'errore nell'approssimazione risultante.

fx = inline('x.^2'); 
[fp,errest] = derivest(fx,3) 

fp = 
      6 
errest = 
    3.6308e-14 

Answer 5

Se la funzione è noto per essere due volte differenziabile, utilizzare

f'(x) = (f(x + h) - f(x - h))/2h 

che è secondo ordine preciso nella h. Se è solo una volta differenziabile, utilizzare

f'(x) = (f(x + h) - f(x))/h  (*) 

che è il primo ordine in h.

Questa è teoria. In pratica, le cose sono abbastanza complicate. Prenderò la seconda formula (primo ordine) poiché l'analisi è più semplice. Fai il secondo ordine come esercizio.

La prima osservazione è che tu devi assicurarti che (x + h) - x = h, altrimenti ottieni enormi errori. Infatti, f (x + h) ep (x) sono vicini tra loro (per esempio 2.0456 e 2.0467), e quando li sottraggono, perdi molte cifre significative (qui è 0.0011, che ha 3 cifre significative in meno di x). Quindi qualsiasi errore su h ha probabilmente un impatto enorme sul risultato.

Quindi, primo passo, fissa un candidato h (ti mostrerò in un minuto come sceglierlo) e prendi come h per il tuo calcolo la quantità h '= (x + h) - x. Se stai usando un linguaggio come C, devi fare attenzione a definire h o x come volatile per quel calcolo che non deve essere ottimizzato.

Successivamente, la scelta di h. L'errore in (*) ha due parti: l'errore di troncamento e l'errore di arrotondamento. L'errore di troncamento è perché la formula non è esatto:

(f(x + h) - f(x))/h = f'(x) + e1(h) 

dove e1(h) = h/2 * sup_{x in [0,h]} |f''(x)|.

L'errore di arrotondamento deriva dal fatto che f (x + h) ef (x) sono vicini l'uno all'altro. Può essere stimato approssimativamente come

e2(h) ~ epsilon_f |f(x)/h| 

dove epsilon_f è la relativa precisione nel calcolo di f (x) (o f (x + h), che è vicino). Questo deve essere valutato dal tuo problema. Per le funzioni semplici, epsilon_f può essere preso come epsilon della macchina. Per quelli più complicati, può essere peggio di quello per ordine di grandezza.

Quindi si desidera h che riduce al minimo e1(h) + e2(h). Collegare tutto insieme e l'ottimizzazione in h produce

h ~ sqrt(2 * epsilon_f * f/f'') 

che deve essere definita sulla base vostra funzione. Puoi fare delle stime approssimative. In caso di dubbio, prendi h ~ sqrt (epsilon) dove epsilon = precisione della macchina. Per la scelta ottimale di h, l'accuratezza relativa a cui è nota la derivata è sqrt (epsilon_f), es. la metà delle cifre significative è corretta.

In breve: troppo piccolo a h => errore di arrotondamento, troppo grande a h => errore di troncamento.

Per la seconda formula ordine, stesse rese calcolo

h ~ (6 * epsilon_f/f''')^(1/3) 

e una precisione frazionaria (epsilon_f)^(2/3) per il derivato (che è tipicamente uno o due cifre significative meglio della prima formula ordine, assumendo una doppia precisione).

Se questo è troppo impreciso, non esitate a chiedere altri metodi, ci sono un sacco di trucchi per ottenere una migliore precisione. L'estrapolazione di Richardson è un buon inizio per le funzioni regolari. Ma quei metodi tipicamente computano f un bel po 'di volte, questo può essere o non essere quello che vuoi se la tua funzione è complessa.

Se si intende utilizzare le derivate numeriche molte volte in punti diversi, diventa interessante costruire un'approssimazione di Chebyshev.