 # On Estimating Machine-Zero Residual

In this paper, we propose two techniques to estimate the magnitude of a machine-zero residual for a given problem, which is the smallest possible residual that can be achieved when we solve a system of discretized equations. We estimate the magnitude of the machine-zero residual by a norm of residuals computed with a randomly-perturbed approximate solution that is considered as close in magnitude to an exactly-converged solution. One method uses free-stream values as the approximate solution, and the other uses a current solution during an iterative solve as the approximate solution via the method of manufactured solutions. Numerical results show that these estimates predict the levels of machine-zero residuals very accurately for all equations of the Euler and Navier-Stokes equations in a transonic flow over an airfoil and viscous flows over a cylinder and a flat plate.

## Authors

##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

A machine-zero residual is defined as a residual (e.g., a finite-volume discretization of the Euler equations) with a converged solution substituted having machine-zero-level perturbation. It is the smallest possible residual value that a solver can achieve; it is difficult to predict especially before a calculation. Also, it varies widely depending on the grid unit and the solution units especially when solving a dimensional form of flow equations. However, if it were possible to estimate a machine-zero residual before or during a calculation, then it would be very useful. For example, if an initial solution happens to be a converged solution (e.g., a restart) or close to it, then a machine-zero residual estimate would help us determine that the residual is well satisfied already and thus we would stop even before calling a solver. On the other hand, if there is no estimate available, we will find that a solver generates very small solution updates but cannot tell if the solution is already converged (it could be stalling). One can look at the corresponding residual level, but again one cannot tell whether it is the smallest possible value (it could be stalling). Note again that the machine-zero residual level depends on many different factors and thus it can be very small or very large: e.g., being e-09 does not necessarily mean a machine-zero residual.

In practical fluid-dynamics solvers, iterative convergence is checked based on the reduction of the residual or iterative-solution-difference from their initial values. This is a practical thing to do, but it can be inefficient. For example, if an initial solution is very close to a converged solution (e.g., a solution at a previous time step in an inner iteration of an implicit time-stepping scheme with a very small time step), then the residual is already very small at the beginning of an iteration. Reducing the residual by two or three orders of magnitude may be much more than necessary; one order of magnitude reduction may suffice. If an accurate estimate of a machine-zero residual is available, then we can determine how well the discrete equations are satisfied (e.g., five orders of magnitude above the machine-zero level) and stop the iteration earlier and could save computing time by a considerable amount over an entire unsteady simulation.

To the best of the author’s knowledge at the time of writing this, no techniques are currently available for estimating a machine-zero residual level for a given flow problem. It is the objective of this paper to propose two such techniques and investigate how well they estimate machine-zero residual levels for inviscid and viscous flow problems in two dimensions.

## 2 Machine-Zero Residual

Consider a nonlinear system of equations, which arise from a finite-volume-type discretization of the dimensional Euler/Navier-Stokes equations,

 R(U)=0, (1)

where and

denote global vectors of residuals and numerical solutions in a given grid, which we solve by an iterative solver:

 Uk+1=Uk+ΔUk, (2)

where is the iteration counter and is a correction suggested by an iterative method (e.g., Newton’s method). Let be an exactly converged solution for a given problem with appropriate boundary conditions:

 R(¯¯¯¯¯U)=0. (3)

In reality, the iteration will not find exactly but converge to with perturbation of the order of machine zero, , where is a machine zero for a value of (e.g., 1.0E-16 in double precision). This is the best we can hope for. Then, the minimum residual that can be attained, which we call the machine-zero residual, would be

 Rmin=R(¯¯¯¯¯U+ϵ¯¯¯¯¯U), (4)

which may be expanded as

 Rmin=R(¯¯¯¯¯U)+¯¯¯¯¯¯¯¯¯¯∂R∂Uϵ¯¯¯¯¯U+O(ϵ2)=¯¯¯¯¯¯¯¯¯¯∂R∂Uϵ¯¯¯¯¯U+O(ϵ2), (5)

where is the residual Jacobian evaluated at . It shows that the machine-zero residual depends on the magnitude of the solution as well as the Jacobian. Hence, there are various factors that affect the magnitude of the machine-zero residual: e.g., grid units and physical units. Therefore, it is not necessarily of the order of 1.0E-16.

At this point, it is clear that the machine-zero residual cannot be computed because the exactly-converged solution is not known.

## 3 Estimates of Machine-Zero Residual

