Iterative solution procedures are said to be converging when the changes in solved-for values become smaller with each iteration and are finally negligible. Convergence is thus what one wishes all flow-simulating computations to exhibit.

The following picture illustrates convergence. They are generated by PHOENICS Input-File-Library Case 249, concerned with laminar flow and heat in a square cavity with a moving lid.

The curves on the left plot the so-called "spot values",

That these curves become horizontal as the abscissa increases proves that indeed the changes have become very small; for the abscissa measures the iteration (called 'sweep') number.

The curves on the right represent the (logarithms of the) sums of the absolute residuals of each variable. These fall satisfactorily to the user-set target value; whereupon the computation terminates.

Not all computations converge so rapidly, as the next picture shows. It refers to the same flow situation, that of library case 249; but CONWIZ, which controls the built-in convergence promoter of PHOENICS, has been set to F rather than T.

Convergence does still occur, but more slowly; and the residuals do not fall far enough to terminate the calculation.

Convergence has occurred because the default values of the relevant relaxation parameters happened to suffice. Inspection of the RESULT file shows these to be:

If now the velocity-related defaults are over-written by inserting in the Q1 file:Group 17. Relaxation RELAX(P1,LINRLX,1.) RELAX(U1,FALSDT,1.) RELAX(V1,FALSDT,1.) RELAX(H1,FALSDT,1.0E+09)

the following picture is generated:RELAX(U1,FALSDT,1.E6) RELAX(V1,FALSDT,1.E6)

In this case the spot-values oscillate continuously; and so do the residuals. Because the values of the latter do not fall sufficiently, the calculation continues until the maximum allowed number of sweeps have been performed.

This run did not converge; but nor did it diverge. Divergence **is** shown
however by the following picture which resulted when the even more unsuitable
settings:

were made.Group 17. Relaxation RELAX(P1,LINRLX,2.0) RELAX(U1,FALSDT,1.E6) RELAX(V1,FALSDT,1.E6)

Evidently, the amplitudes of oscillations have begun to increase without limit; both for spot values and residuals. Such divergence renders the results of the calculations valueless.

The divergence of the last example was contrived deliberately by (as it
is said) **over-** relaxing the pressure variable, P1, through
choosing a value of its linear-relaxation factor in excess of 1.0.

Unfortunately divergence can appear in practice when the relaxation settings are far from perverse; indeed even when they are the best that the user knows how to select.

The rest of this article is devoted to providing advice on what is then to be done.

The following discussion is organized under the headings:

- (a) The choice of relaxation parameters
- (b) The choice of solution method
- (c) The choice of maximum and minimum values
- (d) The choice of initial guesses
- (e) Easily-removed causes of slow convergence
- (f) Solution of one volume fraction equation
- (g) Other measures to be considered

- Linear relaxation is almost always preferable to the false-time-step alternative. The reason is that the appropriate value for this time step is hard to choose; and what has been found suitable for one size of domain may be entirely unsuitable for another, even though other aspects of the flow may be identical.
A common error of inexperienced users is to set very small values of the false time step,

*e.g.*1.0e-3 sec for the velocity. They then find that the oscillations which have alarmed them disappear; and they are comforted by the apparently smooth progress of the solution procedure, not recognising that it is proceeding so slowly that an intolerably large number of sweeps will be needed to attain the**true**solution.For a wind-farm simulation for example, in which a 10 km width and a 10 km/sec wind imply that effects take 1000 sec to be convected across the domain, even a 1.0 sec false time step imposes a severe delay on the attainment of the correct solution.

- To set CONWIZ=T in the Q1 file, and to make no other settings is usually beneficial; for that introduces under-relaxation of the so-called dveldp's,
*i.e.*the rate-of-change-of-velocity-with-pressure coefficients, which play an important role in the so-called 'SIMPLEST' algorithm used by PHOENICS to handle efficiently the links between the continuity and momentum equations.The default value of the dveldp relaxation parameters is 0.5; but users can use change these by inserting into the Q1 statements such as:

