# Rapport de Anim3D (20/11/2020)

# Autheurs

- Audran Doublet
- Sami Issaadi

# Introduction

Le projet présenté dans ce rapport est un générateur d'animations 3D. Il génère une vidéo à partir d'informations sur une scène (les différents objets présents, les mouvements de caméra, ...).

Les différentes techniques de simulation implémentées sont:

- Smoothed Particle Hydrodynamics (SPH) pour la simulation de fluide
- Simulation de tissu par des ressorts
- Simulation de mouvements d'objets rigides

Les différents de la simulation sont alors transformés en un ensemble de mesh.
Le rendu final est ensuite affiché via un algorithme de lancer de rayons (ray tracing).

Ce projet a été entièrement développé en Rust. Nous utilisons un certain nombre de bibliothèques open source, principalement:

- `Kiss3D` pour un affichage des particules (debug)
- `nalgebra` pour les opérations d'algèbre linéaire
- `serde` pour la lecture et écriture de nos configurations et de nos scènes
- `image` pour la lecture d'images
- `rayon` pour le multithreading

Pour un aperçu visuel des différentes possibilités, nous vous invitons à visionner cette vidéo:

**Ancienne demo**:
[](https://www.youtube.com/watch?v=4S8efTqTbFo)

**Nouvelle demo** (pour cette matière):

[](https://www.youtube.com/watch?v=6EsMR0Xt3eE)

**Code source**: https://github.com/AudranDoublet/opr

### Disclaimer

Nous tenons à préciser que nous n'avons pas implémenté tout ce programme dans le cadre du cours d'Animation 3D. Nous avions déjà implémenté:

- la simulation de fluide avec SPH
- les interactions entre les fluides et des objets rigides
- Transformation des particules de SPH en maillage
- Algorithme de ray tracing pour le rendu de SPH et des objets solides
- Systèmes d'animations basiques pour la caméra, certains solides, ...

Nous avons néanmoins implémenté cela il y a quelques mois (uniquement nous deux).
Ce que nous avons implémenté dans le cadre de ce cours est la simulation de tissus.

Aussi, nous ne sommes ni mathématiciens, ni physiciens, et nous sommes entièrement basés sur des papiers scientifiques sur le sujet, mis à part peut-être pour certaines optimisations algorithmiques.

# Simulation de tissu

Comme dit plus haut, une partie non négligeable de ce projet a été réalisé dans un cadre antérieur à ce cours. Notre objectif a été d'intégrer un des thèmes vu en cours (la simulation de tissus) au projet. Cela signifie, pour nous, de développer la méthode et de l'intégrer aux autres éléments simulés (objets solides & fluides).

## Première méthode: mass spring simulation

Nous avons implémenté la méthode vue en cours qui consiste à représenter le tissu comme un ensemble de particules liées par des ressorts.

Pour l’intégration des vitesses et positions, nous utilisions originellement la méthode semi-implicit Euler, puis nous avons implémenté l’intégration RK4 (Runge-Kutta 4ème ordre) pour tenter d’améliorer la stabilité mais n’avons pas observé de grands changements.

Nous avons tenté de gérer les collisions internes en considérant chaque particule comme des sphères pouvant entrer en collision entre elles. Pour accélérer le calcul des collisions, nous avons implémenté un découpage de l’espace en grille (hash-grid). Malheureusement cette méthode naïve ne permettait pas de correctement gérer les collisions internes du tissu.

Nous avons changé de stratégie pour deux raisons:

* la méthode nouus a paru peu stable, on pas suffisement
* l'intégration avec SPH, bien que probablement faisable, n'est ni très pertinente ni bien documentée dans la littérature scientifique

## Deuxième méthode: utilisation de SPH

Après le peu de succès obtenu avec la première méthode, nous avons décidé de partir de l'algorithme SPH. Cela a plusieurs intérêts:

- l'idée de base est la même: on discrétise le tissu en un nombre limité de points
- des self-collisions stables étaient faciles à implémenter grâce aux force de pressions (voir DFSPH plus bas)
- les collisions avec les solides et fluides sont elles aussi gérées par DFSPH

Nous avions déjà implémenté une force permettant de simuler des objets déformables (basés sur [4]), cependant notre implémentation de l'algorithme contenait des erreurs générant de trop grosses instabilités pour simuler correctement un tissu.

Il était par exemple possible de simuler un tatou élastique mais pas un tissu, ou du moins, pas sans tricher sur certains hyperparamètres.

Nous avons donc corrigé notre implémentation du papier.

Pour faire simple, l’idée est de discrétiser un objet en ensemble de particules et d’enregistrer leur configuration initiale: voisinage, distances aux particules voisines. On peut alors calculer les forces d’élasticité en comparant l’état présent d’une particule à son état initial.

> Ceci est un bref résumé, la méthode décrite par le papier est plus compliquée afin de permettre les rotations correctes de l'objet solide élastique (en réponse aux frictions de l'objet avec une surface par exemple).

Ensuite, la pipeline déjà en place permet de déterminer les forces de pressions en chacune des particules et donc de gérer les collisions internes ainsi que les collisions avec les fluides et les solides.

## Animation: changement d'état de la matière

Nous avons implémenté un nouveau système d'animation, permettant de changer l'état des fluides (et objets déformables) au cours de la simulation.

Ce système nous a permis de générer les scènes où des objets déformables se transforment en fluides, puis deviennent de nouveau unis. Cependant, le design de scène est un peu compliqué et la simulation peut demander du temps bien que tout le processus soit parallélisé (cpu), nous n'avons donc pu générer de scène illustrant parfaitement la reconstruction (sans explosion).

## Meshing des tissus

La méthode de meshing que nous avions pour les fluides n'est absolument pas adapté à des tissus (voir meshing des fluides). En effet, cela génère un tissu très épais. De plus, si la densité des particules est faible à un endroit du tissu, celui-ci aurait été troué: autrement dit, un tissu étiré serait un tissu troué.

Nous avons donc implémenté deux nouvelles méthodes de meshing:

- une très simple, où les particules du simulateur sont directement rendu via raytracing
- une plus complexe, spécifique au rendu de tissu. Elle part du principe que le tissu est un rectangle d'une seule couche de particules (comme lorsqu'il est simulé avec des ressorts). Nous créons donc deux triangles entre chaque groupe de 4 particules. Enfin, nous calculons les normales aux points des triangles comme étant la moyenne des normales des triangles dont ce point est un sommet

