[usv001]

Dear all,

I am trying to implement an implicit solver for solving steady-state compressible flows in OpenFOAM. With regards to this, I have coded a matrix-free LU-SGS method [1] with optional preconditioning. I have also coded Venkatakrishan slope limiter [2] which is supposed to provide improved convergence and the SLAU scheme [3] which is capable of handling low Mach number regions.

I tried running the solver with the 1st order upwind scheme for M~0.01 flow over a cylinder. While the solution converged, there was some asymmetry in the solution (see attached image). The journal papers report perfectly symmetric profiles. So, I am not sure where I am going wrong. I increased farfield boundary distance from 10D to 40D but it had no effect. By the way, I am applying fixedValue boundary conditions.

When I switched to using TVD schemes such as van Leer and van Albada, the solution did not converge. I reported this problem in this thread and a user suggested that I use Venkatakrishnan's limiter. However, when I tried to use the 2nd order scheme with Venkatakrishnan's slope limiter, the solution does not converge either. I have tried varying the parameter K from 0.1 to 5 but there has been no improvement in convergence thus far. I have tried running with the HLL flux scheme too but the same thing happens. It converges for 1st order but not for second order. This led me to believe that there are two problems here: one with the flux scheme and the other with the spatial scheme.

To summarise:

1. Running SLAU with 1st order upwind scheme results in asymmetric solution for M~0.01 flow past cylinder.
2. Using Venkatakrishnan's slope limiter does not seem to improve convergence.

Any comments/suggestions would be of great help.

USV

PS: I am using OpenFOAM v2.4.

References

1. H. Luo, J.D. Baum, and R. Löhner. A Fast, Matrix-free Implicit Method for Compressible Flows on Unstructured Grids. J. Comput. Phys., 146:664-690, 1998
2. V. Venkatakrishnan. Convergence to Steady State Solutions of the Euler Equations on Unstructured Grids with Limiters. J. Comput. Phys., 118:120-130, 1995
3. E. Shima and K. Kitamura. Parameter-Free Simple Low-Dissipation AUSM-Family Scheme for All Speeds. AIAA Journal. 49:1693-1709, 2011

[usv001]

Dear all,

I realized that I didn't provide much information to help with the troubleshooting. So, I decided to include parts of my code.

I have implemented Venkatakrishnan limiter as a limitedGradScheme following the same general layout as the cellLimitedGrad scheme. The important function is the limitFace function which is given by:

template<>
inline void VKcellLimitedGrad<scalar>::limitFace
(
    scalar& limiter,
    const scalar& maxDelta,
    const scalar& minDelta,
    const scalar& extrapolate,
    const scalar& epsSqr
)
{
    const scalar& deltaNeg = extrapolate;
    
    scalar deltaPos(0);
    
    if (extrapolate > SMALL)
    {
        deltaPos = maxDelta;
    }
    else if (extrapolate < -SMALL)
    {
        deltaPos = minDelta;
    }
    
    scalar deltaPosSqr(sqr(deltaPos)), 
    	   deltaNegSqr(sqr(deltaNeg)),
    	   deltaPosNeg(deltaPos*deltaNeg);
    	   
    scalar limiter0 
    	= 
    	(
    		  deltaPosSqr 
    		+ epsSqr
    		+ 2*deltaPosNeg
    	)/
    	(
	          deltaPosSqr 
		+ 2*deltaNegSqr
		+ deltaPosNeg
		+ epsSqr
    	);
    	
    limiter = min(limiter, min(limiter0, 1));
}

This is supposed to implement the Venkatakrishnan limiter which is given by:

I have simplified the expression by dividing the numerator by . I don't think this is critical. Please correct me if I am wrong though.

In the main function calcGrad, the important piece of code that is different from the calcGrad function in cellLimitedGrad scheme is coloured red:

// create limiter
scalarField limiter(vsf.internalField().size(), 1.0);
    
const surfaceScalarField& deltaCoeffs = mesh.surfaceInterpolation::deltaCoeffs();

forAll(owner, facei)
{
    label own = owner[facei];
    label nei = neighbour[facei];
    
    // eps^2 = (K*dx)^3
    scalar epsSqr = pow3(K_/deltaCoeffs[facei]);
        
    // owner side
    limitFace
    (
        limiter[own],
        maxVsf[own],
        minVsf[own],
        (Cf[facei] - C[own]) & g[own],
        epsSqr
    );

    // neighbour side
    limitFace
    (
        limiter[nei],
        maxVsf[nei],
        minVsf[nei],
        (Cf[facei] - C[nei]) & g[nei],
        epsSqr
    );
}

The above limiter runs fine (i.e. it is accurate and stable) but it does not produce the convergence it is supposed to regardless of the value of K_. I'd be grateful if someone could tell me where I am going wrong or what I should do differently. I'll be glad to provide any other relevant information.

USV

[naffrancois]

I am not an open foam user, what is behind fixed value boundary condition? Is it a dirichlet bc? Do you prescribe everything at the bc? If so I would rather use a characteristic bc such as a Farfield bc based on Rankine Hugoniot relations or acoustic relations. Given the trials you have already done, a wrong bc set up could explain convergence issues. Dirichlet bc can yield wave reflections, not easily smoothed out using 2nd order reconstruction while characteristic based bc help a lot in allowing pressure waves to exit the domain

As far as asymmetry is concerned, if your solution is converged, it may be caused by the wall bc which generates too much entropy. Do the paper about slau mention if they use a specific bc such as curvature corrected wall bc?

[usv001]

Hello naffrancois,

You're right. fixedValue prescribes everything at the boundary. I modified it to waveTransmissive which should allow the acoustic waves to leave the domain but it didn't have any effect on the asymmetry. Since the farfield is quite some distance away from the cylinder, I don't think that it is the main thing.

On the other hand, you're quite right about the too much entropy being generated near the wall. I was using a simple symmetry technique (ST) at the wall. I went through the papers on curvature corrected symmetry technique (CCST). I was able to get symmetric results for a M=0.3 flow over a cylinder with Venkatakrishnan limiter (using 1st order upwind scheme, generated the asymmetry). For this case, it was crucial for the farfield BC to be waveTransmissive for good convergence. I compared ST and CCST. It turns out that the differences were quite subtle but CCST was more accurate. Thanks for bringing that to my attention! I've attached am image of the results.

By the way, I checked the paper on SLAU and it says the following:

For the solid wall, the density and pressure of ghost cells are given by zeroth-order extrapolation. Point symmetric or mirror symmetric velocity vectors are set for nonslip and slip walls, respectively.

It is not exactly clear as to which method (ST or CCST) they are using.

In summary,

• I think my current implementation of Venkatakrishnan limiter is fine.
• I am guessing that the asymmetry was caused by the high numerical diffusion of the 1st order upwind schemes and, in a small part, due to the wall BC.
• Using waveTransmissive was important for good convergence.

So, the main problem now is that the solver is not able to handle very low Mach numbers but it performs fine at the higher subsonic values. Does anyone has any thoughts on this?

Many thanks,
USV

[dario_z]

Hi USV,


Which low Mach preconditioner are you using (if any)?

Have you been able to progress further?

[AshwaniAssam]

Does it run for case M=0.1, 0.2. SLAU may have some issues at low Mach number. But 0.1 seems fine.



When you run second order scheme. Did you get symmetric solution, even though your simulation might have not converged?


Also, what do you mean by asymmetry? Have you tried to match with Cp or something at the wall.