Algoritmi e Strutture Dati per la Modellazione Interattiva...

UNIVERSITA DEGLI STUDI DI PISA

Facolta di Scienze

Matematiche Fisiche e Naturali

Corso di LaureaSpecialistica in Informatica

Tesi di laurea

Algoritmi e Strutture Dati per la

Modellazione Interattiva di Tagli su

Oggetti Deformabili

Candidato:

Nico Pietroni

Relatori:

Dott. P. Cignoni

Dott. F. Ganovelli

Controrelatore:

Prof. Francesco Romani

Anno Accademico 2002/03

Prova a riconoscere i differenti poliedri presenti nel quadro di M.C Escher.

Ringraziamenti

Il primo ringraziamento va alla mia famiglia ed in particolare ai miei geni-

tori per la fiducia mostrata durante questi anni di studio. Un grazie speciale

a Paola per il grande affetto donato durante questi giorni di scleri. Grazie a

Dinstensao per avermi mostrato la via del magico mondo dell’Universita di

Pisa. Grazie a Fabio, Paolo e Mr.Mendoza (che parla ispano-franco-italiano)

per l’aiuto offerto durante lo svolgimento di questa tesi. Non posso dimen-

ticare infine di includere nei ringraziamenti tutta casa Notari, gli amici di

Sangi e Roncordi per l’energia creativa che mi hanno dato.

Algoritmi e Strutture Dati per la Modellazione Interattiva

di Tagli su Oggetti Deformabili

Nico Pietroni

Dipartimento di Informatica

Universita di Pisa

Sommario

Nel campo della simulazione del comportamento di oggetti deformabili,

uno dei problemi aperti e quello della modellazione e simulazione di tagli.

Questa tesi presenta delle strutture dati ed algoritmi per simulare interatti-

vamente l’applicazione di tagli su oggetti deformabili rappresentati con mesh

volumetriche tetraedrali.

L’algoritmo puo simulare l’effetto dell’azione di taglio sulla superficie con di-

versi livelli di dettaglio, permettendo di bilanciare precisione e performance

in base ad una misura dell’errore percepito dall’utente. La scelta di model-

lare le modifiche sugli oggetti per mezzo di operatori atomici permette di

allocare le risorse di calcolo disponibili con un’alta granularita.

Indice

1 Introduzione 1

1.1 Schema generico di un simulatore . . . . . . . . . . . . . . . . 4

1.2 Obiettivi della tesi . . . . . . . . . . . . . . . . . . . . . . . . 6

1.3 Organizzazione della tesi . . . . . . . . . . . . . . . . . . . . . 7

2 Stato dell’arte 8

2.1 Modellazione del tessuto . . . . . . . . . . . . . . . . . . . . . 8

2.1.1 Concetti di base . . . . . . . . . . . . . . . . . . . . . . 9

2.1.2 Elasticita lineare . . . . . . . . . . . . . . . . . . . . . 10

2.1.3 Modelli continui: FEM . . . . . . . . . . . . . . . . . . 12

2.1.4 Modelli discreti . . . . . . . . . . . . . . . . . . . . . . 14

2.1.5 Particle systems . . . . . . . . . . . . . . . . . . . . . . 14

2.1.6 Il modello Mass Spring . . . . . . . . . . . . . . . . . . 15

2.1.7 Modelli ibridi: Tensor mass . . . . . . . . . . . . . . . 16

2.1.8 Tensor Mass e Mass Spring a confronto . . . . . . . . . 18

2.1.9 Modelli continui e modelli discreti . . . . . . . . . . . . 18

2.2 Definizione di mesh tetraedrali . . . . . . . . . . . . . . . . . . 19

INDICE ii

2.3 Utilizzo di mesh tetraedrali nella rappresentazione dei modelli 20

2.4 Modifiche geometriche e topologiche . . . . . . . . . . . . . . 21

2.4.1 Il taglio di una mesh tetraedrale . . . . . . . . . . . . . 22

2.4.2 Modifiche topologiche nella mesh tetraedrale . . . . . . 23

2.4.3 Approccio con distruzione di tetraedri . . . . . . . . . 24

2.4.4 Approccio con duplicazione di facce . . . . . . . . . . . 24

2.4.5 Approccio con suddivisione di tetraedri . . . . . . . . . 25

2.5 Metodi di integrazione numerica . . . . . . . . . . . . . . . . . 28

2.5.1 Metodi di Runge Kutta . . . . . . . . . . . . . . . . . 30

2.5.2 Metodi impliciti . . . . . . . . . . . . . . . . . . . . . . 31

3 Tecniche di taglio su mesh tetraedrali 33

3.1 Strutture dati e algoritmi per mesh tetraedrali . . . . . . . . . 34

3.1.1 Relazioni di adiacenza . . . . . . . . . . . . . . . . . . 36

3.1.2 Operatori di navigazione sulla mesh . . . . . . . . . . . 42

3.2 Modifiche topologiche: tagli . . . . . . . . . . . . . . . . . . . 50

3.2.1 Gli operatori di taglio . . . . . . . . . . . . . . . . . . 52

3.2.2 Problemi connessi all’uso degli operatori di taglio . . . 53

3.2.3 Vertex-split . . . . . . . . . . . . . . . . . . . . . . . . 56

3.2.4 Edge-split . . . . . . . . . . . . . . . . . . . . . . . . . 64

3.2.5 Riassegnamento dei parametri fisici dopo il taglio . . . 70

3.3 Semplificazione . . . . . . . . . . . . . . . . . . . . . . . . . . 71

3.3.1 Edge Collapse . . . . . . . . . . . . . . . . . . . . . . . 73

3.3.2 Correttezza topologica . . . . . . . . . . . . . . . . . . 74

INDICE iii

3.3.3 Correttezza geometrica . . . . . . . . . . . . . . . . . . 76

3.3.4 Riassegnamento dei parametri fisici dopo la semplifi-

cazione . . . . . . . . . . . . . . . . . . . . . . . . . . 78

3.3.5 Implementazione della decimazione . . . . . . . . . . . 78

4 La modellazione del contatto 81

4.1 Collision Detection . . . . . . . . . . . . . . . . . . . . . . . . 81

4.1.1 Tecniche di Collision Detection . . . . . . . . . . . . . 82

4.1.2 La Collision Detection nella simulazione chirurgica . . 83

4.1.3 La Collision Detection con hardware grafico . . . . . . 84

4.1.4 L’intersezione geometrica . . . . . . . . . . . . . . . . . 89

4.2 Il contatto . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91

4.2.1 Le forze di di deformazione . . . . . . . . . . . . . . . . 92

4.2.2 Approccio con molla . . . . . . . . . . . . . . . . . . . 93

4.2.3 Approccio con Force Grid . . . . . . . . . . . . . . . . 94

4.2.4 Approccio con spostamento delle facce . . . . . . . . . 95

4.3 Il taglio . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98

4.3.1 Modellazione del taglio guidata dall’errore . . . . . . . 98

4.3.2 Vertex-split, calcolo della sezione di taglio . . . . . . . 101

5 Risultati sperimentali e sviluppi futuri 104

5.1 Risultati sperimentali . . . . . . . . . . . . . . . . . . . . . . . 104

5.2 Sviluppi futuri . . . . . . . . . . . . . . . . . . . . . . . . . . . 108

5.2.1 Interazione con un solido . . . . . . . . . . . . . . . . . 109

5.2.2 Lacerazioni . . . . . . . . . . . . . . . . . . . . . . . . 110

INDICE iv

5.2.3 Self Collision . . . . . . . . . . . . . . . . . . . . . . . 111

5.2.4 Sutura . . . . . . . . . . . . . . . . . . . . . . . . . . . 111

5.3 Conclusioni . . . . . . . . . . . . . . . . . . . . . . . . . . . . 113

5.4 Alcune immagini del simulatore . . . . . . . . . . . . . . . . . 115

Capitolo 1

Introduzione

Questa tesi affronta il problema della simulazione realistica dell’interazio-

ne con oggetti deformabili mediante contatto con particolare riferimento al

campo di applicazione della chirurgia virtuale.

Con l’aggettivo deformabile si indica un oggetto la cui forma cambia in

modo sensibile se sottoposto ad uno stress. La difficolta principale nel mo-

dellare questo tipo di oggetti e il costo computazionale di simulare il com-

portamento fisico del materiale, infatti solo nell’ultimo decennio l’evoluzione

delle macchine calcolatrici ha reso possibile l’implementazione delle prime

applicazioni interattive anche su personal computer.

La ricaduta applicativa di maggior rilievo in questo campo e quella della

modellazione di scenari chirurgici a scopo didattico e di addestramento. In

questo contesto si introduce la necessita di modellare non solo la deformazio-

ne ma anche il taglio interattivo del materiale. Questo requisito ha generato

una serie di problemi sia nella modellazione della fisica, sia nel rilevamen-

to del contatto e delle collisioni con altri oggetti sia, in modo minore, nella

visualizzazione, poiche tutte le soluzioni proposte in questo campo si sono

sviluppate nell’ipotesi che la topologia degli oggetti simulati rimanga inva-

riata e che questi possano essere sempre rappresentati con lo stesso modello

geometrico.

Il contributo principale di questa tesi e la definizione e implementazio-

2

ne di strutture dati ed algoritmi per rappresentare oggetti deformabili che

possono subire interattivamente tagli o lacerazioni. I problemi da affrontare

sono i seguenti:

frammentazione Gli oggetti sono rappresentati come un complesso di pri-

mitive geometriche, nella maggior parte dei casi tetraedri, e l’applica-

zione di tagli inevitabilmente porta ad un incremento del loro numero,

con conseguente aumento del costo computazionale del tempo simu-

lazione della fisica, che e proporzionale dal numero di primitive che

compongono l’oggetto;

frame-rate costante L’aggiornamento della vista ad una frequenza non in-

feriore ai 20Hz e un requisito irrinunciabile per poter definire un’opera-

zione interattiva. Mentre il costo di un passo di simulazione della fisica

e costante per tutti i passi, quello di aggiornare la topologia dipende da

quanta parte dell’oggetto e interessata al taglio ed anche alla velocita

con cui tale taglio viene effettuato, quindi si va ad aggiungere come

overhead e in generale puo causare un rallentamento percepibile della

simulazione;

L’approccio seguito si basa su due strategie. La prima e quella di definire

algoritmi basati sulla concatenazione di operazioni atomiche che occupino un

tempo trascurabile e terminino lasciando il sistema in uno stato corretto e

consistente, cosı da ottenere algoritmi interrompibili. La seconda e di avere

piu modi di rappresentare uno stesso taglio, con diversa accuratezza visiva

e costo computazionale, e di poter valutare l’errore percepito dall’utente in

modo da scegliere un modo meno accurato ma piu veloce di rappresentare il

taglio se la differenza in termini di immagine sullo schermo e trascurabile. Al

fine di valutare sperimentalmente le soluzioni proposte si e implementato un

3

prototipo di simulatore che consente all’utente di gestire il bisturi virtuale ed

interagire con gli oggetti che compongono la scena deformandoli e tagliandoli.

1.1 Schema generico di un simulatore 4

Figura 1.1: Schema generico di un simulatore.

1.1 Schema generico di un simulatore

In un simulatore chirurgico l’utente interagisce col modello fisico utilizzando

un tool virtuale (bisturi) . Il contatto tra il tool ed il modello simulato

provochera deformazioni, tagli o lacerazioni di tessuto restituendo sul bisturi

virtuale (tramite periferica) una forza di feedback in modo tale da dare la

sensazione toccare veramente il modello con il tool.

Il simulatore, se analizzato ad alto livello, puo essere composto da una serie

di moduli (Figura 1.1) :

Physical engine Data una descrizione dell’oggetto in termini di geometrie

e di proprieta fisiche del materiale (quali le costanti di elasticita del

materiale cui e composto il modello ) l’engine si occupa di calcolare

quali deformazioni intervengono in conseguenza di sollecitazioni esterne

durante, ad esempio, il contatto con un altro oggetto, e interne, come

la forza di gravita. Il modulo di topological modification si occupa di

modellare le modifiche che intervengono sul modello deformabile (tagli,

1.1 Schema generico di un simulatore 5

lacerazioni. . . ) .

Rendering Questo modulo visualizza gli oggetti che compongono la scena.

Collision detection Il modulo Collision Detection rileva quando esiste un

contatto tra i vari elementi presenti nella scena corrente. Fornisce dati

indispensabili per individuare dove e come sia avvenuto il contatto:

quali siano le entita interessate direttamente, quale sia superficie di

contatto e le forze esercitate tra loro; tutti dati indispensabili all’engine

fisico durante il calcolo e la conseguente deformazione degli oggetti

coinvolti.

Haptic feedback L’interfaccia attraverso cui l’utente interagisce con la sce-

na e costituita di solito da un sistema di force feedback. Questo modulo

deve comunicare col sistema di feedback esterno fornendo al sistema i

dati relativi alla posizione ed all’orientamento del tool nello spazio. Ad

ogni contatto rilevato con il modello deformabile questo modulo utiliz-

zara i dati forniti dall’engine fisico per simulare il contatto tramite una

forza di feedback sull’interfaccia.

Avere una buona frequenza di aggiornamento del sistema (quindi un buon

frame-rate) e l’obiettivo primario da perseguire affinche la simulazione risulti

il piu possibile realistica. I sistemi proposti finora per modellare oggetti

deformabili sono molto efficienti in situazioni in cui il modello non e mai

soggetto a modifiche di tipo topologico (ad esempio i tagli) ; quando pero il

modello viene modificato in tempo reale l’engine fisico deve essere in grado

di reinterpretare e simulare di conseguenza le nuove situazioni.

1.2 Obiettivi della tesi 6

1.2 Obiettivi della tesi

In una simulazione chirurgica l’utente deve essere in grado almeno di tagliare

liberamente il tessuto su cui sta operando. Modifiche topologiche di questo

tipo incrementano la complessita della descrizione del modello aumentan-

do di conseguenza anche il tempo necessario all’engine fisico per simulare il

comportamento dei vari tessuti. La vera difficolta sta fondamentalmente nel

creare sistemi di modifica topologica che non incrementino troppo la descri-

zione del modello fisico, accompagnati tecniche di semplificazione geometrica

che riducano il set di dati su cui dovra lavorare il framework. Le caratteri-

stiche principali del sistema di taglio quı proposto possono essere riassunte

in questi 3 punti fondamentali:

1. Le operazioni non dovranno incrementare inutilmente la complessita del

modello deformabile. Ogni taglio introduce sempre una quantita di in-

formazioni in piu nel modello ; piu informazioni introduce e maggiore e

il realismo fisico dell’azione di taglio; purtroppo pero tutto questo dimi-

nuisce il rate dell’aggiornamento del sistema. Rendere veramente reale

lo scenario per l’utente significa anche avere un sistema che reagisca in

modo istantaneo alle sollecitazioni esterne (il tempo di computazione

non deve essere percepito) . Il giusto bilanciamento tra realta e sem-

plicita sara quindi la filosofia che sta alla base della costruzione degli

operatori di taglio.

2. Le operazioni di modifica saranno accompagnate da operazioni di sem-

plificazione del modello. Questa caratteristica puo far tendere la com-

plessita del modello ad una situazione costante.

3. Le operazioni saranno, nella loro semplicita abbastanza versatili da

consentirne in futuro l’utilizzo nella simulazione altre operazioni (ad

esempio le lacerazioni del tessuto) .

1.3 Organizzazione della tesi 7

1.3 Organizzazione della tesi

Nel capitolo 2 di questa tesi verra illustrato lo stato dell’arte relativo alla

modellazione di oggetti deformabili. In particolar modo verranno illustrati i

possibili modelli matematici utilizzati nella simulazione di oggetti deformabili

e le principali tecniche di rappresentazione di taglio proposte finora.

Nel capitolo 3 verranno presentate inizialmente le strutture dati utilizzate per

la rappresentazione dell’oggetto deformabile oltre alle strutture necessarie ad

effettuare interrogazioni sulle sue proprieta geometriche. Sulla base di queste

strutture verranno poi definite le operazioni di modifica e semplificazione

necessarie all’esecuzione del taglio.

Nel capitolo 4 verra definita in dettaglio l’interazione tra il bisturi e l’oggetto

deformabile mettendo in particolare rilievo le differenze con gli altri metodi

finora utilizzati.

Infine il capitolo 5 sara dedicato all’analisi dei risultati ottenuti e ai possibili

sviluppi futuri del simulatore.

Capitolo 2

Stato dell’arte

La realizzazione di ognuno dei componenti di un simulatore propone problemi

che spaziano in vari campi della computer graphics e della robotica. Questo

capitolo presenta un breve stato dell’arte sul tema, con particolare rilievo sui

problemi studiati in questa tesi.

2.1 Modellazione del tessuto

Un modello deformabile deve reagire in maniera realistica agli stress impo-

sti dagli eventi esterni relativamente alle sue proprieta elastiche e geometri-

che. Nel simulatore che ci proponiamo di realizzare la tecnica adottata per

modellare un tessuto deve soddisfare 2 caratteristiche:

• L’interattivita dell’ applicazione comporta che la deformazione del tes-

suto sia elaborata in tempo reale.

• La possibilita di effettuare tagli impone che il modello sia facilmente

adattabile alle modifiche geometriche.

2.1 Modellazione del tessuto 9

2.1.1 Concetti di base

Ogni modello per una simulazione del comportamento fisico di un oggetto si

basa sulla seconda legge del moto di Newton:

mx + γx +∂Erg(x)

∂u= f (2.1)

dove x indica la posizione del punto di massa m e γ il damping factor, ov-

vero la resistenza al movimento dovuta all’attrito. Il termine ∂Erg(x)∂u

indica

la derivata dell’energia potenziale rispetto alla posizione del punto, infine

f e la forza esterna esercitata sulla massa. Se, ad esempio, il punto fosse

legato all’origine del sistema di coordinate da una molla di coefficiente k e

si stesse muovendo lungo l’asse x, allora E(x) =∫ x0 kudu = 1

2ku2 e quindi

∂Erg(x)∂u

= kx. Supponiamo che un oggetto sia rappresentato da un numero

di punti x0, . . . , xn di massa m0, . . . ,mn, distribuiti in uno spazio continuo.

Se riusciamo a trovare l’equazione dell’energia associata alla posizione di cia-

scuno di questi punti potremmo predire e quindi calcolare il comportamento

dell’oggetto risolvendo le equazioni differenziali risultanti.

Modelli discreti e modelli continui Anche se adesso fossimo in grado

di descrivere il comportamento di questo oggetto per ogni punto attraverso

l’equazione 2.1, bisogna stabilire quale sia l’insieme di punti scelto in modo

tale da descrivere comportamento del modello in maniera adeguata.

Le strategie di rappresentazione di un modello possono essere raggruppate

in due grandi classi: I modelli discreti ed i modelli continui. Nei modelli

discreti la rappresentazione del modello viene ridotta ad un sottoinsieme dei

suoi punti tali da caratterizzare in modo piu realistico possibile le proprieta

meccaniche dell’oggetto. La complessita fisica dell’oggetto richiedera una

distribuzione di punti adeguata: Se ad esempio l’oggetto da rappresentare

avesse la forma di una barra di gomma; il sistema potrebbe essere forma-

to da un insieme di punti uniformemente distribuiti lungo l’estensione della


barra; ognuno di essi sarebbe legato agli adiacenti tramite una molla. Piu

complesso e l’oggetto da descrivere e piu complessa sara la definizione dell’e-

nergia. I modelli continui si basano invece su una rappresentazione continua

dell’oggetto formata da infiniti punti distribuiti nello spazio del modello. In

questo caso l’equazione del moto dovra essere riscritta come :

µ(u)x + λ(u)x +∂Erg(x)

∂u= f(x) (2.2)

dove x = x(u, t) e la posizione del punto u al tempo t. Siccome questa equa-

zione e risolvibile per via analitica solo in casi ideali (non e possibile descrive-

re analiticamente la forma di molti oggetti), il problema verra discretizzato

attraverso l’analisi ad elementi finiti.

2.1.2 Elasticita lineare

Quasi tutti i materiali possiedono una certa quantita di di elasticita. Se forze

esterne producono una deformazione che non eccede un certo limite, la de-

formazione scomparira appena le forze saranno rimosse. Lo stress σ, definito

come forza per unita di area e lo strain ε (tensione) , definito come lo spo-

stamento per unita di lunghezza, soddisfano la seguente relazione conosciuta

come la legge di Hook [TG85]:

εx =σx

E(2.3)

dove E e il modulo di elasticita. L’estensione dell’elemento lungo la di-

rezione delle x e accompagnata da componenti laterali di strain (l’oggetto

allungato si contrae e si assottiglia) :

εy = −νυσx

E(2.4)

εz = −νυσx

E(2.5)


dove ν e una costante chiamata coefficiente di Poisson. Le equazioni 2.3,

2.4 e 2.5 possono essere utilizzate anche per descrivere forze di compressione

sul modello (il modulo di elasticita ed il coefficiente di Poisson sono gli stessi

che in caso di tensione) . Le due costanti E e ν definiscono completamente

le proprieta elastiche di un corpo.

L’Energia elastica Prima di tutto definiamo il piu semplice modello bio-