# Smoothed Particle Hydrodynamics

SPH (Smoothed Particule Hydrodynamics) est une méthode pour simuler la mécanique de milieux continue. Pour ce faire, le milieu est décomposé en un nombre donné de particule (l’augmentation de ce nombre rendant la simulation plus réaliste), on peut alors estimer les propriétés du milieu en chacune des particules à l’aide d’un noyau de convolution (smoothing kernel) appliqué en son voisinage.


En utilisant les propriétés ainsi estimées (e.g. densité, volume), on peut approcher itérativement une solution aux équations de Navier-Stokes qui décrivent les mouvements d'un fluide.

On calculera alors les différentes forces s'appliquant sur chacune de ces particules ainsi que leurs vitesses, et donc, leurs positions.

Les principales forces calculées sont:

- les forces de pressions entre les particules ou avec d'autres types d'objets
- des forces reliées aux propriétés physiques des fluides (tension de surface, viscosité, vorticité, ...)
- la gravité (et oui!)

L'algorithme ce déroule comme suis:

1. on calcule, en fonction de la vitesse maximale, la `time step` entre la step précédente et celle-ci. Cela permet à l'algorithme d'être stable avec des vitesses relativement élevées tout en étant relativement rapide avec des vitesses lentes
2. on calcule la nouvelle vitesse après application des forces de pressions
3. on calcule les autres forces (forces externes)
4. on applique ces forces, et obtient donc une nouvelle vitesse et une nouvelle position
5. on recommence à la première étape

## Calcul de pression: Divergence Free SPH

Il existe plusieurs variation du calcul des forces de pressions pour la simulation de fluides incompressibles, mais celle que nous avons retenu aura été Divergence-Free SPH [2].

