---

Geometrie.JPG

Hallo zusammen!  

Ich bin neu hier im Forum, allerdings ist es mir schon länger bekannt und der ein oder andere Beitrag konnte mir in der Vergangenheit schon weiterhelfen.

Diesmal habe ich allerdings ein Problem, dass ich momentan leider nicht so richtig gebacken bekomme...

Mein Vorhaben

Ich möchte einen Teil eines Luftleitungssystem  mit OpenFoam (simpleFoam) berechnen.

Es handelt sich um Rohrleitungen mit Krümmern und Durchmessern von 250 und 350 mm die in einen großen Rechteckkanal münden. Die Geschwindigkeiten liegen in der Größenordnung von maximal ~ 24 m/s bis auf 2 m/s

Das Netz erstelle ich mittels SnappyHexMesh und einem reinen Hexaeder-Grundgitter aus gmsh. Hier habe ich bereits mehrere Netzvarianten ausprobiert. Grob / fein / mit Layer oder ohne.. CheckMesh ergab immer das mein Netz in Ordnung ist.

Vorgegeben wird die Geschwindigkeit am Einlass und es wird am Ende des Rechtecks frei "ausgeblasen", hier setze ich den Druck auf "0".

Ich verwende das k-Epsilon Turbulenzmodell mit Wandfunktionen, habe allerdings auch schon kOmegaSST ausprobiert.

Ich initialisiere mittels Potentialfoam.

Mein Problem

Ich habe nun allerdings das Problem, dass ich mit meinem Residuum von p kaum unter den Wert 0,005 - 0,001 oder etwas in der Größenordnung komme. Auch die anderen Größen erreichen meistens keine kleineren Werte.

Das sollte doch eigentlich noch weiter sinken, das ist doch eigentlich nicht akzeptabel so?

Die Ergebnisse sehen eigentlich im groben und ganzen plausibel aus. Geschwindigkeiten sind im korrekten Bereich, beim Druck sind auch keine großen Ausreißer dabei.

Ich mir einfach unsicher und unschlüssig ob ich hier wirklich alles korrekt mache und alle Randbedingungen etc. passen  

Ich wäre wirklich sehr dankbar, falls jemand einfach mal über meinen Case drübersieht und mir eventuell noch die ein oder andere Idee mit auf den Weg gibt!     

Beste Grüße!

Hier ein beispielhafter Auszug meines Rechenoutputs. Ich habe auch schon über 1000 Iterationen gerechnet ohne Erfolg.

Code:

---

Time = 150

smoothSolver:  Solving for Ux, Initial residual = 0.00378163, Final residual = 3.17913e-005, No Iterations 19
smoothSolver:  Solving for Uy, Initial residual = 0.00795081, Final residual = 6.85687e-005, No Iterations 23
smoothSolver:  Solving for Uz, Initial residual = 0.00467368, Final residual = 4.40989e-005, No Iterations 18
GAMG:  Solving for p, Initial residual = 0.0132319, Final residual = 8.3468e-009, No Iterations 23
GAMG:  Solving for p, Initial residual = 0.00198833, Final residual = 8.91338e-009, No Iterations 14
GAMG:  Solving for p, Initial residual = 0.000139475, Final residual = 5.24627e-009, No Iterations 10
GAMG:  Solving for p, Initial residual = 1.69987e-005, Final residual = 5.30083e-009, No Iterations 7
time step continuity errors : sum local = 2.43281e-007, global = 6.73135e-009, cumulative = 1.22496e-008
smoothSolver:  Solving for epsilon, Initial residual = 0.0105538, Final residual = 0.000100227, No Iterations 14
smoothSolver:  Solving for k, Initial residual = 0.0144642, Final residual = 0.000144077, No Iterations 18
ExecutionTime = 2162.26 s  ClockTime = 2162 s

---

Der Vollständigkeit halber auch alle anderen Codes meines Case's.

fvSchemes

Code:

---

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  5                                     |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
    default         steadyState;
}

gradSchemes
{
    default         Gauss linear;
    grad(U)         cellLimited Gauss linear 1;
}