meccanico: il modello elastico lineare. Denotiamo come Minitial il modello

allo stato iniziale.

Sotto eventi esterni il modello si deforma. Rappresentiamo la deformazione

di un modello dalla sua forma iniziale con un vettore di spostamento U(x, y, z)

con (x, y, z) ∈ Minitial; quindi:

Mdeformed = Minitial + U(x, y, z) (2.6)

L’energia elastica Welastic misura la quantita di deformazione di Mdeformed

e puo essere espressa come[CA99]:

Welastic = λ2(trE)2 + µ(trE)2 (2.7)

dove:

• la matrice 3 x 3 simmetrica E, detta anche green tensor e definita come:

E =1

2(∇U +∇UT ) (2.8)

• λ e µ sono i coefficienti di Lame che caratterizzano la rigidezza di un

materiale, definiti come:

λ =νE

(1 + ν)(1− 2ν)(2.9)


µ =E

2(1 + ν)(2.10)

dove E e ν rappresentano il modulo di elasticita ed il coefficiente di

Poisson.

2.1.3 Modelli continui: FEM

Questa sezione e dedicata principalmente ad illustrare la soluzione proposta

dal metodo ad elementi finiti (FEM)[Zie67]. Il modello sui cui e stata definita

l’equazione 2.1 e composto da masse puntiformi. Questa equazione puo essere

estesa in un dominio continuo: sia Ω il dominio dei punti del volume occupato

dal modello, estendendo l’equazione 2.1 nel continuo otteniamo:∫Ω

∂Erg(x)

∂udu−

∫Ω

fdu = 0 (2.11)

L’energia potenziale dovuta alla deformazione equivale al lavoro esercitato

dalla forza che agisce sul corpo. Consideriamo tale corpo sotto l’azione di

carichi distribuiti, ad esempio la forza di gravita, e sia f la somma di tali

forze. Se il sistema e in equilibrio, quando uno spostamento viene imposto

dall’esterno, il lavoro esterno dato dal carico deve equivalere al cambio di

energia di strain interna del corpo:

1

2

∫Ω[σ]T εdΩ−

∫Ω[σ]T pdΩ−

∫§[σ]T q = 0 (2.12)

dove σ e ε sono il vettore di strain e di stress, δ e lo spostamento, p e

la densita della forza nel corpo, q e la densita della forza sulla superficie S

(boundary stress) . Questa situazione non e risolvibile analiticamente, eccetto

per casi ideali. L’idea principale del FEM consiste nel decomporre il dominio

Ω in un insieme di elementi finiti di cui si conosca la descrizione analitica

del dominio di integrazione (e quindi del volume stesso occupato da questi

elementi) . L’integrazione su tutto il dominio, cosı come avevamo descritto


precedentemente, sara sostituita dalla sommatoria di tutti gli integrali che

costituiscono gli elementi in cui viene decomposto il modello:

1

2

∑e

∫Ωe

[σ]T εdΩ−∑e

∫Ωe

[σ]T pdΩ−∑e

∫Se

[σ]T qdS = 0 (2.13)

dove Ωe e la parte di dominio locale dell’elemento finito, che sara definito

da un insieme di punti (ad esempio un elemento esaedrale o tetraedrale) .

L’integrazione dell’intero dominio si riduce cosı ad una serie di integrali sui

singoli elementi. Possiamo esprimere lo spostamento δ come interpolazione

lineare di degli spostamenti dei punti che definiscono l’elemento cella. Sia

quindi definita la funzione:

δ = Nδe (2.14)

dove δe e un vettore di spostamenti di nodi con 3 componenti per ogni

vertice. Sostituendo δ in 2.13 possiamo ottenere la formulazione discreta:

Kδ = F (2.15)

Naturalmente i nodi di ciascun elemento sono condivisi dagli elementi

adiacenti, per questo le equazioni non sono indipendenti. La matrice K di

3n × 3n e chiamata matrice di rigidita (Stiffness matrix ) . La soluzione di

questo sistema definisce la posizione dei vertici nello spazio, e quindi la defor-

mazione dell’oggetto. Una volta definita l’energia, diventa di fondamentale

importanza definire la funzione di interpolazione δ. Tale funzione risultera

definita sul tipo di elemento primitivo utilizzato per suddividere il dominio.

L’idea e di approssimare la funzione δ tramite una somma pesata dei vertici

che definiscono la forma della primitiva:

δ(x, y, z) =l−1∑i=0

δihi(x, y, z) (2.16)


dove l e il numero di elementi che definiscono la forma (4 per un tetraedro)

. Una funzione di interpolazione deve essere almeno C1 continua lungo gli

spigoli dell’elemento (la tensione non puo essere infinita) . Aumentando la

granularita della suddivisione il metodo tendera verso la soluzione corretta.

2.1.4 Modelli discreti

I modelli discreti rappresentano il comportamento di un materiale attraverso

un insieme di parti meccaniche elementari e un insieme di relazioni tra essi. In

questo tipo si modelli i piu utilizzati sono quelli che si basano su una struttura

di tipo particle system[Wita]; la loro popolarita e dovuta soprattutto alla

semplicita di realizzazione oltre che alla loro versatilita.

2.1.5 Particle systems

I Particle sono oggetti che hanno massa, posizione, velocita, rispondono alle

forze pur non avendo una dimensione spaziale. Nonostante la loro semplicita

i particles possono presentare una grande quantita di caratteristiche inte-

ressanti. Come tutti i corpi il moto delle particle e governato dalla legge

f = m · a. Le particles possono essere messe in relazione tra loro in modo

tale che alcune forze interne vengano trasmesse da un insieme di particles

ad un altro; questa caratteristica puo essere utilizzata per simulare proprieta

meccaniche di un fluido ad esempio, oppure (nel caso non ci sia alcun legame

tra le particelle) di un gas. La semplicita di queste strutture diventa effi-

cienza quando si tratta di simulare il comportamento di strutture non rigide

(deformabili ed interattive nel nostro caso) . Le forze interne si trasmettono

tra le varie particle secondo un modello numerico che simuli l’elasticita del

modello. L’imposizione di una qualsiasi forza esterna al particle system cau-

sera la deformazione del modello ed il conseguente aggiornamento dello stato

delle particelle.


1. La situazione iniziale di equilibrio viene perturbata da forze esterne

cambiando la posizione iniziale di un insieme di particles.

2. La forza scaturita si trasmette alle altre particle secondo un modello

elastico lineare.

3. Il solver di equazioni differenziali calcola l’accelerazione la velocita e la

prossima posizione della particle.

4. Il sistema perturbato continua a mutare fino a raggiungere una posi-

zione di equilibrio.

2.1.6 Il modello Mass Spring

Una particolare istanza del Particle system per la modellazione degli ogget-

ti deformabili e il modello Mass spring. In questo modello le particle che

costituiscono il modello sono interconnesse tra di loro tramite molle. Le

connessioni tra particles possono seguire ad esempio una rete spigoli di celle

tetraedrali che descrivano l’intero volume dell’oggetto deformabile. Il valore

della forza interna agente sulla particella Pi assume la forma:

Fi =∑

j∈N(Pi)

kij(‖PiPj‖ − l0ij)PiPj

‖PiPj‖(2.17)

Dove N(Pi) e l’insieme dei vertici Pj adiacenti a Pi, kij e il coefficiente di

rigidita tra i vertici Pi e Pj, l0ij e la lunghezza iniziale della molla tra Pi e Pj.

Il modello Mass spring costituisce nel complesso un metodo abbastanza sem-

plice da implementare; purtroppo la resa dal punto di vista del realismo fisico

e piu scarsa rispetto ai metodi continui come il FEM. Nel FEM la variazio-

ne di energia interna dell’oggetto viene calcolata basandosi sul volume del

dominio di integrazione prescelto all’interno del quale la massa e distribuita


in maniera uniforme. Questo tipo di modello fisico costituisce una rappre-

sentazione piu accurata della realta; le differenze col modello Mass Spring

si notano poi in termini di realismo durante la simulazione. Un metodo di-

screto che simuli il comportamento dell’oggetto in maniera piu efficiente e il

modello Tensor Mass.

2.1.7 Modelli ibridi: Tensor mass

Un altro modello discreto (ma che puo essere considerato un modello ibrido)

e il modello Tensor Mass[CA99]. Supponiamo che il volume dell’oggetto

deformabile sia suddiviso in tetraedri (questa rappresentazione sara trattata

piu in dettaglio piu nella prossima sezione) . Denotiamo come Minitial il

modello allo stato iniziale a con P 0i la posizione del vertice i-esimo all’istante

zero. Procediamo in 3 passi:

1. definiamo le equazioni di interpolazione che definiscono il vettore di

spostamento di un punto (x, y, z) dentro al tetraedro Ti come funzione

dei 4 vettori di spostamento dei 4 vertici di Ti.

2. Scriviamo l’energia elastica del tetraedro come funzione dei 4 vettori di

spostamento dei 4 vertici di Ti.

3. Calcoliamo la forza elastica Fi che agisce sul vertice Pi.

Vettore di spostamento Dato Ti, il tetraedro definito dai quattro vertici

Ti(j), j=0,. . . ,3, allora il vettore di spostamento sara definito come:

UTi(x, y, z) =

3∑j=0

bj(x, y, z)P 0Ti(j)PTi(j) (2.18)

Dove bj(x, y, z) sono le coordinate baricentriche del punto (x,y,z) dentro

il tetraedro Ti. Poiche le coordinate baricentriche sono lineari rispetto alle


coordinate (x,y,z), la matrice E(x, y, x) definita in 2.8 rimane costante per

ogni punto del tetraedro Ti.

Energia elastica Se λi e µi sono i coefficienti di Lame del tetraedro Ti,

allora possiamo esprimere l’energia elastica Welastic(Ti) del tetraedro Ti come:

WElastic(Ti) = V (Ti)(λi

2(div(U))2 + µitr(∇UT∇U)− µi

2‖curlU‖2) (2.19)

Dove div e la funzione di divergenza, mentre curl e la funzione roto-

re. L’energia totale Welastic(Mdeformed) sara la somma delle energie di ogni

tetraedro.

Forza lineare elastica Data l’espressione dell’energia elastica deriviamo

la forza Fi applicata su un vertice Pi tramite la relazione:

Fi = −∂Welastic(Mdeformed)

∂Pi=

∑Tj∈L(i)−

∂Welastic(Tj)

∂Pi

dove L(i) e l’insieme dei tetraedri adiacenti al vertice Pi. Mentre all’inter-

no del tetraedro la forza FTj(j) applicata sul vertice PTj(j) assume la seguente

forma:

FTj(j) =∑3

k=0[KTijk ]P 0

Ti(k)PTi(k)

dove KTijk e una matrice 3 x 3 chiamata matrice di rigidita definita come:

dato un tetraedro Ti ed i suoi vertici P 0Ti(0)

, P 0Ti(1), P

0Ti(2), P

0Ti(3)

, i sei tensori

sono definiti come:

[KTijk ] = 1

36V (Ti)(λiMkM

Tj + µiMjM

Tk + µi(MjMk)[Id3x3])

dove [Id3x3] e una matrice identita di 3 righe e 3 colonne.

E importante sottolineare come tale matrice dipenda solamente dalle ca-

ratteristiche del materiale (λi, µi) e dalla geometria dell’ oggetto allo stato


iniziale (tale matrice non deve essere quindi aggiornata ad ogni passo della

simulazione) .

Per un dato vertice Pi, la forza elastica Fi sara allora la somma dei

contributi FTj(j) dei tetraedri adiacenti Ti:

Fi = [KiiP0i Pi] +

∑j∈N(Pi)[KijP

0j Pj]

dove Kij e la somma dei tensori [KTj

ii ] associati ai tetraedri adiacenti a Pj,

Kij e la somma dei tensori KTj

ij associati ai tetraedri adiacenti attraverso

l’edge (i, j) e N(Pi) e la lista dei vertici vicini a Pj.

2.1.8 Tensor Mass e Mass Spring a confronto

I due metodi hanno la stessa complessita, che risulta lineare rispetto al nume-

ro degli edge (il metodo del tensor mass necessita di qualche calcolo iniziale in

piu nella valutazione della matrice K) . Il modello mass spring si basa su una

rappresentazione discreta e il suo realismo e fortemente influenzato dalla di-

sposizione della rete di molle all’interno volume dell’oggetto. In questo caso,

effettuare un taglio puo significare rimuovere una o piu molle; cio potrebbe

stravolgere notevolmente il comportamento del modello stesso. D’altro lato

il modello tensor mass offre una simulazione piu realistica e, essendo basato

su una descrizione continua del modello, tale metodo viene influenzato in

maniera minore dalla topologia interna della mesh.

2.1.9 Modelli continui e modelli discreti

Il Metodo ad Elementi Finiti, come gia esposto, scompone il volume dell’og-

getto in un insieme di elementi di cui si riesca a descrivere analiticamente

il dominio di integrazione; di fatto di compie quindi una discretizzazione di

un volume; ma la rappresentazione dell’oggetto rimane di fatto sempre con-

tinue e costituita da infiniti punti. Nei metodi discreti invece il modello a

2.2 Definizione di mesh tetraedrali 19

cui facciamo riferimento durante la simulazione e costituito da un insieme di

masse puntiformi che interagiscono tra loro. Partendo da questo tipo di ap-

prossimazione della realta e ovvio che la resa dal punto di vista del realismo

e peggiore di quella dei modelli continui. Per quanto riguarda la complessita

dei calcoli coinvolti nella simulazione risultano piu efficienti i metodi discreti.

In questi metodi la massa e concentrata nelle particelle di dimensione infini-

tesime costituenti il modello, di conseguenza ogni calcolo di forze o energie e

basato sulle proprieta di questo insieme finito di punti. Nei metodi continui

invece la massa e distribuita nei volumi delle celle, il calcolo viene quindi

ripartito ad un insieme di integrazioni sulle singole celle. I modelli discreti

hanno inoltre la importante caratteristica di essere piu flessibili e piu facil-

mente adattabili a modifiche geometriche del modello. Questo deriva dal

fatto che molti dei metodi continui necessitano spesso una significativa pre-

computazione di strutture dati ausiliarie; una eventuale modifica del modello

(ad esempio un taglio) produrrebbe una incoerenza tra i dati e quello che e

la reale condizione dell’oggetto (a meno che tali dati non vengano anch’essi

ricalcolati) .

2.2 Definizione di mesh tetraedrali

Esponiamo alcune definizioni formali utilizzate in seguito.

Definizione 1 Un simplesso σ a d dimensioni e l’insieme dei punti

dati dalla combinazione convessa di d+1 punti linearmente indipendenti (un

segmento per d=1, un triangolo per d=2. . . ) .

Definizione 2 La faccia di un simplesso σ e un simplesso i cui punti

sono un sottoinsieme di σ (un punto e una faccia di un segmento e di un

triangolo, un segmento e una faccia di un triangolo. . . )

2.3 Utilizzo di mesh tetraedrali nella rappresentazione dei modelli 20

Definizione 3 Dato il simplesso σ di dimensione d, le facce proprie di σ

sono l’insieme di facce di dimensione d-1 appartenenti a σ.

Definizione 4 Dati due simplessi σ di dimensione n e σ′ di dimensione

m con (n > m) allora si dice che σ′ ∈ σ se esiste una sequenza di facce

S = F0, F1, . . . , Fn−m tali che F0 = σ′, Fn−m = σ e per ogni Fi ∈ S − Fn−m

Fi e faccia propria di Fi+1;

Definizione 5 Un complesso simpliciale Γ a d dimensioni (o d-

complex) e un insieme di simplessi tali che:

• per ogni simplesso σ appartenente a Γ , tutte le facce di σ appartengono

a Γ.

• per ogni coppia di simplessi σ, σ’ appartenenti a Γ, tali che σ∩σ′ 6= ∅,σ ∩ σ′ e una faccia di entrambi (σ e σ′).

• d e il massimo ordine di simplessi appartenenti a Γ.

Una mesh tetraedrale e un complesso simpliciale di dimensione 3. Ogni faccia

di tetraedro (triangolo) appartene alla mesh. L’intersezione tra 2 tetraedri

adiacenti se non e nulla e una faccia del tetraedro o un edge (segmento)

oppure un punto, in tutti i casi sono facce del tetraedro stesso. In una mesh

di tetraedri non esistono simplessi di ordine maggiore di tre.

2.3 Utilizzo di mesh tetraedrali nella rappre-

sentazione dei modelli

Si possono pensare diversi modi per modellare organi con i quali interagire.

La strategia piu utilizzata e quella di rappresentare gli organi come mesh

tetraedrali irregolari. Tali mesh, sebbene difficili da generare, offrono notevoli

vantaggi:

2.4 Modifiche geometriche e topologiche 21

• Le rappresentazioni basate su dati volumetrici sono superiori in termini

di qualita rispetto agli approcci basati sulla modellazione di superfici

in quanto la deformazione calcolata risulta piu realistica.

• Rispetto ai sistemi basati sui voxel utilizzando mesh tetraedrali la

rappresentazione della superficie e piu accurata.

• Il tetraedro e l’oggetto volumetrico con la topologia piu semplice. In

caso di tagli sul tetraedro stesso la topologia generata risulta la piu

semplice da gestire.

2.4 Modifiche geometriche e topologiche

Il taglio di tessuto e una delle attivita principali effettuate nella chirurgia, il

risultato delle modifiche geometriche conseguenti ad un taglio puo incidere

in modo rilevante sul realismo del processo di taglio stesso e sulla comples-

sita del modello deformabile. Ci sono alcune caratteristiche fondamentali da

rispettare durante la simulazione di un taglio affinche esso risulti realistico:

Accuratezza I tagli devono essere piu fedeli possibile alle forme tracciate

con l’haptic-device dall’utente. Non bisogna costringere mai l’uten-

te a tagliare lungo posizioni prestabilite (come se si muovesse lungo

una griglia) , perche il senso di realismo ne risulterebbe gravemente

compromesso.

Interattivita L’intervallo di tempo che intercorre tra quando l’utente esegue

il taglio e quando questo viene rappresentato non deve essere percepi-

bile.

Frammentazione Le operazioni di taglio tendono ad aumentare la com-

plessita della descrizione del modello, per mantenere il rate di aggior-


Figura 2.1: a) La superficie di taglio compresa tra due istanti successivi t0 e t1. b)

La triangolazione della superficie di taglio.

namento della simulazione costante, il numero dei nuovi elementi creati

deve essere minimo.

2.4.1 Il taglio di una mesh tetraedrale

Prima di definire come si possono effettuare tagli dobbiamo stabilire una

rappresentazione dei bisturi nello spazio; la letteratura esistente a riguardo

descrive il movimento del bisturi come:

Definizione 6 La posizione del bisturi nello spazio e definita da un

segmento (tail,head) che definisce la linea della lama.

Definizione 7 La superficie di taglio St′,t′′ e la superficie (tail, head)t, t ∈ [t′, t′′];

ovvero la parte di spazio spazzata dalla lama nell’intervallo di tempo.

Durante la simulazione la posizione del bisturi e conosciuta ad intervalli di-

screti ∆ t. La superficie di taglio tra due istanti consecutivi e ricavata come


Figura 2.2: a) Approccio con distruzione di tetraedri. b) Approccio con duplicazione

delle facce. c) Approccio con suddivisione di tetraedri.

interpolazione lineare dei due estremi del bisturi. Nel caso pratico, essen-

do l’intervallo temporale piccolo, questo piano puo essere descritto con due

triangoli quasi complanari. La superficie di taglio St′,t′′ sara cosı descritta da

una serie di triangoli (Figura 2.1) .

2.4.2 Modifiche topologiche nella mesh tetraedrale

Gli approcci proposti finora si possono raggruppare in 3 tipologie:

• Approccio con distruzione di tetraedri (Figura 2.2.a) .

• Approccio con duplicazione di facce (Figura 2.2.b) .

• Approccio con suddivisione di tetraedri (Figura 2.2.c) .


2.4.3 Approccio con distruzione di tetraedri

Il taglio viene effettuato eliminando i tetraedri che intersecano con il tool[FDA02].

Cancellando i tetraedri il tempo di calcolo della simulazione naturalmente

diminuisce; come conseguenza pero tale metodo viola i principi fisici di con-

servazione della massa.

Se si utilizza l’approccio con distruzione dei tetraedri per ottenere un effetto

realistico la mesh deve essere molto dettagliata. Alcuni approcci propongono

di dettagliare maggiormente la mesh nei punti dove si effettuano i tagli (per

fare cio dovremmo conoscere in anticipo dove i tagli avverranno) .

2.4.4 Approccio con duplicazione di facce

Il taglio viene effettuato duplicando i nodi, e di conseguenza le facce, che sono

vicini al bisturi[Nie03, MLBdC01], rimuovendo poi le relazioni di adiacenza

lungo il taglio:

• I nodi della mesh piu vicini alla traiettoria di taglio sono spostati in

modo tale che le loro facce interne siano lungo la traiettoria di taglio.

• I vertici che sono stati spostati vengono duplicati ed associati alle

relative facce interne.

Diversamente dall’approccio precedente non c’e perdita di volume, e il nu-

mero di nodi rimane ancora pressoche costante.

La duplicazione dei nodi comporta la modifica della matrice di rigidita so-

