[hjasak]

You can't help it: in order to satisfy mass conservation, you have to solve the pressure equation implicitly. Whatever else you do will give you checkerboarding.

As for the convection term in the momentum equation, thing about what it means in terms of cost to make the convection implicit. Also, using upwind differencing introduces a lot of discretisation error and will kill your LES.

Hrv

[eugene]

I suggest you try and reformulate your convection term into implicit part with an explicit correction.

i.e.
Ustar = c*U_implicit + f(U_explicit)
phi_star = fvc::surfaceInterpolate(Ustar) <- can remain fully explicit

thus:

fvm::div(c*phi_star, U_implicit) +
fvc::div(phi_star, f(U_explicit)

(I offer no guarantee that this analysis is correct!)

Of course this would mean redoing your gridFilter function to output the components of c and f(U_explicit), instead of just doing an explicit filter operation.