in order to slow down the sweep-to-sweep changes in values.**SPEDAT(RLXFAC,RLXDU1,R,0.25) SPEDAT(RLXFAC,RLXDV1,R,0.25)** - The greater the number of interacting dependent variables, the more likely is it that convergence may be difficult fo achieve. Consequently it is sometimes useful to apply the relaxation not to the dependent variables themselves but to the auxiliary variables which effect the interactons. Two examples are:
- In
**natural-convection heat-transfer**problems, especially those with heat sources, the fluid flow influences the temperature distribution, which affects the density, which affects the gravitational force, which then influence the flow.

In these circumstances under-relaxation of the**density**may promote convergence. - When
**turbulent flows**are being simulated by way of, say, the k-epsilon turbulence model, the fluid flow influences the distributions of k and epsilon, which affects the turbulent-viscosity distribution which in turn affects the flow.

In these circumstances it is under-relaxation of the**turbulent viscosity**may promote convergence.

- In
- Two final points to make about the choice of relaxation parameter are:
- Because the variety of flows which users of PHOENICS may expect the code to simulate is vast, it is not possible to prescribe infallibly successful rules for choosing relaxation parameters.
- In particular, should users introduce,
*via*In-Form, PLANT or their own GROUND coding, representations of unusual forces acting on the fluid, they should themselves take steps to prevent divergence.

For example, If the In-Form source formulation is being employed the 'with line' (*i.e.*linearisation) option should be employed whenever the source formula contains the dependent variable itself.

- Because the variety of flows which users of PHOENICS may expect the code to simulate is vast, it is not possible to prescribe infallibly successful rules for choosing relaxation parameters.

- In three-dimensional problems, or even two-dimensional ones for which NZ exceeds 1, switching from the 'slab-wise' to the 'whole-field' solution procedure usually promotes convergence.

This is something which happens automatically if CONWIZ is set 'true'. - For body-fitted-co-ordinate problems, I is always wise to choose the
GCV=T option; for this takes especially accurate account of the departures from orthogonality which may sometimes occur.
It is especially important for flow-over-natural-terrain problems; for steep gradients of just one part of the ground surface may suffice to cause divergence if GCV is left at its default value, namely F.

An example is this:

- It is desired to calculate the steady-state flow and tenperature within a living room.
- There is a heat source within the room; and there is an open window. But air circulation results from natural convection alone.
- The user relies on PHOENICS to find out what the velocity distribution and therefore leaves the initial values at their defaults.
- At the very first sweep, PHOENICS finds that there are cells with the finite heat sources which the user prescribed, but no flow to carry the heat away. It therefore computes enormous temperatures for those cells; and large buoyancy forces ensue, especially if the otherwise reasonable Boussinesq approximation is in use.
- Resulting next-sweep high velocities trigger divergence.

The setting of the PIL variables VARMIN and VARMIX can prevent this from occurring.

Of the two modes of using these variables, that which sets limits to the **increments** of variables is the preferable. Thus, if VARMIN(temperature) is set to -1.e-11 (to trigger that mode) and VARMAX(temperature) to 10.0 deg Celsius, the above-mentioned second-sweep divergence will simply be unable to occur.

Less drastic but still troublesome failure to converge can usually be controllled in this way.

If one knew the solutions to the equations at the start, and could insert the corresponding variable-field values as initial guesses, convergence would be immediate. Of course, the solutions never are known; however, if initial guesses can be made which are close to these solutions, convergence is likely to be rapid.

It is therefore usually worthwhile to set values of the initial fields, by way of
FIINIT, which are as close as possible to the expected average values of the solutions. It
is possible to do even better, if the PATCH and INIT commands are used; for these permit
linearly-varying values (and indeed more complicated piece-wise linear distributions) to
be inserted rather easily. Even more powerful is
In-Form's (initial of *variable_name* is *formula*) command.

How far to go in this direction should be determined by reference to the man-time cost of determining and introducing suitable initial values and of the computer-time cost that is thereby saved.

A means of saving computer time that requires very little man-time, once the system has been set up, is to employ the SPINTO facility. This launches a series of runs, on increasingly finer grids. Each run uses as its initial field values obtained by interpolation from the output of the previous run.

It often occurs that convergence proceeds smoothly, but not fast enough. In order to establish whether acceleration may be possible, the cause of slow convergence should first be sought among the following:-