lamente dei nodi duplicati; la localita delle modifiche rende questo metodo

facilmente integrabile con i metodi di analisi ad elementi finiti. Muovere i

nodi nella mesh significa cambiare la posizione di energia minima, invalidan-

do cosı la simulazione per tali nodi. Piu grandi saranno le distanze tra le

particle e maggiore sara lo spostamento da effettuare; tale spostamento puo

provocare in alcuni casi comportamenti anomali.


2.4.5 Approccio con suddivisione di tetraedri

Il taglio viene effettuato suddividendo i tetraedri che sono intersecati dal

tool[GCMS00, BMG99].

Definizione 8 Dato un tetraedro σ ∈ Γ e una superficie di taglio S tale che

S ∩ σ 6= ∅, chiamiamo sostituzione valida una sostituzione di σ con un

insieme di elementi T tali che:

• ∀σ′ ∈ T, S ∪ σ′ ⊂ δσ′: non c’e intersezione propria tra la superficie dei

taglio ed ognuno dei tetraedri di T.

• Ω(σ) =⋃

σ′∈T Ω(σ′): Il volume occupato dall’insieme dei tetraedri in T

e lo stesso volume che occupava σ.

• Γ \ σ ∪ T e un complesso simpliciale.

Definizione 9 L’operatore di taglio e una funzione tale che, utilizzando

una superficie di taglio S ed una mesh Γ, applica una azione di sostituzione

valida a tutti i tetraedri σ ∈ Γ intersecati da S

In questo metodo assumiamo innanzitutto che che la superficie di taglio

interseca un tetraedro se e solo se taglia uno dei suoi edge, inoltre se appros-

simiamo S∩σ con una superficie, la suddivisione del tetraedro sara piu facile

da trovare che nel caso generale.

Suddivisione di un tetraedro Date le precedenti assunzioni ci sono solo

5 modi in cui gli edge di un tetraedro possono essere intersecati da una

superficie di taglio, queste sono chiamate configurazioni di taglio ( Figura

2.3 ).

Nelle prime due configurazioni S interseca rispettivamente 3 e 4 spigoli,

dividendo il tetraedro in due parti, mentre negli altri casi il tetraedro e inter-

secato senza essere diviso completamente. Questa soluzione mira a cercare


Figura 2.3: Suddivisioni possibili di un tetraedro

la suddivisione valida che crei il minor numero possibile di tetraedri.

Per fare questo si utilizza una tabella, chiamata look up table che viene in-

dirizzata dagli edge intersecati dalla superficie di taglio. Per ogni possibile

configurazione di taglio scriviamo nella lookup table quali vertici andranno

inseriti e quali facce del tetraedro andranno divise, preservando per costru-

zione le proporieta del complesso simpliciale. La terza condizione per una

sostituzione valida puo essere riscritta come:

• T e un complesso simpliciale.

• T ∩ (Γ \ σ) e un complesso simpliciale.

Le facce del tetraedro diviso dal piano di taglio devono essere suddivise

in triangoli (si puo notare che nelle divisioni mostrate nell’esempio di Figura

2.3 si crea sempre almeno faccia che non e piu un triangolo) . In pratica, il

problema resta nei bordi della nuova triangolazione che sono condivisi con

il resto della mesh. Affinche le proprieta del complesso simpliciale vengano


Figura 2.4: a) La superficie di taglio raffigurata in puo essere triangolarizzata nelle

due maniere riportate in b) . c) Una faccia che non e completamente

divisa d) Faccia in cui due spigoli vengono divisi

preservate i tetraedri che hanno facce in comune con il tetraedro suddiviso

andranno a loro volta suddivisi in modo tale che l’intersezione tra di essi sia

ancora una volta una faccia. Le Figure 2.4.a 2.4.b e mostrano un caso prati-

co: supponiamo che la superficie di taglio intersechi due tetraedri adiacenti

σ1 e σ2 nei punti a,b,c,d. Dopo due sostituzioni valide i due tetraedri sono

sostituiti rispettivamente dagli insiemi T1 e T2. Siccome la faccia (p1, p2, p3)

e condivisa da σ1 e σ2 le nuove facce (a, b, p2) e (p1, a, b, p3) andranno trian-

golarizzate nella giusta maniera. Nell’esempio illustrato due triangolazioni

sono possibili per la faccia (p1, a, b, p3), introducendo lo spigolo (p1, b) oppure

(a, p3).

In generale il problema viene risolto garantendo una triangolazione fissa

in questi casi. Osserviamo che la faccia condivisa f puo trovarsi solo in due

2.5 Metodi di integrazione numerica 28

situazioni:

• Solo uno spigolo di f viene diviso.

• Due spigoli di f vengono divisi.

Nel primo caso viene generato un poligono a 6 vertici, esso viene triango-

larizzato utilizzando le diagonali che connettono il nuovo vertice a con i tre

vertici p1, p2, p3. Nel secondo caso viene generata una faccia triangolare ed

un quadrilatero. Utilizzando un ordine totale dei vertici sceglieremo in modo

univoco la stessa diagonale. A seconda dell’ordine in cui sono messi i vertici

possiamo quindi avere differenti divisioni per la stessa configurazione di ta-

glio. Alla fine otterremo una lookup table di 33 entrate. Chiaramente ogni

caso va considerato per ogni possibile rotazione del tetraedro ( ad esempio il

caso A genera 4 differenti tipi di taglio, ognuno per ogni possibile rotazione

del tetraedro) .

2.5 Metodi di integrazione numerica

I vertici che compongono l’oggetto deformabile si muovono in relazione alle

forze interne scaturite durante la deformazione. Utilizzando le equazioni

del moto possiamo essere in grado di stabilire, conoscendo le forze interne

al modello, la posizione di ciascun vertice v in un istante di tempo t. Per

risolvere le equazioni del moto dobbiamo risolvere delle equazioni differenziali

di secondo grado.

Ogni equazione differenziali con grado maggiore di uno puo essere riscritto

in un sistema di equazioni di primo ordine, il problema si riconduce quindi

alla soluzione del sistema[Witb]:

x(t) = φ(x, t) (2.20)


Figura 2.5: In alto a sinistra viene raffigurato un esempio di taglio di tipo A: i tre

spigoli intersecati dalla superficie di taglio producono il numero 42 utiliz-

zato per accedere alla lookup table che definisce l’insieme dei 4 tetraedri

utilizzati per la sostituzione. In basso a destra viene raffigurato lo schema

necessario all’interpretazione del numero dei vertici restituiti dalla tabella;

infine in basso a sinistra la decomposizione nei 4 tetraedri.


I metodi piu utilizzati per la risoluzione numerica di questi sistemi id

equazioni sono i metodi di Runge Kutta.

2.5.1 Metodi di Runge Kutta

L’idea base e quella di esprimere x(t + h) utilizzando la serie di Taylor :

x(t + h) = x(t) + hx′(t) +h2

2!x′′(t) +

h3

3!x′′′(t) +

hn

n!

dnx

dtn(2.21)

Assumendo che x(t) sia C∞, possiamo troncare questa serie ad un dato

termine, ad esempio il secondo:

x(t + h) = x(t) + hx′(t) + O(h2) (2.22)

Otterremmo cosı una approssimazione di x(t + h) con un errore O(h2).

Chiamiamo φ la funzione che calcola x′ possiamo scrivere la formula prece-

dente come:

x(t + h) = x(t) + hφ(x, t) + O(h2) (2.23)

Questa approssimazione del metodi di Runge-Kutta al primo ordine vie-

ne chiamato anche metodo di Eulero. Il metodo di Runge-Kutta al secondo

ordine viene ottenuto troncando la serie al terzo termine della serie. L’appros-

simazione di x(t + h) generera un errore O(h2) . Chiaramente maggiore sara

l’ordine e piu accurato sara il risultato, d’altra parte il numero di valutazioni

di φ(x, t) e un punto cruciale nella complessita computazionale. Questi me-

todi vengono detti espliciti in quanto cio che non e conosciuto viene espresso

come soluzione di un qualcosa che sappiamo gia.


2.5.2 Metodi impliciti

Diversamente dai metodi espliciti nei metodi impliciti esprimiamo cio che

non e conosciuto come soluzione di un sistema di coppie di equazioni[Bar].

Per le equazioni differenziali questo significa esprimere il valore x(t + h)

utilizzando φ(x(t + h), t + h), in particolare il passo di Eulero puo essere

riscritto come:

x(t + h) = x(t) + hφ(x(t + h), t + h) + O(h2) (2.24)

dove il termine φ(x(t + h), t + h) viene approssimato attraverso la serie

di Taylor troncata al primo termine:

φ(x(t + h), t + h) = φ(x(t), t) + hφ′(x(t), t) (2.25)

Consideriamo adesso che x e una variabile singola appartenente all’insie-

me dei numeri reali, x(t) e φ(x, t) sono vettori, mentre φ(x, t) e una matrice.

Sostituendo 2.25 in 2.24 e risolvendo per x(t + h):

x(t + h) = x(t) + (1

hI − φ′(x(t), t))−1φ(x(t), t) (2.26)

Questo metodo richiede la computazione di una matrice inversa (o equi-

valentemente risolvere un sistema lineare) ad ogni passo della simulazione.

Questo rende i metodi impliciti piu complessi in termini di numero di opera-

zioni rispetto ai metodi espliciti per un singolo passo di integrazione. D’altra

parte sono metodi che possono essere utilizzati in contesti pratici in cui il

valore di h e elevato. Entrambi i metodi (espliciti ed impliciti) fanno parte

della classe di metodi chiamati single step, ovvero i metodi che esprimono

il valore di x(t + h) usando solamente lo stato dello stato precedente (ovve-

ro x(t)) . Al contrario i metodi multistep utilizzano una interpolazione dei

precedenti n punti:


x(t + h) = x(t) + hφ(β0x(t + h) + β1x(t) + . . . + βnx(t− (n− 1)h)) (2.27)

Possiamo notare come questa formulazione esprime anche i metodi single

step (in questo caso i = 1) sia espliciti (β0 = 0) che impliciti (β0 6= 0).

Capitolo 3

Tecniche di taglio su mesh tetraedrali

Le tecniche di taglio su mesh tetraedrali presentate nel capitolo precedente

curavano ognuna un singolo aspetto del problema:

• La tecnica di cancellazione di tetraedri riduce il numero di tetraedri

della mesh ma non preserva il volume e non rappresenta il taglio in

maniera accurata.

• Nell’approccio con duplicazione delle facce e dei nodi viene modifica-

ta la posizione di equilibrio della mesh e di conseguenza la matrice di

rigidita. In alcuni casi grandi spostamenti di vertici possono rivelare al-

l’utente cambiamenti nella forma della mesh, degradando notevolmente

il fattore di realismo.

• La tecnica di rimpiazzamento dei tetraedri, nella sua precisione di rap-

presentazione dell’azione di taglio, ha lo svantaggio di creare un numero

elevato di entita incrementando la complessita di calcolo.

L’aumento di complessita nella descrizione della mesh puo causare un

ritardo percepibile, dovuto alle operazioni di aggiornamento, tra l’azione di

taglio da parte dell’utente e la divisione del tessuto.

Nella prima sezione di questo capitolo illustreremo la struttura della mesh

tetraedrale e definiremo le funzioni necessarie alla sua navigazione. Nella

3.1 Strutture dati e algoritmi per mesh tetraedrali 34

Figura 3.1: Convenzione sulla numerazione locale dei vertici degli edge e delle facce:

a) Orientamento del tetraedro: La numerazione locale dei vertici e tale

che: [(T (V1)− T (V0))× (T (V3)− T (V0))] · (T (V2)− T (V0) > 0. b) Le

facce associate alle triple di vertici. c) Gli edge associati alle coppie di

vertici.

seconda sezione invece verranno definite le operazioni modifica della mesh

utilizzate nella modellazione del taglio. Infine, nella terza sezione viene illu-

strato il metodo di semplificazione della mesh utilizzato nello svolgimento di

questa tesi.

3.1 Strutture dati e algoritmi per mesh te-

traedrali

Le entita principali che compongono una mesh di tetraedri sono i vertici, gli

edge, le facce ed i tetraedri. Una struttura dati costruita per rappresentare


una mesh di tetraedri deve rappresentare queste entita e le relazioni di appar-

tenenza che intercorrono tra esse. Questa struttura dati deve essere costruita

in modo tale da mantenere la consistenza della connettivita delle entita della

mesh (questa caratteristica puo permettere ad esempio di conoscere tutti i

tetraedri che incidono in un edge). Piu in generale possiamo dire che una

caratteristica fondamentale della struttura dati e la possibilita di effettuare

query sui dati rappresentanti la mesh. Il contesto del problema impone inol-

tre che questa struttura dati possa essere aggiornata in tempo reale a causa

delle modifiche geometriche della stessa (deformazioni, tagli, semplificazioni)

.

L’efficienza di questa struttura dati si puo quindi esprimere in termini di:

• Utilizzo di memoria.

• Efficienza nelle esecuzione delle query.

• Efficienza nelle modifiche della connettivita.

Una mesh di tetraedri viene descritta da due entita esplicite fondamen-

tali: Il vertice ed il tetraedro. Ogni tetraedro viene descritto dalla posizione

dei quattro vertici nello spazio. Tali vertici devono essere condivisi in modo

appropriato dai tetraedri adiacenti. Questa condizione e di vitale importan-

za in quanto la modifica della posizione di un vertice della mesh durante la

sua deformazione si deve riflettere su tutti i tetraedri che incidono sul vertice

stesso.

In dettaglio una mesh contiene un insieme di vertici e di tetraedri; ogni te-

traedro e dotato di quattro puntatori, ciascuno per ogni vertice; in questo

modo l’istanza di ogni vertice risulta unica e condivisa dagli opportuni te-

traedri. Per la completezza della struttura stessa non esisteranno vertici che

appartengano a nessun tetraedro. I puntatori a vertice in ogni tetraedro sa-

ranno sistemati in modo tale che il vertice 0 veda i vertici 1 2 e 3 in senso


antiorario. In particolare varra il seguente criterio: Sia T (Vi) il vertice i-

esimo del tetraedro T , allora affinche l’ordinamento dei vertici sia corretto

dovra valere il criterio (Figura 3.1.a) :

[(T (V1)− T (V0))× (T (V3)− T (V0))] · (T (V2)− T (V0) > 0 (3.1)

Se questo criterio non viene rispettato il vertice 2 ed il vertice 0 verranno

scambiati di posizione.

Le facce e gli spigoli del tetraedro esistono in quanto facce del tetraedro

descritte dai vertici. Risulta superfluo e complica la coerenza nell’aggiorna-

mento dei dati creare oggetti che modellino esplicitamente queste entita. Ai

fini della nostra applicazione una faccia o un edge sara solo l’insieme dei suoi

vertici. Le tabelle di Figura 3.1 riportano l’ordine dei vertici per ogni faccia

e l’ordine dei vertici per ogni spigolo.

Seguendo queste tabelle possiamo dedurre ad esempio che la faccia 1 di un

tetraedro e quella descritta dall’insieme dei vertici T (V0), T (V1), T (V3), op-

pure che l’edge 1 e quello che va dal vertice T (V0) a T (V2). Naturalmente

l’ordine in cui sono specificati i vertici sulle facce per ragioni di correttezza

e antiorario, mentre su uno spigolo l’ordine risulta superfluo. Da queste due

tabelle fondamentali si possono costruire altre tabelle in modo da fare un

numero minore di accessi possibile per ogni ricerca. La tabella 3.1.a riporta

le facce incidenti su ogni vertice, la 3.1.b gli spigoli incidenti su ogni vertice,

la 3.1.c gli spigoli di ogni faccia, la 3.1.d le facce per ogni spigolo. Queste

tabelle vengono spesso utilizzate nella navigazione della mesh.

3.1.1 Relazioni di adiacenza

Adiacenza tetraedro-tetraedro (via faccia) Se il tetraedro Ti condivi-

de tre vertici col tetraedro Tj, ci sara una faccia F , descritta da questi tre


Vertici F0 F1 F2

0 0 1 2

1 0 1 3

2 0 2 3

3 1 2 3

Vertici E0 E1 E2

0 0 1 2

1 0 3 4

2 1 3 5

3 2 4 5

a b

Facce E0 E1 E2

0 0 3 1

1 2 4 0

2 1 5 2

3 4 5 3

Edge F0 F1

0 0 1

1 0 2

2 1 2

3 0 3

4 1 3

5 2 3

c d

Tabella 3.1: a) Tabella vertice-faccia. b) Tabella vertice-edge. c) Tabella faccia-edge.

d) Tabella edge-faccia.

Figura 3.2: Un primo esempio di mesh tetraedrale: a sinistra vengono messi in risalto

i tetraedri che compongono la mesh


vertici che apparterra ad entrambi i tetraedri (Ti,Tj) . Per navigare all’inter-

no della mesh diventa fondamentale definire una funzione di adiacenza tra

tetraedri che condividono la stessa faccia.

Definizione 10 Dato un tetraedro Ti ed una faccia Fj ∈ Ti, la funzione di

adiacenza tetraedro-tetraedro TT (Ti, Fj) restituisce il tetraedro Tk 6= Ti,

tale che Fj ∈ T (k), se T (k) non esiste allora TT (Ti, Fj) = Ti.

L’ultimo caso della definizione precedente e relativo alle facce esterne del-

la mesh. Le facce esterne (quindi quelle che dovranno essere disegnate nella

fase di rendering) saranno quindi caratterizzate dal fatto che la funzione di

adiacenza tetraedro tetraedro restituira lo stesso tetraedro che figura come

paramento della funzione, cioe se TT (Ti, Fj) = Ti allora la faccia Fj del te-

traedro Ti e esterna.

Affinche la funzione possa essere invertibile dobbiamo definire un’altra fun-

zione:

Definizione 11 La funzione ZT (Ti, Fj) restituisce la faccia Fh tale che esi-

ste un tetraedro Tk 6= Ti tale che TT (Tk, Fh) = Ti

In altre parole la funzione ZT restituisce l’indice della faccia Fj del te-

traedro Ti vista dal tetraedro adiacente attraverso tale faccia. Nel caso di una

faccia esterna tale funzione avra valore uguale all’indice della faccia stessa.

Date queste due funzioni possiamo definire la proprieta invariante:

Proprieta 1 Data la mesh tetraedrale TM composta da tn tetraedri, allora:

∀(Ti, Fj)(i = 0 . . . tn), (j = 0 . . . 3), TT (TT (Ti, Fj), ZT (Ti, Fj)) = Ti


Figura 3.3: a) I tetraedri T0 e T1 condividono la stessa faccia b) TT (T0, F2) = T1;

TT (T1, F3) = T0; ZT (T0, F2) = F3; ZT (T1, F3) = F2

Implementazione Per implementare tali funzioni e necessario utilizzare

strutture dati che mantengano nel tempo tale connettivita. Queste strutture

dati verranno inizializzate al momento della creazione della mesh stessa e

modificate localmente in caso di modifiche topologiche (tagli o collassi) .

Ogni istanza di tetraedro conterra un array di quattro puntatori a tetrae-

dro; ogni puntatore specifica l’adiacenza attraverso la faccia indicata della

posizione dell’array; naturalmente una faccia esterna e caratterizzata da un

puntatore allo stesso tetraedro. La funzione ZT verra anch’essa rappresen-

tata attraverso l’uso di un array di quattro interi, la faccia esterna sara

caratterizzata in questo caso dall’indice della faccia stessa.

Adiacenza vertice-tetraedro In alcune operazioni di modifica della mesh

risulta necessario conoscere tutti i tetraedri che condividono un vertice. De-

finire una funzione del genere sulla base della adiacenza tetraedro-tetraedro

e poco efficiente. Data una mesh tetraedrale TM di tn tetraedri e vn vertici,

allora:

Definizione 12 La funzione di adiacenza vertice-tetraedro V T (Vi), (i =


0 . . . vn) restituisce l’insieme dei tetraedri Tj(0 ≤ j < tn) tale che Vi ∈ Tj

Implementazione della topologia vertice-tetraedro La funzione di

adiacenza vertice-tetraedro mette in relazione una entita esplicita (il vertice)

con in insieme di entita (i tetraedri che lo condividono) , questo dettaglio

rende l’implementazione di questa funzione piu problematica rispetto all’im-

plementazione della topologia tetraedro-tetraedro dove la relazione coinvolge

al massimo due tetraedri. In questo caso risulta necessario costruire strutture

dati ausiliarie che mantengano, per ogni vertice, il sottoinsieme di tetraedri

che incidono su di esso. Tale sottoinsieme non richiede un ordinamento par-

ticolare; le operazioni piu frequenti su questa struttura dati saranno la visita

completa e l’inserimento di nuovi elementi.

La struttura dati scelta per modellare questo tipo di relazione e la lista. Ogni

vertice delle mesh contiene il puntatore al primo tetraedro che compone la

lista e un intero che specifica l’indice di quel vertice all’interno del tetraedro.

Ogni tetraedro conterra sua volta due array di quattro elementi, ciascuno per

ogni vertice:

L’array Tp Per ogni posizione j (j = 0 . . . 4), un puntatore indichera il

prossimo tetraedro che condivide il vertice j-esimo.

