Detecting Topological Changes in Dynamic Delaunay Triangulations - - PowerPoint PPT Presentation

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA

GTC 2017 Speaker: Dr. Matteo Lulli

• Prof. M. Bernaschi and Prof. M. Sbragaglia

May the 10th, 2017

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 1 / 26

Outlook

1

Motivation and Definitions

2

Physical system and Simulations

3

Algorithm Description

4

Physics Results

5

Conclusions and Perspectives

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 2 / 26

Outline for section 1

1

Motivation and Definitions

2

Physical system and Simulations

3

Algorithm Description

4

Physics Results

5

Conclusions and Perspectives

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 3 / 26

Motivation...

Many systems are composed of a number of entities which dynamically interact with each other in real continuous space, i.e. not on a lattice...

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 4 / 26

Motivation...

Many systems are composed of a number of entities which dynamically interact with each other in real continuous space, i.e. not on a lattice... Pedestrians on the street Birds in a flock Droplets in an emulsion Many others...

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 4 / 26

Motivation...

Many systems are composed of a number of entities which dynamically interact with each other in real continuous space, i.e. not on a lattice... Pedestrians on the street Birds in a flock Droplets in an emulsion Many others...

How to associate a topology?

There are different ways to associate a topology to a given set of points: triangulations, quadrangulations... In the literature it is common to use Voronoi tesselation (dual to Delaunay) Delaunay triangulation (dual to Voronoi)

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 4 / 26

Motivation...

Many systems are composed of a number of entities which dynamically interact with each other in real continuous space, i.e. not on a lattice... Pedestrians on the street Birds in a flock Droplets in an emulsion Many others...

How to associate a topology?

There are different ways to associate a topology to a given set of points: triangulations, quadrangulations... In the literature it is common to use Voronoi tesselation (dual to Delaunay) Delaunay triangulation (dual to Voronoi)

Topology changes

Many important properties on the dynamics of such sets of points can be studied by analyzing changes occurring in the topology during the time evolution of the system, e.g., topological changes are associated to stress redistribution in emulsions...

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 4 / 26

...and Definitions

Commonly used Voronoi tesselation: given the set of points {zi}, the Voronoi cell associated to zi is a convex region of the plane defined as Ci = {x ∈ R2 | |x − zi| < |x − zk| ∀ k = i} The simplest topology change is the T1 neighbours switching:

• ne shared edge of two adjacent Voronoi cells degenerates in a

point

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 5 / 26

...and Definitions

Commonly used Voronoi tesselation: given the set of points {zi}, the Voronoi cell associated to zi is a convex region of the plane defined as Ci = {x ∈ R2 | |x − zi| < |x − zk| ∀ k = i} The simplest topology change is the T1 neighbours switching:

• ne shared edge of two adjacent Voronoi cells degenerates in a

point We follow an alternative approach based on the dual of the Voronoi tesselation: the Delaunay triangulaton in which topology changes correspond to changes in the triangulation There can be more complex events in bulk and in the presence

• f boundaries
• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 5 / 26

...and Definitions

Commonly used Voronoi tesselation: given the set of points {zi}, the Voronoi cell associated to zi is a convex region of the plane defined as Ci = {x ∈ R2 | |x − zi| < |x − zk| ∀ k = i} The simplest topology change is the T1 neighbours switching:

• ne shared edge of two adjacent Voronoi cells degenerates in a

point We follow an alternative approach based on the dual of the Voronoi tesselation: the Delaunay triangulaton in which topology changes correspond to changes in the triangulation There can be more complex events in bulk and in the presence

• f boundaries
• A. Corbetta PhD thesis, “Multiscale crowd

dynamics: physical analysis, modeling and applications”

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 5 / 26

Our Approach

We propose a unified approach for tracking topological changes occurring in the Delaunay triangulations of the points composing 2D systems We built up a fast CUDA library which can be integrated with already existing simulations or data analysis programs

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 6 / 26

Our Approach

We propose a unified approach for tracking topological changes occurring in the Delaunay triangulations of the points composing 2D systems We built up a fast CUDA library which can be integrated with already existing simulations or data analysis programs

Features

Lattice based algorithm Library interface only needs a i) density field, ii) a threshold, iii) (optional) a map for irregular geometries Arbitrary obstacles geometry handling Periodic and non-periodic boundary conditions handled in both directions Negligible overhead on our Lattice Boltzmann simulations: check for topology changes every ∆t = 102 steps on