divSchemes
{
    div(phi,U)      Gauss linearUpwind grad(U);
    div(phi,k)      Gauss upwind;
        div(phi,epsilon)      Gauss upwind;
    div(phi,omega)  Gauss upwind;
    div((nuEff*dev2(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected;
}

fluxRequired
{
    default         no;
    p               ;
}

// ************************************************************************* //

---

fvSolution

Hier habe ich auch mit den Relaxationen gespielt - Standard-Werte p 0.3 Rest 0.7 wie im Tutorial, das hat allerdings auch alles keinen großen Unterschied gemacht.

Code:

---

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  5                                     |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
    p
    {
solver GAMG;
tolerance 1e-8;
relTol 0;
smoother GaussSeidel;
nPreSweeps 0;
nPostSweeps 2;
cacheAgglomeration on;
agglomerator faceAreaPair;
nCellsInCoarsestLevel 120;
mergeLevels 1;
maxIter 500;
    }

    Phi
    {
        $p;
    }

    U
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-8;
        relTol          0.01;
        nSweeps         1;
    }

    k
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-8;
        relTol          0.01;
        nSweeps         1;
    }

    omega
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-8;
        relTol          0.01;
        nSweeps         1;
    }
   
    epsilon
    {
        solver          smoothSolver;
        smoother        GaussSeidel;
        tolerance       1e-8;
        relTol          0.01;
        nSweeps         1;
    }
}

SIMPLE
{
    nCorrectors 1;
    nNonOrthogonalCorrectors 3;
    consistent               yes;

}

potentialFlow
{
    nNonOrthogonalCorrectors 10;
}

relaxationFactors
{
    equations
    {

        U   0.9;
        k   0.9;
        epsilon 0.9;
      
    }
}

cache
{
    grad(U);
}

// ************************************************************************* //

---

Und meine Randbedingungen:

U

Code:

---

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 1 -1 0 0 0 0];

internalField   uniform (0 0 0);

boundaryField
{

    outlet
{
        type pressureInletOutletVelocity;
        value uniform (0 0 0);
}

    wall
    {
        type            fixedValue;
        value           $internalField;
    }

inlet
{
        type                fixedValue;
        value               uniform (-12 0 0);
}
}

// ************************************************************************* //

---

p

Code:

---

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{

    inlet
    {
        type            zeroGradient;
    }

    outlet
    {
            type            fixedValue;
        value           uniform 0;

    }

    wall
    {
        type            zeroGradient;
    }

}

// ************************************************************************* //

---

k

Code:

---

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0.02;

boundaryField
{
    inlet
    {
        type            fixedValue;
        value           uniform 0.02;
    }

    outlet
    {
        type            zeroGradient;
    }

    wall
    {
        type            kqRWallFunction;
        value           $internalField;
    }

}

// ************************************************************************* //

---

epsilon

Code:

---

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -3 0 0 0 0];

internalField   uniform 0.038;

boundaryField
{

   wall
    {
        type            epsilonWallFunction;
          value $internalField;
    }

    inlet
    {
        type            fixedValue;
        value           $internalField;
    }

    outlet
    {
        type            zeroGradient;
    }

}

// ************************************************************************* //

---

nut

Code:

---

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -1 0 0 0 0];

internalField   uniform 0;

boundaryField
{

    inlet
    {
        type            calculated;
        value           uniform 0;
    }

    outlet
    {
        type            calculated;
        value           uniform 0;
    }

    wall
    {
        type            nutkWallFunction;
        value           uniform 0;
    }

}

// ************************************************************************* //

---

Und zuletzt mein letztes CheckMesh-File

Code:

---

Create time

Create polyMesh for time = 0

Time = 0

Mesh stats
    points:           1768526
    faces:            4478326
    internal faces:   4058490
    cells:            1368541
    faces per cell:   6.2379
    boundary patches: 4
    point zones:      0
    face zones:       0
    cell zones:       1