L’array Zp Per ogni posizione j (j = 0 . . . 4), un intero indichera l’indice

del vertice j-esimo sul prossimo tetraedro della lista (Tp[j]) . tetraedro

L’inserimento di un nuovo tetraedro alla lista comportera quindi solo una

modifica alla testa della lista stessa, cioe l’istanza del vertice stesso (oltre alla

modifica degli array Tp e Zp del tetraedro da inserire) . Va precisato inoltre

che la fine della lista viene contrassegnata da un puntatore nullo nell’array

Tp ed un indice -1 nello Zp.

template < class VERT_TYPE > class TETRA_TYPE

public:


a) b)

Figura 3.4: a) I tetraedri T0, T1 e T2 condividono lo stesso vertice V; V T (V ) =

T0,T1,T2 b) Implementazione della lista per il caso raffigurato

...

/// The vertex type

typedef typename VERT_TYPE MVTYPE;

typedef typename VERT_TYPE::scalar_type scalar_type;

MVTYPE *v[4];

...

///pointers to tetra for vertex-tetra topology

TETRA_TYPE *tv[4];

short int zv[4];

TETRA_TYPE *t[4];

int z[4];

...

scalar_type _volume;

;//end class


3.1.2 Operatori di navigazione sulla mesh

Utilizzando le strutture di connettivita appena specificate sono stati defini-

ti alcuni operatori di navigazione sulla mesh. Gli operatori che definiremo

permettono di eseguire delle query (ad esempio se si volesse conoscere tutti i

tetraedri che condividono un vertice, oppure un edge. . . ) in modo efficiente

e astratto dalla struttura dati sottostante.

Data la frequenza e l’importanza d’utilizzo di tali oggetti nell’implementa-

zione delle operazioni di modifica che seguono, risulta utile descrivere la loro

struttura in dettaglio.

Iteratore vertice-tetraedro: NextV Risulta utile all’utente poter na-

vigare attraverso la connettivita di vertice-tetraedro in modo trasparente

rispetto alla sua rappresentazione in termini di strutture dati. Creare un og-

getto dedito a tale funzione, significa creare l’iteratore vertice-tetraedro che

scorre la lista dei tetraedri che condividono uno stesso vertice.

Un oggetto di questo tipo conterra un puntatore al tetraedro corrente ( t ) ed

l’indice del vertice sul tetraedro stesso. Inizializzando questo iteratore con i

valori dell’indice e del tetraedro presenti sul vertice stesso (relativamente Zp e

Tp del vertice) possiamo scorrere tutta la catena utilizzando una operazione

NextV cosı definita:

tmp=t;

t=tmp->Tp[z];//il puntatore al tetraedro viene aggiornato

al successivo nella lista

z=tmp->Zp[z];//l’indice z viene aggiornato coll’indice del

vertice sul nuovo tetraedro

Una importante proprieta invariante caratterizza questo operatore:

Proprieta 2 Siano t0 e z0 i valori corrispondenti a z e t dell’iteratore vertice

tetraedro e sia V il vertice di indice z0 del tetraedro t0; siano poi t1 e z1 i


valori di z e t dopo una operazione di Next e sia V 1 il vertice di indice z1 del

tetraedro t1, allora:

• o z1 = −1 e t1 = NULL,

• oppure V 1 = V

In pratica il vertice indicato dall’iteratore sul tetraedro deve essere sem-

pre lo stesso, oppure diventa nullo nel caso in cui non ci siano piu tetraedri

da visitare. Questa proprieta e garantita naturalmente solo nel caso in cui sia

corretta la costruzione della topologia vertice-tetraedro. Nel caso di modifi-

che della mesh un test che sfrutti tale proprieta invariante, iterato su tutte

le catene, si e rivelato molto utile nella individuazione di alcuni errori.

Half-edge L’half-edge e un oggetto utilizzato nella navigazione di una me-

sh. Come spiega la parola stessa (meta edge) l’half-edge rappresenta un edge

orientato. Come rappresentato in Figura 3.5 in un complesso simpliciale a 1

dimensione (spezzata) rappresentare un edge orientato significa arricchire la

descrizione di un edge con qualcosa in piu che ne identifichi la direzione. I

due half-edge rappresentati, pur riferendosi allo stesso edge, sono differenti.

Potremmo cosı esprimere l’half-edge come una coppia (edge, vertice) nel no-

stro esempio di complesso simpliciale a una dimensione. Naturalmente gli

half-edge sono rappresentati solo da un sottoinsieme di tutte le possibili cop-

pie edge, vertice; cioe quelle in cui e il vertice e incluso nell’edge specificato

nella coppia. Attraverso le funzioni di adiacenza, possiamo definire adesso

delle funzioni di iterazione basate su questo tipo di concetto. Nella naviga-

zione della mesh potremmo cosı conoscere per ogni edge i vertici (e viceversa)

; inoltre operatori di iterazione permetterebbero la transizione dell’half-edge

da un edge all’altro oppure da un vertice all’altro.

Nella Figura 3.5 il concetto di half-edge viene anche espresso in un comples-

so simpliciale a due dimensioni (mesh di triangoli). In questo caso l’half-edge


Figura 3.5: L’half-edge viene rappresentato graficamente come un triangolo rettan-

golo: a) I due half-edge He0 e He1 pur essendo definiti sullo stesso edge

hanno vertici differenti: He0 = (b, Edge1) mentre He01 = (c, Edge1)

b) L’half-edge definito sulla mesh di triangoli He2 = (l, Edge11, F5)

definito precedentemente viene esteso ad una tripla faccia, edge, vertice;nell’esempio viene mostrato in particolare l’half-edge l, Edge11, F5. Come

nel caso precedente solo un sottoinsieme di queste triple assume un signifi-

cato coerente con la mesh.

In generale un half-edge in un complesso simpliciale a n dimensioni puo essere

rappresentato come una di (n+1)-tupla dove:

Definizione 13 1. Hedge = Fn, Fn−1, . . . , F0, Fi e faccia del complesso

simpliciale a n-dimensioni.

2. F0 ∈ F1 ∈ F2 . . . Fn−1 ∈ Fn. Dove Fi ∈ Fi+1significa che Fi e una

faccia di Fi+1

Nel nostro caso, essendo la mesh tetraedrale (quindi un complesso sim-

pliciale a 3 dimensioni) , l’half-edge sara rappresentato da una quadrupla


Figura 3.6: L’half-edge viene rappresentato graficamente come un triangolo rettan-

golo con un vettore indicante il tetraedro su cui e definito. In questo caso

l’half-edge e costituito dalla quadrupla (o,Edge3, F3, T1)


tetraedro, faccia, edge, vertice. La condizione 2 implica che, per ogni pos-

sibile half-edge, il vertice specificato dovra appartenere all’edge, il quale a sua

volta sara parte della faccia anch’essa appartenente al tetraedro specificato.

L’half-edge rivela la sua vera importanza se utilizzato come base per co-

struire una serie di operatori di navigazione all’interno della mesh. Nell’im-

plementazione della quadrupla half-edge nella mesh tetraedrale, la faccia e

l’edge saranno due interi, mentre il vertice ed il tetraedro saranno puntatori

all’istanza del vertice e del tetraedro stesso.

Operazioni base sugli half-edge: flip Dato un half-edge Ti, Fj, Eh, Vk,possiamo notare che per ogni elemento della quadrupla esiste una sola pos-

sibile variazione che generi un half-edge come definito in 2. Ciascuna di

queste variazioni comporta la modifica di un solo elemento della quadrupla.

Ad esempio, come mostrato in Figura 3.7, modificando il vertice passiamo

dall’half-edge 0 all’half-edge 1: la successiva modifica del vertice riportera

l’half-edge nello stato 0. Le due configurazioni ottenute (0 e 1) sono le uni-

che due che possono essere ottenute variando solo il vertice nella quadrupla

dell’esempio. Una operazione di questo tipo viene chiamata flip dell’half-edge

sul vertice.

Effettuando la stesso tipo di operazione cambiando l’edge della quadrupla 0

(flip dell’half-edge sull’edge) passiamo dallo stato 0 allo stato 2, mentre il flip

sulla faccia dallo stato 0 al 3, infine il flip sul tetraedro dallo stato 0 al 4.

Se la mesh tetraedrale e una sola componente connessa, allora si puo provare

come, tramite successioni di questi operatori si possa di fatto visitare tutti

i possibili half-edge. La quantita di informazione presente dell’half-edge e

le operazioni di flip costituiscono il punto di partenza nella costruzione di

iteratori di navigazione della mesh. Le operazioni di flip appena esposte si

basano esclusivamente sulla topologia tetraedro-tetraedro.


Figura 3.7: Possibili flip su una mesh di tetraedri: Partendo dallo stato 0 arriviamo

all’half-edge dello stato1 tramite flip sul vertice e viceversa. Ogni transi-

zione di stati viene raggiunta attraverso l’operazione indicata sulla freccia

(ad ogni operazione corrisponde anche l’inversa).


Figura 3.8: Possibili Next su una mesh di tetraedri: a) Il next su una faccia corri-

sponde ad un flip sul tetraedro. b) Il next sull’edge restituisce in modo

ordinato tutti i tetraedri che condividono un edge.

Operazioni derivate sugli half-edge: Next Combinando tra loro le

operazioni di flip possiamo ottenere iteratori che permettono la visitare tutti

i possibili tetraedri della mesh. Operazioni di questo tipo sono chiamate

Next. In particolare possiamo implementare due tipi di Next significativi su

una mesh tetraedrale:

Definizione 14 (NextFace) L’half-edge viene ribaltato sul tetraedro che

condivide la stessa faccia. Come mostrato in Figura 3.8.a, si passa dallo

stato 0 allo stato 1 e viceversa; in questo caso, siccome esistono solo due te-

traedri che condividono la stessa faccia questa operazione si riduce ad un flip

sul tetraedro. Se nessun tetraedro condivide la faccia l’operazione di NextFace

non produrra nessun risultato.


Figura 3.9: L’operazione di Next edge equivale ad effettuare una flip sul tetraedro

seguita da un flip sulla faccia.

NextFace = FlipTetra

Definizione 15 (Next edge) Questa operazione restituisce il prossimo half-

edge con tetraedro differente che condivide lo stesso edge.

Come mostrato in Figura 3.8.b, iterando questa operazione si possono ot-

tenere, ad ogni passo, l’insieme ordinato dei tetraedri che condividono tale

edge. L’ordinamento viene garantito da fatto che questo tipo di Next viene

effettuato combinando un flip sul tetraedro ed un flip sulla faccia. Come mo-

strato in Figura 3.9 la prima di queste due operazioni sposta l’half-edge sul

tetraedro che condivide la stessa faccia (e quindi lo stesso edge), mentre il

secondo flip (quello sulla faccia) muove l’half-edge sull’altra faccia dell stesso

tetraedro che condivide l’edge: In sostanza la prima operazione consente di

ottenere il prossimo tetraedro, mentre la seconda crea le condizioni per arri-

vare sul prossimo tetraedro (senza che si torni allo stato iniziale come nelle

operazioni di flip) .

Nextedge = FlipFace(FlipTetra)

Ciascuna di queste due operazioni di iterazione, essendo basate sulle ope-

razioni di flip, utilizzano esclusivamente le funzioni di adiacenza tetraedro-

etraedro; questo dettaglio si rivelera molto importante in un contesto di

3.2 Modifiche topologiche: tagli 50

modifica: dobbiamo essere in in grado di visitare un gruppo di tetraedri

modificandoli, senza che tale modifica influisca sulla visita stessa.

3.2 Modifiche topologiche: tagli

Prima di introdurre la strategia di taglio utilizzata conviene soffermarsi e

spiegare piu nel dettaglio i punti deboli delle altre strategie utilizzate.

Il punto di forza della strategia con duplicazione di facce consiste nel poter

effettuare tagli senza introdurre nuovi tetraedri all’interno della mesh. Per

contro, e richiesto che alcuni vertici vengano artificiosamente spostati sull’in-

tersezione tra il bisturi e la superficie, questa operazione puo dar luogo ad

anomalie. Consideriamo l’esempio illustrato in Figura 3.10: secondo l’imple-

mentazione dovremmo proiettare il vertice 0 (che e il piu vicino alla lama)

sulla superficie di taglio; in questo caso il taglio non potrebbe essere effettua-

to in quanto la superficie di taglio non interseca piu la mesh e sdoppiare il

vertice 0 non avrebbe quindi piu senso. Si potrebbe pensare a questo punto

di avvicinare il vertice 1 sulla superficie di taglio; questo lungo spostamento

puo risultare pero una modifica troppo grande e percepibile.

Nell’implementazione con rimpiazzamento di tetraedri si assume di poter

suddividere un tetraedro (ed eventualmente quelli a lui adiacenti) basandosi

sul movimento del bisturi all’interno del tetraedro stesso. Come mostrato

nel capitolo 2, il gruppo di tetraedri introdotto per rimpiazzare uno tagliato

dipende da quali edges di quest’ultimo sono intersecati dal tool. Questa

dipendenza crea due problemi:

1. Occorre gestire la mesh nell’intervallo di tempo in cui il tool ha attraver-

sato alcuni degli edges del tetraedro senza attraversarlo completamente

(e quindi non e stata individuata ancora una sostituzione valida) .


Figura 3.10: a) Il taglio illustrato in non puo essere effettuato con l’approccio con

duplicazione di facce. b) Se si seleziona il vertice piu vicino , ovvero il

vertice 0, avremmo un vertice e una faccia corrispondente da duplicare

che sono esterni e quindi non condivisi da nessuna altra entita dalla quale

sia possibile separare. c) Lo spostamento di un altro vertice, ad esempio

il vertice 1, puo causare modifiche al modello percepibili dall’utente.

2. Il numero e la qualita dei tetraedri possono degenerare velocemente (ad

esempio effettuando ripetutamente piu tagli nella stessa regione).

A nostro avviso, l’origine dei problemi sta nel fatto che in questo modello

il tool non interagisca con l’oggetto deformabile in maniera fisica, cioe tra-

mite il contatto con la superficie esterna, bensı il suo movimento determina

direttamente le modifiche della struttura geometrica (mesh) che rappresenta

l’oggetto. Per fare cio viene consentito che il tool si trovi all’interno dell’og-

getto, cioe che essi occupino lo stesso spazio fisico.

L’approccio proposto in questa tesi si basa invece su una modellazione del

taglio in funzione dello stress. Come avviene nella realta, un tessuto non si

lacera a causa del contatto con il bisturi, bensı a causa del fatto che la pres-

sione, o la tensione, esercitata da un oggetto esterno supera un certo limite.

Nel metodo proposto il volume della mesh e del bisturi non dovranno mai


intersecarsi, la deformazione ed il taglio avvengono proprio affinche questa

proprieta venga preservata.

Modellare il taglio attraverso questi accorgimenti fa si che il metodo propo-

sto (oltre che essere piu realistico del metodo con separazione di tetraedri)

sia adattabile a qualsiasi tool per il quale si conosca lo stress esercitato sulla

mesh durante il contatto, o piu in generale a una qualsiasi causa di stress,

come ad esempio quello indotto dalla forza di gravita. I due operatori de-

scritti nella sezione che segue saranno utilizzati per modellare il taglio in

corrispondenza dei punti di stress piu alto.

3.2.1 Gli operatori di taglio

Ci sono due funzioni per implementare il taglio che sono stati chiamati vertex-

split e edge-split:

Vertex-split: Separazione di vertici Come mostrato in Figura 3.11.a que-

sto operatore prende in input una sezione di taglio formata da un in-

sieme connesso di facce Fs0, Fs1, . . . , Fsn (evidenziate nell’esempio) e

produce come risultato la mesh stessa dove il vertice t e stato duplicato

e le facce di sezione sono diventate esterne.

Edge-split: Divisione di un edge Questo operatore prende in input un

punto su un edge della mesh Et. Come mostrato in Figura 3.11.b

l’applicazione di questo operatore sostituisce l’insieme dei tetraedri che

condividono l’edge Et con due nuovi sottoinsiemi di tetraedri ciascuno

dei quali condivide o l’edge E ′t oppure l’edge E ′′

t . Le facce interne

condivise dai due sottoinsiemi diventano come nel caso precedente facce

esterne ed il vertice condiviso dai due nuovi edge viene duplicato.

La mesh tetraedrale modella il volume di un oggetto, questo complesso

simpliciale (ovvero la mesh stessa) deve essere (e rimanere) manifold , cioe


Figura 3.11: a) Vertex-split duplica il vertice t utilizzando la superficie di taglio evi-

denziata. b) Edge-split suddivide in due l’edge Et, di conseguenza ogni

tetraedro che condivide l’edge genera due tetraedri.

omeomorfa ad una sfera in ogni punto. I due operatori di taglio appena

presentati devono preservare queste proprieta pur modificando la mesh.

3.2.2 Problemi connessi all’uso degli operatori di taglio

Riprendiamo gli operatori di taglio introdotti precedentemente, come possia-

mo notare in Figura 3.12, l’applicazione del vertex-split alla sezione di taglio

evidenziata produrra una mesh non piu manifold.

La mesh in Figura 3.12 rappresenta un caso non realistico nella modellazione

di un oggetto deformabile: i due tetraedri della mesh non possono fisicamen-

te essere attaccati l’uno all’altro attraverso due edge. Oltre al fatto che la

mesh di Figura 3.12 rappresenti un modello non realistico va precisato che

la funzione di adiacenza T-T non puo collegare i due tetraedri tra loro; di

conseguenze nessun half-edge definito sulla mesh permetterebbe la visita di


Figura 3.12: L’applicazione del vertex-split alla mesh tetraedrale raffigurata produce

una situazione in cui la mesh non e piu omeomorfa ad una sfera lungo

gli edge (V0, V2) e (V1, V2)

tutti i tetraedri, andando cosı a violare le nostre ipotesi.

Cerchiamo allora di definire sotto quali condizioni possa essere effettuato il

taglio effettuato dal vertex-split e come questo possa essere applicato in mo-

do che dopo la modifica la mesh sia ancora manifold.

Sia S una superficie di taglio composta da facce interne della mesh te-

traedrale T e sia σ(S) = v0, . . . , vh l’insieme dei vertici di bordo di S.

Definiamo una nuova mesh e(S) ottenuta da S coma segue:

Per ogni sequenza di vertici appartenenti a σ(S) che sono anche sul bordo

della mesh tetraedrale si aggiunga un vertice fittizio connesso a tutti i vertici

della sequenza e le relative facce (Figura 3.13) .

Diciamo che separiamo T attraverso S, in breve V Split(T, S), per indicare

l’operazione duplicazione di tutti i vertici interni (Figura 3.15) .

Lemma 1 Sia T una mesh 3-manifold e S una superficie di taglio 2-manifold

in T . Se e(S) e 2-manifold e non esiste un edge interno di e(S) tra due

vertici in σ(S), allora V Split(T, S) e 3-manifold.


Figura 3.13: Esempio di estensione di una sezione di taglio:I vertici fittizi vengono

connessi con le sequenze di vertici esterni nella mesh tetraedrale.

Figura 3.14: Una sezione di mesh tetraedrale che evidenzia la superficie di taglio S

(evidenziata): La superficie S viola le condizione esposte nel Lemma 1

in quanto non e 2 -manifold (a causa dell’edge (e, d0)) e gli edge interni

(a, b) e (c, d) connettono due vertici di bordo.


Figura 3.15: Esempio di V Split di una sezione di taglio: I vertici interni i relativi

edge e facce vengono duplicati.

Dim: Se la mesh tetraedrale e 3-manifold, gli unici punti in cui l’operazione

V Split(T, S) puo creare dei punti non omeomorfi ad una sfera sono gli edge

ed in vertici della superficie S. Data la definizione di V Split i vertici interni

ed i relativi edge verranno duplicati.

Un edge e non manifold quando l’operazione di V Split duplica le facce che

lo condividono senza duplicare l’edge stesso. Cio avviene solo nel caso in cui

un edge interno ad S colleghi due vertici di bordo, caso escluso per ipotesi.

I vertici appartenenti a σ(S) saranno tutti omeomorfi ad una semi-sfera dopo

l’operazione di V Split, in quanto bordo della nuova mesh tetraedrale.

Se, al contrario, c’e almeno un vertice non omeomorfo ad una semi-sfera

significa che lo era precedentemente (e quindi esterno nella mesh tetraedrale)

e unico elemento di una sequenza di vertici di bordo ma in questo caso la

mesh e(S) non e manifold come richiesto nelle ipotesi (Figura 3.14) .

3.2.3 Vertex-split

I due operatori di taglio qui descritti modificano la mesh utilizzando so-

lamente le informazioni di adiacenza fornite dalle funzioni TT e VT pre-

cedentemente descritte. L’azione del vertex-split puo essere suddivisa in 4

fasi:


Verifica delle Condizioni Questa fase controllera la sezione di taglio della

mesh S in accordo con le condizioni individuate nel Lemma 1. Una

sezione di taglio viene chiamata valida se rispetta tali condizioni.

Individuazione dei vertici da duplicare Se la sezione di taglio e valida

allora dobbiamo individuare i vertici da duplicare.

Individuazione dei tetraedri coinvolti nel taglio Per ogni vertice che

intendiamo duplicare dobbiamo individuare i tetraedri che condividono