• verall simulation time ∼ 107 steps

Variable number of generating points... The implementation is composed by several elements that can be individually replaced according to the actual needs Our code can be downloaded from http://twin.iac.rm.cnr.it/dynvorcuda.tgz

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 6 / 26

Our Approach

We propose a unified approach for tracking topological changes occurring in the Delaunay triangulations of the points composing 2D systems We built up a fast CUDA library which can be integrated with already existing simulations or data analysis programs

Features

Lattice based algorithm Library interface only needs a i) density field, ii) a threshold, iii) (optional) a map for irregular geometries Arbitrary obstacles geometry handling Periodic and non-periodic boundary conditions handled in both directions Negligible overhead on our Lattice Boltzmann simulations: check for topology changes every ∆t = 102 steps on

• verall simulation time ∼ 107 steps

Variable number of generating points... The implementation is composed by several elements that can be individually replaced according to the actual needs Our code can be downloaded from http://twin.iac.rm.cnr.it/dynvorcuda.tgz

Current limitations

Lattice based Triangulation relies on a non-compressed adjacency matrix Order ∼ 104 points

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 6 / 26

Outline for section 2

1

Motivation and Definitions

2

Physical system and Simulations

3

Algorithm Description

4

Physics Results

5

Conclusions and Perspectives

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 7 / 26

What are emulsions?

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 8 / 26

What are emulsions?

Ensembles of Droplets that do not coalesce - surfactants...

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 8 / 26

What are emulsions?

Ensembles of Droplets that do not coalesce - surfactants... Display both solid and liquid features: yield stress σY

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 8 / 26

What are emulsions?

Ensembles of Droplets that do not coalesce - surfactants... Display both solid and liquid features: yield stress σY Solid-like behaviour - elastic response: σ < σY Liquid-like behaviour - flow is non-Newtonian: σ > σY , σ ∝ ˙ γn, n = 1

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 8 / 26

What are emulsions?

Ensembles of Droplets that do not coalesce - surfactants... Display both solid and liquid features: yield stress σY Solid-like behaviour - elastic response: σ < σY Liquid-like behaviour - flow is non-Newtonian: σ > σY , σ ∝ ˙ γn, n = 1 Emulsions flow when the droplets configuration undergoes irreversible topological changes: plastic events

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 8 / 26

Lattice Boltzmann Model for Emulsions

Starting point: Shan-Chen multi-component forced Lattice Boltzmann fσ,i(x + ci, ci, t + 1) = − 1 τσ

• (1 − τσ)fσ,i(x, ci, t) − feq

σ,i(x, ci, t)

+ Fσ,i(x, ci, t) with the usual equilibrium distribution and intercomponent forcing F a

µ (r) = −ψµ(r)

• ν,i

Gµν(r)wiψν(r + ci)ca

i

• 1R. Benzi, M. Bernaschi, M. Sbragaglia, S. Succi and S. Chibbaro, Journal of Chemical Physics 131, 104903 (2009).
• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 9 / 26

Lattice Boltzmann Model for Emulsions

Starting point: Shan-Chen multi-component forced Lattice Boltzmann fσ,i(x + ci, ci, t + 1) = − 1 τσ

• (1 − τσ)fσ,i(x, ci, t) − feq

σ,i(x, ci, t)

+ Fσ,i(x, ci, t) with the usual equilibrium distribution and intercomponent forcing F a

µ (r) = −ψµ(r)

• ν,i

Gµν(r)wiψν(r + ci)ca

i

How to avoid droplets coalescence?

• 1R. Benzi, M. Bernaschi, M. Sbragaglia, S. Succi and S. Chibbaro, Journal of Chemical Physics 131, 104903 (2009).
• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 9 / 26

Lattice Boltzmann Model for Emulsions

Starting point: Shan-Chen multi-component forced Lattice Boltzmann fσ,i(x + ci, ci, t + 1) = − 1 τσ

• (1 − τσ)fσ,i(x, ci, t) − feq

σ,i(x, ci, t)

+ Fσ,i(x, ci, t) with the usual equilibrium distribution and intercomponent forcing F a

µ (r) = −ψµ(r)

• ν,i

Gµν(r)wiψν(r + ci)ca

i

How to avoid droplets coalescence? Augment the set of velocities for the discretization of the equilibrium velocity distribution, double belt D2Q251, and impose frustrated self-interactions

1

Short-range attractive interactions Gµµ < 0

2

