Blog : Mettre une sphère dans un tube
Cet article fait suite à la série de simulations sur l’allée de tourbillons de
Von Kármán. Cette fois-ci, nous étendons le domaine en 3D et générons le
maillage de calcul à l’aide de snappyHexMesh. La géométrie étant relativement
simple, nous utilisons les surfaces recherchables (searchable surfaces)
intégrées à OpenFOAM au lieu de fichiers STL externes.
Géométrie et maillage
Cette simulation est une extension en 3D de la classique allée de tourbillons de Von Kármán. Au lieu d’une géométrie plane en 2D, nous utilisons un tube cylindrique d’un diamètre de 4D et y plaçons une sphère de diamètre D. Nous appliquons une vitesse uniforme en entrée et observons le développement tridimensionnel de l’écoulement dans le sillage.
Le maillage est généré avec snappyHexMesh. Puisque notre domaine est constitué
de formes géométriques de base, nous utilisons des objets géométriques intégrés
(searchableCylinder et searchableSphere) plutôt que d’importer des fichiers
STL externes :
geometry
{
tube
{
type searchableCylinder;
point1 (0 0 -100);
point2 (0 0 100);
radius #eval "2*$D";
}
sphere
{
type searchableSphere;
origin (0 0 0);
radius #eval "0.5*$D";
}
}
Pour l’étape de maillage créneaux (castellated mesh), nous ne raffinons que les limites des surfaces. Les formes géométriques étant lisses et continues, l’extraction des arêtes n’est pas nécessaire (même l’entrée et la sortie sont définies par la boîte englobante du maillage de fond) :
refinementSurfaces
{
sphere
{
level (1 1);
}
tube
{
level (1 1);
}
}
Pour résoudre les gradients de la couche limite le long des surfaces lisses, nous ajoutons des couches d’inflation prismatiques :
layers
{
sphere
{
nSurfaceLayers 3;
}
tube
{
nSurfaceLayers 3;
}
}
Le reste du fichier snappyHexMeshDict utilise des paramètres standards.
L’exécution de snappyHexMesh produit un maillage propre avec une
non-orthogonalité maximale d’environ 50° :
Note : Pour éviter que ParaView n’affiche de manière erronée les cellules limites coupées, l’option de découpe Crinkle clip est activée.
Checking geometry...
Overall domain bounding box (-2 -2 -4) (2 2 12)
Mesh has 3 geometric (non-empty/wedge) directions (1 1 1)
Mesh has 3 solution (non-empty) directions (1 1 1)
Boundary openness (-3.15769e-16 2.9073e-16 7.81195e-17) OK.
Max cell openness = 4.07388e-16 OK.
Max aspect ratio = 4.04672 OK.
Minimum face area = 0.000805246. Maximum face area = 0.013633. Face area magnitudes OK.
Min volume = 4.0276e-05. Max volume = 0.00125246. Total volume = 200.51. Cell volumes OK.
Mesh non-orthogonality Max: 50.0132 average: 7.10234
Non-orthogonality check OK.
Face pyramids OK.
Max skewness = 0.377543 OK.
Coupled point location match (average 0) OK.
Conditions d’écoulement et propriétés du fluide
Comme pour le cas 2D, nous spécifions une vitesse constante en entrée (Uin) et
faisons varier la viscosité cinématique (nu) pour atteindre le nombre de
Reynolds ciblé (Re). Ceci est configuré de manière dynamique via l’entrée
#eval dans constant/transportProperties :
#include "<case>/settings"
transportModel Newtonian;
nu #eval "$Uin*$D/$Re";
Les variables Uin, D et Re sont définies dans le dictionnaire settings à
l’intérieur du dossier du cas.
Visualisation de l’écoulement
Dans les simulations 2D, tracer les coefficients de force ou les cartes de contours de vitesse sur le plan médian fournit une image claire du champ d’écoulement. Cependant, en 3D, les coupes dans le plan médian ne parviennent souvent pas à capturer le détachement tourbillonnaire complexe autour de la sphère.
Une approche de visualisation bien supérieure consiste à tracer les lignes de courant en 3D.
Génération des lignes de courant
Approche sous ParaView
ParaView propose deux méthodes principales pour l’ensemencement (seeding) des lignes de courant :
- Stream Tracer : Place les points de départ en utilisant des formes géométriques de base (ligne ou sphère).
- Stream Tracer With Custom Source : Utilise une géométrie source personnalisée pour définir l’emplacement des points de départ.
Pour étudier les structures d’écoulement prenant naissance juste devant l’obstacle, il est utile d’ensemencer les lignes de courant à l’intérieur d’une zone circulaire sur le plan d’entrée correspondant à la projection de la sphère. Nous pouvons créer cette source d’ensemencement personnalisée dans ParaView grâce au pipeline suivant :
- Extract Block: Isole le patch de la limite d’entrée.
- Delaunay 2D: Triangule la géométrie du patch pour la préparer au raffinement.
- Loop Subdivision: Raffine le maillage de surface triangulé.
- Clip (Cylinder): Découpe le patch raffiné pour créer un disque d’ensemencement circulaire du rayon souhaité.
En réduisant le rayon de cette source personnalisée, nous pouvons visualiser clairement les motifs tourbillonnants complexes qui se forment dans le sillage de la sphère :
Enfin, l’application du filtre Stream Tracer permet de produire des animations fluides de l’écoulement :
Objet de fonction OpenFOAM streamLine
Les lignes de courant peuvent également être calculées directement pendant
l’exécution ou en post-traitement à l’aide de l’objet de fonction streamLine
d’OpenFOAM. Cet objet de fonction génère des particules lagrangiennes, suit
leurs trajectoires à travers le champ de vitesse et enregistre les parcours
résultants au format VTK.
Cependant, les options d’ensemencement intégrées sont assez limitées. L’option
native la plus proche est patchSeed, qui distribue un nombre fixe de points sur
l’ensemble d’un patch limite. Si on l’applique directement à notre patch
d’entrée, patchSeed génère de nombreuses lignes de courant rectilignes sans
grand intérêt près des parois extérieures du tube.
Pour restreindre l’ensemencement à une région spécifique, il nous faut une
méthode du type faceZoneSeed. En combinant cela avec topoSet, nous pouvons
définir une petite sous-zone sur l’entrée.
Implémentation de faceZoneSeedSet
L’examen de la classe patchSeedSet d’OpenFOAM montre qu’elle s’appuie sur
triangulatedPatch pour distribuer les points d’ensemencement sur les faces d’un
patch. Puisqu’une zone de faces (faceZone) peut être évaluée de manière
similaire, nous pouvons facilement implémenter une classe personnalisée
faceZoneSeedSet.
Nous adaptons patchSeedSet en remplaçant patchSet_ par faceZoneName_ et en
supprimant selectedLocations_ :
class faceZoneSeedSet
:
public sampledSet
{
// Private Data
//- Face zone name to sample
const word faceZoneName_;
//- Maximum number of zone faces to seed
const label maxPoints_;
//- Random number generator
autoPtr<Random> rndGenPtr_;
// Private Member Functions
//- Calculate distribution points on the face zone
void calcZoneSamples
(
const label nAvailable,
const label nZonePoints,
DynamicList<point>& samplingPts,
DynamicList<label>& samplingFaces,
DynamicList<label>& samplingSegments,
DynamicList<scalar>& samplingCurveDist
);
public:
//- Runtime type information
TypeName("faceZoneSeed");
// Constructors
faceZoneSeedSet
(
const word& name,
const polyMesh& mesh,
const meshSearch& searchEngine,
const dictionary& dict
);
};
Le positionnement des points d’ensemencement s’effectue dans calcZoneSamples,
qui simule calcPatchSamples de la classe patchSeedSet :
void Foam::faceZoneSeedSet::calcZoneSamples
(
const label nAvailable,
const label nZonePoints,
DynamicList<point>& samplingPts,
DynamicList<label>& samplingFaces,
DynamicList<label>& samplingSegments,
DynamicList<scalar>& samplingCurveDist
)
{
if (nAvailable < 1) return;
Random& rndGen = *rndGenPtr_;
globalIndex globalSampleNumbers(nAvailable);
label nGlobalPatchPoints = returnReduce(nZonePoints, sumOp<label>());
point pt;
label facei;
const faceZone& fz = mesh().faceZones()[faceZoneName_];
triangulatedFaceZone tfz(fz);
const label np = nAvailable * fz.size() / scalar(nGlobalPatchPoints);
for (label i = 0; i < np; ++i)
{
tfz.randomLocalPoint(rndGen, pt, facei);
samplingPts.append(pt);
samplingFaces.append(fz[facei]);
samplingSegments.append(0);
samplingCurveDist.append(globalSampleNumbers.toGlobal(i));
}
}
L’algorithme appelle randomLocalPoint via triangulatedFaceZone. Puisque les
particules de lignes de courant résident dans le volume des cellules, nous
associons les faces aux cellules intérieures à l’aide de faceOwner():
labelList samplingCells(samplingFaces.size(), -1);
forAll(samplingFaces, i)
{
samplingCells[i] = mesh().faceOwner()[samplingFaces[i]];
}
Le code source complet des classes est disponible dans le dépôt d’exemples.
Configuration et Résultats
Nous configurons topoSetDict pour définir la faceZone cible (smallInlet) :
actions
(
{
name f0;
type faceSet;
action new;
source patchToFace;
sourceInfo
{
patch inlet;
}
}
{
name f0;
type faceSet;
action subset;
source cylinderToFace;
sourceInfo
{
point1 (0 0 -100);
point2 (0 0 100);
radius 0.2;
}
}
{
name smallInlet;
type faceZoneSet;
action new;
source setToFaceZone;
sourceInfo
{
faceSet f0;
}
}
);
Et la configuration suivante pour l’objet de fonction streamLine :
type streamLine;
libs (fieldFunctionObjects libfaceZoneSeed);
setFormat vtk;
direction forward;
fields (p U);
lifeTime 10000;
nSubCycle 5;
cloud particleTracks;
seedSampleSet
{
type faceZoneSeed;
axis xyz;
zone smallInlet;
maxPoints 200;
}
L’exécution de l’objet de fonction génère les trajectoires de lignes de courant suivantes :
Note : Le maillage de fond au niveau de l’entrée étant relativement grossier, la
limite de la zone d’ensemencement (faceZone) apparaît crantée/carrée.
Dans ParaView, le fait de raffiner et de trianguler le patch extrait avant la
découpe produit un profil d’ensemencement circulaire nettement plus lisse.
Les fichiers du cas de simulation sont disponibles dans le dépôt d’exemples.