Overall number of cells of each type:
    hexahedra:     1146266
    prisms:        48425
    wedges:        0
    pyramids:      0
    tet wedges:    93
    tetrahedra:    0
    polyhedra:     173757
    Breakdown of polyhedra by number of faces:
        faces   number of cells
            4   26247
            5   14342
            6   20343
            7   2
            8   3
            9   86303
           10   12
           12   18976
           13   3
           15   7379
           18   147

Checking topology...
    Boundary definition OK.
    Cell to face addressing OK.
    Point usage OK.
    Upper triangular ordering OK.
    Face vertices OK.
    Number of regions: 1 (OK).

Checking basic patch addressing...
                   Patch    Faces   Points
                    wall   350661   378310
                   inlet     4354     4505
                  outlet    51327    51795

Checking geometry...
    Overall domain bounding box (-3.88806 -3.63696 -0.632848) (2 -0.04195 1.917)
    Mesh has 3 geometric (non-empty/wedge) directions (1 1 1)
    Mesh has 3 solution (non-empty) directions (1 1 1)
    Boundary openness (-1.95926e-016 -1.8415e-015 1.707e-016) OK.
    Max cell openness = 3.4469e-016 OK.
    Max aspect ratio = 7.21486 OK.
    Minimum face area = 2.51639e-006. Maximum face area = 0.00302785.  Face area magnitudes OK.
    Min volume = 4.5284e-009. Max volume = 0.000163703.  Total volume = 2.06708.  Cell volumes OK.
    Mesh non-orthogonality Max: 62.0921 average: 10.7761
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 1.80749 OK.
    Coupled point location match (average 0) OK.

Mesh OK.

---

[Diese Nachricht wurde von Absynthe23 am 28. Feb. 2019 editiert.]

Eine Antwort auf diesen Beitrag verfassen (mit Zitat/Zitat des Beitrags) IP

---

Hallo Absynthe (das Parfüm?), ...

ich geh nur kurz auf deine Aussagen ein.

• Schön gegliedert, und code Tags verwendet. Find ich sehr lobenswert
  
• Diese Aussage ist nicht korrekt
  

Zitat:

---

Ich habe nun allerdings das Problem, dass ich mit meinem Residuum von p kaum unter den Wert 0,005 - 0,001 oder etwas in der Größenordnung komme. Auch die anderen Größen erreichen meistens keine kleineren Werte.

Das sollte doch eigentlich noch weiter sinken, das ist doch eigentlich nicht akzeptabel so?

---

Der SIMPLE Algorithmus (du verwendest ja den SIMPLEC) kann transiente Erscheinungen nicht berechnen. Wenn du also alle 5 oder 10 Iterationen rausschreiben lässt und dir einen Schnitt durch deinen Rechteck-Kanal legst, gehe ich stark davon aus, dass du hier irgendwelche transienten Erscheinungen vorfinden wirst. Zumindest irgendwo in deiner Geometrie. Entsprechend können deine Residuuen nicht weiter sinken auch wenn die Lösung eigentlich schon gut approximiert[i] ist.

Allgemein rechnest du ja sehr stabil (1 Ordnung - Gauss Upwind). Beim Outlet würde ich beim Geschwindigkeitsfeld ein inlet/Outlet verwenden und ein Rückströmen verhindern (du weißt ja nicht wie und was eigentlich rückströmen würde). Ferner hast du hier einen fixierten Druck, wodurch diese [i]pressureXY Randbedingung wenig Sinn macht. Siehe Quellcode (sofern ich mich jetzt richtig erinnere).

• Essentiell ist natürlich auch der Fehler, dass du den Druck nicht unterrelaxierst (auch wenn du SIMPLEC) verwendest solltest du dein Druckfeld unterrelaxieren
  
• Stabilisierung kannst du auch damit erreichen, dass du deine Gleichungen nicht mit 0.9 sondern mit 0.7 unterrelaxierst, damit fallen deine Residuuen im Allgemeinen auch
  