Long- and mid-range repulsive interactions Gµµ > 0 Disjoining pressure between interfaces: coalescence is hampered

• 1R. Benzi, M. Bernaschi, M. Sbragaglia, S. Succi and S. Chibbaro, Journal of Chemical Physics 131, 104903 (2009).
• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 9 / 26

Outline for section 3

1

Motivation and Definitions

2

Physical system and Simulations

3

Algorithm Description

4

Physics Results

5

Conclusions and Perspectives

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 10 / 26

Algorithm Description

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 11 / 26

Flagging

We associate to each mesh node, a flag taking on a limited number of values: EMPTY, BUBBLE, UW, DW, LW, RW, OBS

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 12 / 26

Flagging

We associate to each mesh node, a flag taking on a limited number of values: EMPTY, BUBBLE, UW, DW, LW, RW, OBS The assignment of the EMPTY or BUBBLE flag depends on the value of the density with respect to a pre-set threshold ρσ(xi) ≶ ρ∗

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 12 / 26

Flagging

We associate to each mesh node, a flag taking on a limited number of values: EMPTY, BUBBLE, UW, DW, LW, RW, OBS The assignment of the EMPTY or BUBBLE flag depends on the value of the density with respect to a pre-set threshold ρσ(xi) ≶ ρ∗ This is the entry point of the library, from now on all the analysis is based on this data structure only

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 12 / 26

Labelling

The points in the bulk of the domain assumes binary values (EMPTY, BUBBLE): Ising model-like

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 13 / 26

Labelling

The points in the bulk of the domain assumes binary values (EMPTY, BUBBLE): Ising model-like We can use GPU algorithms for Swendsen-Wang cluster-flip dynamics in order to label the droplets in our Lattice Boltzmann simulations [Komura Y, Okube Y, Comput. Phys. Comm, 185 (2014).]

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 13 / 26

Labelling

The points in the bulk of the domain assumes binary values (EMPTY, BUBBLE): Ising model-like We can use GPU algorithms for Swendsen-Wang cluster-flip dynamics in order to label the droplets in our Lattice Boltzmann simulations [Komura Y, Okube Y, Comput. Phys. Comm, 185 (2014).] We had to modify them in order to move along all 2D-LB velocities

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 13 / 26

Centers of Mass

Once we have the labels we can compute the centers of mass by using reduction operations

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 14 / 26

Centers of Mass

Once we have the labels we can compute the centers of mass by using reduction operations We resorted to the Thrust library to compute the size of each labelled domain (as if it were a histogram)

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 14 / 26

Centers of Mass

Once we have the labels we can compute the centers of mass by using reduction operations We resorted to the Thrust library to compute the size of each labelled domain (as if it were a histogram) We mark each domain size with a binary mask if it sits across a periodic boundary x-Periodic atomicOr(d labcnt + d label[k], 0x80000000); y-Periodic atomicOr(d labcnt + d label[k], 0x40000000); We unmark the domain size when computing the center of mass cond = -((labcnt & 0x80000000) >> 31); labcnt = ((cond & (labcnt ˆ 0x80000000)) | ((˜ cond) & labcnt));

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 14 / 26

Digital Voronoi

Once we have the centers of mass we can compute the Digital Voronoi tesselation

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 15 / 26

Digital Voronoi

Once we have the centers of mass we can compute the Digital Voronoi tesselation We naively find the nearest center of mass for each point of the flag lattice

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 15 / 26

Digital Voronoi

Once we have the centers of mass we can compute the Digital Voronoi tesselation We naively find the nearest center of mass for each point of the flag lattice Voronoi cells inherit droplets labels

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 15 / 26

Vertices Detection

Once we have the Digital Voronoi tesselation we can find the vertices of the Voronoi cells

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 16 / 26

Vertices Detection

Once we have the Digital Voronoi tesselation we can find the vertices of the Voronoi cells In our approach we follow a different route, with respect to the literature, that reduces the burden

• f managing Euclidean colouring pathologies on

the lattice [G. Rong, et al., I3D ’08, ACM, 2008,

• pp. 89-97; M. Qi, et al., I3D ’12, ACM, 2012, pp.

39-46.]

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 16 / 26

Vertices Detection

Once we have the Digital Voronoi tesselation we can find the vertices of the Voronoi cells In our approach we follow a different route, with respect to the literature, that reduces the burden

• f managing Euclidean colouring pathologies on