tale vertice e, per ognuno di essi stabilire a quale dei due nuovi vertici

dovra riferirsi.

Sdoppiamento dei vertici e ripristino delle topologie Alla fine prov-

vederemo a duplicare i vertice stesso e ad aggiornare di conseguenza

adiacenza in modo tale da modellare il taglio.

Condizioni sulla sezione di taglio Il parametro su cui verranno poste

condizioni sull’ applicabilita del taglio sara la sezione di taglio vista come

mesh triangolare. Per l’efficienza dell’implementazione la sezione di taglio

estesa S ′ = φ(S) non viene esplicitamente costruita, bensı le condizioni per

cui S ′ e 2-manifold vengono espresse in termini di S con dei test che com-

prendono tutti i casi.

In seguito utilizzeremo la seguente notazione: se un vertice V o un edge E

sono di bordo nella mesh triangolare allora diremo che tali entita apparten-

gono all’insieme 2dB, rispettivamente 3dB se tali entita sono di bordo nella

mesh tetraedrale.

Data la mesh triangolare S indicante la sezione di taglio, allora l’opera-

zione di vertex-split puo essere applicata su S se:

1. La mesh S e una superficie composta da una sola componente connessa

omeomorfa ad un disco.


Figura 3.16: Condizioni sulla sezione di taglio del vertex-split (i vertici e gli edge

esterni sono evidenziati) . a) La sezione di taglio non e manifold. b)

Esiste un edge interno che connette due vertici esterni. c) Esiste un

vertice di bordo sulla mesh tetraedrale che sta tra due edge non di

bordo. d) Esiste una sequenza di edge di bordo composta da un solo

edge.


Figura 3.17: a) Tra i vertici interni alla sezione di taglio esiste un vertice v di bordo

nella mesh tetraedrale. b) Risultato dello split applicato alla mesh estesa:

la mesh M1 non e piu omeomorfa da una sfera in ogni punto.

2. Non esistono edge interni alla mesh S che uniscano 2 vertici apparte-

nenti all’insieme 2dB.

3. Non esiste tra due edge 2dB ma non 3dB (lungo cioe il bordo di S) un

vertice 3dB .

4. Non esiste una sequenza di edge 2dB e 3dB composta da un solo edge

3dB.

5. Tra i vertici interni di S non esistono vertici che sono 3dB.

Queste condizioni derivano dalle ipotesi del Lemma 1. In particolare i casi

esposti nelle condizioni 1 e 2 sono illustrati in Figura 3.16.a e 3.16.b, mentre le

condizioni 3 e 4 corrispondono alle mesh triangolari illustrate rispettivamente

in Figura 3.16.c e 3.16.d.

La condizione 5 risulta necessaria affinche il vertex-split non crei la situazione

mostrata in Figura 3.17.b.


Figura 3.18: a) Solamente i vertici b e c sono da duplicare in quanto si trovano

entrambi tra due edge esterni. b) La mesh prodotta dall’esecuzione del

vertex-split.

Individuazione dei vertici da duplicare Data una sezione di taglio S

che rispetta le condizioni 1. . . 5 definiamo i vertici V ∈ S che vanno duplicati.

Sia S’ l’estensione della sezione di taglio S (S ′ = ϕ(S)), vanno duplicati

(secondo la definizione dello split) tutti i vertici interni ad S’. Questa ultima

affermazione si traduce sulla mesh non estesa S come:

• Vanno duplicati tutti i vertici interni a S.

• Vanno duplicati tutti i vertici 2dB di S che sono anche 3dB e si trovano

tra due edge di S che sono sia 2dB che 3dB.

Nell’esempio rappresentato in Figura 3.18 i vertici di bordo a, b, c e d sono

di bordo sia nella mesh S, sia nella mesh tetraedrale; solamente i vertici b e

c sono da duplicare in quanto sono rispettivamente inclusi tra gli edge E0 ed

E1 (nel caso di b) e tra E0 ed E1 (nel caso di c).

Da notare l’importanza dell’half-edge nella navigazione sia della mesh di


triangoli sia della mesh di tetraedri: un esempio di utilizzo riguarda l’in-

dividuazione dell’edge esterno. In questa fase, attraverso l’operazione di

Next edge, possiamo verificare se esiste almeno una faccia esterna che incide

sull’edge.

Individuazione dei tetraedri coinvolti nel taglio Data la sezione di

taglio S sulla mesh tetraedrale, esiste sicuramente un insieme di tetraedri

T 6= ∅ che condividono alcune entita con essa. Ci possono essere tetraedri

che hanno una faccia sulla sezione di taglio, altri che condividono un edge,

altri ancora solamente un vertice. L’insieme T puo essere suddiviso in due

sottoinsiemi disgiunti Ta e Tb, rispetto alla sezione di taglio stessa. I due

sottoinsiemi devono essere tali che le uniche entita da essi condivise sono

quelle presenti sulla sezione di taglio, piu formalmente:

Sia S una sezione di taglio valida e sia T l’insieme dei tetraedri che

condividono almeno un’entita con S allora i due sottoinsiemi Ta e Tb devono

godere delle seguenti proprieta:

1. Ta ∩ Tb = ∅

2. Ta ∪ Tb = T

3. ∀Ti ∈ Ta, Tj ∈ Tb : Ta ∩ Tb 6= ∅ =⇒ Ta ∩ Tb ∈ S.

Un esempio di scomposizione in questi due sottoinsiemi e riportato in Fi-

gura 3.19. I due sottoinsiemi di tetraedri Ta e Tb vengano costruiti in modo

tale che ognuno di essi possa successivamente sostituire le entita che condi-

vide con S con una delle due che si creano dopo la duplicazione.

La proprieta 3 risulta fondamentale nella suddivisione dei tetraedri, non ri-

spettarla significherebbe creare un taglio non corretto rispetto alla sezione

definita (ci potremmo ritrovare ad esempio un tetraedro che attraversa il ta-

glio). La creazione dei due sottoinsiemi Ta e Tb viene effettuata attraverso

tre passi successivi:


Figura 3.19: a) Tutti i tetraedri che hanno una entita condivisa sulla sezione di

taglio (evidenziata) sono coinvolti nella divisione. b) Le due mesh

rappresentano i sottoinsiemi Ta e Tb.

Figura 3.20: a) La sezione di taglio S nella mesh tetraedrale con tutti i tetraedri

che condividono almeno una entita con la sezione. b) Suddivisione

dei tetraedri che condividono una faccia con la sezione. c) I tetrae-

dri che condividono vertici o edge con S vengono suddivisi sfruttando la

condivisione dei vertici.


1. Vengono inseriti in un insieme Ta0 tutti i tetraedri che hanno una fac-

cia sulla superficie di taglio; attraverso l’adiacenza tetraedro-tetraedro

l’altro tetraedro che condivide la faccia sara a lui opposto ed apparterra

all’insieme Tb0 (Figura 3.20.b) .

2. Per ogni tetraedro appartenente a T controlliamo se questo condivi-

de almeno un vertice V 6∈ S con un tetraedro Tx che appartiene a Ta

oppure Tb (una opzione esclude l’atra visto che un tetraedro non puo

avere contemporaneamente due vertici che stanno su due lati differen-

ti rispetto ad S) . Se tale vertice esiste, allora inseriamo il tetraedro

nell’insieme corrispondente Tx (Figura 3.20.b) .

3. Ripetiamo il passo 2 fino che tutti i tetraedri Ti ∈ T apparterranno a

Ta oppure a Tb.

Da notare che la divisione dei due sottoinsiemi puo essere descritta come:

Ta = t|∃v 6 ∈S, v ∈ t, v ∈ t′, t′ ∈ Ta ∨ t′ ∈ Ta0Tb = t|∃v 6 ∈S, v ∈ t, v ∈ t′, t′ ∈ Tb ∨ t′ ∈ Tb0

Sdoppiamento del vertice e ripristino delle topologie Quando arri-

viamo a questa fase sappiamo gia quale sia l’insieme dei vertici da duplicare

ed i due sottoinsiemi Ta e Tb divisi dalla sezione di taglio S.

A questo punto come prima cosa dobbiamo duplicare i vertici stabiliti; per

ogni tetraedro di Ta sostituiamo il vecchio vertice che stava sulla sezione as-

sociando una delle due copie create, di conseguenza nei tetraedri di Tb questa

sostituzione avviene utilizzando l’altra copia del vertice. Appena effettuate

questa sostituzioni dobbiamo modificare le funzioni di adiacenza:

• Le facce che si trovano sulla sezione di taglio vengono poste a fac-

ce esterne; modificando quindi la funzione di adiacenza tetraedro -

tetraedro.


• La funzione di adiacenza vertice-tetraedro viene inizializzata per ognu-

no dei vertici duplicati con il tetraedro cui e stato assegnato.

Queste due operazioni sulla topologia garantiscono la correttezza degli

iteratori che verranno riutilizzati all’interno della mesh.

3.2.4 Edge-split

L’operazione di vertex-split consente si modellare solo tagli che passano at-

traverso i vertici e gli edge della mesh originale. Cio puo non essere sufficiente

qualora il punto dove deve comparire il taglio e lontano dai vertici selezio-

nati, poiche risulterebbe visivamente inaccurato. Inoltre, come illustrato in

Figura 3.10, esistono casi in cui il taglio non sarebbe possibile.

L’edge-split viene introdotto per ovviare a queste limitazioni, modificando

la mesh in modo che l’applicazione del vertex-split risulti accurata. Tale

modifica consiste nell’inserimento di un nuovo vertice all’interno dell’edge

da dividere. Questa operazione portera ad un riassetto locale della mesh

suddividendo tutte le entita che condividevano l’edge suddiviso (tetraedri e

facce) oltre che ad una modifica delle funzioni di adiacenza. Piu in dettaglio

l’edge-split puo essere suddiviso in 2 fasi:

Edge Split La divisione in due di un edge coincide con la divisione dei te-

traedri ai quali appartiene. Cosı come l’edge viene sostituito da una

coppia di nuovi edge, anche i tetraedri che lo contengono saranno so-

stituiti da una coppia di tetraedri ciascuno contenente uno dei due

edge creati.Tutti i nuovi tetraedri dovranno essere connessi con gli altri

tramite le funzioni di adiacenza.

Modellazione del taglio Definiamo una sezione di taglio che comprenda

tutte le facce interne che stanno a meta del vecchio edge da dividere


Figura 3.21: a) Mesh allo stato iniziale. b) I tetraedri originali vengono sostituiti da

un nuovo insieme di tetraedri in modo tale che l’edge venga diviso in due

parti. c) Le funzioni di adiacenza vengono ripristinate. d) Utilizzando il

vertex-split si effettua la separazione dei tetraedri sulla sezione di taglio

appena definita.

(tutte queste facce apparterranno ai nuovi tetraedri) ed utilizziamo

l’operazione di vertex-split per modellare il taglio.

La prima fase non modifica l’apparenza del modello. Effettuano queste

operazioni sostituiamo un insieme di tetraedri con un altro che occupa lo

stesso spazio, passando situazione come mostrata in Figura 3.21.a quella

mostrata in Figura 3.21.b e poi in 3.21.c. la prima fase prepara la mesh

affinche si possa utilizzare , nella fase successiva, l’operazione di vertex-split

su una sezione di facce ben definita. Dopo quest’ultima fase si arriva alla

situazione mostrata in Figura 3.21.d: Il taglio ottenuto divide quindi l’edge


Figura 3.22: a) Divisione dell’edge 1: I due nuovi tetraedri saranno formati dalle

quadruple di vertici (0, 1, Vn, 3) e (V n,1, 2, 3). b) Divisione dell’edge 2:

I tetraedri creati sarebbero formati dai vertici (0, 2, Vn, 1) e (Vn, 2, 3, 1).

nelle due parti richieste.

Divisione di un edge Dato un edge Ei i tetraedri che dovranno essere

sostituiti sono tutti quelli che condividono l’edge Ei stesso. Lo strumento

adatto a scorrere in modo ordinato, rispetto ad un senso di rotazione, tutti

i tetraedri che condividono l’edge e l’half-edge (attraverso next edge) .

Supponiamo inizialmente di dover dividere l’edge 1; per suddividere tale

edge in due come prima cosa dobbiamo inserire un vertice nel punto in cui

l’edge va suddiviso. Come mostrato in Figura 3.22.a il vertice Vn e inserito

lungo edge 1 in modo tale da suddividerlo in due; il primo segmento va dal

vertice 0 a Vn, mentre il secondo da Vn al vertice 2. Per ogni tetraedro che

condivide l’edge dovremo creare due tetraedri come illustrato nell’esempio in

Figura 3.22.a; come possiamo notare i due tetraedri creati condivideranno i

vertici 1,3, e Vn, i primi due sono condivisi in quanto non fanno parte dell’edge

da separare. Il tetraedro Ta sara formato quindi dai vertici (0, 1, Vn, 3) mentre

il tetraedro Tb dai vertici (V n,1, 2, 3).

Se invece volessimo separare l’edge 2 i tetraedri creati sarebbero formati dai


Figura 3.23: Affinche la divisione di un edge sia effettuata in modo topologicamente

corretto e necessario suddividere tutti i tetraedri che lo condividono.

vertici (0, 2, Vn, 1) e (Vn, 2, 3, 1) Figura 3.22.b. Il caso ottenuto tagliando

l’edge 2 puo essere ricondotto al caso dell’edge 1 tramite una rotazione del

tetraedro, rotazione che deve essere effettuata in modo tale da far coincidere

l’edge 2 nella posizione dell’edge 1.

In pratica ogni caso di taglio di edge possibile puo essere ottenuto attraverso

rotazioni e mirroring del caso illustrato in Figura 3.22.a. Otterremo quindi

12 possibili casi di taglio (uno per ogni edge in entrambe le direzioni) , per

ognuno dei quali dovranno essere specificate le quadruple ordinate di vertici

che definiscono i tetraedri Ta e Tb.

Dopo che tutti i tetraedri che condividevano l’edge da dividere sono stati

sostituiti dai tetraedri Ta, i e Tb, j dobbiamo ripristinare le funzioni di adia-

cenza in modo corretto. Siano Told l’insieme dei tetraedri sostituiti e TA e

TB rispettivamente l’insieme dei tetraedri Ta e Tb definiti come illustrato nel

paragrafo precedente.


Per quanto riguarda la funzione di adiacenza vertice-tetraedro le operazioni

da effettuare sono 3:

1. Inserire tutti i tetraedri Ti ∈ TA ∪ TB nell’insieme di tetraedri asso-

ciato al nuovo vertice Vn (e ovvio che ogni nuovo tetraedro condividera

questo vertice).

2. Togliere ogni tetraedro Ti ∈ Told dagli insiemi associati ai suoi vertici.

3. Inserire ogni nuovo tetraedro Ti ∈ TA ∪ TB negli insiemi associati ai

suoi vertici.

Per quanto riguarda la funzione di adiacenza tetraedro-tetraedro, per ogni

nuovo tetraedro Ti ∈ TA ∪ TB:

1. Come illustrato nelle Figure 3.24.a e 3.24.b, le facce F1 ed F3 coinci-

deranno con le facce segnate come F ′1 ed F ′

3. La funzione di adiacenza

TT per queste facce connettera i nuovi tetraedri Ta e Tb con i tetraedri

cui era connesso T; se F1 ed F3 erano facce esterne allora anche F ′1 ed

F ′3 lo saranno.

2. Le faccia Fnew di Figura 3.24.c, sta tra Ta e Tb; la funzione TT per

questa faccia dovra quindi connettere i due tetraedri tra loro.

3. La funzione TT sulle facce F0 ed F2 (Figura 3.24.a) connetteva il tetrae-

dro T con gli altri che condividevano l’edge da dividere, oppure erano

facce esterne. Dopo la sostituzione ognuna di queste due facce viene

divisa in F0a, F0b e F2a, F2b. La funzione di adiacenza TT connettera

ognuna di queste due nuove facce con quelle che vengono create dalla

suddivisione del tetraedro successivo; oppure sara una faccia esterna.


Figura 3.24: a) Le adiacenze originali del tetraedro prima della modifica. b) Le adia-

cenze delle facce F1 ed F3 sono equivalenti a quelle delle facce F ′1 ed

F ′3. c) I due tetraedri creati vanno connessi tra loro tramite la faccia

condivisa. d) Tra le coppie di nuovi tetraedri creati vanno connessi in

senso ordinato rispetta all’edge diviso.


Modellazione del taglio La mesh si presenta in questa fase formalmente

corretta e dotata di tutte le funzioni di adiacenza. L’edge e gia stato sud-

diviso ma non e stato effettuato il taglio vero e proprio. In virtu di queste

ipotesi possiamo adesso utilizzare l’operazione di vertex-split definendo co-

me superficie di taglio le facce condivise tra ogni coppia di tetraedri di ogni

sostituzione (le facce Fnew di Figura 3.24.c) . Il vertex-split provvedera a

realizzare il resto del taglio restituendo una mesh omeomorfa ad una sfera in

ogni punto.

3.2.5 Riassegnamento dei parametri fisici dopo il taglio

La specifica degli operatori di taglio appena descritta astrae dal modello fi-

sico. In questa fase dobbiamo essere in grado di determinare i valori dei

parametri fisici da associare ai nuovi tetraedri utilizzando i parametri geo-

metrici ottenuti della divisione dei tetraedri.

Dobbiamo preservare il contributo di massa del dominio tetraedrale qualun-

que sia la sua decomposizione, assegnamo cioe il contributo di massa dei

nuovi tetraedri in maniera tale che la loro somma sia equivalente alla somma

dei contributi di massa dei tetraedri originali.

Supponiamo che M0 sia il valore di massa di un vertice della mesh prima

della sua divisione. I valori delle masse dei due vertici creati dalla sua divi-

sione m0 ed m1 devono essere tali che la loro somma sia equivalente ad M .

La ripartizione dei contributi di massa sui due nuovi vertici viene calcolata

in base alla somma dei volumi dei tetraedri che incidono sui due differenti

vertici (come se la massa fosse quindi concentrata all’interno del volume oc-

cupato dai tetraedri). Siano quindi T0 e T1 i due insiemi di tetraedri che

incidono sui vertici v0 e v1 creati dopo l’applicazione del vertex-split, allora

definiamo V ol0 e V ol1 come:

V ol0 =∑

ti∈T0V ol(ti)

3.3 Semplificazione 71

V ol1 =∑

ti∈T0V ol(ti)

Dove V ol(ti) e una funzione che restituisce il volume di un tetraedro. Dato

che T0 e T1 sono due insiemi disgiunti il volume V occupato dall’insieme di

tutti i tetraedri che incidevano sul vertice originale equivale alla somma dei

volumi ripartiti nei due sottoinsiemi creati dopo il taglio.I contributi di massa

dei due vertici saranno allora calcolati come:

m0 = (V ol0/V )M

m1 = (V ol1/V )M

3.3 Semplificazione

Il metodo di taglio proposto in questa tesi puo creare alcuni problemi relativi

alla creazione di nuovi tetraedri:

Problema della frammentazione Nell’utilizzo dell’edge-split, per ogni di-

visione di edge viene raddoppiato il numero dei tetraedri che condivi-

devano quell’edge. Anche se questi operatori di taglio sono stati creati

in modo tale da ridurre al minimo l’incremento di tetraedri in caso di

taglio (rispetto ad esempio all’approccio con suddivisione di tetraedri

gia illustrato nel secondo capitolo) , tale incremento, anche se minimo,

puo degradare col tempo le prestazioni della simulazione.

Tetraedri degeneri L’edge-split puo creare tetraedri degeneri. Un tetrae-

dro degenere e un tetraedro che ha una aspetto non regolare, nell’esem-

pio di Figura 3.25.a il tetraedro T2 puo essere considerato degenere (il

criterio per stabilire l’aspetto di un tetraedro sara esposto in seguito) .

Questi tetraedri possono generare comportamenti anomali nel calcolo

della deformazione del modello.

Il degrado delle prestazioni della simulazione e legato soprattutto al nu-

mero di tetraedri presenti sulla mesh: piu tetraedri ci sono e maggiore sara


a) b)

Figura 3.25: a) Le operazioni di taglio possono creare tetraedri dall’aspetto non re-

golare (degeneri), in questo caso il tetraedro T2 non e un tetraedro di

buona qualita. b) La semplificazione diminuisce la quantita di tetraedri

eliminando i tetraedri di bassa qualita , migliorando di conseguenza la

qualita della mesh stessa. Le operazioni di taglio tendono invece ad

incrementare la complessita e peggiorare la qualita della mesh.


la complessita nel calcolo delle forze interne al modello. Affinche il nume-

ro di tetraedri rimanga in un ordine costante e vengano eliminati i tetraedri

degeneri, le operazioni di taglio sulla mesh devono essere compensate da ope-

razioni di semplificazione in modo tale che il numero di tetraedri della mesh

tenda a stabilizzarsi (Figura 3.25.b) .

Un algoritmo di semplificazione costruisce una rappresentazione semplifica-

ta della mesh con un numero minore di entita cercando di preservare il piu

possibile il volume e la forma del modello iniziale. L’approccio adottato in

questa tesi e basato su una strategia di semplificazione di tipo incrementa-

le: La semplificazione procede attraverso una sequenza di modifiche locali

