Premessa

Il corso "Analisi matematica attraverso esercizi" recentemente inserito nel sito di Electroportal, illustra molto bene e con appropriati esempi, ciò che sull'argomento viene insegnato nelle scuole italiane e, naturalmente, non fa alcun cenno alle possibilità offerte dall'impiego di calcolatori elettronici per superare le difficoltà che normalmente si incontrano nella pratica applicazione della matematica ai reali problemi industriali. E' questo uno degli argomenti che da tempo cerco di proporre (finora senza successo) agli studenti desiderosi non solo di superare gli esami, ma anche di formarsi una base di conoscenza delle tecniche di calcolo fondamentali per l'impiego dei calcolatori  nella progettazione (Computer Aided Design). D'altra parte, da un'indagine svolta in Internet, mi rendo conto che mancano pubblicazioni adatte ad illustrare elementarmente questi principi, di per sé semplici e che risalgono alle origini del calcolo infinitesimale. Le basi del calcolo numerico furono infatti trattate da Newton, Eulero, Boole e molti altri negli ultimi secoli in modo puramente speculativo, ed applicate in pratica solo da pochi decenni, dopo l'avvento del calcolatore elettronico che permette il rapido svolgimento della miriade di calcoli necessari alla loro implementazione. Scopo di questo articolo e' dunque un'introduzione al metodo delle differenze finite, esemplificandone l'applicazione ai casi principali, in alternativa alle soluzioni "simboliche" dell'analisi matematica. Segnalo che gli esempi qui forniti si basano sull'impiego di MathCad, ma che possono essere svolti in qualsiasi software matematico in grado di elaborare variabili indicizzate.

Il concetto delle differenze finite

Nell'analisi matematica il concetto di derivata viene introdotto come limite della differenza fra due punti della funzione, divisa per il loro intervallo, quando questo tende a zero:

Data infatti una qualsiasi funzione di x, calcolata in un punto (x₀) ed in secondo punto (x₀ + l'intervallo h), è possibile ricavare la pendenza della retta che passa per i due punti. E' chiaro che se h tende a 0, i due punti tendono a coincidere e la retta a diventare la tangente della funzione stessa nel punto x₀. E' quindi questo "passaggio al limite" l'essenza del calcolo infinitesimale, così fecondo per tutti gli sviluppi che ne sono seguiti. Ma se ci fermiamo un momento prima dello zero, cioè diamo ad h un valore piccolo ma finito, diciamo Dx, e sostituiamo y ad f(x), potremo applicare lo stesso criterio:

Cioè la derivata prima di y nel punto di ascissa x₀, è approssimativamente uguale alla differenza dei due valori di y, divisa per Dx (e sarà più approssimato, più  Dx tende a zero). Questo implica tuttavia la necessità di calcolare con grande precisione i valori di y, la cui differenza diventa ovviamente evanescente al diminuire di Dx. Per i moderni mezzi ci calcolo, questo non è però un problema: anche il più modesto PC offre precisioni tali da soddisfare le normali esigenze pratiche. Diventano così superflui i vari metodi più o meno elaborati che sono citati nella letteratura matematica per aumentare la precisione di calcolo.

Derivata in un punto

Per passare subito alla più semplice applicazione del metodo visto, possiamo applicarlo al calcolo della derivata prima in un punto dato di una funzione. Scegliamo per semplicità una parabola e valutiamone la pendenza nel punto x=1:

La fig. mostra la procedura classica,  ricavando prima la derivata della funzione (y'= -4x² +5) e ponendo in questa x=1. Ed ecco la soluzione con le differenze finite:

Un primo calcolo, con Dx  =  1/1000, dà un risultato approssimato al  2 per 1000, il secondo dà addirittura il risultato esatto (dy sta per y' e naturalmente ya è l'ordinata della funzione per x=1, mentre yb è quella per x= 1+Dx ). Non credo ci sia bisogno di altri commenti su questo esempio, salvo sottolineare che se per una parabola il metodo può apparire del tutto inutile (data la semplicità della soluzione convenzionale), non così sarebbe se la funzione fosse non immediatamente derivabile. La procedura alle differenze finite rimane infatti così semplice anche se le funzioni trattate sono estremamente complesse.

Grafico della derivata di una funzione

Applicando la procedura vista si può ottenere direttamente l'andamento in forma grafica della derivata prima di una funzione data. Si tratta di utilizzare variabili indicizzate (variabili con indice in apice) che consentono di accumulare serie di valori che le esprimono

Con N si fissa il numero di punti dell'indice n e con Dx  l'intervallo di "campionamento" della funzione. In questo modo si può far variare la variabile x fra due limiti stabiliti (in questo caso da -5 a +5). Si calcolano poi i valori di y corrispondenti alla funzione data (in questo caso una cubica), che può essere tradotta in grafico (linea blu). Per ottenere l'andamento della derivata, si applica la forma delle differenze finite (con l'unica avvertenza di considerare il punto attuale e quello precedente), ottenendo la funzione dy (cioè y', linea rossa). Anche in questo caso non ritengo si debba sottolineare che sono immediatamente rilevabili massimo, minimo e flesso, semplicemente osservando i grafici (e nelle applicazioni pratiche, questo significa un'immediatezza nel rilievo dei risultati). E' poi ovvia l'estensione del metodo alle derivate di ordine superiore: se al posto di y si utilizzano i campionamento di dy si ottiene l'andamento della derivata seconda (y''), e così via

Calcolo di aree

Il campionamento delle funzioni permette anche l'applicazione al calcolo di integrali definiti, geometricamente rappresentanti aree racchiuse da  queste. Come esempio si veda il corso citato (esercizio n.3):

L'integrazione dell'area (A) racchiusa fra le due parabole y1 ed y2, fra i limiti 0 e 8/3 (che rappresentano le ascisse di intersezione) avviene semplicemente "sommando" le ordinate moltiplicate per Dx.  Questo rappresenta infatti la somma di 10000 (N)  rettangoli di altezza y e larghezza  Dx  per entrambe le parabole. Come si può vedere il risultato è esatto alla terza cifra decimale rispetto al calcolo dell'integrale ricavato nell'esempio del corso.

Soluzioni grafiche di equazioni differenziali

L'applicazione forse più interessante delle differenza finite è quella della soluzione delle equazioni differenziali. Supponiamo infatti di conoscere la funzione della derivata prima (dy) e di voler ricavare da questa l'andamento grafico della funzione originale (y): con un riarrangiamento algebrico dell'espressione delle differenze finite si ottiene immediatamente la serie dei valori di y. Se infatti   

                    dy_n = (y_n - y_{n -1})/Dx     si ha:      y_n  =  y_{n -1} + Dx* dy_n

Vediamone subito un esempio, in cui x è sostituito dal tempo (t) (e quindi Dx da DT). La funzione derivata è una parabola, quindi la funzione originale è una cubica, osservata nell'intervallo 0-10  (N*DT), con 10000 (N) campionamenti, distanziati di 1/1000 (DT).

Si osservi che si è dovuto fissare l'origine (y₀=0, cioè la costante di integrazione) e che il confronto fra i valori della funzione calcolata dopo 10000 iterazioni e quella "esatta" dall'integrale, differiscono solo dell'1.5  per mille.

Possiamo quindi affrontare esempi più significativi, quale quello trattato in questo articolo (vedi: "esempio di equazione differenziale", con svolgimento simbolico). L'equivalente soluzione numerica è riportata nella figura seguente, in cui si può notare la "coincidenza" dei grafici ottenuti dalla procedura numerica (in rosso) e da quella simbolica (a punti blu). Si osservi che per aggirare la condizione iniziale di x=0 (che al denominatore porterebbe ad un valore infinito), si è posto x = ad 1 miliardesimo.

Come si vede, la procedura è rimasta elementare anche in presenza di un'equazione tutt'altro che semplice da risolvere.

Calcoli dei circuiti elettronici

Chiudo questa breve panoramica di applicazioni con un cenno alle possibilità offerte da questo metodo nei calcoli che illustrano i comportamenti di circuiti elettronici. Nella maggioranza dei casi, questi comportamenti sono determinati da equazioni differenziali difficilmente risolvibili con la matematica classica. Anche l'applicazione delle trasformate di Laplace, non sempre conduce a risultati pratici e comunque presenta, se non nei casi più elementari, notevoli complicazioni. La pratica più diffusa è oggi quella della simulazione con software dedicati ormai talmente perfezionati da rendere obsoleti i vari mezzi di elaborazione fin qui utilizzati. Metteremo quindi qui a confronto una simulazione eseguita con Micro-Cap 9 di un circuito con induttanza, resistenza e capacità, alimentato in onda quadra, con una procedura di soluzione alle differenze finite di un sistema di equazioni differenziali che definiscono appunto il comportamento di questo circuito. Iniziamo dalla simulazione:

V1 è un generatore di onda quadra a 10V, con periodo di 100ms, e si vuole conoscere la tensione di uscita V2, dati i valori dei 3 componenti. Micro-Cap è già così in grado di tracciare sia V1 (traccia blu), sia il risultato cercato V2 (traccia rossa), mostrando l'andamento oscillatorio dovuto agli effetti di L e C. L'approccio di calcolo, pur con l'impiego delle differenze finite, è naturalmente più complicato. Chiamando il la corrente che percorre l'induttanza, ic quella che percorre la capacità, ed ovviamente V1 e V2 le tensioni ai rispettivi punti, possiamo scrivere:

Che possono facilmente essere ricondotte nella forma di differenze finite (l'equazione contenente l'integrale e' stata derivata per questo). Occorre però osservare che le 3 equazioni, per essere calcolate simultaneamente, devono essere scritte in forma di matrice. Ecco l'implementazione in MathCad, con i necessari chiarimenti. N è il numero di campionamenti, ed n l'indice di questi. L'espressione di V1 tiene conto della possibilità di espandere il grafico a più onde quadre (struttura ricorsiva). Vengono poi impostati i valori di R,L,C, tenendo conto delle loro unità di misura. Tutti gli altri parametri del calcolo vengono poi inizializzati, e viene preparata la scale dei tempi (in ms): 10000 (N) campionamenti per 10us (DT) corrisponde infatti a 100ms. Solo alla fine viene scritta la matrice che riproduce il sistema di equazioni differenziali

Non resta ora che rappresentare in grafico il risultato del calcolo di V2 (insieme all'andamento della tensione di ingresso V1):

Che dimostrano la perfetta corrispondenza con il grafico ottenuto con Micro-Cap.