CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM > OpenFOAM Programming & Development

All Mach number implicit solver with Kurganov-Tadmore scheme - pisoCentralFoam

Register Blogs Members List Search Today's Posts Mark Forums Read

Like Tree24Likes

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   February 28, 2022, 12:49
Default
  #101
Senior Member
 
mkraposhin's Avatar
 
Matvey Kraposhin
Join Date: Mar 2009
Location: Moscow, Russian Federation
Posts: 355
Rep Power: 21
mkraposhin is on a distinguished road
slicedSurfaceScalarField works as pointer to data of array (field) that was specified during construction.

Quote:
Originally Posted by strakey View Post
Matvey - Thanks for the reply. I think I understand this better now. I downloaded the OF6 version and compiled it for OF7 with some minor changes. I also downloaded the Yeqn.H and heqn.H from the 2112 branch and this also compiled fine under OF7. The MULES modifications you made here did have a big impact on diffusion of temperature. I also have it outputting the hLambdaCoeff as a volScalarField with an fvc::average() procedure and the coefficients are near 1 now as expected. I still, however see too much diffusion in the species, much more than energy. From the Yeqn.H code, it looks like all species equations are being modified by the same set of lambdaCoeffs, which are 0 regardless of what I set my inertSpecies to. I'm not sure how a slicedSurfaceScalarField works though. I'm going to insert some more output fields to see what exactly is being computed for lambdaCoeffs for each of the species (4 in my case) but it seems as though lambdaCoeffs should be different for each species field?


I also tried backward, vanLeer and other things to minimize diffusion but these all have a small effect.


Lastly, I don't know if you have seen this work (https://www.appliedccm.com/wp-conten...M-2019-DWS.pdf), but this seems similar to what you are doing here with your hybrid solver.


Thanks
Pete
mkraposhin is offline   Reply With Quote

Old   March 1, 2022, 10:54
Default
  #102
Senior Member
 
mkraposhin's Avatar
 
Matvey Kraposhin
Join Date: Mar 2009
Location: Moscow, Russian Federation
Posts: 355
Rep Power: 21
mkraposhin is on a distinguished road
Quote:
Originally Posted by strakey View Post
Matvey - Thanks for the reply. I think I understand this better now. I downloaded the OF6 version and compiled it for OF7 with some minor changes. I also downloaded the Yeqn.H and heqn.H from the 2112 branch and this also compiled fine under OF7. The MULES modifications you made here did have a big impact on diffusion of temperature. I also have it outputting the hLambdaCoeff as a volScalarField with an fvc::average() procedure and the coefficients are near 1 now as expected. I still, however see too much diffusion in the species, much more than energy. From the Yeqn.H code, it looks like all species equations are being modified by the same set of lambdaCoeffs, which are 0 regardless of what I set my inertSpecies to. I'm not sure how a slicedSurfaceScalarField works though. I'm going to insert some more output fields to see what exactly is being computed for lambdaCoeffs for each of the species (4 in my case) but it seems as though lambdaCoeffs should be different for each species field?



Update: It looks like the lambdaCoeffs (allFacesLambda) being used to modify the transport equations in Yeqn.H is that for the last species in the loop calling the MULES limiter,
Code:
forAll(Y, iCmpt)
    {
        volScalarField& Yi = Y[iCmpt];
        if ( Yi.name() != inertSpecie )
        {
            Info << "Limiting for: " << Yi.name() << endl;
            surfaceScalarField& rhoPhiYCorr = phiY[iCmpt];
        
            mulesWithDiffusionImplicitLimiter
            (
                rho,
                Yi,
                phi_own,
                phi_nei,
                allFacesLambda,
                rhoPhiYCorr,
                diffusiveFlux[iCmpt],
                mDCf[iCmpt],
                SuSp[iCmpt]
            );
        
        if (Yi.name() == "CH4")  //pete 2-28-2022
        {
        lambdavolH2O=fvc::average(lambdaCoeffs);
        }
        }
    }
it just keeps getting overwritten.



For my case, the lambdaCoeff for the last species is 0 everywhere, so 0 is being used for all of the species fields everywhere. I think hLambdaCoeff is just left at its initial value (1.0).



I also tried backward, vanLeer and other things to minimize diffusion but these all have a small effect.


Lastly, I don't know if you have seen this work (https://www.appliedccm.com/wp-conten...M-2019-DWS.pdf), but this seems similar to what you are doing here with your hybrid solver.


Thanks
Pete

Hi, Pete!


This is how MULES algorithm for multicomponent flow was designed. To make approximation of transport equations consistent between different phases, it takes minimum value of lambdaCoeff (the same approach is used for V-version of limited schemes in OpenFOAM, such as vanLeerV, MinmodV, etc):

\phi_f = \phi_f^{Upwind} + \lambda \times (\phi_f^{TVD} - \phi_f^{Upwind}) .

However, if \phi = 0, then MULES sets lambda to 0 and all other fields are transported with upwind. This is a very diffusive, but at the same time robust approach.





I think, this is rather technical issue, but it's solution requires some efforts.