La condition d'incompressibilité d'un fluide (pas de variation de densité) implique qu'il n'y a [pas de divergence de vélocité](https://fr.wikipedia.org/wiki/Champ_sol%C3%A9no%C3%AFdal) dans le fluide quand la masse volumique est constante. DFSPH se base sur cette contrainte pour mieux forcer la condition d'incompressibilité du fluide.

Nous avons choisi cette méthode car elle est plutôt récente (2015), stable et permet d'utiliser de plus grands pas de temps que d'autres méthodes plus anciennes, rendant ainsi les simulations plus rapides.

## Autres forces liées au fluide

Dans un souci de synthèse, nous ne rentrerons pas dans le détail des différentes forces et de leur implémentation. Néanmoins, les voici:

- une force de viscosité stable uniquement pour les fluides peu visqueux (comme l'eau), mais rapide à calculer
- une force de viscosité stable pour les fluides hautement visqueux (basé sur [3])
- une force d'élasticité, permettant de simuler des objets déformables (basé sur [4])
- une force de traînée, pour simuler les interactions avec des particules d'air sans simuler des millions de particules; utilisé principalement pour simuler un tatou qui danse (basé sur [5])
- une force de tension (basé sur [6])
- une force vorticité (basé sur [7])
- des interactions entre différents types de fluides ayant des propriétés différentes

## Animations

Nous avons également mis en place des fonctionnalités d'animation
Au-delà de simuler en soit des évènements physiques, notre objectif était de créer des scènes. Pour cela, nous avons mis en place un système permettant de configurer des scènes complexes et de les animer.

Ces différents objets animés ont des propriétés (position, vitesse, accélération, rotation, vitesse angulaire, accélération angulaire). Ces propriétés peuvent être modifiées de différentes manières:

- la valeur de ces variables peut être fixée pour un temps donnée
- la valeur de ces variables peut suivre une progression linéaire ou en courbe de Béziers
- pour les mouvements, un effet de "fondu" peut être ajouté, afin de rendre l'effet plus lisse
- la rotation peut être définie afin que l'objet "regarde" dans une direction
- ces différentes modifications peuvent se suivre dans le temps et être exécutées en boucle

Ce système générique nous a permis d'implémenter un certain nombre d'objets animés, notamment:

- une caméra animée
- des objets solides animés, qui permettent par exemple de faire un moteur
- un vent dont la force change au cours du temps (toujours pour un certain tatou)
- des émetteurs de fluides animés

Les particules de fluide peuvent aussi avoir des propriétés précises:

- possibilité d'ajouter du fluide dans une forme prédéfinie (soit en cube, soit suivant une mesh)
- possibilité d'avoir des particules fixées (encore et toujours pour un tatou dansant)

## Maillage du fluide

Bien que l'affichage d'un fluide soit possible simplement en affichant les particules en tant que sphère avec un certain rayon, cela ne permettrait pas un rendu réalistique.

Pour obtenir un rendu réalistique du fluide, il existe plusieurs méthodes. Dans notre cas, nous convertissons le fluide en maillage polygonal (mesh) à partir des différentes informations qu'on a (e.g. densité du fluide à une position donnée). *Cela n’a malheureusement pas permis d’obtenir le résultat escompté.*

### Marching Cube

Pour générer un maillage depuis particules de fluide obtenus à partir de la simulation, nous utilisons l'algorithme de marching cube.
Il s'agit d'abord d'estimer la surface du fluide en parcourant son domaine à l'aide d'un cube et en estimant la densité du fluide en chacun des sommets de ce cube.
Si un ou plusieurs sommets ont une densité trop faible, on considère alors qu'ils n'appartiennent pas au fluide et on génère une frontière composée de triangles.
On obtient alors ainsi un maillage pour notre fluide.


### Interpolation linéaire

L'algorithme ne nous donne pas les positions des sommets des triangles mais uniquement les segments auxquels ils appartiennent. Pour améliorer le maillage, on interpole linéairement la position des sommets des triangles dans les segments auxquels ils appartiennent.

### Noyeau anisotropique

L'estimation de la densité au niveau des sommets du cube est faite avec un noyau de convolution comme pour la simulation, cependant, ce noyau est originellement isotropique: il a les mêmes propriétés en tous points de l'espace.


Nous remplaçons ce noyau par un noyau anisotropique [10] qui tient compte de la déformation du fluide et permet alors une meilleure estimation du mesh.

## Interactions avec SPH

# Rigid Body Physics

Nous avons implémenté le mouvement des objets physiques en utilisant les équations explicites de la mécanique newtonienne. Malheureusement, celles-ci n’ont pas toujours été suffisantes lorsque nous développions et que nos objets ont atteint des vitesses relativistes !

## Rotation de l'objet

Pour simuler la rotation d'un objet non déformable, il faut connaître des propriétés plus complexes que la masse:

1. **le tenseur d'inertie** de l'objet, qui, de notre compréhension, mesure la résistance qu'oppose le corps à toute accélération angulaire.
2. **le centre de masse** de l'objet, qui est donc le point autour duquel cet objet tournera.

Nous récupérons ces informations à partir de la description de l'objet en mesh et de la densité de l'objet. Cela permet aussi d'obtenir la masse. Pour calculer ces propriétés, nous avons implémenté l'algorithme proposé par GeometricTools dans [1].

Enfin, une dernière information nécessaire est le point auquel une force s’exerce. Cela permet de calculer le moment de cette force, et donc l’impact sur l’accélération l’objet, et l’impact sur son accélération angulaire.

## Récupération rapide des propriétés du corps en un points de l'espace

Afin d’interagir au mieux avec les autres composants du simulateur, nous avons besoin d’accéder à plusieurs propriétés physiques de l’objet:

- le vecteur entre un point de l'espace et le point le plus proche de l'objet
- savoir si un point de l'espace est à l'intérieur ou à l'extérieur de l'objet
- connaître le volume en un point de l'espace de l'objet (voir interaction avec SPH)

Ces informations sont très coûteuses à calculer pour chaque particule et à chaque étape de simulation.

Pour cela nous avons implémenté [8]. Ce papier propose de récupérer ces informations avant la simulation pour une grille de points autour de l'objet, puis de faire une interpolation pour le point qui nous intéresse.

### Champs de distance signée (Signed Distance Field)

Ainsi, nous avons commencé par implémenter cette grille afin qu'elle contienne un champ de distance signée à l'objet.

Une distance signée est une distance qui est positive si le point est à l'extérieur de l'objet, et négative si le point est dans l'objet.

#### Récupération de la distance à l'objet

Pour récupérer la distance à un objet, nous parcours tous les triangles de celui-ci et nous récupérons les distances entre ces triangles et notre point. La distance minimale obtenue sera la distance gardée. Pour calculer la distance minimale entre un point et un triangle (et donc obtenir la coordonnée barycentrique du triangle la plus proche), nous avons implémenté l'algorithme proposé par GeometricTools dans [9].

En réalité, nous ne parcourons pas tous les triangles de l'objet, puisque nous utilisons une structure de recherche optimisée: des Bounding Sphere Hierarchy. L'idée est de construire un arbre ou chaque niveau sépare les triangles en deux groupes, centrés autour d'un point. On peut ensuite calculer la distance minimale entre cette sphère et notre point (la distance au centre moins le rayon), et ainsi beaucoup optimiser la recherche.

#### Récupération de la distance signée

L'objectif est de vérifier si le point est à l'intérieur ou à l'extérieur de l'objet.

**La première approche** que nous avons utilisée est de calculer le produit scalaire entre la normale du triangle de l'objet le plus proche, et la direction entre ce triangle et notre point. Si ce produit scalaire est positif, alors le point est à l'extérieur.
Cette approche par du postulat que les normales des triangles de l'objet devraient être tournées vers l'extérieur afin que celui-ci soit visible.

Nous avons abandonné cette approche car elle était trop sensible aux erreurs de calcul sur les flottants. N'ayant pas trouvé de source sur cette méthode, nous en avons choisi une autre.

**La deuxième approche** que nous avons utilisée est de faire un lancé de rayon à partir d'un point à l'extérieur de l'objet vers notre point. Si le nombre d'intersections est pair, on sera à l'extérieur de l'objet. Sinon, on sera à l'intérieur.


Cette méthode est plus coûteuse, mais fonctionne mieux. De plus, elle ne nécessite pas que les normales de l'objet soient correctes.

#### Résultats

Une fois la distance signée en tous les points récupérés, on peut construire notre grille et la sauvegarder pour la simulation.

Ainsi, le calcul de l'interpolation donne la distance à l'objet. Le calcul de la dérivée de l'interpolation donne la direction vers le point le plus proche. Cela permet d'obtenir le vecteur que nous recherchions.

Ci-dessous, une reconstruction d'un slice de l'objet à partir du SDF.

### Champs de volume

La deuxième donnée dont nous avons besoin est le volume de l'objet en toute zone de l'espace. Cette méthode est proposée dans [8] afin de gérer les interactions avec les fluides de SPH.

Pour calculer ce volume, nous appliquons la formule analytique du calcul de volume, à savoir une triple intégrale (par x, y et z). Nous calculons donc, pour tous les points négatifs de notre grille, la triple intégrale de la distance, via la méthode de Gauss-Legendre.

Ci-dessous, une reconstruction du volume d'un slice de l'objet à partir du SDF:

# Interactions des objets solides avec les fluides

Pour implémenter l'interaction des objets solides avec les fluides, nous nous sommes basé sur [8].

Celui-ci nécessite de connaître à la position de la particule, la distance signée à l'objet, la direction de l'objet, et le volume de l'objet en ce point. Nous expliquons plus haut comment nous récupérons ces données.

L'implémentation est alors relativement simple, car il suffit de prendre en compte l'objet comme une particule lors du calcul des forces.

Calculer le mouvement de l'objet dû au fluide est aussi extrêmement simple, puisqu'il suffit d'appliquer la troisième loi de Newton: si l'on calcule une force exercée par le solide sur une particule, on applique l'opposée de cette force sur le solide.

# Raytracer

Après les différentes étapes que nous avons décrites ci-dessus, nous avons:

- les informations de déplacement de la caméra
- les différents objets de la scène aux différents temps

Nous avons donc souhaité pouvoir faire le rendu final, et créer une image pour chaque sortie du programme. Pour cela, nous avons choisi un algorithme de rendu photoréaliste, à savoir le lancé de rayon.

S'il nous est évident qu'utiliser un moteur de rendu dédié comme Blender ou Maya aurait permi d'obtenir de bien meilleurs résultats, nous avons tout de même décidé d'implémenter le raytracing, pour deux raisons:

- l'objectif initial de ce projet était un projet de synthèse d'image, et il fallait bien trouver un lien
- pourquoi pas ? (ce projet n'ayant aucun autre objectif de production)

## Configuration

Pour pouvoir rendre complètement les scènes, le raytracer nécessite d'avoir quelques informations supplémentaires:

- la liste des lumières de la scène (lumière d'ambiance & directionnelles)
- des objets non présents dans la simulation venant ajouter du décor: un ciel, un sol, ou autres
- des informations sur comment rendre les différents tissus et fluides; nous avons utilisé pour cela le format de fichier .mtl

## Implémentation

Notre ray tracer fonctionne avec des formes en triangle. Nous avons donc utilisé un BVH (Bounding Volume Hierarchy) pour optimiser les lancer de rayons.

Nous avons implémenté différents effets de lumière:

- Réflection, refraction, spéculaire
- la loi de Fresnel
- la loi de Beer-Lambert

Nous avons aussi implémenté la gestion des textures et les déformations de normales.

# Références

[1] https://www.geometrictools.com/Documentation/PolyhedralMassProperties.pdf
[2] Divergence-Free Smoothed Particle Hydrodynamics; 2016; Bender and al.
[3] A Physically Consistent Implicit Viscosity Solver for SPH Fluids; 2018; Weiler and al
[4] An Implicit SPH Formulation for Incompressible Linearly Elastic Solids; 2018; Peer and al.
[5] Approximate Air-Fluid Interactions for SPH; 2017; Glisser and al.
[6] Versatile Surface Tension and Adhesion for SPH Fluids; 2013; Akinci and al.
[7] Turbulent Micropolar SPH Fluids with Foam; 2017; Bender and al.
[8] Volume Maps: An Implicit Boundary Representation for SPH; 2019; Bender and al.
[9] https://www.geometrictools.com/Documentation/DistancePoint3Triangle3.pdf
[10] Reconstructing surfaces of particle-based fluids using anisotropic kernels; 2013; Yu, J. & Turk, G.