To estimate a machine-zero residual, we need to approximate the exactly-converged solution . As our interest is only in estimating the magnitude of the machine-zero residual in some norm, the approximation just has to be close in magnitude to the exactly-converged solution . There are two possibilities.

### 3.1 Free-stream estimate Rc

One is a constant free-stream solution:

 Rc=||R(U∞+ϵrU∞)||, (6)

where is a random number () and denotes a constant solution, which is an exact solution with a free-stream condition applied at all boundaries: . Note that in practice we apply different random numbers to solutions at different nodes (or cells) as we will show below. It is reasonable to set a free-stream condition at all boundaries, but we did not do so for the test cases shown in this paper. This estimate can be computed before a calculation and also separately for different equations in a target system.

In our dimensional Euler/Navier-Stokes solver, we store the primitive variables at nodes, where and are velocity components, is the temperature, and is the gauge pressure , where is a free-stream pressure. Hence, the residual is a function of the primitive variables stored at nodes: , where denotes the number of nodes in a grid. We define the perturbed free-stream solution at a node as

 (7)

where is a random number defined at the node and the maxmod function is defined by

 maxmod(a,b)={a,if|a|≥|b|,b,otherwise, (8)

Note that the free-stream value of the gauge pressure is zero. The velocity components and are computed for a given set of a free-stream Mach number and an angle of attack. It is emphasized again that the perturbed free-stream solution is different at different nodes.

### 3.2 Current solution as a manufactured solution Rm

The other possibility is to use a current solution as an exact solution, which is possible by the method of manufactured solutions. Let be a source term defined by

 S=R(Uk), (9)

Then, is an exact solution to the problem,

 R(U)=S, (10)

which leads to the following estimate:

 Rm=||R(Uk+ϵrUk)−S||. (11)

This is a machine-zero residual for a slightly different problem, but expected to serve well as an estimate, at least in magnitude, for the original problem. In fact, if we expand it as

 Rm≈∣∣∣∣∣∣∂R∂UϵrUk∣∣∣∣∣∣, (12)

we find that the estimate is an approximation to the machine-zero residual in the form (5). It would be an accurate estimate if the current solution is close to a converged solution in magnitude. In practice, we again apply different random numbers to solutions at different nodes. As in the previous case, this estimate can also be computed for each equation in a target system. It may be updated at every iteration if needed. Note that this technique is more general than the free-stream estimate ; it can be used with a free-stream solution by setting (then there is no need to specify a free-stream condition at all boundaries).

In our dimensional solver, we define the perturbed solution in the primitive variables, at a node , as

 wkj+ϵrwkj=⎡⎢ ⎢ ⎢ ⎢ ⎢⎣maxmod(p′j(1+ϵrj),ϵrj)maxmod(uj(1+ϵrj),ϵrj)maxmod(vj(1+ϵrj),ϵrj)Tj(1+ϵrj)⎤⎥ ⎥ ⎥ ⎥ ⎥⎦, (13)

where the variables are values at the -th iteration.

## 4 Results

In this section, we will present numerical results for inviscid and viscous test cases. The Euler and Navier-Stokes equations are solved in the dimensional form. The discretization method is the edge-based method  with the Roe inviscid flux  and the alpha-damping viscous flux : the residual is defined at a node by

 Resj=∑k∈{kj}Φjk(njk)+sjVj, (14)

where is the measure of a dual control volume around the node , is a set of neighbor nodes of the node , is the sum of inviscid and viscous numerical fluxes, and is the directed area vector along the edge. The system of residual equations are then solved by an implicit defect-correction solver. See Ref. for further details. We store the primitive variables at nodes, reconstruct them to achieve second-order accuracy, and directly update them in a solver. For all problems, the free-stream pressure and temperature are defined by

 p∞=101325.0 [Pa] ,T∞=288.15 [K] , (15)

and the velocity will be determined for a given Mach number and an angle of attack.

The machine-zero residual estimate is computed only once before a calculation and is computed at every iteration to see how it changes. The estimates and actual residual norms will be computed in the norm for all cases (i.e., the arithmetic average of the absolute value of a residual over all nodes in a grid). The machine zero is set to be e-16 unless otherwise stated. We will show also the norm of the iterative solution error defined for each primitive variable , , at the -th nonlinear iteration by

 dw(i)=1N∑j∈{nodes}|wj(i)k−wj(i)k−1|~w(i), (16)

where

 ~w(i)=⎧⎨⎩maxj∈{nodes}|wj(i)|,if maxj∈{nodes}|wj(i)|≥% 1.0e-05,1,otherwise. (17)

This quantity is computed here only to confirm machine-zero convergence of the solution variables. We would not determine convergence based on it since a small iterative solution error does not necessarily mean convergence (e.g., it may be stalling). We wish to determine convergence by the residual norm, which tells us how well the discrete equations are satisfied.

### 4.1 Inviscid transonic flow over a Joukowsky airfoil: M∞=0.85 at an angle of attack 1.25∘

We solve the Euler equations for a transonic flow over a Joukowsky airfoil at an angle of attack with . Results are shown in Figure 1. As can be seen, the machine-zero residuals are accurately predicted by both estimates. It is seen also that the estimate does not change very much during the iteration. The iterative solution difference is also plotted. The results confirm that the solution is indeed converged to machine zero.

### 4.2 Viscous flow over a cylinder: M∞=0.001 and Re∞=20

This is a low-Mach viscous flow test case. The free-stream velocity is determined for the Mach number 0.001 and a zero angle of attack. The viscosity is a constant and determined from . To overcome well-known low-Mach problems, we employed here a preconditioned technique of Weiss and Smith  and solved an preconditioned Navier-Stokes system. Results are shown in Figure 2. As can be seen, the estimates are very accurately predicted and the iterative solution difference has reached machine-zero convergence.

### 4.3 Viscous flow over a flat plate: M∞=0.15 and Re∞=104

We consider a viscous flow over a flat plate. We again solve the preconditioned Navier-Stokes equations in the dimensional form. The velocity and the viscosity are determined from (zero angle of attack) and . The grid is a triangular grid with 274194 nodes. Results are shown in Figure 3. It can be seen that the estimates very accurately predict the4 machine-zero residuals and the iterative solution difference has reached machine-zero convergence.

For this test case, we explored also the use of larger values of to see if these estimates can predict the residual level for partially converged solutions. Here, the solver is terminated when the iterative solution difference is less than a specified value of . Results are shown in Figure 4. As expected, the solver is terminated earlier for a larger value of . The solutions are very similar as shown in Figure 5, which are sampled along a vertical grid line at . Therefore, the residuals have been reduced sufficiently for producing accurate solutions. However, the residual estimates are not very accurate as can be seen in Figures 3(b) and 3(c). In fact, we performed the same test for both the airfoil and cylinder cases and found similar results (not shown). One possible way to use them is to check first the iterative solution difference whether the solution is converged to a desired level and then use these estimates to confirm that the current residual is small enough (smaller than these estimates); if the residual is larger than the estimates, the solver may be stalling.

## 5 Conclusions

We have proposed and demonstrated two techniques to estimate the norm of machine-zero residuals for a given problem. Numerical results show that these estimates can predict machine-zero residual levels very accurately, at least for the inviscid and viscous flow problems considered. For practical applications, these estimates may be employed to check the convergence of an iterative solver. Focusing on how close the magnitude of a current residual is to a machine-zero-residual estimate (instead of how many orders of magnitude it is reduced from an initial residual), we may terminate the iteration when the residual norm is, say, five orders of magnitude larger than the machine-zero-residual estimate.

## Acknowledgments

The second author gratefully acknowledges support from Software CRADLE (part of Hexagon).

## References

•  Liu, Y. and Nishikawa, H., “Third-Order Inviscid and Second-Order Hyperbolic Navier-Stokes Solvers for Three-Dimensional Inviscid and Viscous Flows,” 46th AIAA Fluid Dynamics Conference, AIAA Paper 2016-3969, Washington, D.C., 2016.
•  Roe, P. L., “Approximate Riemann Solvers, Parameter Vectors, and Difference Schemes,” J. Comput. Phys., Vol. 43, 1981, pp. 357–372.
•  Nishikawa, H., “Beyond Interface Gradient: A General Principle for Constructing Diffusion Schemes,” Proc. of 40th AIAA Fluid Dynamics Conference and Exhibit, AIAA Paper 2010-5093, Chicago, 2010.
•  Weiss, J. M., Maruszeski, J. P., and Smith, W. A., “Implicit Solution of Preconditioned Navier-Stokes Equations Using Algebraic Multigrid,” AIAA J., Vol. 37, No. 1, 1999, pp. 29–36.