- Possibly too much under-relaxation is being applied. This is often the case when the user has had a struggle to procure convergence at all, has achieved it in a rather unsystematic manner by applying under-relaxation to every variable he could think of, and, being thankful that his calculation now is converging, is fearful of changing any of the settings.
- Possibly too much time is being spent in obtaining intermediate solutions of excessively
high accuracy for variables which have little influence on the number of sweeps which are
needed for attainment of the final solution.

For example, there is no point in solving for a scalar variable at all, if it does not directly or indirectly influence the velocity field, until after the velocity field has converged; and one should certainly not set low values of ENDIT and high values of LITER for such variables while velocity-field adjustments are in progress. - Even for variables which do participate significantly in the velocity-field solution, such as the pressure, or one of the velocity components, excessive demands may be being made at the linear-equation-solution level.
- Similarly, LSWEEP may have been set to a very large number, and at least some of the RESREF's to excessively low values, so that the computer is either producing solutions of a higher accuracy than are needed, or is continuing to perform adjustment cycles without (because of round-off error) increasing the accuracy of the solution at all. It is best indeed not to set RESREFs at all, but to allow EARTH to calculate then by settting SELREF=T.
- If the flow being simulated is a time-dependent one, the time step chosen may be too small. This will often be the case if the time step was selected so as to ensure convergence, rather than because the phenomenon is a rapidly-varying one.

How can one determine whether one of the above is contributing to excessive computer charges? By changing the relevant settings and observing what happens. Every reduction in the amount of under- relaxation, even an enlargement in time step, is to be welcomed, provided that it does not on the one hand lead to divergence or on the other entail a significant falling off in numerical accuracy.

It is recommended that only the second-phase volume-fraction equation is solved when
the two phases coexist everywhere in finite proportions, *e.g.* the volume fraction of one of
the phases never falls below 0.001. The following PIL instruction:

SOLUTN(R1,Y,N,...); SOLUTN(R2,Y,Y,...),

causes EARTH to solve for R2 and to calculate R1 as 1-R2. In addition to being more economical than solving for R1 and R2, this practice has proved to give faster convergence in several 2-phase simulations.

When the phase volume fractions do not coexist everywhere in finite proportions, it is necessary to solve R1 and R2, by way of the settings:

SOLUTN(R1,Y,Y,...); SOLUTN(R2,Y,Y,...).

Were R2 alone solved in this situation, loss of accuracy would occur in the calculation of R1 from 1-R2 in those location in which R2 is close to unity.

There are other measures which can be taken at the SATELLITE level in order to reduce computer time. One important one is to change the choice of solution procedure in the light of the considerations discussed in Section 5.4. The point-by-point procedure should be used for the velocity components and volume fractions when the diffusion terms are of little importance; but it should definitely not be used otherwise.

Different velocity components may need different solution procedures. For example, a boundary-layer flow in which z is the predominant flow direction is best handled with point- by-point solution for u and v, but slab-wise solution for w. The reason is that w is strongly influenced by viscosity in these circumstances,whereas u and v are not.

In some fluid-flow simulations, it is best to use the TERMS command in order to switch off the viscous terms entirely. More surprisingly, perhaps, it is sometimes permissible to switch off the convection terms as well. For example, when the flow takes place in a porous medium such as EARTH or sand, the only important terms are the pressure gradient and the resistance exerted by the medium; so there is no point in computing terms which are negligible. The Encyclopaedia entry, Darcy, explains how to activate this feature.

A number of other devices *e.g.* ADDDIF, GALA, ZDIFAC, NEWRH1 etc. are available for
solution control and economy. These are in Group 8; they are described in the Encyclopaedia.

Finally, attention is drawn once more to the possibility of varying the starting-time and the frequency of updating of the various dependent variables, and for that matter of the auxiliary variables as well. PHOENICS possesses many devices for effecting such variations, which you are advised to explore.

The guiding principle should be this: quantities which both influence and are influenced by the velocity (or other dominant-variable) field should be updated frequently throughout the computation; those which influence but are not influenced should be set at the beginning and never updated; and those which are influenced but do not influence should be updated only at the end.