>diaw: But, surely we have then set about to solve the wrong governing equation?

Exactly. It's good to realize that, because that's what discretization does: You end up solving a slightly different equation.

>This is my contention with the addition of artificial viscosity (upwinding etc) - leading to 'false diffusion' in the resuting flow field.

Not exactly. Don't confuse numerical dissipation (the discretization error of you base scheme) with artificial dissipation. Artificial dissipation is meant as a fix to the above problem. You are trying to eliminate the error terms in your equation to get back to the original ones, which would be a possible solution to the problem mentioned above. However, the error term is not known exactly, so to be on the safe side we usually add a little too much diffusion. If you can find a way to exactly eliminate the discretization error without excessive dissipation (and do so efficiently) you'll be famous. (good luck!) I am not convinced this is possible. You'll always have to pay a price for discretization.

>In the limit as 'dx->0', a correct numerical scheme must surely arrive back at the original pde?

Absolutely. This condition is referred to as "consistency". Consistency and stability will guarantee convergence (according to Lax, I think). Otherwise your scheme is no good.

>Adding terms to the original equation is not equivalent to transform-solve-inv_transform process, which is where I hope to eventually end up.

I disagree. Whatever you do to modify your discrete equations, you could always recast it into an addition of certain terms. Your transformation might make it easier or more efficient, but in the end it could be done otherwise.

>The manipulation of viscosity is what upwinding seems to do.

That's over-simplified (at least in the case of higher order upwinding). I see upwinding as a smart way of interpolation. Interpolation is necessary because of discretization. For example, with finite volumes, if you describe flow variables at a cell center, you'll need some way of propagating those cell values to the cell faces in order to evaluate fluxes (equivalent problems you have with FEM or FDM). Upwinding takes the physical wave propagation (direction and/or information) into account in order to progagate information in a physically correct way.

One point on the side: If you do find that your instability is caused by singularities of the N-S equations, what are you going to do about it? If you're aiming for an exact solution of the N-S equations, shouldn't you want your scheme to blow up in case of a singularity (because that would be the correct solution)? Or would you rather have a robust scheme that gives you a solution even when there is none in reality.

I am interested in you cylinder experiment from the above point of view, but maybe it comes all back to your research. You say that some codes will give you a steady state even beyond the Reynolds number at which such state exists in reality, while some other codes refuse to converge beyond that point. I would like to understand (just like you), what the reason for this is. What are the characteristics of those two classes of schemes: a) the ones that force a steady state, even when it's unrealistic, b) the ones that seem to mimic real instabilities by refusing to converge in an unstable case? If your reasearch can answer this question in a general way (not just for FEM, FVM or FDM, but for any numerical treatment) it will be a great contribution.

Thanks Mani for an excellent thought process...

--------- Sub-section 1:

>diaw: But, surely we have then set about to solve the wrong governing equation?

>Mani: Exactly. It's good to realize that, because that's what discretization does: You end up solving a slightly different equation.

>diaw: This is my contention with the addition of artificial viscosity (upwinding etc) - leading to 'false diffusion' in the resuting flow field.

>Mani: Not exactly. Don't confuse numerical dissipation (the discretization error of you base scheme) with artificial dissipation. Artificial dissipation is meant as a fix to the above problem. You are trying to eliminate the error terms in your equation to get back to the original ones, which would be a possible solution to the problem mentioned above. However, the error term is not known exactly, so to be on the safe side we usually add a little too much diffusion. If you can find a way to exactly eliminate the discretization error without excessive dissipation (and do so efficiently) you'll be famous. (good luck!) I am not convinced this is possible. You'll always have to pay a price for discretization.

