Un’introduzione al metodo degli elementi finiti

In questo e nei prossimi quaderni di progettazione andremo ad introdurre il metodo degli elementi finiti. Questo è probabilmente l’approccio numerico più utilizzato per la risoluzione per via numerica dei problemi strutturali in ambito ingegneristico. Attraverso un’introduzione teorica e degli esempi pratici, nei prossimi mesi andremo ad analizzare nel dettaglio la formulazione di varie tipologie di elemento finito. Mostreremo, per ciascuno, le limitazioni ed i benefici attraverso esempi applicativi.

di Franco Concli

Il metodo agli elementi finiti (FEM) è un approccio numerico utilizzato per risolvere le equazioni differenziali che governano pressoché ogni fenomeno fisico. Per questo, negli anni, grazie anche ad una costante crescita della potenza di calcolo, il metodo ha trovato una sempre più larga diffusione in molte discipline ingegneristiche. Le applicazioni più comuni includono l’analisi strutturale, lo scambio termico, la meccanica dei fluidi ed i campi elettromagnetici.

Dai sistemi complessi agli elementi finiti

Il metodo agli elementi finiti altro non è che un approccio sistematico per scomporre problemi complessi in “sottosistemi” più piccoli e gestibili. Questi sono appunto i cosiddetti elementi finiti. Questa scomposizione/semplificazione è possibile attraverso un processo noto come discretizzazione spaziale. Il problema o dominio di interesse viene suddiviso in elementi (segmementi, superfici o volumi elementari) che costituiscono la così detta mesh (o griglia di calcolo). In questo modo, un problema anche complicato viene scomposto in elementi geometricamente semplici come mattoncini di Lego. La soluzione viene quindi calcolata per ogni elementino.

La formulazione di un problema mediante il metodo agli elementi finiti si traduce quindi in un sistema di equazioni algebriche. L’unione delle soluzioni calcolate per ogni elemento finito, nel suo insieme, altro non è che la soluzione approssimata del problema ingegneristico complesso di partenza.

In alcuni casi, strutture o sistemi possono essere risolti per via analitica, soprattutto quando isostatici. Tuttavia, tali scenari sono rari e spesso lo l’onere richiesto per risolvere i problemi per via analitica è particolarmente gravoso. Si pensi ad una struttura come la Torre Eiffel. Qualcuno, prima dell’avvento dei computer, ha progettato tale struttura – con le sue 18.038 parti unite da circa 2.500.000 rivetti – completamente per via analitica!

Ma queste sono imprese di altri tempi che oggi non hanno più senso. Il lavoro che, pressoché 15 decadi orsono, ha richiesto anni e squadre di ingegneri, oggi può essere gestito da un piccolo team in tempi molto più ristretti.

Il metodo agli elementi finiti

Il metodo agli elementi finiti è una tecnica di approssimazione e soluzione numerica di funzioni incognite. I campi di applicazione più comuni sono quello strutturale e quelle termico anche se altri domini come i fluidi o i campi elettromagnetici possono anch’essi essere risolti mediante un approccio FEM. Esistono varie tipologie di elementi finiti, ma alla base di ognuno troviamo semplici equazioni algebriche. Quando un problema complesso viene discretizzato in elementi finiti, le equazioni algebriche di ogni singolo elemento vengono “assemblate” in un sistema di equazioni più ampio rappresentativo dell’intero problema ingegneristico. Nello studio o analisi di un fenomeno mediante approccio agli elementi finiti, la dimensione della griglia di calcolo gioca un ruolo fondamentale:

• Se questa è troppo grossolana, la geometria del sistema non viene rappresentata correttamente. Inoltre, il problema viene risolto in un numero di punti insufficiente a rappresentare correttamente le variabili di campo.
• Se la mesh è troppo fine, la soluzione del problema richiederà invece troppo tempo.

Come vedremo in seguito, una prima classificazione degli elementi finiti può essere fatta a secondo dalla loro natura spaziale: elementi 1D, elementi 2D ed elementi 3D. In molti casi, è possibile arrivare ad una buona rappresentazione del problema ingegneristico con elementi finiti “semplici” che portino a risultati accurati in tempi di calcolo relativamente ridotti.