the lattice [G. Rong, et al., I3D ’08, ACM, 2008,

• pp. 89-97; M. Qi, et al., I3D ’12, ACM, 2012, pp.

39-46.] We identify Voronoi vertices by counting how many gaps among colours are present around each point

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 16 / 26

Vertices Detection

k k k k Voronoi cell edge: 2 gaps Voronoi cell vertex: 3 gaps k Voronoi cell bulk: no gaps

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 16 / 26

Vertices Detection

Isolated point: 2 gaps 3 gaps 4 gaps: prolongued vertex 3 gaps 4 gaps: prolongued vertex

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 16 / 26

Vertices Detection

k spx spy smypx k spx spy spxpy spxpy ak ak a) b) smy spymx smypx

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 16 / 26

Delaunay Triangulation

We look at the triangulation as a graph represented by its adjacency list

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 17 / 26

Delaunay Triangulation

We look at the triangulation as a graph represented by its adjacency list When we find a vertex we store the topological information in a structure similar to an adjacency matrix: each non-zero entry of a row takes on the value of the column position, each zero entry takes on a very high value MAXLINK = 0x77777777

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 17 / 26

Delaunay Triangulation

We look at the triangulation as a graph represented by its adjacency list When we find a vertex we store the topological information in a structure similar to an adjacency matrix: each non-zero entry of a row takes on the value of the column position, each zero entry takes on a very high value MAXLINK = 0x77777777 This allows us to compress the matrix contents via a simple sort which we perform by using the moderngpu library

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 17 / 26

Extraction of the Topological Change

Once we have two subsequent Delaunay triangulations (labelled as Now and Past) we can compare them

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 18 / 26

Extraction of the Topological Change

Once we have two subsequent Delaunay triangulations (labelled as Now and Past) we can compare them We need to find an isomorphism ΦP →N (very hard problem in the general case [L´ aszl´

• Babai,
• Proc. STOC ’16 Proceedings 48th ACM symp.

Theory of Computing Pag. 684-697]) between the labels: we naively check which center of mass in the Past configuration is closest to another one in the Now configuration

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 18 / 26

Extraction of the Topological Change

Once we have two subsequent Delaunay triangulations (labelled as Now and Past) we can compare them We need to find an isomorphism ΦP →N (very hard problem in the general case [L´ aszl´

• Babai,
• Proc. STOC ’16 Proceedings 48th ACM symp.

Theory of Computing Pag. 684-697]) between the labels: we naively check which center of mass in the Past configuration is closest to another one in the Now configuration We translate the Past adjacency matrix using ΦP →N so that after sorting we can compare them

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 18 / 26

Extraction of the Topological Change

Once we have two subsequent Delaunay triangulations (labelled as Now and Past) we can compare them We need to find an isomorphism ΦP →N (very hard problem in the general case [L´ aszl´

• Babai,
• Proc. STOC ’16 Proceedings 48th ACM symp.

Theory of Computing Pag. 684-697]) between the labels: we naively check which center of mass in the Past configuration is closest to another one in the Now configuration We translate the Past adjacency matrix using ΦP →N so that after sorting we can compare them If they differ we extract the plastic event details

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 18 / 26

Extraction of the Topological Change

Once we have two subsequent Delaunay triangulations (labelled as Now and Past) we can compare them We need to find an isomorphism ΦP →N (very hard problem in the general case [L´ aszl´

• Babai,
• Proc. STOC ’16 Proceedings 48th ACM symp.

Theory of Computing Pag. 684-697]) between the labels: we naively check which center of mass in the Past configuration is closest to another one in the Now configuration We translate the Past adjacency matrix using ΦP →N so that after sorting we can compare them If they differ we extract the plastic event details

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 18 / 26

Extraction of the Topological Change

b c d DP DN ΦP →N(a) a ΦP →N(b) ΦP →N(c) ΦP →N(d) a b c ΦP →N(a) ΦP →N(b) ΦP →N(c) Bulk Boundary DP DN

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 18 / 26

Outline for section 4

1

Motivation and Definitions

2

Physical system and Simulations

3

Algorithm Description

4

Physics Results

5

Conclusions and Perspectives

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 19 / 26

Experiments and numerical comparison

In collaboration with Prof. Matteo Pierno and Dr. Laci Derszi from Padua University Poiseuille flow in rough microchannels2