----- diaw: Do you have a few references to the ideas used ot get a handle on the size, or estimate, of the discretization & numerical dissipation errors? (I have Prof. Jasak's thesis as a starter - a few additonal references would be very useful).

Part of my research work is the estimation of the discretization residuals. I have been deriving them major term group at a time. I still have a long way to go.

---------

Sub-section 2:

>diaw: In the limit as 'dx->0', a correct numerical scheme must surely arrive back at the original pde?

>Mani: Absolutely. This condition is referred to as "consistency". Consistency and stability will guarantee convergence (according to Lax, I think). Otherwise your scheme is no good.

------ diaw: Thanks. The discretization scheme - in the limit as 'dx->0' must give the pde, otherwise we will be introducing some level of permanent offset, I would imagine - leading to eventual numeric overflow.

----------

Sub-section 3:

>diaw: Adding terms to the original equation is not equivalent to transform-solve-inv_transform process, which is where I hope to eventually end up.

>Mani: I disagree. Whatever you do to modify your discrete equations, you could always recast it into an addition of certain terms. Your transformation might make it easier or more efficient, but in the end it could be done otherwise.

---------- diaw: Ok, fair comment - I see where you are going with this logic. As long as one is keeping track of the relevant terms. Good.

----------

Sub-section 4:

>diaw: The manipulation of viscosity is what upwinding seems to do.

>Mani: That's over-simplified (at least in the case of higher order upwinding). I see upwinding as a smart way of interpolation. Interpolation is necessary because of discretization. For example, with finite volumes, if you describe flow variables at a cell center, you'll need some way of propagating those cell values to the cell faces in order to evaluate fluxes (equivalent problems you have with FEM or FDM). Upwinding takes the physical wave propagation (direction and/or information) into account in order to progagate information in a physically correct way.

--------- diaw: I was refering (perhaps a little simplistically) to how 1D simple upwinding affects the singularity point - in the FEM method. I have observed that the 'alfa' value basically squeezes the effect of the convection term to zero (at 'alfa'=1), & hence pretty-much guarantees stability. Granted, this is for one particular scheme. I do not know of any current 2d, 3d upwinding schemes for triangular elements in FEM.

I plan to thoroughly investigate the carry-through of simple-upwinding & higher-order schemes for FVM, & FDM.

----------

Sub-section 5:

>Mani :One point on the side: If you do find that your instability is caused by singularities of the N-S equations, what are you going to do about it? If you're aiming for an exact solution of the N-S equations, shouldn't you want your scheme to blow up in case of a singularity (because that would be the correct solution)? Or would you rather have a robust scheme that gives you a solution even when there is none in reality.

---------- diaw: 1) If you know where the predicted 'dominant' singularity is located - you aim to avoid it. In a practical sense this translates into setting bounds on the element dimension/s used. From what I have observed in 2d, one singularity is dominant. You aim to avoid him...

and/or,

2) You use a suitable transformation scheme to set up a suitable calculation space in which dominant singularities are pushed far-enough out of reach to allow stable-accurate computation. (This is my 'Holy Grail'

-----------

Sub-section 6:

>Mani: I am interested in you cylinder experiment from the above point of view, but maybe it comes all back to your research. You say that some codes will give you a steady state even beyond the Reynolds number at which such state exists in reality, while some other codes refuse to converge beyond that point. I would like to understand (just like you), what the reason for this is. What are the characteristics of those two classes of schemes: a) the ones that force a steady state, even when it's unrealistic, b) the ones that seem to mimic real instabilities by refusing to converge in an unstable case? If your reasearch can answer this question in a general way (not just for FEM, FVM or FDM, but for any numerical treatment) it will be a great contribution.

------

diaw: I have observed the following, at this stage:

a) The schemes that seem to force a steady state - FEM = CBS scheme (splitting as I understand it), for FVM, in the main so far the Star-CD-type solvers seem to be rather solid, although eventually they also succumb.

b) The schemes that seem to mimic the real instabilities are : FEM - fully coupled schemes. I have tested an internal 6-node triangular scheme & quad scheme. The internal code will not under any circumstances converge beyond the instability, at least not into any sensible field - there are some simple problems where it will appear to settle - but one look at the flow field will show total garbage. In most cases, they simply yoyo back & forth.

For the split-scheme, I have noted that using a 'pseudo-time-step' that is too small often ends up with a garbage converged flow-field.

These are really more observations at this point.

----------

Thanks again for your really constructive input.

diaw...

Some feedback on N-S singularities & the complex plane.

After Tom's comment regarding the presence of singularities being present in the 'complex plane' I performed a mathematical analysis & transformed an ODE (1D N-S) into the complex plane... There they are - to be sure, but... They are actually 'unstable' as they appear in a simplistic sense...

Now this is getting interesting...

---------- Previous thread: Q: diaw: (1) "Is the reason why typical mathematics approaches cannot find the singularities is that the 'differential' 'transformation' process simply moves all the singularities to the origin as we decrease delta'x' to zero."

A: Tom: I suspect that the singularities shoot off into the complex plane as the grid spacing is refined. The appearance of them on the real line is an artefact of the discretization. A good example of a discrete problem behaving very differently to the continuous one is to discretize the logistic equation using a forward Euler scheme. The resulting discrete equation possesses a period doubling cascade and chaos as the free parameter (timestep) is varied while the continous problem has no such behaviour. In the limit of the timestep going to zero the discrete problem approaches the continuous soultion but outside this range there is no reason for them to share the same behaviour.

diaw's reply: I would tend to agree with the singularities being squeezed somewhere as dx is reduced to zero - perhaps they then squeeze onto the complex plan to be reflected there - nice point.