Inoltre, è interessante notare come quando si crea un modello FEM, all’atto pratico è possibile visualizzare “gli elementi finiti” sullo schermo. Tuttavia, va sottolineato come tali elementi non siano “reali”: il vero fondamento del metodo risiede non tanto negli elementi ma nei nodi (vertici degli elementi), che rappresentano i punti in cui la soluzione viene di fatto calcolata.

Gli elementi 1D

Nell’analisi agli elementi finiti si possono distinguere due tipi di elementi monodimensionali: gli elementi asta e gli elementi trave. La differenza principale riguarda il modo in cui le forze e le deformazioni vengono modellate. Gli elementi asta possono venire utilizzati per modellare componenti capaci di resistere solamente a carichi assiali (trazione o compressione). Sono, ad esempio, utilizzati per modellare le strutture reticolari.

Oltre agli elementi asta esistono i cosiddetti elementi trave che, a loro volta, possono essere classificati in elementi trave piani (capaci di resistere solo a flessione nel piano e a taglio trasversale), elementi trave 2D (capaci di resistere anche a forze assiali) ed elementi trave 3D (capaci di resistere a tutte le componenti di forza e momenti nodali).

Elementi asta

Consideriamo una barra elastica a sezione uniforme di lunghezza L, con modulo elastico E e sezione A. Questo avrà una rigidezza assiale pari a EA/L.

Un elemento finito che possa modellare tale semplice sistema può essere, nella sua versione più semplice, un segmento che congiunga i due nodi alle estremità. I gradi di libertà (GDL) di questo elemento finito “asta” sono solo due: gli spostamenti assiali dei due nodi di estremità (ovvero u₁ e u₂).

Per imporre uno spostamento d ad uno dei nodi e deformare la trave devono essere applicate delle forze su entrambi i nodi in modo da garantire l’equilibrio dell’elemento. In particolare, dovrebbe essere applicata una forza su uno dei due nodi in modo da indurne lo spostamento mentre sull’altro, la cui traslazione è vincolata, si dovrebbe andare a generare una reazione vincolare. Le forze corrispondenti sui due nodi dovranno essere contrapposte in modo da garantire l’equilibrio dell’elemento. Tali forze possono essere calcolate partendo dalla definizione di rigidezza di una asta lineare.

in cui F_ij è la forza associata al nodo i-esimo che si ha a seguito di uno spostamento del nodo j-esimo. Per garnatire l’equilibrio, le forze esterne F₁ ed F₂ dovranno risultare pari a F₁ = F₁₁ – F₁₂ e F₂ = F₂₁ – F₂₂.

Questi due equilibri, possono essere scritti in forma matricale (tenendo conto della rigidezza dell’elemento) come:

k prende il nome di matrice di rigidezza, d è il vettore degli spostamenti nodali mentre r rappresenta il vettore dei carichi esterni.

Una colonna di k rappresenta i carichi nodali che devono essere applicati all’elemento per sostenere uno stato di deformazione in cui il GDL corrispondente assuma un valore unitario mentre tutti gli altri GDL nodali risultino nulli.

In questo caso la soluzione è semplice poiché la relazione tra carichi, spostamenti e rigidezza di un asta è nota dalla teoria di de Saint Venant. Questo approccio per la scrittura delle matrici di rigidezza dell’elemento è detto Metodo Diretto. Tuttavia, per la maggior parte degli elementi finiti è necessario utilizzare un approccio più generale che non richieda la conoscenza delle proprietà dell’elemento (Procedura Formale).

Ciò è possibile alla luce della considerazione per cui il lavoro svolto dai carichi nodali sia immagazzinato nell’elemento come energia di deformazione elastica.

in cui B è la matrice che lega deformazione e spostamento, E è la matrice delle proprietà del materiale (detta anche matrice costitutiva) e ???? è il volume dell’elemento. Per ottenere B per l’elemento asta, lo spostamento assiale ???? di un punto arbitrario all’interndo dell’asta è ottenuto mediante un’interpolazione lineare tra i valori nodali ????₁ e ????₂. Ciò può essere scritto come:

in cui N è detta matrice delle funzioni di forma. Le funzioni di forma N_i permettono di descrivere la variazione di u con x all’interno dell’elemento in funzione degli spostamenti nodali u_i.

In questo caso le funzioni interpolanti sono lineari. Definito il campo di spostamento u per ogni punto all’interno dell’elemento, è possibile andare a calcolare la deformazione assiale ε come:

Nel caso dell’elemento asta, la matrice B risulta:

I campi di spostamento e deformazione evidenziano la limitazione di un elemento asta a due nodi: può rappresentare solo uno stato di deformazione costante. Inoltre, tale elemento non risulta adatto a rappresentare sistemi in cui le forze (assiali) siano distribuite (anziché applicate direttamente ai nodi), o aste a sezione variabile.

Un approccio pratico consiste nel modellare l’asta reale mediante l’impiego di più elementi finiti e applicare il carico distribuito ai vari nodi. In questo modo è possibile arrivare ad una soluzione approssimata la cui accuratezza cresce all’aumentare del numero di elementi finiti.

Elementi trave

Per determinare la matrice di rigidezza di un elemento trave 2D, si consideri innanzitutto l’elemento trave piano i cui gradi di libertà sono costituiti dai soli spostamenti verticali dei due nodi (v₁ e v₂) e dalle rotazioni attorno all’asse z ( q_z1 e q_z2 ). La Figura mostra l’applicazione di uno spostamento o rotazione unitario e le rigidezze k associate.

La matrice di rigidezza può essere costruita colonna per colonna secondo

la regola generale vista nel paragrafo precedente. Per ottenere i termini di una colonna, è necessario risolvere un problema staticamente indeterminato (struttura iperstatica). Ciò può essere facilmente fatto sfruttando il Principio dei Lavori Virtuali.

In alternativa, è possibile procedere, come nel caso precedente, con un procedimento formale.

Lo spostamento ???? di un elemento trave può essere espresso come un polinomio del terzo ordine in ???? (l’elemento ha 4 gradi di libertà, quindi 4 costanti)_.

La derivata rispetto a ???? rappresenta la rotazione

Ogni costante ????_???? può essere espressa in termini di gradi di libertà nodali effettuando le opportune sostituzioni.

Ad esempio, con riferimento al primo grado di libertà ????₁ = 1 è possibile scrivere

da cui:

Ripentendo l’operazione per tutti i gradi di libertà si arriva a scrivere

La curvatura della trave è la derivata seconda di ???? rispetto ad x.

Il momento flettente e le sollecitazioni flessionali sono calcolati dalla curvatura, che a sua volta è funzione dei gradi di libertà nodale:

Il momento flessionale varia, dunque, linearmente con x in ciascun elemento. Un modello agli elementi finiti con elementi trave fornisce una soluzione esatta quando la forza e/o il momento sono applicati ai nodi. Un carico uniformemente distribuito, invece, produrrebbe una deflessione di quarto grado in x, mentre la funzione di forma è solamente di terzo ordine. Nel caso di elementi trave piani, un carico distribuito può essere applicato, come forze equivalenti applicate ai nodi.

La combinazione delle matrici di rigidezza dell’elemento asta (GDL u) con quella dell’elemento trave piana (GDL v e θ_z) permette di ottenere la matrice di rigidezza dell’elemento trave 2D.

Conclusioni

In questo primo Quaderno di Progettazione dedicato al FEM abbiamo presentato il metodo degli elementi finiti applicato a problemi monodimensionali, mettendo in evidenza le fasi fondamentali del processo.
L’approccio 1D, sebbene concettualmente semplice, rappresenta un passaggio essenziale per comprendere i principi alla base del metodo e per sviluppare familiarità con gli strumenti matematici e computazionali utilizzati.

L’uso di elementi lineari, funzioni di forma e formulazione debole forniscono una base solida per affrontare sistemi più complessi in due o tre dimensioni. Inoltre, l’analisi 1D è spesso sufficiente per risolvere problemi reali in cui una sola dimensione domina il comportamento fisico (come barre, travi o condotti), semplificando molto la vita degli ingegneri…

Il prossimo passo consisterà nell’estendere questi concetti a geometrie e carichi più complessi, affrontando il metodo degli elementi finiti in dimensioni superiori.