della mesh dove, ad ogni passo, si riduce la descrizione della mesh. La sem-

plificazione incrementale di una mesh e puo essere effettuata seguendo due

strategie:

Strategie di Raffinamento Partendo da una rappresentazione di massima

della regione da semplificare si aggiungono dettagli in modo da ottenere

un modello che si avvicini il piu possibile a quello reale.

Strategie di Decimazione Si eliminano dalla regione da semplificare en-

tita scelte secondo una strategia in modo da generare un modello ap-

prossimato che sia piu fedele possibile a quello reale.

L’algoritmo di semplificazione implementato in questa tesi si basa su una

strategia di decimazione.

3.3.1 Edge Collapse

Semplificare localmente una mesh significa togliere entita in modo da ottenere

una rappresentazione della mesh composta da un minor numero di tetraedri.

L’approccio adottato in questa tesi e basato sull’edge collapse. Ogni passo di

decimazione corrisponde col collasso di un edge; collassare un edge significa


unificare i vertici dell’edge modificando di conseguenza la mesh, i tetraedri

che contenevano quell’edge saranno cancellati.

La Figura 3.26.a mostra il collasso di un edge su una mesh di triangoli: i

vertici Ve0 e Ve1 sono sostituiti dal solo vertice V col risultato che i triango-

li Fe0 e Fe1, che condividevano l’edge, vengono eliminati, mentre i triangoli

F0 . . . F7 vengono modificati in modo da condividere il vertice V. Non sem-

pre e possibile pero collassare un edge; la Figura 3.26.b mostra un caso in

cui non e possibile collassare l’edge evidenziato; un collasso di questo tipo

genererebbe una situazione topologicamente scorretta.

3.3.2 Correttezza topologica

Dato un edge da collassare c’e un insieme di condizioni necessarie e sufficien-

ti che preservano la correttezza topologica della mesh tetraedrale. Queste

condizioni vengono introdotte in [DEGN99] come link conditions.

Definiamo come Σi la mesh tetraedrale in cui ogni vertice esterno e connesso

ad un vertice fittizio VDummy; secondo questa rappresentazione ogni vertice

esterno genera un edge con VDummy, ogni edge esterno una faccia ed ogni fac-

cia esterna un tetraedro. Sia Ei l’edge da collassare e Σi+1 la mesh ottenuta

dopo il collasso. Sia σ una entita della mesh Σi, definiamo:

St(σ) = τ ∈ Σi|σ ∈ τ

Se σ ad esempio fosse un vertice allora St(σ) sarebbe l’insieme di tutti gli

edge, le facce ed i tetraedri a cui σ appartiene.

Definiamo adesso Lk(σ) come l’insieme di tutte gli elementi appartenenti a

St(σ) ma non incidenti in σ. Siano V0 V1 i vertici che delimitano l’edge da

collassare Ei, allora le seguenti affermazioni sono equivalenti:

• Lk(V0) ∩ Lk(V1) = Lk(Ei)


a)

b)

Figura 3.26: a) Il collasso dell’edge evidenziato produce una nuova situazione in cui

il numero di triangoli e ridotto di due (ovvero coloro che condividevano

l’edge. b) Il collasso dell’edge evidenziato produrrebbe una situazione

in cui le due facce F0 ed F1 finirebbero schiacciate nel collasso.


Figura 3.27: a) Lk(a) ∩ Lk(b) = x, y = Lk(Eab). b) Lk(a) ∩ Lk(b) =

x, y, z, Ezx 6= Lk(Eab) = x, y, z.

• Se Σi e un poliedro, allora anche Σi+1 e un poliedro.

Nella Figura 3.27 viene riportato un esempio di collasso su mesh trian-

golare. In 3.27.a Lk(a) ∩ Lk(b) = x, y che e uguale a Lk(Eab), quindi sa-

rebbe possibile effettuare il collasso dell’edge (Eab), la situazione e diversa

per l’esempio mostrato in 3.27.b: Lk(a) ∩ Lk(b) = x, y, z, Ezx, mentre

Lk(Eab) = x, y, z; non e possibile quindi collassare l’edge (Eab) .

Un semplice esempio di applicazione delle link conditions su una mesh

tetraedrale viene rappresentato in Figura 3.28: il vertice fittizio (Dummy)

viene posto al centro del tetraedro per semplicita di rappresentazione. In que-

sta situazione si cerca di collassare un edge di una mesh tetraedrale formata

da un solo tetraedro (e ovvio che tale collasso non debba essere effettuato

affinche la mesh rimanga poliedrica) ; Lk(a) ∩ Lk(b) 6= Lk(Eab) in quanto

Lk(a) ∩ Lk(b)− Lk(Eab) = Ex,Dummy, Ey,Dummy, Fx,y,Dummy.

3.3.3 Correttezza geometrica

La correttezza geometrica del modello risulta piu difficile da garantire. In

ordine di importanza e facilita di implementazione i seguenti test possono


Figura 3.28: Le Link Conditions non vengono rispettate in quanto Lk(a) ∩ Lk(b)−Lk(Eab) = Ex,Dummy, Ey,Dummy, Fx,y,Dummy

essere effettuati sui tetraedri modificati [CCM+00]:

Flip In alcuni casi il collasso di un edge puo provocare l’inversione di un

tetraedro.

Shape Alcuni collassi possono creare tetraedri deformi (ci sono diversi me-

todi per valutare la qualita della forma di un tetraedro).

Auto Intersezioni In rari casi il collasso di un edge puo provocare una auto

intersezione della mesh.

Nel framework proposto questa tesi sono stati risolti i primi due problemi:

Per quanto riguarda il flipping del tetraedro e sufficiente effettuare il test

sull’ordinamento dei vertici definito in 3.1.

Il problema relativo ai tetraedri degeneri viene risolto attraverso l’ordine in

cui vengono effettuati i collassi durante la decimazione; i collassi cui vengono

eliminati tetraedri degeneri sono quelli con priorita piu alta.

Il problema delle auto intersezioni e un problema molto complesso in quan-

to e l’unico caso i cui effetti possono essere non locali rispetto al collasso (si


possono creare cioe auto-intersezioni tra qualsiasi tetraedro appartenente alla

mesh anche se questi non e coinvolto direttamente nel collasso) . Fortuna-

tamente questi tipi di situazioni sono piuttosto rare e possono percio essere

trascurate.

3.3.4 Riassegnamento dei parametri fisici dopo la sem-

plificazione

Come avviene per le operazioni di taglio, la massa, alla fine di una operazione

di semplificazione, deve essere redistribuita sui vertici in modo proporzionale

al volume dei tetraedri vi incidono

L’oggetto deformabile deve essere trattato come se la sua massa fosse unifor-

memente distribuita all’interno del volume dei tetraedri. Ciascun tetraedro

contribuisce, in maniera proporzionale al suo volume ad alla densita di massa

del materiale cui e composto, alla massa presente su ciascun vertice. Sup-

ponendo che densita di massa sia la stessa in ciascun tetraedro della mesh,

dopo ciascuna operazione di edge collapse dobbiamo effettuare le seguenti

operazioni :

• Eliminare il contributo di massa dei tetraedri eliminati o modificati

durante il collasso.

• Aggiungere il contributo di massa dei tetraedri modificati dopo il col-

lasso.

3.3.5 Implementazione della decimazione

L’edge collapse e uno strumento che viene utilizzato per migliorare localmente

la mesh nei punti in cui sono stati effettuati dei tagli. Il dominio di questi

collassi sono tutti i tetraedri creati durante il taglio, ogni collasso interessa

un edge di questi tetraedri. L’approccio utilizzato consiste nel trovare quali


siano i tetraedri di qualita peggiore e quindi quali siano gli edge da collassare

in modo da rimuovere quei tetraedri.

Come gia presentato precedentemente quando avviene un collasso l’edge viene

sostituito da un vertice; questo vertice puo essere posto sulla posizione di

uno dei due vertici, a meta dell’edge collassato, oppure possiamo cercare il

punto per il quale otteniamo come risultato la mesh migliore. Quest’ultima

soluzione risulta la piu efficace ma anche quella di complessita maggiore;

del resto gli edge che andremo a collassare (in quanto condivisi da tetraedri

degeneri) saranno di solito edge molto corti, percio la posizione su cui porre il

vertice durante il collasso non rappresenta un fattore importante nell’efficacia

della decimazione.

Aspect Ratio di un tetraedro Una misura della qualita di un tetraedro

puo essere ottenuta attraverso l’aspect ratio. Ci sono diversi metodi per

esprimere l’aspect ratio di un tetraedro, i piu utilizzati sono due:

• Misurare l’aspect ratio attraverso il rapporto tra la sfera inscritta e

quella circoscritta in un tetraedro.

• Misurare l’aspect ratio come il minimo angolo solido di un tetraedro.

Ogni collasso candidato avra quindi un indice di priorita calcolato in base

alla media aspect ratio dei tetraedri che lo condividono. L’algoritmo di de-

cimazione sceglie tra tutti i possibili collassi quello che elimina il gruppo di

tetraedri di peggiore qualita, reinserendo nell’insieme di tutti i candidati i

collassi relativi ai tetraedri appena modificati.

Il meccanismo delle marche temporali I collassi candidati devono es-

sere eseguiti nell’ordine giusto, prediligendo nel nostro caso quelli che coin-

volgono il gruppo di tetraedri che ha un aspect ratio medio peggiore.

Di solito si rappresenta l’insieme dei possibili collassi attraverso una coda con


priorita (di solito un heap) contenente tutti i possibili collassi con la relativa

stima di priorita. Questa coda con priorita deve essere aggiornata dopo ogni

collasso modificando la priorita di ogni collasso adiacente a quello appena

eseguito (si possono modificare notevolmente le proprieta geometriche dei

tetraedri coinvolti e quindi anche il loro aspect ratio) .

Un buon metodo per evitare l’aggiornamento esplicito degli elementi della

coda con priorita e quello di utilizzare delle marche temporali sui vertici per

stabilire se un collasso sia aggiornato oppure no. Piu precisamente:

• Inizialmente ogni collasso candidato ed ogni vertice hanno la stessa

marca temporale.

• Per ogni collasso tolto dalla coda con priorita stabiliamo se e up to date.

Un collasso e up to date se la marca temporale del collasso e maggiore

o uguale a quella dei suoi vertici.

• Se il collasso appena estratto e up to date allora possiamo effettuare il

collasso marcando i vertici dell’edge coinvolto con la marca temporale

corrente.

Capitolo 4

La modellazione del contatto

Nel capitolo precedente sono state definite le operazioni di modifica e sempli-

ficazione della mesh tetraedrale. Ciascuna di queste operazioni (vertex-split,

edge-split e collasso) modifica la mesh tetraedrale preservando le sue pro-

prieta. Queste operazioni sono pertanto operazioni atomiche.

I due operatori di modifica sono stati realizzati con lo scopo di modellare

dei tagli sull’oggetto deformabile. Per simulare il taglio effettuato dal bistu-

ri sull’oggetto dobbiamo innanzitutto modellare il contatto tra il bisturi e la

mesh. Attraverso l’analisi degli stress imposti dal bisturi sulla mesh possiamo

stabilire dove utilizzare gli operatori di modifica simulando di conseguenza

il taglio. L’utente finale utilizzara una periferica di tipo force-feedback per

muovere il bisturi nello spazio e provare, quando avviene un contatto con

l’oggetto deformabile, la sensazione del contatto.

4.1 Collision Detection

Il termine collision detection indica il problema di rilevare le collisioni tra gli

oggetti presenti nella scena e di determinare le informazioni sul contatto quali

la posizione dove questo occorre e la profondita di penetrazione tra gli oggetti

coinvolti. La difficolta di questo problema varia con il comportamento degli

oggetti: dal caso piu facile, quello di oggetti rigidi, a quello piu complesso,

4.1 Collision Detection 82

quando gli oggetti sono fortemente deformabili. Nel caso preso in esame in

questa tesi gli oggetti che compongono la scena sono due: il bisturi (o un

qualsiasi altro attrezzo utilizzato in chirurgia) e il modello deformabile.

4.1.1 Tecniche di Collision Detection

La collision detection tra modelli generici e stata studiata negli anni in molti

campi: CAD/CAM, robotica, simulazioni e computer animation. Le soluzio-

ni proposte variano a seconda della rappresentazione geometrica degli oggetti

compongono la scena e dei dati che l’algoritmo di collision detection dovra

fornire.

Per esempio, gli algoritmi che mantengono la minima distanza euclidea tra

due oggetti sono spesso utilizzati nelle applicazioni di progettazione. Nel

nostro caso siamo interessati ai metodi che rilevano la collisione tra modelli

poligonali in tempo reale, non abbiamo quindi bisogno della distanza Eucli-

dea tra due oggetti che non collidono.

Molto del lavoro svolto nella collisione tra modelli poligonali e orientato verso

gli algoritmi che operano su poliedri convessi. Questi algoritmi si basano su

specifiche strutture dati create per la ricerca dei punti di contatto tra polie-

dri mantenendo la loro coerenza con la geometria durante l’animazione. Le

possibili ottimizzazioni su questi tipi di algoritmi si basano su organizzazioni

gerarchiche di bounding volumes: le soluzioni variano da alberi di box alli-

neati con gli assi,alberi di sfere, BSPTree.

Nonostante la loro efficienza questo tipo di algoritmi non puo essere utilizzato

nel caso di un simulatore chirurgico in quanto gli organi sono generalmente

non convessi e sottoposti in tempo reale a forti deformazioni.

L’algoritmo utilizzato nello svolgimento di questa tesi si basa sull’utilizzo

dell’hardware grafico.


Figura 4.1: a) Esempio di collisione statica: nel time step di controllo t0 il bisturi

collide con la mesh.b) Esempio di collisione dinamica: la collisione avviene

tra due time step successivi (t0, t1) ma istantaneamente nessuno di essi

si trova in collisione con la mesh. ).

4.1.2 La Collision Detection nella simulazione chirur-

gica

Nel contesto della simulazione chirurgica, non possono essere applicate tec-

niche di collision detection basate sulla precomputazione di strutture dati

dipendenti dalla forma degli oggetti, in quanto sono presenti nella scena og-

getti deformabili, ovvero la cui forma non e costante.

Contatti tra il bisturi ed il modello deformabile possono avvenire ad ogni

istante (in una situazione reale il chirurgo operando provoca contatti molto


di frequente) . Per questa ragione il controllo sulle intersezioni deve essere

effettuato ad ogni passo discreto della computazione; la collision detection

deve quindi essere effettuata in maniera molto efficiente. Non e sufficiente

pero testare la collisione in ogni passo discreto della computazione: il nostro

algoritmo di collision detection deve essere in grado di rilevare anche le colli-

sioni che avvengono tra due passi discreti consecutivi: oggetti che si muovono

possono anche entrare in collisione tra due time steps consecutivi senza che

in ognuno dei due istanti gli oggetti siano siano realmente in contatto tra

loro (Figura4.1) .

Le collisioni che vengono rilevante in un istante discreto sono dette stati-

che mentre quelle che avvengono tra due istanti successivi sono dinamiche.

Determinare l’intersezione geometrica tra le primitive che compongono gli

oggetti presenti nella scena e una operazione di elevato costo computazionale

e puo costituire un vero e proprio collo di bottiglia se tale controllo viene

effettuato ad ogni time step della simulazione.

Utilizzando un metodo basato sull’hardware grafico possiamo ridurre il nu-

mero di poligoni coinvolti nel test geometrico. Questo ultimo test infine

rivelera l’intersezione tra un sottoinsieme delle primitive che compongono la

scena, fornendo i dati necessari all’implementazione del contatto tra il bisturi

e l’oggetto deformabile.

4.1.3 La Collision Detection con hardware grafico

Rilevare una collisione tra due oggetti consiste nel verificare che il volume

del primo (ad esempio il tool che ha una forma regolare) intersechi quello del

secondo. Questo processo ha delle similitudini col processo di visualizzazione

di una scena.

In questo caso l’utente specifica un volume di vista, il view frustum, caratte-

rizzato dalla locazione, orientamento e proiezione della camera; l’algoritmo


Figura 4.2: a) La camera ortografica in OpenGL. b) La camera Prospettica. In en-

trambi i casi il volume di vista viene limitato dai parametri Top, bottom,

left, right, near, far ).


di frustum culling scarta i poligoni che non intersecano il volume di vista.

Il metodo che impieghiamo e stato proposto in [LCN99] e consiste nell’ap-

prossimare la forma del bisturi (che avendo una forma abbastanza semplice

puo del resto essere schematizzato come un parallelepipedo) con un view

frustum e sfruttare il vantaggio che l’algoritmo di frustum culling viene ese-

guito dall’hardware.In questo modo selezioniamo il sottoinsieme di primitive

componenti il modello che sulle quali effettuare il test geometrico.

Implementazione con OpenGL Il processo di visualizzazione puo esse-

re diviso in fase geometrica e fase di rasterizzazione. La parte geometrica

converte tutte i poligoni nel sistema di coordinate di camera, seleziona tutte

le facce che intersecano il volume di vista e calcola la proiezione dei vertici

sullo schermo. Successivamente, la fase di rastering determina quali pixel

devono essere accesi (utilizzando lo Z-Buffer ed il color buffer) . OpenGL

offre la possibilita di cambiare modalita di rendering in select mode. In que-

sto tipo di modalita di rendering, diversamente da quella classica non viene

effettuata alcuna fase di raster e vengono eliminate inoltre altre operazioni di

complessita notevole (ad esempio il calcolo dell’illuminazione) . Il select mode

(che si ottiene tramite la funzione glRenderMode(GL_SELECT)) restituisce le

informazioni relative ai gruppi di facce visibili. Ogni gruppo viene definito

utilizzando la funzione glPushName(), chiamata prima che ciascun gruppo

di facce sia disegnato. OpenGL riempie un array durante la fase geometrica

scrivendo tutti i gruppi che appaiono nel volume di vista. Associando un

nome per ogni tetraedro si possono cosı conoscere i tetraedri che appaiono

sul volume di vista.

La Collision Detection statica Chiamiamo adesso Head e Tail rispet-

tivamente la punta e la base del tool. La collisione statica tra il bisturi e

a mesh poligonale rappresentante l’organo puo essere effettuata associando


Figura 4.3: a) I punti head e tail del bisturi b) La camera posizionata col volume di

vista corrispondente alla forma del bisturi

una camera ortografica al bisturi (Figura 4.3) . La camera e posizionata sul

punto Head ed e direzionata lungo l’asse (Head,Tail) (utilizzando al funzione

gluLookAt()) . I parametri Near e Far (che definiscono gli estremi del volu-

me di vista) sono rispettivamente 0 e ‖Tail −Head‖. La sezione del bisturi

e compresa utilizzando i parametri Left, Right, Top e Bottom della camera

in accordo con la forma del bisturi.

Il movimento del bisturi ed i dati restituiti La soluzione presentata

nel paragrafo precedente verifica la collisione tra il bisturi in una determinata

posizione e il modello deformabile in un preciso time step. Se l’utente muove

il bisturi in maniera veloce parte della superficie di taglio definita dal movi-

mento potrebbe non venire considerata nella verifica della collisione (Figura

4.4) .

Per risolvere questo di problema approssimiamo il volume descritto dal bi-


Figura 4.4: a) Il tetraedro evidenziato viene individuato nel caso di collisione statica

b) Caso in cui un tetraedro che entra in collisione tra due time step

consecutivi non viene rilevato.

sturi durante il suo movimento. Piu in dettaglio creiamo un volume di vista

utilizzando la proiezione ortografica una posizione del bisturi sul quella del

time step successivo. Se in due time step successivi ti e t(i+1) abbiamo le

due coppie (Headti , Tailti) e (Headt(i+1), Tailt(i+1)

) il volume di vista sara co-

struito in modo tale da comprendere i volumi definiti nei due singoli istanti

(Figura 4.5.a) . Con questa estensione vengono riscontrate anche le collisioni

avvenute tra due time step successivi.

Ci possono essere pero situazioni in cui questo tipo di algoritmo individua

collisioni anche se in realta esse non accadono (Figura 4.5.b) . Il successi-

vo passo dell’algoritmo individuera tra i tetraedri restituiti da questo primo

passo di collision detection quelli che realmente sono stati toccati dal bisturi.


Figura 4.5: a ) La regione di spazio individuata e la proiezione di un volume sul

successivo. b ) Caso in cui l’estensione rivela collisioni che in realta non

