Capitolo 1
___________________________________________________ __________________________
1.1 Introduzione
In questo capitolo si introduce il metodo agli elementi finiti nell’ipotesi di materiali a
comportamento lineare.
La gran parte dei problemi di ingegneria civile non possiede
forma chiusa perchØ, a causa della loro complessità in termini di condizioni al contorno
e di geometria, risulta impossibile risolvere le equazioni di equilibrio
Per questo motivo, si è pensato di adottare dei met
come il Metodo agli Elementi Finiti (FEM). Quest’ultimo permette
problemi ingegneristici, anche complessi, mediante
equazioni differenziali (alle derivate parziali o tota
discretizzazione della geometria del continuo attribuendo, ad esso, solo un numero
finito di gradi di libertà. Conseguentemente, è necessario discretizzare anche le suddette
equazioni differenziali (che, in questo modo, divengo
L’aspetto evidente è che, piø fitta è la discretizzazione operata sul continuo (e quindi
sulle equazioni che reggono il problema) maggiore è la precisione della soluzione
Figura 1 . geotecnico
___________________________________________________ __________________________
si introduce il metodo agli elementi finiti nell’ipotesi di materiali a
La gran parte dei problemi di ingegneria civile non possiede una soluzione analitica in
forma chiusa perchØ, a causa della loro complessità in termini di condizioni al contorno
e di geometria, risulta impossibile risolvere le equazioni di equilibrio
Per questo motivo, si è pensato di adottare dei metodi approssimati di tipo numerico,
come il Metodo agli Elementi Finiti (FEM). Quest’ultimo permette
, anche complessi, mediante la risoluzione numerica di sistemi di
equazioni differenziali (alle derivate parziali o totali). Il metodo si basa sulla
discretizzazione della geometria del continuo attribuendo, ad esso, solo un numero
finito di gradi di libertà. Conseguentemente, è necessario discretizzare anche le suddette
equazioni differenziali (che, in questo modo, divengono equazioni algebriche).
L’aspetto evidente è che, piø fitta è la discretizzazione operata sul continuo (e quindi
sulle equazioni che reggono il problema) maggiore è la precisione della soluzione
. 1. Condizioni al contorno in un problema
geotecnico
2
___________________________________________________ __________________________
si introduce il metodo agli elementi finiti nell’ipotesi di materiali a
soluzione analitica in
forma chiusa perchØ, a causa della loro complessità in termini di condizioni al contorno
e di geometria, risulta impossibile risolvere le equazioni di equilibrio (vedi figura 1.1).
odi approssimati di tipo numerico,
come il Metodo agli Elementi Finiti (FEM). Quest’ultimo permette lo studio di
la risoluzione numerica di sistemi di
Il metodo si basa sulla
discretizzazione della geometria del continuo attribuendo, ad esso, solo un numero
finito di gradi di libertà. Conseguentemente, è necessario discretizzare anche le suddette
no equazioni algebriche).
L’aspetto evidente è che, piø fitta è la discretizzazione operata sul continuo (e quindi
sulle equazioni che reggono il problema) maggiore è la precisione della soluzione
Condizioni al contorno in un problema
Capitolo 1 3
___________________________________________________ __________________________
numerica; è ovvio che, a questo, corrisponde una maggiore onere in termini di tempo
impiegato per la risoluzione del problema.
1.2 Variabili primarie del problema
In un qualsiasi problema numerico, è necessario stabilire quali siano le “variabili
primarie” rispetto alle quali risolvere il problema (ad esempio, nei problemi di
ingegneria Geotecnica, in genere vengono assunte come variabili primarie gli
spostamento, le pressioni interstiziali, … vedi figura 1.2). In un secondo momento è,
poi, necessario definire in che modo queste variabili si evolvono all’interno del nostro
dominio geometrico; questa necessità equivale ad assegnare una legge di interpolazione
che definisca la loro variazione a partire dal valore noto in corrispondenza dei nodi.
Tale legge di interpolazione viene denominata Funzione di forma.
La trattazione matematica del metodo FEM nel caso di un generico problema
geotecnico è complessa, in quanto è necessario risolvere un sistema di equazioni
differenziali alle derivate parziali le cui incognite sono gli spostamenti nei nodi “solidi”
e le pressioni interstiziali nei nodi “fluidi” della “mesh”. Il numero di variabili per ogni
nodo definisce il numero di gradi di libertà per nodo (i nodi comuni a due elementi
devono avere gli stessi gradi di libertà e le loro funzioni di forma devono essere
Figura 1.2. Variabili primarie per un
problema geotecnico
Capitolo 1
___________________________________________________ __________________________
caratterizzate da polinomi dello stesso ordine
libertà (i due spostamenti) e per la parte fluida h
interstiziale).
1.3 Discretizzazione in elementi finiti
La discretizzazione dello spazio consiste nel dividere il dominio (continuo) in un
numero finito di elementi connessi tra loro attraverso dei nodi disposti in
corrispondenza degli angoli e lungo i lati degli elementi in funzione del tipo di elemento
che si utilizza. La geometria viene schematizzata attraverso una
composta da un numero finito di elementi
del valore finale delle incognite nodali. In geotecnica si usano due
una per la parte solida e una per la parte fluida. La forma dei vari elementi può essere
mono-, bi- e tridimensionale ed il numero di nodi solidi
dell’elemento (vedi figura
forma triangolare o quadrangolare.
Figura 1.3. Geometria degli
elementi per la discretizzazione
___________________________________________________ __________________________
caratterizzate da polinomi dello stesso ordine); il nodo per la parte solida ha due gradi di
libertà (i due spostamenti) e per la parte fluida ha un grado di libertà (la pressione
retizzazione in elementi finiti
La discretizzazione dello spazio consiste nel dividere il dominio (continuo) in un
numero finito di elementi connessi tra loro attraverso dei nodi disposti in
corrispondenza degli angoli e lungo i lati degli elementi in funzione del tipo di elemento
utilizza. La geometria viene schematizzata attraverso una
composta da un numero finito di elementi, tale da assicurare una sufficiente precisione
del valore finale delle incognite nodali. In geotecnica si usano due
la parte solida e una per la parte fluida. La forma dei vari elementi può essere
e tridimensionale ed il numero di nodi solidi dipende dall’ordine
(vedi figura 1.3). Per problemi bidimensionali usualmente l’elemento è di
forma triangolare o quadrangolare.
Geometria degli
elementi per la discretizzazione
del modello 4
___________________________________________________ __________________________
; il nodo per la parte solida ha due gradi di
a un grado di libertà (la pressione
La discretizzazione dello spazio consiste nel dividere il dominio (continuo) in un
numero finito di elementi connessi tra loro attraverso dei nodi disposti in
corrispondenza degli angoli e lungo i lati degli elementi in funzione del tipo di elemento
utilizza. La geometria viene schematizzata attraverso una mesh (maglia)
tale da assicurare una sufficiente precisione
del valore finale delle incognite nodali. In geotecnica si usano due mesh sovrapposte:
la parte solida e una per la parte fluida. La forma dei vari elementi può essere
dipende dall’ordine
Per problemi bidimensionali usualmente l’elemento è di
Capitolo 1 5
___________________________________________________ __________________________
Il passo successivo è quello di attribuire, ai nodi dei vari elementi, un numero in
maniera sequenziale. Nel momento in cui si deve definire la maglia di elementi finiti,
spesso, un aiuto viene fornito dalla presenza di discontinuità, naturali o artificiali. ¨
evidente che il tipo di mesh, così come le condizioni al contorno, influenzano il
risultato.
La forma ed il numero di elementi dipende dal comportamento del materiale. Quando il
materiale ha un comportamento lineare, la soluzione è immediata perchØ la matrice
costitutiva non varia durante l’analisi. Piø complessa è, invece, la situazione quando il
materiale ha un comportamento non lineare e la soluzione dipende dalla storia del
carico, in questi casi la matrice di rigidezza viene aggiornata durante l’analisi. Un
aspetto importante, che bisogna considerare per migliorare la precisione dell’analisi, è
relativo alla dimensione relativa tra i vari elementi della mesh che deve variare
progressivamente in modo da garantire una transizione graduale dalle zone con elementi
piø piccoli a quelle con elementi piø grandi. ¨ opportuno non usare elementi distorti
(ovvero con rapporto tra la dimensione maggiore e quella minore superiore a 2) in modo
da non rendere singolare la matrice di rigidezza dell’elemento e causare instabilità
numerica.
In zone nelle quali sono presenti discontinuità geometriche e di materiale, in zone di
applicazione del carico (come anche in quelle in cui di prevedono notevoli gradienti di
stato tensionale o di spostamento), la mesh deve essere “raffinata” ovvero costituita da
elementi piø piccoli (vedi figura 1.4).
Figura 1.4. Schematizzazione della geometria della mesh [1]
Capitolo 1 6
___________________________________________________ __________________________
1.4 Equazioni di equilibrio, locale e globale
L’equilibrio del singolo elemento si ottiene risolvendo numericamente un sistema
matriciale del tipo:
nullnull null null null ∆null null null null null ∆null null null (1.1)
dove:
• nullnull null null rappresenta la matrice di rigidezza del singolo elemento, funzione della legge
costitutiva scelta;
• null ∆null null null è il vettore delle incognite nodali incrementali del singolo elemento (es.:
spostamenti o pressioni interstiziali);
• null ∆null null null è il vettore delle forze nodali incrementali del singolo elemento.
Le matrici di rigidezza dei singoli elementi vengono espresse in termini di tutte le
variabili nodali dell’intera mesh e assemblate nella matrice di rigidezza globale del
sistema. Allo stesso modo, il vettore delle forze nodali del sistema è ottenuto sommando
le singole forze agenti in ogni nodo della mesh. Combinando le equazioni dei singoli
elementi si ottiene un sistema matriciale del tipo:
nullnull null null null ∆null null null null null ∆null null null (1.2)
dove:
• nullnull null null rappresenta la matrice di rigidezza del sistema globale, ottenuta assemblando le
matrici di rigidezze dei singoli elementi;
• null ∆null null null è il vettore delle incognite nodali incrementali del sistema globale (es.:
spostamenti o pressioni interstiziali);
• null ∆null null null è il vettore delle forze nodali incrementali del sistema globale.
Risolvendo il quale, si ottiene il vettore delle incognite globali null ∆null null null . I termini della
matrice di rigidezza si ottengono considerando i nodi in comune dei singoli elementi. I
termini del vettore delle forze sono ottenuti considerando i singoli carichi applicati in
corrispondenza dei nodi. Nel momento in cui si va a scrivere la matrice di rigidezza, per
i nodi in comune le singole rigidezze locali vengono sommate. Se la matrice di
Capitolo 1 7
___________________________________________________ __________________________
rigidezza del singolo elemento è simmetrica anche la matrice globale è simmetrica (vedi
figura 1.5).
Indicato con m il numero totale di gradi di libertà del sistema, la matrice di rigidezza
globale avrà dimensioni mxm. Se la mesh è costituita per es. da 5000 nodi, se ogni nodo
possiede due gradi di libertà, la matrice globale sarà costituita da 10000x10000
Figura 1.5. Procedura di assemblaggio della matrice di rigidezz a
globale per una mesh semplice a due elementi [1]
Capitolo 1 8
___________________________________________________ __________________________
coefficienti (10
8
). Soltanto pochi coefficienti sono non nulli e, di solito, sono disposti
attorno alla diagonale principale della matrice. Dunque, i termini della matrice di
rigidezza crescono al crescere del numero di elementi della mesh.
1.5 Valutazione degli spostamenti
Di solito nel metodo degli elementi finiti, lo spostamento è la quantità incognita. Forze
e tensioni sono le quantità secondarie che possono essere ricavate una volta ottenuti gli
spostamenti. La deformata dell’elemento è legata al tipo di funzione di forma scelta.
Nello spazio delle due dimensioni (condizioni Piane), individuato da un sistema di
coordinate x e y, il comportamento deformativo è caratterizzato da due spostamenti
globali u e v. Negli elementi lineari le funzioni di forma sono polinomi del primo
ordine, ovvero il valore delle incognite primarie varia in maniera lineare tra due nodi
consecutivi. Se, per esempio, si considera un elemento con tre nodi:
noto il valore degli spostamenti in corrispondenza dei nodi, null null, null null , nullnull , nullnull, nullnull null nullnull ,
ottengo un sistema di sei equazioni in sei incognite ed è possibile determinare le
costanti null null , null null , nullnull , nullnull , nullnull null nullnull :
nullnullnull null nullnull null nullnullnull null null nullnullnull null nullnull null nullnullnull null null (1.3)
Figura 1.6. Elemento a tre
nodi [1]
Capitolo 1 9
___________________________________________________ __________________________
Quindi note le coordinate di un punto all’interno dell’elemento, è possibile determinare
gli spostamenti di un generico punto:
null null null nullnull null nullnull null null null nullnull null null null null null nullnull null nullnull null null null nullnull null null null null null nullnull null nullnull null null null nullnull null null null null (1.4)
null null null nullnull null nullnull null null null nullnull null null null null null nullnull null nullnull null null null nullnull null null null null null nullnull null nullnull null null null nullnull null null null null (1.5)
In generale è possibile scrivere le incognite in un punto dell’elemento in funzione del
valore in corrispondenza dei nodi:
null null null nullnull nullnull nullnull null null null null null null
null null null null null null null null null nullnull null null null null null nullnullnullnull (1.6)
dove nullnull null rappresenta la matrice delle funzioni di forma; le componenti dello
spostamento u e v sono espresse in termini del loro valore in corrispondenza dei nodi.
Negli elementi quadratici le funzioni di forma sono polinomi del secondo ordine,
mentre negli elementi cubici le funzioni di forma sono polinomi del terzo ordine.
L’accuratezza di un’analisi agli elementi finiti dipende dalla dimensione degli elementi.
Se la geometria del problema presenta contorni curvilinei, è bene usare elementi con
nodi lungo i lati (elementi di ordine superiore al primo).
1.6 Valutazione della matrice di rigidezza del singolo elemento
Le equazioni che governano il comportamento deformativo del singolo elemento sono
le equazioni di equilibrio, di congruenza e la legge costitutiva del materiale. In un
problema piano, in cui le variabili primarie siano soltanto gli spostamenti, è possibile
scrivere:
null ∆null null nullnull ∆null ∆null nullnull nullnull null null ∆null ∆null null null null nullnull nullnull ∆null null null (1.7)
Capitolo 1 10
___________________________________________________ __________________________
Noti gli spostamenti è possibile determinare le deformazioni, ottenute come derivate
degli spostamenti. Per un problema di deformazione piana in cui ∆null null null0 ; ∆null null null0 :
null null ∆null null null null null ∆null null ∆null null ∆null nullnull ∆null null null null null ∆null null null null null ∆null null ∆null null ∆null null ∆null nullnull null null null (1.8)
dove ∆null null null∆null nullnull null∆null nullnull null0 .
Il precedente sistema di equazioni, riscritto in forma compatta fornisce la seguente
relazione matriciale:
null ∆null null null nullnull nullnull ∆null null null (1.9)
In cui nullnull null è la matrice che contiene le derivate parziali degli spostamenti.
Se si scelgono, come variabili principali, gli spostamenti incrementali, per semplice
derivazione del vettore degli spostamenti nodali globale, è possibile ottenere le
deformazioni incrementali e, per integrazione della legge costitutiva, il vettore delle
tensioni incrementali. Tale operazione è eseguita in corrispondenza dei punti di Gauss
di ogni elemento. Il problema è a controllo di deformazione, in quanto data la
deformazione si vuole conoscere la tensione:
null ∆null null null nullnull nullnull ∆null null (1.10)
• In cui nullnull null è la matrice costitutiva del materiale;
• null ∆null null è l’incremento di tensione;
• null ∆null null è l’incremento di deformazione.
1.7 Integrazione numerica
Per materiali a comportamento lineare, invocando il principio del Minimo dell’Energia
Potenziale, si ottiene che alla scala locale:
nullnull null null
null null nullnull null null nullnull nullnull null nullnullnullnullnull
nullnullnull (1.11)
null ∆null null null null null nullnull null null null ∆null null nullnullnullnull
nullnullnull null null nullnull null null null ∆null null nullnullnullnull
nullnullnull (1.12)
Capitolo 1 11
___________________________________________________ __________________________
dove:
• null ∆null null è il vettore delle forze di volume applicate al contorno;
• null ∆null null è il vettore delle forze di trazione superficiali applicate al contorno.
La matrice di rigidezza del singolo elemento ed il vettore delle forze nodali incrementali
non possono essere determinati analiticamente, ma si deve ricorrere ad un metodo di
integrazione numerica. Lo scema di integrazione numerica piø comunemente utilizzato,
nel FEM, è quello dell’integrazione Gaussiana. Gli integrali si valutano soltanto in
specifici punti di integrazione, detti punti di Gauss dell’elemento. Il numero di punti di
Gauss dell’elemento definisce l’ordine di integrazione (ridotta o completa).
Per valutare la matrice di rigidezza ed il vettore delle forze è necessario calcolare gli
integrali relativi alle relazioni (1.11) e (1.12). Essi, tuttavia, non sono facilmente
valutabili, di conseguenza è necessario avvalersi di metodi di calcolo numerici. Ad
esempio, discretizzando il dominio di integrazione in intervalli di ampiezza regolare, è
possibile valutare l’area, sottesa ad una qualsiasi funzione continua nel dominio di
integrazione, attraverso la formula dei trapezi (vedi figura 1.7) conoscendo il valore
della null null null null nei punti intermedi dei predetti intervalli. Al crescere del numero di intervalli
(quindi al ridursi della loro ampiezza), aumenta l’onere in termini di tempo di
computazione, ma allo stesso tempo migliora la precisione del risultato.
Figura 1.7. Integrazione con la
formula dei trapezi [1]
Capitolo 1 12
___________________________________________________ __________________________
1.8 Soluzione delle equazioni globali
Per determinare le incognite globali bisogna risolvere il sistema:
nullnull null null null ∆null null null null null ∆null null null (1.13)
Intuitivamente, saremmo portati a risolvere il sistema invertendo la matrice di rigidezza
globale del sistema, nullnull null null , ma viste le sue dimensioni tale operazione non è facile. La
soluzione adottata è, pertanto piø semplice; essa si articola nell’utilizzo del metodo di
decomposizione triangolare di nullnull null null e del forward elimination and back-substitution
method. Se la matrice di rigidezza globale del sistema nullnull null null non è singolare (ovvero il
suo determinante non è nullo), essa può essere sempre scritta come prodotto di due
matrici:
nullnull null null null nullnull nullnull null null (1.14)
dove:
• nullnull null è una matrice triangolare bassa con coefficienti unitari lungo la diagonale
principale;
• nullnull null è una matrice triangolare alta.
L’equazione (1.13) si può riscrivere, dunque, come:
nullnull nullnull null nullnull ∆null null null null null ∆null null null (1.15)
Definiamo,ora, la matrice intermedia:
nullnull nullnull ∆null null null null nullnull null (1.16)
Sostituendo nella relazione (3.15) ottiene:
nullnull nullnull null null null null ∆null null null null nullnull null null null ∆null null nullnull null null nullnull (1.17)
Calcolato nullnull null lo si va a porre nella relazione (1.16) e si riesce a valutare l’incremento
del vettore degli spostamenti incogniti:
nullnull nullnull ∆null null null null nullnull null null null ∆null null null null nullnull nullnull null null nullnull (1.18)
Capitolo 1 13
___________________________________________________ __________________________
Da questa operazione se sono stati scelti, come variabili principali, gli spostamenti
incrementali null ∆null null null , allora per semplice derivazione del vettore degli spostamenti nodali
globale, è possibile ottenere le deformazioni incrementali e per integrazione della legge
costitutiva, il vettore delle tensioni incrementali.
CAPITOLO 2
TEORIA DEL METODO AGLI ELEMENTI FINITI IN
CONDIZIONI STATICHE:
Materiali a comportamento non lineare
Capitolo 2 14
___________________________________________________ __________________________
2.1 Introduzione
Questo capitolo descrive la teoria del metodo degli elementi finiti nell’ipotesi piø
realistica di materiali a comportamento non lineare. In questo caso esistono due modi di
soluzione:
• Metodo tangent stiffness;
• Metodo di Newton-Raphson.
Quando il materiale ha un comportamento elastico non lineare o elasto-plastico, la
matrice costitutiva non è piø costante ma varia durante l’analisi.
2.2 Analisi non lineare con il metodo agli elementi finiti
L’abbandono della linearità trova spiegazione nel fatto che, in realtà, il terreno esibisce
un comportamento meccanico caratterizzato da: forte non-linearità del legame tensioni-
deformazioni, irreversibilità, influenza rispetto alla storia tensionale subita. La non-
linearità viene introdotta scrivendo l’equazione di equilibrio, l’equazione costitutiva e di
congruenza in forma incrementale. L’equazione di equilibrio globale scritta in forma
incrementale è la seguente:
nullnull null null null null ∆null null null null null null ∆null null null null (2.1)
dove:
• i è il numero di incrementi o steps;
• null null null null null è la matrice di rigidezza globale all’i-esimo step, che dipende dallo stato tenso-
deformativo;
• null ∆null null null null è il vettore di incremento degli spostamenti dei nodi;
• null ∆null null null null è il vettore di incremento delle forze nodali.
E’ necessario adottare particolari strategie risolutive che permettano di tener conto della
non-linearità dovuta alla legge costitutiva.
Capitolo 2 15
___________________________________________________ __________________________
2.3 Tangent stiffness method
Il carico totale applicato al contorno viene suddiviso in un numero i di incrementi.
Applicato il primo incremento di carico null ∆null null null null , usando come prima matrice di
rigidezza del sistema quella elastica null null null null null , è possibile ricavare il primo incremento di
spostamento nodale null ∆null null null null , applicando il primo step di carico. Si impone nota la curva
reale nel diagramma null ∆null null null null null ∆null null null del sistema, che rappresenta la soluzione vera (vedi
figura 2.1).
Per la realizzazione dell’analisi si utilizza, come prima matrice per risolvere il problema
quella elastica, cioè viene considerata come prima matrice la tangente all’origine null null null null null ,
che è proprio la tangente della curva del comportamento reale nella parte iniziale a
piccolissime deformazioni. Noto l’incremento di carico, e la matrice di rigidezza
elastica, da queste due informazioni si ricava il primo spostamento globale null ∆null null null null . Si
osserva che lo spostamento ottenuto non è quello che segue la curva reale, questo
perchØ si è applicato l’incremento di carico utilizzando una rigidezza costante. Dunque
Figura 2.1. Tangent Stiffness Method [1]