Ansonsten, sagen die Residuuen nicht immer das aus, ob eine Lösung konvergiert ist oder nicht. Siehe hierzu die Berechnung der Residuuen in OpenFOAM + Skalierung. Irgendwo auf cfd-online.com hab ich das mal beschrieben.

Zuletzt:

Zitat:

---

Mesh non-orthogonality Max: 62.0921 average: 10.7761

---

Das ist mit Abstand das Wichtigste. Ein Wert von 62 ist eigentlich noch gar nicht so hoch, daher würde ich auch deine 3x NonOrthoCorrection auf 0 setzen und eher die nCorrectors auf 2. Auswirkungen sind in der Druck-Korrekturgleichung vorzufinden.

Ach ja, für eine bessere Approximation der Turbulenzgrößen würde ich dir die turbulentKineticEnergy... für k empfehlen und die [i]mixingLengthDissipationRateInlet...[i] für epsilon und entsprechend die MixingLength bestimmen (=d_hyd * 0.037).

------------------
Viele Grüße,
Tobi

OpenFOAM® Tutorials | Training | Publikationen | Für Anfänger wiki.openfoam.com

Eine Antwort auf diesen Beitrag verfassen (mit Zitat/Zitat des Beitrags) IP

---

Hallo Tobias,

vielen lieben Dank erstmal dafür, dass Du Dir kurz Zeit genommen und meinen Case angesehen hast!   

Zitat:

---

Der SIMPLE Algorithmus (du verwendest ja den SIMPLEC) kann transiente Erscheinungen nicht berechnen. Wenn du also alle 5 oder 10 Iterationen rausschreiben lässt und dir einen Schnitt durch deinen Rechteck-Kanal legst, gehe ich stark davon aus, dass du hier irgendwelche transienten Erscheinungen vorfinden wirst. Zumindest irgendwo in deiner Geometrie. Entsprechend können deine Residuuen nicht weiter sinken auch wenn die Lösung eigentlich schon gut approximiert[i] ist.

---

Okay, da gebe ich dir natürlich Recht! Vermutlich werden nach der plötzlichen Aufweitung des Rohres in den Rechteck-Kanal große Ablösegebiete mit instationären Verwirbelungen/Verhalten auftreten..und diese kann der SIMPLEC Algorithmus nicht lösen. Das würde dann eventuell auch die doch deutlich schwankenden Residuen erklären?

Ich habe deine Ratschläge beachtet und meine Simulation dahingehend angepasst. Es ändert sich kaum etwas, nur ein klein wenig weiter sinken die Residuen. Ich werde das Ergebnis (für verschiedene Iterationsanzahlen) auch noch untereinander vergleichen.

Jetzt hätte ich dann allerdings noch eine Frage. Um die transienten Effekte abbilden zu können, würde ich dann wohl an dieser Stelle als Solver pimpleFoam verwenden, richtig?

Das würde ich mir dann die nächsten Tage mal genauer ansehen/probieren.

Gibt es hier etwas zu beachten? Wie groß ist der Aufwand? Oder einen Art "Ratschlag" den Du mir vorab geben könntest? 

Beste Grüße!

Timo

P.S:

Absynthe ist u.a. ein Lied 

Eine Antwort auf diesen Beitrag verfassen (mit Zitat/Zitat des Beitrags) IP

---

Hi,

• Wie sieht deine Unterrelaxation vom Druck aus (Eintrag in fvSolutions)
  
• pimpleFoam für transiente Berechnungen, ja. Wie man die Algorithmen richtig verwendet (PISO / PIMPLE / SIMPLE), findest du in meinem Buch Mathematics, Derivations, Numerics and OpenFOAM
  
• Aufwand ist meines Erachtens nicht sehr groß (also um den Case umzubauen)
  
• Du kannst auch transiente Tutorials von meiner Seite holen (holzmann-cfd.de)
  

Gutes Gelingen. Grüße Tobi

------------------
Viele Grüße,
Tobi

OpenFOAM® Tutorials | Training | Publikationen | Für Anfänger wiki.openfoam.com

Eine Antwort auf diesen Beitrag verfassen (mit Zitat/Zitat des Beitrags) IP