usnistgov / fipy

FiPy is a Finite Volume PDE solver written in Python
http://pages.nist.gov/fipy/en/latest
Other
504 stars 148 forks source link

Stokes example Rhie-Chow is correct? #689

Open evstigneevnm opened 4 years ago

evstigneevnm commented 4 years ago

Using stokesCavity.py added verification of div(V) as:

velocity_d = FaceVariable(mesh=mesh, rank=1)
velocity_d[0] = xVelocity.arithmeticFaceValue;
velocity_d[1] = yVelocity.arithmeticFaceValue;
velocity_d[..., mesh.exteriorFaces.value] = 0.
velocity_d[0, mesh.facesTop.value] = U
div_value = max(abs(velocity_d.divergence.value))

and it turns out to be constant (near the top corners of the cavity) with value about 3.713 for 40X40 grid with almost no decrease on sweeps. Is that valid?

wd15 commented 4 years ago

I think that the example does point out that there is error near the boundaries due to a number of missing features in FiPy so this isn't that surprising.

evstigneevnm commented 4 years ago

I checked OpenFOAM cavity test example and there's the same divergence error, about BigO(1) and it doesn't depend on the mesh resolution. I guess it's the flaw of the Rhie-Chow method for this particular problem due to singularities in corners. If div(V)!=0 it means that there are sources and drains locally available that may corrupt correct solution of the momentum eq-s. Maybe add another test case (and remove this one) where Rhie-Chow performs better? If you don't agree the issue can be closed.

wd15 commented 4 years ago

I'm curious. Did you get the same issues in OpenFOAM even with the features switched on that FiPy is missing?

I probably won't be adding any more flow examples unless I get back into FiPy development, which won't be happening anytime soon. Thanks for raising this issue though.

evstigneevnm commented 4 years ago

Dear wd15, Yes, it's the same behaviour for a test example that is provided in OpenFOAM in: $FOAMTUTORIALS/incompressible/icoFoam/cavity/cavity The point is that Rhie-Chow method recovers continuity only up to the order of approximation for periodic boundary conditions. However, at the presence of sharp corners (singularity) such convergence is no longer true. The behaviour is as follows: ||div(V)||{C} grows, but the problem localizes itself near singularity, hence with the decrease of grid spacing ||div(V)||_{L_2} remains bounded from above by a constant. So the convergence of the whole method in terms of honest div(V) is as BigO(1) in numerical L2 norm. I will put a test project with some matlab examples, divergence estimation for OpenFOAM and div(V) norms later on my github and put a link here. Best regards, Nick

evstigneevnm commented 4 years ago

Dear wd15, sorry for the delay. Here is the link to the test using OpenFOAM. divergence test in 2D lid driven cavity flow Best regards, Nick