0.5 1 1.5 2 2.5 3 3.5 0.2 0.4 0.6 0.8 1 1.2 v z 2/2 Comparison ∆P = 1bar H = 9, ∆P = 0.6 · 10−7lbu ∆P = 2bar H = 9, ∆P = 1.2 · 10−7lbu ∆P = 3bar H = 9, ∆P = 1.6 · 10−7lbu

2Ladislav Derzsi, Massimo Bernaschi, ML, Mauro Sbragaglia et al. Phys. Rev. E 95, 052602 (2017)

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 20 / 26

Experiments and numerical comparison

In collaboration with Prof. Matteo Pierno and Dr. Laci Derszi from Padua University Poiseuille flow in rough microchannels2

0.5 1 1.5 2 2.5 3 3.5 0.2 0.4 0.6 0.8 1 1.2 v z 2/2 Comparison ∆P = 1bar H = 9, ∆P = 0.6 · 10−7lbu ∆P = 2bar H = 9, ∆P = 1.2 · 10−7lbu ∆P = 3bar H = 9, ∆P = 1.6 · 10−7lbu

Roughness induced fluidization measure ∆vslip = vflat − vrough Particularly for Couette flows, constant stress σ ∆vslip =

H/2

−H/2

dz ˙ γ(z) ∝

H/2

−H/2

dz f(z) ∝

H/2

−H/2

dz Γ(z)

2Ladislav Derzsi, Massimo Bernaschi, ML, Mauro Sbragaglia et al. Phys. Rev. E 95, 052602 (2017)

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 20 / 26

Fluidization scaling: Experiments and Numerics

Experimentally and numerically one finds the scaling relation3 ∆vslip ∼ λ−1

3Ladislav Derzsi, Massimo Bernaschi, ML, Mauro Sbragaglia et al. Phys. Rev. E 95, 052602 (2017)

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 21 / 26

Fluidization scaling: Experiments and Numerics

Experimentally and numerically one finds the scaling relation3 ∆vslip ∼ λ−1 Poiseuille Flow (Experiments)

10 20 40 60 80 100 0.1 1 10 1 10 0.1 1

1 2 3 1 2 3 1 1

, , g = 15 µm , , g = 20 µm , , g = 37.5 µm

(A) ∆vslip (mm s

• 1)

λ (µm)

1 1 , , g/d = 3.15 , , g/d = 4.21 , , g/d = 7.89

∆vslip / vplug λ/d (B)

1 bar 2 bar 3 bar 1 bar 2 bar 3 bar vplug (mm s-1) ∆p (bar)

3Ladislav Derzsi, Massimo Bernaschi, ML, Mauro Sbragaglia et al. Phys. Rev. E 95, 052602 (2017)

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 21 / 26

Fluidization scaling: Experiments and Numerics

Experimentally and numerically one finds the scaling relation3 ∆vslip ∼ λ−1 Poiseuille Flow (Experiments)

10 20 40 60 80 100 0.1 1 10 1 10 0.1 1

1 2 3 1 2 3 1 1

, , g = 15 µm , , g = 20 µm , , g = 37.5 µm

(A) ∆vslip (mm s

• 1)

λ (µm)

1 1 , , g/d = 3.15 , , g/d = 4.21 , , g/d = 7.89

∆vslip / vplug λ/d (B)

1 bar 2 bar 3 bar 1 bar 2 bar 3 bar vplug (mm s-1) ∆p (bar)

Couette Flow (Numerics)

1 10 1 10

2 4 6 8 10 16 0.1 1 10

∆vslip (x 10

• 4 lbu)

g/d=1 g/d=1.5 g/d=2 g/d=4 g/d=6

ΓΤΟΤ (x 10

• 5 lbu)

1 1 ∆vslip(x 10-4 lbu) λ/d 1 1

3Ladislav Derzsi, Massimo Bernaschi, ML, Mauro Sbragaglia et al. Phys. Rev. E 95, 052602 (2017)

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 21 / 26

Two Scenarios: trapped vs. get-in-get-out

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 22 / 26

Two Scenarios: trapped vs. get-in-get-out

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 23 / 26

Outline for section 5

1

Motivation and Definitions

2

Physical system and Simulations

3

Algorithm Description

4

Physics Results

5

Conclusions and Perspectives

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 24 / 26

Conclusions and Perspectives

Overview & Results

We built a unified implementation for detecting topological changes in dynamic 2D Delaunay triangulations The method is general and can be applied in different contexts (emulsions, pedestrian tracking1...) Complex boundaries are automatically handled Negligible overhead for Lattice Boltzmann Simulations We obtained the highest statistics so far for emulsion plastic events in complex confined channels The code is composed by replaceable components Our code can be downloaded at http://twin.iac.rm.cnr.it/dynvorcuda.tgz