Maybe you will find time to solve it
mkraposhin is offline   Reply With Quote

Old   March 8, 2022, 12:46
Default
  #103
New Member
 
Pete
Join Date: Feb 2016
Posts: 8
Rep Power: 10
strakey is on a distinguished road
Matvey - I'm working on this but haven't found a good way to fix this problem without going into the MULESTemplates.C code and I was hoping to not have to mess with that.


Also, I think I found a bug in your Y.eqn (dev-2112) code on line 187:
Code:
owner = mesh.owner()[iFace];
neighbour = mesh.owner()[iFace];
shouldn't the second line be:
Code:
neighbour = mesh.neighbour()[iFace];
I made that change but then the code threw an "index out of range error" on line 232. Not sure what is happening there.


Thanks
Pete
strakey is offline   Reply With Quote

Old   March 8, 2022, 14:14
Default
  #104
Senior Member
 
mkraposhin's Avatar
 
Matvey Kraposhin
Join Date: Mar 2009
Location: Moscow, Russian Federation
Posts: 355
Rep Power: 21
mkraposhin is on a distinguished road
Yes, you are right. This is a bug, thank you.

To circumvent indexing problem, wrap next code into if (mesh.isInternalFace(iFace)):

Code:
    if (mesh.isInternalFace(iFace))
    {
            owner = mesh.owner()[iFace];
            neighbour = mesh.neighbour()[iFace];
            deltaY = mag(Y[iCmpt][owner] - Y[iCmpt][neighbour]);
            if (deltaY > maxDeltaY.primitiveField()[iFace])
            {
                 maxDeltaY.primitiveFieldRef()[iFace] = deltaY;
            }
    }
Quote:
Originally Posted by strakey View Post
Matvey - I'm working on this but haven't found a good way to fix this problem without going into the MULESTemplates.C code and I was hoping to not have to mess with that.


Also, I think I found a bug in your Y.eqn (dev-2112) code on line 187:
Code:
owner = mesh.owner()[iFace];
neighbour = mesh.owner()[iFace];
shouldn't the second line be:
Code:
neighbour = mesh.neighbour()[iFace];
I made that change but then the code threw an "index out of range error" on line 232. Not sure what is happening there.


Thanks
Pete
mkraposhin is offline   Reply With Quote

Old   March 8, 2022, 15:21
Default
  #105
New Member
 
Pete
Join Date: Feb 2016
Posts: 8
Rep Power: 10
strakey is on a distinguished road
Thanks, that did the trick with the indexing error. Now I am seeing realistic values of hLambdaCoeffs and they are much better behaved than the lambdaCoeffs (only less than 1 in regions of steep phi gradients).


What I am testing now is moving the section of code: "Solve for Components" below the calculation of hLambdaCoeffs and using hLabmdaCoeffs on the Y equations (instead of lamdaCoeffs) as well as h. That way everything is using the same limiter and it seems to get around the phi=0 issue with MULES. I'll keep you posted.
Pete
strakey is offline   Reply With Quote

Old   March 16, 2022, 02:56
Default
  #106
Senior Member
 
mkraposhin's Avatar
 
Matvey Kraposhin
Join Date: Mar 2009
Location: Moscow, Russian Federation
Posts: 355
Rep Power: 21
mkraposhin is on a distinguished road
Please write if you get any progress
mkraposhin is offline   Reply With Quote

Old   February 27, 2024, 17:38
Default
  #107
New Member
 
Join Date: Aug 2022
Location: Mexico
Posts: 16
Rep Power: 3
Leckie is on a distinguished road
Quote:
Originally Posted by mkraposhin View Post
C = sqrt(Cp/Cv/(drho/dp)) becomes sqrt(Cp/Cv/psi) = sqrt(Cp/Cv*p/rho) in case of perfect gas EOS. Where psi = (drho / dp) = 1/(R*T) and rho = psi*p
This formulation is true for all compressible implicit OF solvers. If you will change EOS, for example - rho = A*(p/p0)^k, then you have to change equation for pressure. for example - see sonicLiquidFoam

Hello, I know it's been some years since the last post on this thread. I'm still learning OpenFoam and I'm actually doing some work at Mach 6-7 based on the rhoCentralFoam/forwardStep tutorial and I'm having issues when applying different initial conditions aside of P=1Pa and T=1K... What changes do I need to make in order to simulate flow at Mach 6-7 with P=200Pa & T=220K?
Leckie is offline   Reply With Quote

Reply

Thread Tools Search this Thread
Search this Thread:

Advanced Search
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
DPMFoam - Serious Error --particle-laden flow in simple geometric config benz25 OpenFOAM Running, Solving & CFD 27 December 19, 2017 20:47
foam-extend_3.1 decompose and pyfoam warning shipman OpenFOAM 3 July 24, 2014 08:14
Solver is finishing with huge Mach number Fonzie CFX 1 March 12, 2007 14:15
High Mach number solver error Luke CFX 3 January 31, 2007 22:26
TVD scheme at low Mach number Axel Rohde Main CFD Forum 5 August 6, 1999 02:01


All times are GMT -4. The time now is 16:42.