esistono (la superficie di taglio evidenziata non interseca nessun tetraedro

4.1.4 L’intersezione geometrica

I tetraedri che sono stati individuati nel passo precedente non son necessa-

riamente tetraedri che sono stati toccati dal bisturi. Per individuare questo

sottoinsieme e necessario basarsi sulle proprieta geometriche della superficie

di taglio e del modello deformabile. Come gia accennato nel capitolo 2 il mo-

vimento del bisturi, e quindi la superficie di taglio, puo essere schematizzata

come una serie di triangoli connessi.

Chiamiamo Tetra0 l’insieme dei tetraedri individuati dalla collision detection

al primo passo di selezione e siano (Headti , Tailti) e (Headt(i+1), Tailt(i+1)

)

i punti head e tail del bisturi agli istanti ti e t(i+1). Tra due time step

successivi si puo dividere la superficie di taglio nei due triangoli Tr0 =

(Headti , Tailti , Headti) e Tr1 = (Headt(i+1), Tailt(i+1)

, Tailti). L’insieme dei

tetraedri Tetra1 che sono intersecati dalla superficie di taglio sara:


Figura 4.6: In questo schema sono sintetizzati i due passi della collision detection: Un

primo passo utilizza l’hardware grafico e produce un sovrainsieme di tutti

i tetraedri coinvolti nella collisione; mentre nel secondo passo vengono

determinati, tramite il test geometrico, quali tra essi sono effettivamente

entrati in contatto con la superficie di taglio

Tetra1 = Ti ∈ Tetra0 : (Ti ∩ Tr0) ∪ (Ti ∩ Tr1) 6= ∅ (4.1)

Ovvero l’insieme dei tetraedri selezionati al primo passo di collision de-

tection che si intersecano con uno dei due triangoli della superficie di taglio.

In caso tale collisione avvenisse l’algoritmo deve produrre i dati necessari

alla modellazione della deformazione.

Per modellare una deformazione (o successivamente un taglio) l’algoritmo di

collision detection, per ogni tetraedro coinvolto nel contatto, dovra restituire:

• Su quale faccia esterna del tetraedro sia avvenuto il contatto.

• Quale e il segmento di intersezione tra questa faccia e il triangolo della

superficie di taglio

• Quantificare la penetrazione del bisturi nella mesh.

Una buona quantificazione della penetrazione potrebbe essere la minima

distanza che occorre tra il segmento di taglio (head, tail) ed il segmento di

intersezione sulla faccia.

4.2 Il contatto 91

Figura 4.7: Il segmento (S1, S2) indica il punto in cui la superficie di taglio e entrata

in contatto con il tetraedro T, mentre il vettore d indica la penetrazione

all’interno del tetraedro stesso

Supponiamo di avere la situazione illustrata in Figura 4.7: la faccia da resti-

tuire e la faccia Fa del tetraedro T , il segmento di intersezione e costituito dal-

la coppia S1 e S2, e il vettore d rappresenta la penetrazione. Saranno questi

tutti dati indispensabili al calcolo della deformazione e della corrispondente

forza di feedback da restituire all’utente.

4.2 Il contatto

Simulare la deformazione di un tessuto significa calcolare la sua forma in ogni

istante di tempo. Per una realistica modellazione dell’intersezione con un og-

getto, la sola visualizzazione della deformazione del modello non e sufficiente;

e fondamentale fornire una forza di feedback all’haptic device in modo che

l’utente possa provare la sensazione di toccare il tessuto.

4.2 Il contatto 92

Affinche la forza di feedback della simulazione sia percepita in modo continuo

e realistico dall’utente, dovremmo ottenere una frequenza di aggiornamento

del sistema nell’ordine di 1KHz. Purtroppo la frequenza di aggiornamento

del framework puo variare nell’ordine delle decine di Hz. Questa differenza

tra le frequenza di aggiornamento cui necessita l’interfaccia haptic feedback

e la frequenza di aggiornamento del sistema, puo essere causa di instabilita

nell’haptic device (queste oscillazioni possono creare notevoli problemi all’u-

tente durante la simulazione) .

Il problema della deformazione del modello a contatto con il tool e connesso

al calcolo dell’approssimazione della forza di feedback restituita all’utente.

Durante lo svolgimento di questa tesi non e stato possibile purtroppo dispor-

re di un haptic device, mi sono limitato quindi, in questa fase del lavoro,

alla parte grafica della deformazione. Nelle sezioni che seguono illustrero le

varie tipologie di approcci alla deformazione ed al calcolo della forza di feed-

back, con particolare rilievo all’approccio utilizzato durante lo svolgimento

di questa tesi.

4.2.1 Le forze di di deformazione

Nel sistema di scambio di forze tra utente e oggetto deformabile che stia-

mo descrivendo, l’utente imprime una certa forza Fext al tool che produce

alcune deformazioni sul tessuto. Il tessuto resistera a tale deformazione con

una forza variabile dipendente dai suoi parametri fisici, restituendo cosı una

forza di reazione Ftissue. Maggiore sara la forza impressa dall’utente Fext e

maggiore sara la forza di reazione del tessuto Ftissue che continuera a defor-

marsi. Quando le due forze Fext e Ftissue risultano opposte ma con la stessa

intensita la deformazione avra raggiunto il suo massimo. Mentre nel caso in

cui |Fext| < |Ftissue| il tessuto tornera verso la situazione iniziale di equilibrio.

La forza che viene restituita all’utente equivale a Fext.

4.2 Il contatto 93

Figura 4.8: Un esempio di haptic device: Il phantom 1.5/6 DOF e un force feedback

device che fornisce una forza di feedback con 6 gradi di liberta. Questo

modello, oltre a restituire una forza all’utente offre la possibilita di definire

anche forze di torsione.

Su questo tipo di modello ideale sono state create differenti strategie per il

calcolo della deformazione e della forza Ftissue. Queste strategie approssima-

no il contatto tra il tool ed il tessuto modellando la superficie di contatto

come se fosse un punto (non una superficie come sarebbe nella realta) .

4.2.2 Approccio con molla

Questo metodo propone di modellare l’iterazione tra il tool ed il tessuto

ponendo una molla che connetta il punto in cui si trova l’estremita del tool

Ptool con il punto in cui il tool ha penetrato la superficie di contatto Pproxy.

La forza elastica della molla viene adesso utilizzata per calcolare la forza Fext:

Fext = −Fspring = c · (Ppoint − Pproxy) (4.2)

La forza esterna trovata viene poi applicata su un nodo della mesh producen-

do la deformazione desiderata.La costante c rappresenta la rigidita dinamica

4.2 Il contatto 94

Figura 4.9: Le fasi della deformazione del Tessuto: a) Deformazione b) Stato di

equilibrio c) Rilassamento

del materiale.

Mentre il tessuto e deformato la sue forza elastica Ftissue tendera a riportare il

tool fuori dall’intersezione con il tessuto. Il sistema raggiunge lo stato di equi-

librio quando Fext = −Ftissue (Figura 4.10) . Utilizzando questo approccio,

se non viene utilizzata una costante di elasticita della molla adeguato, pos-

siamo trovarci in situazioni in cui, durante la deformazione, la parte frontale

del tool svanisce all’interno della mesh.

4.2.3 Approccio con Force Grid

Consideriamo un haptic device con 3 gradi di liberta che produca una for-

za 3D in un punto Pproxy. L’oggetto deformabile e immerso in una griglia

regolare W . Ogni vertice di W contiene le forze di reazione che dovrebbe

essere restituita se la superficie dell’oggetto passasse par tale vertice. Tali

valori sono ottenuto proiettando la forza di reazione dei vertici dell’oggetto

sui vertici di W con una funzione basata sulla distanza. Ne segue che i valori

di W vengono aggiornati alla stessa frequenza del rendering visuale.

Disponendo della griglia W , la forza di reazione in un punto generico P vie-

4.2 Il contatto 95

Figura 4.10: Fasi del calcolo dell deformazione utilizzando una molla:La forza

Fext = −Ftissue, l’allungamento della molla da proxypoint a toolpoint

restituisce al tool la forza Fspring

ne ottenuta da quella memorizzata nei vertici della cella che contiene P per

interpolazione trilineare.

4.2.4 Approccio con spostamento delle facce

L’approccio utilizzato in questa tesi si concentra maggiormente sulla mo-

dellazione della deformazione piuttosto che sulla modellazione del contatto

stesso tramite l’haptic feedback device.

In questo metodo il modello viene modificato geometricamente in modo che

il tool non sia mai dentro il volume del modello deformabile; come reazio-

ne a questa deformazione geometrica il tessuto tendera a tornare alla forma

iniziale (di equilibrio) . Il principio su cui si basa questo metodo puo essere

visto come opposto rispetto all’approccio con la molla. Nell’approccio con la

molla viene posta una resistenza artificiale alle forze esterne, la deformazione

avviene quindi attraverso l’allungamento o la contrazione della molla stessa

e la forza di feedback viene calcolata di conseguenza. Il metodo dello sposta-

mento delle facce e invece un metodo di fatto adattivo: lo spostamento delle

4.2 Il contatto 96

Figura 4.11: Una parte di Force Grid

faccia esterna produce una perturbazione nello stato iniziale della mesh; le

proprieta elastiche intrinseche al modello deformabile creano una forza che

tende a riportare la mesh allo stato iniziale, le componenti normali alla su-

perficie costituiranno la forza di feedback.

Come illustrato in Figura 4.12 le facce che entrano in collisione con il tool

vengono spostate affinche i due volumi non si intersechino piu. Questo mo-

vimento crea delle forze interne al modello deformabile Fi0 ed Fi1 che ten-

derebbero a riportare il modello nella posizione di equilibrio. Annullando la

componente di tali forze lungo la normale delle faccia N impediamo che tool

e modello deformabile entrino di nuovo in contatto tra loro.

La forza di feedback sul tool puo essere calcolata tramite un’interpolazione

delle forze presenti sui tre vertici della faccia rispetto alla posizione attuale

del tool.

4.2 Il contatto 97

Figura 4.12: Le tre fasi dell’approccio con traslazione della faccia: a) Il tool interseca

la mesh. b) La faccia viene traslata affinche la mesh non intersechi piu il

tool, i vertici traslati (grazie alle proprieta elastiche del modello) reagi-

scono allo spostamento con una forza (Fi0 e Fi1) che tende a riportarli

nella posizione iniziale c) Le componenti di Fi0 e Fi1 lungo la direzione

della normale N, ovvero F ′i0 e F ′

i1, vengono annullate tramite due forze

opposte che agiscono sugli stessi vertici Fext0 e Fext1

4.3 Il taglio 98

4.3 Il taglio

Le soluzioni proposte finora nello stato dell’arte sono orientate verso una

concezione di taglio guidata dal movimento del tool all’interno del volume

della mesh.

Nell’approccio con suddivisione dei tetraedri, ad esempio, la suddivisione

di un insieme di tetraedri viene guidata dal movimento compiuto dal tool

all’interno del volume di un tetraedro. Questo tipo di rappresentazione del

taglio si basa su una rappresentazione errata della fisica del contatto tool-

tessuto. In realta i due volumi non si intersecano mai; il taglio avviene dove lo

stress esercitato dal tool sulla superficie del modello supera una certa soglia;

il tipo di taglio proposto in questa tesi si basa su questo principio. Quando

una faccia della mesh subisce uno stress troppo alto allora avviene il taglio.

Lo stress sulla faccia e una forza opposta alla forza di feedback ma di uguale

intensita.

Una volta individuato il vettore che definisce lo stress sulla superficie della

mesh dobbiamo scegliere quale dei due operatori di taglio utilizzare ed infine

come utilizzarli in maniera tale che il taglio sia una rappresentazione piu

fedele possibile degli stress a cui e sottoposta la mesh.

4.3.1 Modellazione del taglio guidata dall’errore

Supponiamo di conoscere il punto in cui deve avvenire il taglio e la dire-

zione dello stress imposto dal tool sulla superficie del modello deformabile.

Applicando l’edge-split potremmo modellare il taglio nel modo piu accurato

possibile: La faccia toccata dal tool verrebbe suddivisa nel punto di maggiore

stress (Figura 4.13.b) .

L’utilizzo dell’edge-split incrementa pero il numero di tetraedri della mesh. Il

vertex-split potrebbe essere invece utilizzato definendo una sezione di taglio

in modo tale da dividere i vertici di maggior stress. L’utilizzo del vertex-

4.3 Il taglio 99

Figura 4.13: differente utilizzo degli operatori di taglio: a) Il tool interseca la mesh,

creando i vettori di stress evidenziati b) Utilizzando l’edge split dividia-

mo la mesh esattamente nel punto di maggiore stress c) Utilizzando

l’operazione di vertex-split il taglio viene effettuato sul vertice che ha lo

stress maggiore

4.3 Il taglio 100

split produrrebbe quindi una rappresentazione piu approssimata del taglio;

il punto in cui il tool entra in contatto con superficie (che e anche il punto

di maggior stress) sarebbe cosı approssimato con il vertice, appartenente alla

stessa faccia, che subisce lo stress esterno maggiore (Figura 4.13.b ) . Pur

fornendo una soluzione al taglio meno accurata dell’edge-split, l’operazione

di vertex-split ha il grande vantaggio di non creare nuove entita nella me-

sh. Risulta quindi necessario formulare una strategia di scelta tra questi due

operatori durante la simulazione.

Ogni qualvolta che utilizziamo il vertex-split approssimiamo il taglio al ver-

tice della stessa faccia dove lo stress risulta maggiore. Questo errore puo es-

sere quantificato ed utilizzato per scegliere quale dei due operatori utilizzare.

Possiamo quindi esprimere l’errore come:

Err = |Pstress − Vstress| (4.3)

dove Pstress indica il punto di maggiore stress, mentre Vstress indica il vertice

nel caso in cui utilizzassimo il vertex-split. Questa misura dell’errore, anche

se corretta, e incompleta: dobbiamo considerare l’errore rispetto al punto di

vista dell’utente durante la simulazione.

Se ad esempio, durante il rendering, la mesh e lontana dal punto di vista

dell’utente l’errore commesso sara visibilmente meno evidente dello stesso

errore nel caso in cui avessimo un punto di vista ravvicinato. Possiamo quindi

considerare l’errore commesso come la proiezione dell’errore precedentemente

trovato sul punto di vista dell’utente, ottenendo cosı una stima dell’errore in

pixel. La funzione di OPENGL gluProject permette di ottenere , date

le coordinate 3d di un punto nel sistema di riferimento dell’ambiente, le

coordinate 2d della sua proiezione sulla finestra visualizzata. L’errore in

pixel potra essere quindi espresso come:

Errpixel = |gluProject(Pstress)− gluProject(Vstress)| (4.4)

4.3 Il taglio 101

Se l’errore stimato in pixel supera la tolleranza allora utilizzeremo l’edge-split,

altrimenti il vertex-split.

bool Tolerance(Point3d p0,Point3d p1)

return (ComputeError(p0,p1)<pixel);

void Docut(MESH_TYPE::MTTYPE *t0,int face)

.....// cut \‘e il punto di maggiore stress

//vcut \‘e il vertice nel caso in cui utilizzassimo

//vertex-split

if (!Tolerance(vcut->P(),cut))

Operator1(t0,face);

else

Operator0(vcut);

4.3.2 Vertex-split, calcolo della sezione di taglio

Supponiamo di dover modellare il taglio sul tessuto deformabile tramite l’ap-

plicazione del vertex-split. Come esposto nel precedente paragrafo questo

operatore va applicato in modo da duplicare il vertice di stress maggiore

Vstress.

Per applicare l’operazione di vertex-split dobbiamo innanzitutto definire una

sezione di taglio valida costruita in maniera tale da dividere il vertice Vstress.

Potremmo definirla utilizzando solamente facce che condividono tale vertice.

Siccome Vstress e vertice esterno questa sezione di taglio conterra sicuramen-

te due edge esterni che lo condividono, questo comporta la duplicazione di

Vstress durante l’applicazione del vertex-split. Possiamo costruire piu sezioni

4.3 Il taglio 102

da taglio di questo tipo utilizzando uno stesso vertice Vstress.

La deformazione del tessuto (ovvero lo spostamento della faccia a contatto

col tool) , date le proprieta elastiche del modello, crea delle forze di stress sui

vertici interni della mesh. La sezione di taglio prescelta deve essere costruita

in modo da separare le facce che contengono i vertici dove lo stress concen-

trato e massimo. La costruzione della sezione di taglio puo essere suddivisa

in due fasi:

Scelta della prima faccia Tra tutte le facce interne incidenti sul vertice

scegliamo quella i cui due vertici diversi da Vstress hanno in media uno

stress interno maggiore. Lo stress interno viene valutato come proie-

zione della forza interna sulla forza esterna impressa sul vertice Vstress

(Figura 4.14.a) .

Completamento della sezione di taglio Supponiamo che Finit sia la fac-

cia iniziale selezionata. Come illustrato in Figura 4.14.b , a partire dai

due edge E0 ed E1 vengono selezionate le facce che condividono questi

due edge e che sono piu complanari possibile rispetto alla faccia ini-

ziale. Iterando questo passo i nuovi edge selezionato saranno E2 ed

E3. Questa operazione di selezione viene ripetuta fino a che verranno

selezionati due edge esterni (Figura 4.14.c ) .

4.3 Il taglio 103

Figura 4.14: Creazione di una superficie di taglio: a) La faccia Finit e quella di

maggiore stress interno. b) A partire dai due edge E0 ed E1 vengono

selezionate le facce che condividono questi due edge e che sono piu

complanari possibile rispetto alla faccia iniziale. c) Nell’ultima iterazione

gli edge E4 ed E5 sono esterni.

Capitolo 5

Risultati sperimentali e sviluppi futuri

5.1 Risultati sperimentali

Le Tabelle che seguono sono i risultati dei test eseguiti su un PC Intel Pen-

tium 4 2.53 GHz. Le Tabelle 5.1 e 5.2 riportano rispettivamente i tempi medi

di esecuzione dell’edge-split e del vertex-split espressi in millisecondi ottenuti

utilizzando due mesh tetraedrali di diversa grandezza e forma.

Le tabella 5.3 illustra il tempo medio necessario all’esecuzione di un pas-

so di decimazione. Il test e stato ottenuto effettuando collassi su un numero

molto elevato di edge della mesh senza che intervenissero operazioni esterne

ad incrementare la dimensione in entita della mesh. I collassi effettuati mo-

tetraedri vertici edge-split Tempo Media

768 243 384 140 0.3645

6000 1331 3000 1375 0.4583

Tabella 5.1: Tabella dei tempi di esecuzione dell’edge-split: prime due colonne della

tabella riportano il numero di tetraedri e di vertici della mesh allo stato

iniziale, mentre la terza indica il numero di volte che l’edge-split e stato

eseguito. Seguono il tempo totale di esecuzione e la media di tempo per

operazione espressi in millisecondi.

5.1 Risultati sperimentali 105

tetraedri vertici Split effettuati Split non effettuati Tempo Media

768 243 127 84 15 0.0712

6000 1331 980 142 63 0.05614

Tabella 5.2: Tabella dei tempi di esecuzione del vertex-split: Le prime due colonne

della tabella riportano il numero di tetraedri e di vertici della mesh allo

stato iniziale, mentre la terza e la quarta indicano il numero di volte che

il vertex-split e stato eseguito ed il numero di volte in cui la sezione di

taglio definita non era valida. Seguono il tempo totale di esecuzione e

la media di tempo per operazione espressi in millisecondi.

tetraedri vertici collassi effettuati collassi non effettuati Tempo Media

1168 288 82 101 500 2.73

1168 288 206 494 953 1.36

1168 288 241 814 1047 0.99

Tabella 5.3: Tabella dei tempi della decimazione con tolleranza massima di

cambiamento di volume pari a 1/50 del volume originale

dificano la mesh per un volume minore di 1/50 del volume totale.

L’elevato numero di collassi non effettuati di questo test deriva dalle link

conditions che si fanno sempre piu difficili da soddisfare all’aumentare dei

colassi effettuati. Nel simulatore la decimazione coinvolge solo localmente la

mesh nelle zone in cui sono stati effettuati divisioni di tetraedri, pertanto

non presenta una percentuale di fallimento cosı alta.

La tabella 5.4 riporta invece i tempi necessari al ripristino della proprieta

fisiche della mesh dopo un collasso (redistribuzione delle masse sui vertici,

assegnazione delle proprieta elastiche dei tetraedri modificati) .

Nel grafico riportato in Figura 5.1 vengono riportate le differenti distribu-

zioni relative all’aspect Ratio dei tetraedri dell’intera mesh dopo le operazioni

di taglio riportate nella tabella 5.5.


tetraedri vertici collassi effettuati Tempo Media

768 243 3840 31 0.00807

1168 288 5840 47 0.00807

6000 1331 30000 266 0.0088

Tabella 5.4: Tabella dei tempi di ripristino delle proprieta fisiche della mesh dopo un

collasso

Figura 5.1: Il grafico di sinistra riporta la distribuzione dell’aspect ratio dei tetraedri

creati dopo il test riportato in tabella 5.5. Il grafico a destra illustra il

numero di tetraedri totali dopo l’esperimento.


ratio media tetraedri vertici e-v e-s collassi

0.201192 660 201 0 0 0

0.159835 1227 303 16 83 0

0.177829 724 239 16 83 64

Tabella 5.5: Risultati finali del test illustrato nel grafico di distribuzione dell’Aspect

Ratio di Figura 5.1. La prima riga della tabella mostra i dati della me-

sh allo stato iniziale, senza che nessuna operazione di taglio sia stata

effettuata. Nella seconda riga invece sono stati effettuati tagli ma non

collassi, mentre la terza riga mostra l’effetto delle stesse modifiche ap-

portate alla mesh con l’utilizzo dei collassi (tolleranza sulla variazione di

volume sui collassi pari a 1/1000 del volume originale) .

Emerge dal grafico che la distribuzione dell’aspect ratio senza l’utilizzo della