• 1A. Corbetta PhD thesis, “Multiscale crowd dynamics: physical analysis, modeling and applications”
• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 25 / 26

Conclusions and Perspectives

Overview & Results

We built a unified implementation for detecting topological changes in dynamic 2D Delaunay triangulations The method is general and can be applied in different contexts (emulsions, pedestrian tracking1...) Complex boundaries are automatically handled Negligible overhead for Lattice Boltzmann Simulations We obtained the highest statistics so far for emulsion plastic events in complex confined channels The code is composed by replaceable components Our code can be downloaded at http://twin.iac.rm.cnr.it/dynvorcuda.tgz

Current Limitations

Lattice based Single-GPU implementation only Full adjacency matrix is needed: ∼ 104 points can be analyzed

• 1A. Corbetta PhD thesis, “Multiscale crowd dynamics: physical analysis, modeling and applications”
• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 25 / 26

Conclusions and Perspectives

Overview & Results

We built a unified implementation for detecting topological changes in dynamic 2D Delaunay triangulations The method is general and can be applied in different contexts (emulsions, pedestrian tracking1...) Complex boundaries are automatically handled Negligible overhead for Lattice Boltzmann Simulations We obtained the highest statistics so far for emulsion plastic events in complex confined channels The code is composed by replaceable components Our code can be downloaded at http://twin.iac.rm.cnr.it/dynvorcuda.tgz

Current Limitations

Lattice based Single-GPU implementation only Full adjacency matrix is needed: ∼ 104 points can be analyzed

Perspectives

Change data structures in order to handle more points Multi-GPU implementation Application in other domains Three-dimensional case

• 1A. Corbetta PhD thesis, “Multiscale crowd dynamics: physical analysis, modeling and applications”
• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 25 / 26

Massimo Bernaschi, ML, Mauro Sbragaglia, “GPU based detection of topological changes in Voronoi diagrams”

• Comput. Phys. Commun. 213, 19 (2017)

Ladislav Derzsi, Massimo Bernaschi, ML, Mauro Sbragaglia et al. Phys. Rev. E 95, 052602 (2017) Videos on Youtube channel: GPUFoam

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 26 / 26

Massimo Bernaschi, ML, Mauro Sbragaglia, “GPU based detection of topological changes in Voronoi diagrams”

• Comput. Phys. Commun. 213, 19 (2017)

Ladislav Derzsi, Massimo Bernaschi, ML, Mauro Sbragaglia et al. Phys. Rev. E 95, 052602 (2017) Videos on Youtube channel: GPUFoam

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 26 / 26

Adjacency matrix comparison

MAXLINK p1n1 MAXLINK MAXLINK MAXLINK MAXLINK MAXLINK p2n2 pNnN p21 p11 pN1 p1 p2 pN MAXLINK MAXLINK MAXLINK MAXLINK MAXLINK MAXLINK ΦP →N(p11) p1 p2 pN ΦP →N(p21) ΦP →N(pN1) ΦP →N(p2n2) ΦP →N(pNnN ) ΦP →N(p1n1) translation via isomorphism (EP)P P (EP)P N Past Adjacency matrix MAXLINK ΦP →N(p11) p1 p2 pN ΦP →N(p21) ΦP →N(pN1) ΦP →N(p2k2) sorting and compressing ΦP →N(p1k1) MAXLINK ΦP →N(pNkN ) MAXLINK MAXLINK n2n2 n21 n11 nN1 n1 n2 nN Now Compressed Adjacency matrix n1n1 nNnN MAXLINK MAXLINK

Let us assume that ΦP →N(p11), . . . , ΦP →N(pN1) are still the smallest labels For example, if n11 = ΦP →N(p11) then a topology change occured

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 27 / 26

Euclidean Colouring Pathologies on lattice

Let us consider the worst possible case: all centers of mass along a circumference very close to each other Configuration very close to degeneracy, and grid colouring pathologies

Isolated point: 2 gaps 3 gaps 4 gaps: prolongued vertex 3 gaps 4 gaps: prolongued vertex k spx spy smypx k spx spy spxpy spxpy ak ak a) b) smy spymx smypx

• Comput. Phys. Commun. 213, 19 (2017)

Detecting Topological Changes in Dynamic Delaunay Triangulations Using CUDA 28 / 26