decimazione e peggiore di quella ottenuta dalle stesse operazioni di taglio ac-

compagnate dalla decimazione. La distribuzione dell’aspect ratio nel caso del

test effettuato con l’utilizzo della decimazione si avvicina alla distribuzione

dell’aspect ratio della mesh allo stato iniziale.

Nella parte destra del grafico viene riportato il numero dei tetraedri compo-

nenti la mesh dopo le operazioni di taglio. Anche in questo caso la decima-

zione riduce in modo evidente il numero totale dei tetraedri.

Da queste tabelle emerge il vantaggio apportato alla simulazione dalla de-

cimazione in termini di riduzione e miglioramento della qualita dei tetrae-

dri della mesh.Come espresso nella tabella 5.3 il passo di decimazione viene

eseguito ad una velocita media di circa 1msec.

La tabella 5.7 riporta il numero di entita create per ciascun utilizzo del

vertex-split. Questa tabella mostra innanzitutto che all’aumentare della tol-

leranza in pixel viene utilizzato in maniera minore l’edge-split. Anche i col-

lassi aumentano al diminuire della tolleranza in pixel. Infatti piu tetraedri

vengono suddivisi e maggiore sara il numero di edge da collassare nell’heap

5.2 Sviluppi futuri 108

vertex-split edge-split collassi tetraedri vertici pixel

19.3 43 19.3 75 33.6 3

23.3 19.6 10.3 30 46.6 4

Tabella 5.6: Tabella delle entita create e del numero di operazioni coinvolte in una

serie di tagli con differenti approssimazioni.

pixel edge-split collassi tetraedri vertici

2 3.222 2.444 0.666 2.222

3 2.375 1.625 2.125 2

4 2.14 1.45 0.71 0.71

Tabella 5.7: Tabella delle medie di entita create e del numero di operazioni coinvolte

per ogni split effettuato con differenti approssimazioni.

di decimazione.

Da questa tabella possiamo dedurre che ciascuna divisione di un vertice com-

porta un aumento del numero di tetraedri di circa 0.666 nel caso di una tol-

leranza di 2 pixel, di 2.125 nel caso di una tolleranza di 3 pixel e di 0.71 nel

caso di una tolleranza di 4 pixel. Il minore aumento di tetraedri nel caso del

taglio con tolleranza pari a 2 pixel e dovuto all’elevato numero di collassi che

ne consegue.

L’efficienza del metodo di taglio proposto in questa tesi risiede proprio nel

basso numero di entita create durante il taglio: I metodi basati sulla Look-up

Table creano in media un numero di tetraedri per taglio di circa 4.5 [GCMS00]

contro un risultato ottenuto in questa Tesi che varia dagli 0.6 ai 2.2.

5.2 Sviluppi futuri

Come illustrato nei capitoli precedenti il taglio del modello viene guidato

dall’interazione delle forze elastiche interne e da quelle esterne impresse dal


Figura 5.2: a) Il solido S entra in contatto con la mesh attraverso la superficie C

esercitando un pressione. b) Il bisturi T avendo una superficie di contatto

minore del solido S, pur imprimendo la stessa forza, esercita una pressione

maggiore.

bisturi. Alla base di questa metodologia le operazioni di modifica geometrica

sulla mesh sono create in maniera tale da rappresentare situazioni reali di

modifica della struttura del tessuto.

Questo tipo di rappresentazione oltre a rendere la simulazione del taglio piu

realistica offre una vasta gamma di possibilita di sviluppi futuri. I paragrafi

che seguono illustrano alcuni degli sviluppi futuri basati sulla struttura dati

e sulle operazioni di modifica su di essa definite.

5.2.1 Interazione con un solido

In una operazione chirurgica esistono altri strumenti oltre che al bisturi che

un chirurgo puo utilizzare. Modellando l’interazione di forze tra un solido

qualsiasi e la mesh deformabile avremmo la possibilita di rappresentare nuovi


Figura 5.3: a) Le forze F0 ed F1 propagano forze di stress di direzione opposta sui

vertici delle due facce evidenziate. b) Tramite vertex-split viene realizzata

la lacerazione.

strumenti utilizzabili dall’utente durante la simulazione.

Come rappresentato in Figura 5.2.a, il solido S modifica la posizione delle

facce con cui entra in contatto, annullando le componenti delle forze interne

che tendono a riportare la faccia in posizione di rilassamento.

La pressione che il solido imprime su una faccia della mesh e il rapporto tra

questa forza e la superficie di contatto tra i due oggetti C. In proporzione alla

pressione esercitata possiamo cosı stabilire quando applicare le operazioni di

modifica. Nel caso di un bisturi che imprime una forza utilizzando la sua

lama, la pressione esercitata sara sicuramente maggiore della stessa forza

impressa dal solido S.

5.2.2 Lacerazioni

Quando due forze di direzione opposta vengono esercitate su un tessuto de-

formabile si potrebbero creare delle lacerazioni.

Generalizzando la definizione della superficie di taglio possiamo essere in


grado di separare due insiemi qualsiasi di tetraedri adiacenti. La lacerazione

potrebbe avvenire in corrispondenza delle facce i cui vertici sono sottoposti

a due diverse forze esterne di direzione opposta. Su queste facce potrebbe

essere definita una sezione di taglio che porterebbe poi, tramite l’utilizzo del

vertex-split, alla separazione dei vertici (Figura 5.3) .

5.2.3 Self Collision

Il termine self collision indica l’autointersezione della superficie di un ogget-

to. Nel caso in cui un oggetto viene tagliato e probabile che le parti separate,

che si trovano in prossimita, si tocchino. Segue quindi la necessita di intro-

durre tra i lavori futuri lo sviluppo di tecniche per la modellazione della self

collision, che e ancora un problema aperto nella simulazione interattiva di

oggetti deformabili[JTT01].

Inoltre questo problema assume particolare rilievo nel caso in cui si volesse

inserire nel simulazione la possibilita di effettuare operazioni di sutura.

5.2.4 Sutura

Un’altra possibilita importante da introdurre in un simulatore di operazioni

chirurgiche sarebbe quella di poter effettuare suture. La simulazione di una

sutura introduce due problemi:

Rappresentazione dell’ago L’inserimento dell’ago nel tessuto non puo es-

sere modellato attraverso l’interazione di un solido con la mesh te-

traedrale. Il foro provocato della pressione di un ago nel tessuto non

puo essere il risultato degli operatori di modifica della mesh definiti e

comporta quindi un trattamento a parte.

Triangolazione delle regioni di tessuto rimarginate Se volessimo simu-

lare la rimarginazione della ferita suturata il problema maggiore sa-


Figura 5.4: a) Le due mesh devono essere attaccate. b) Il vero problema consiste nel

ricreare una triangolazione unica della superficie di contatto.

rebbe quello di costruire una triangolazione unica delle due regioni di

superficie esterna che vengono attaccate tra loro (Figura 5.4).

Modellazione dell’ago Seppure la forza di resistenza offerta dal tessuto

alla penetrazione dell’ago risulta minima l’inserimento di una ago nel tessuto

comporta una maggiore resistenza alle forze di movimento dell’ago che non

siano nella direzione del foro provocato (ad esempio forze di torsione). Non

possiamo quindi utilizzare uno degli operatori definiti in quanto portereb-

be alla creazione di una superficie di taglio lungo la quale l’ago potrebbe

muoversi liberamente. Un’idea di possibile implementazione potrebbe essere

quella di inserire un vertice nel punto della faccia esterna toccata dall’ago.

Questa operazione comportera la creazione di nuovi tetraedro ed in partico-

lare di un edge interno incidente sul vertice appena creato. La traiettoria

dell’ago, e quindi del filo della sutura puo essere approssimato con uno degli

edge interni alla mesh che incidono su questo vertice (Figura 5.5) .

5.3 Conclusioni 113

Figura 5.5: a) L’ago preme con una forza F sulla faccia esterna nel punto evidenziato.

b) Una operazione di remeshing locale crea un nuovo vertice sulla faccia

esterna selezionata, l’ago si muovera lungo l’edge appena creato.

5.3 Conclusioni

I risultati ottenuti mostrano che l’approccio proposto in questa tesi, basato

su tecniche di raffinamento e semplificazione implementati in modo interrom-

pibile e guidato dall’errore percepito dall’ utente permettono di conseguire

l’obiettivo che ci eravamo proposti, cioe di simulare in maniera interattica

deformazioni e tagli.

Infatti le operazioni di taglio, accompagnate dalle operazioni di semplifi-

cazione, limitano il fenomeno della frammentazione nel sistema. Per questo

motivo le prestazioni del sistema non vengono di fatto alterate dalle diverse

operazioni di taglio effettuate dall’utente. Uno stesso taglio puo essere ese-

guito secondo piu livelli di approssimazione: un taglio piu accurato genera

piu frammentazione di un taglio eseguito in maniera piu approssimata. Il

livello di accuratezza viene scelto in maniera tale che l’approssimazione non

sia percepita dall’utente.

I ridotti tempi di esecuzione delle operazioni di modifica e di semplificazione

mostrano inoltre l’efficacia di aggiornamento della struttura dati che mantie-

5.3 Conclusioni 114

ne la geometria e la topologia della mesh tetraedrale.

I metodi di taglio e semplificazione sono interrompibili e restituiscono

una struttura dati consistente. L’alta granularita delle operazioni di aggior-

namento permette di ottimizzare la ripartizione del tempo di esecuzione tra

le varie componenti del sistema (raffinamento, semplificazione, simulazione

. . . ).

Le modifiche dell’oggetto dipendono solo dallo stress cui e sottoposto.

Potremmo quindi utilizzare questo modello per implementare l’interazione

con qualsiasi altro tipo di strumento per il quale siamo in grado di modellarne

il contatto fisico con l’oggetto deformabile.

Quanto esposto in questa tesi e stato oggetto di un articolo proposto per

pubblicazione alla conferenza Medical Simulation 2004 [NCFP04].

5.4 Alcune immagini del simulatore 115

5.4 Alcune immagini del simulatore

Figura 5.6: Un esempio di deformazione ottenuta premendo il bisturi di taglio.


Figura 5.7: Altri esempi di deformazione: A destra le stesse deformazioni sono

rappresentate mettendo in evidenza i tetraedri che compongono la mesh.


a) d)

b) e)

c) f)

Figura 5.8: Fasi di utilizzo dell’edge-split per remeshing con due differenti tipi di

risoluzione: 3 pixel a sinistra e 4 pixel a destra. a) e d) La deformazione

vista dall’utente. b) e e) Vengono suddivisi gli edges dove lo stress e

maggiore. c) e f) Le operazioni di collasso riducono il numero di tetraedri

della mesh.


Figura 5.9: La superficie di taglio (costituita dai triangoli con i bordi evidenziati)

suddivide l’insieme dei tetraedri colorati in due sottoinsiemi.


a)

b)

Figura 5.10: Esempio di taglio con tolleranza di 4 pixel: a) Con una visione della

mesh piu lontana il taglio viene approssimato al vertice piu vicino. b)

Avvicinando l’inquadratura si puo notare l’approssimazione.


Figura 5.11: Due esempi di taglio:In alto viene asportata una parte della mesh, in

basso viene tagliato un angolo


Figura 5.12: Un esempio di taglio su una mesh di forma sferica

Elenco delle figure

1.1 Schema generico di un simulatore . . . . . . . . . . . . . . . . 4

2.1 Superficie di taglio . . . . . . . . . . . . . . . . . . . . . . . . 22

2.2 Differenti tecniche di taglio esistenti . . . . . . . . . . . . . . . 23

2.3 Suddivisione di un tetraedro . . . . . . . . . . . . . . . . . . . 26

2.4 Possibili suddivisioni di una faccia . . . . . . . . . . . . . . . . 27

2.5 Esempio di Lookup Table . . . . . . . . . . . . . . . . . . . . 29

3.1 Topologia di un tetraedro . . . . . . . . . . . . . . . . . . . . 34

3.2 Un primo esempio di mesh tetraedrale . . . . . . . . . . . . . 37

3.3 Adiacenza tetraedro-tetraedro . . . . . . . . . . . . . . . . . . 39

3.4 Adiacenza vertice-tetraedro . . . . . . . . . . . . . . . . . . . 41

3.5 Half-edge in una spezzata ed in una mesh di triangoli . . . . . 44

3.6 Half-edge in a mesh tetraedrale . . . . . . . . . . . . . . . . . 45

3.7 Possibili flip su una mesh di tetraedri . . . . . . . . . . . . . . 47

3.8 Operazioni di next su un half-edge . . . . . . . . . . . . . . . 48

3.9 Next edge come composizione di flip . . . . . . . . . . . . . . 49

3.10 problematica dell’approccio con duplicazione di facce . . . . . 51

ELENCO DELLE FIGURE 123

3.11 Operatori di taglio . . . . . . . . . . . . . . . . . . . . . . . . 53

3.12 Problemi sul vertex-split . . . . . . . . . . . . . . . . . . . . . 54

3.13 Esempio di estensione di una sezione di taglio . . . . . . . . . 55

3.14 Dimostrazione del Lemma sulla correttezza della mesh dopo il

taglio . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55

3.15 Esempio di V Split di una sezione di taglio . . . . . . . . . . . 56

3.16 Condizioni sulla sezione di taglio del vertex-split . . . . . . . . 58

3.17 Condizione di applicabilita del vertex-split . . . . . . . . . . . 59

3.18 Individuazione dei vertici da duplicare del vertex-split . . . . . 60

3.19 Individuazione dei sottoinsiemi Ta e Tb nel esecuzione di un

vertex-split . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62

3.20 Passi per l’individuazione dei sottoinsiemi Ta e Tb nel vertex-split 62

3.21 Passi per l’esecuzione dell’edge-split . . . . . . . . . . . . . . . 65

3.22 Suddivisione di un tetraedro lungo un edge . . . . . . . . . . . 66

3.23 Suddivisione di un edge . . . . . . . . . . . . . . . . . . . . . . 67

3.24 Edge-split: Ripristino delle adiacenze faccia tetraedro . . . . . 69

3.25 Semplificazione schema generale . . . . . . . . . . . . . . . . . 72

3.26 Esempio di collasso di edge su una mesh di triangoli . . . . . . 75

3.27 Esempio di Link Conditions . . . . . . . . . . . . . . . . . . . 76

3.28 Esempio di Link Conditions su una mesh tetraedrale . . . . . 77

4.1 Collisione Statica e Collisione Dinamica . . . . . . . . . . . . . 83

4.2 La Camera in OpenGL . . . . . . . . . . . . . . . . . . . . . . 85

4.3 Collision detection con Camera ortogonale . . . . . . . . . . . 87


4.4 Esempio di Collision detection con Camera ortogonale . . . . . 88

4.5 Estensione della regione di spazio tra due time step successivi

nella collision detection. . . . . . . . . . . . . . . . . . . . . . 89

4.6 Fasi di raffinamento della collision detection . . . . . . . . . . 90

4.7 Quantificazione della penetrazione all’interno di un tetraedro . 91

4.8 Un esempio di haptic device:Il phantom 1.5/6 DOF . . . . . . 93

4.9 Le fasi della deformazione del Tessuto . . . . . . . . . . . . . . 94

4.10 Approccio della deformazione con molla . . . . . . . . . . . . . 95

4.11 Force Grid . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96

4.12 Approccio con spostamento della faccia . . . . . . . . . . . . . 97

4.13 differente utilizzo degli operatori di taglio . . . . . . . . . . . . 99

4.14 Creazione di una superficie di taglio . . . . . . . . . . . . . . . 103

5.1 Distribuzione dell’aspect ratio di una mesh tetraedrale con e

senza collassi dopo operazioni di modifica . . . . . . . . . . . . 106

5.2 Interazione tra un solido e la mesh deformabile . . . . . . . . . 109

5.3 Esempio di lacerazione . . . . . . . . . . . . . . . . . . . . . . 110

5.4 Problema della rimarginazione delle ferite . . . . . . . . . . . . 112

5.5 Esempio di inserimento di un ago . . . . . . . . . . . . . . . . 113

5.6 Esempio di deformazione . . . . . . . . . . . . . . . . . . . . . 115

5.7 Esempi di deformazione . . . . . . . . . . . . . . . . . . . . . 116

5.8 Utilizzo dell’edge-split per remeshing con tolleranza di 3 pixel

e 4 pixel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117

5.9 Schematizzazione del vertex-split . . . . . . . . . . . . . . . . 118


5.10 Esempio di taglio con tolleranza di 4 pixel . . . . . . . . . . . 119

5.11 Due esempi i taglio . . . . . . . . . . . . . . . . . . . . . . . . 120

5.12 Un esempio di taglio su una mesh di forma sferica . . . . . . . 121

Elenco delle tabelle

3.1 Tabelle corrispondenze entita in un tetraedro . . . . . . . . . . 37

5.1 Tabella dei tempi di esecuzione dell’edge-split . . . . . . . . . 104

5.2 Tabella dei tempi di esecuzione del vertex-split . . . . . . . . . 105

5.3 Tabella dei tempi di decimazione . . . . . . . . . . . . . . . . 105

5.4 Tabella dei tempi di ripristino delle proprieta fisiche della mesh

dopo un collasso . . . . . . . . . . . . . . . . . . . . . . . . . . 106

5.5 Risultati finali di operazioni di modifica con l’utilizzo o meno

dei collassi. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 107

5.6 numero di entita create in alcune operazioni di taglio . . . . . 108

5.7 Tabella delle medie di entita create e del numero di operazioni

coinvolte per ogni split . . . . . . . . . . . . . . . . . . . . . . 108

Bibliografia

[Bar] David Baraff. Implicit method for differential equation.

[BK] Andrew B.Mor and Takeo Kanade. Modifyng soft tissue mo-

dels:progressive cutting with minimal element creation. Techni-

cal report, Center for Medical Robotics and Computer Assisted

Surgery,Carnegie Mellon University Pittsbourgh.

[BMG99] D. Bielser, V.A. Maiwald, and M.H. Gross. Interactive cuts

through 3-dimensional soft tissue. Computer Graphics Forum

(Eurographics’99 Proc.), 18(3):C31–C38, September 1999.

[CA99] H. Delingette S. Cotin and N. Ayache. A hybrid elastic model

allowing real-time cutting, deformations and force-feedback for

surgery training and simulation. In CAS99 Proceedings, pages

70–81, May 1999.

[CCM+00] P. Cignoni, D. Costanza, C. Montani, C. Rocchini, and R. Sco-

pigno. Simplification of tetrahedral volume with accurate error

evaluation. Technical Report B4-05, I.E.I. – C.N.R., Pisa, Italy,

March 2000.

[DEGN99] T.K. Dey, H. Edelsbrunner, S. Guha, and D.V. Nekhayev.

Topology preserving edge contraction. Technical Report

RGI-Tech-99, RainDrop Geomagic Inc. Champaign IL., 1999.

BIBLIOGRAFIA 128

[FDA02] C. Forest, H. Delingette, and N. Ayache. Cutting simulation of

manifold volumetric meshes, 2002.

[GCMS00] Fabio Ganovelli, Paolo Cignoni, Claudio Montani, and Roberto

Scopigno. A multiresolution model for soft objects supporting

interactive cuts and lacerations. Computer Graphics Forum,

19(3):??–??, 2000.

[GO01] F. Ganovelli and C. O’Sullivan. Animating cuts with on-the-

fly re-meshing. EuroGraphics Short Presentations, 2001. (J. C.

Roberts, editor), 2001.

[JTT01] P. Jimnez, F. Thomas, and C. Torras. 3D Collision Detection:

A Survey. Computers and Graphics, 25(2):269–285, April 2001.

[LCN99] Jean-Christophe Lombardo, Marie-Paule Cani, and Fabrice

Neyret. Real-time collision detection for virtual surgery. In

CA, pages 82–, 1999.

[MLBdC01] C. Mendoza, C. Laugier, and F. Boux-de Casson. Virtual rea-

lity cutting phenomena using force feedback for surgery simula-

tions. In Interactive Medical Image Visualization and Analisys,

MICCAI, Holland, 2001.

[NCFP04] N.Pietroni, C.Mendoza, F.Ganovelli, and P.Cignoni. A robust

adaptive local remesher (alr) for cutting deformable objects ba-

sed on tetrahedral meshes. Medical Simulation 2004 (submitted),

2004.

[Nie03] H.W. Nienhuys. Cutting in Deformable Objects. Ph.d. thesis,

Utrecht University, 2003.

BIBLIOGRAFIA 129

[SC98] Nicholas Ayache Stephane Cotin, Herve Delingette. Efficient

linear elastic models of soft tissues for real-time surgery si-

mulation. Technical report, Institut National de Recherche en

Informoatique et en Automatique, 1998.

[TG85] S.P. Timoshenko and J.N. Goodier. Theory of elasticity.

McGRAW-HILL International Editions, 1985.

[Wita] Andrew Witkin. Particle system dynamics.

[Witb] David Baraff Andrew Witkin. Differential equation basics.

[Zie67] O. C. Zienkiewicz. The Finite Element Method in Structural

and Continuum Mechanics. McGraw-Hill, 1967.

Algoritmi e Strutture Dati per la Modellazione Interattiva...

Documents

Transcript of Algoritmi e Strutture Dati per la Modellazione Interattiva...