Is variable scaling essential when solving some PDE problems numerically?

In semiconductor simulation, it is common that the equations are scaled so they have normalised values. For example, in extreme cases electron density in semiconductors can vary over 18 order of magnitude, and electric field can change shapely, over 6 (or more) orders of magnitude.

However, the papers never really give a reason for doing this. Personally I am happy dealing with equations in real units, is there any numerical advantage to do this, is it impossible otherwise? I thought with double precision there would be enough digits to cope with these fluctuations.

Both answers are very useful, thanks very much!

• "can vary over 18 orders of magnitude" - and if you consider how many digits are retained in double precision, you'll see if "with double precision there would be enough digits to cope with these fluctuations" is indeed true... – J. M. May 25 '13 at 16:59
• And the real problem starts when you feed these numbers into a numerical algorithm: Take the square, and suddenly you have 36 orders of magnitude difference... – Christian Clason May 25 '13 at 17:07

Solving a (linear) PDE consists in discretizing the equation to yield a linear system, which is then solved by a linear solver whose convergence (rate) depends on the condition number of the matrix. Scaling the variables often reduces this condition number, thus improving convergence. (This basically amounts to applying a diagonal preconditioner, see Nicholas Higham's Accuracy and Stability of Numerical Algorithms.)

Solving nonlinear PDEs in addition requires a method for solving nonlinear equations such as Newton's method, where the scaling can also influence convergence.

Since normalizing everything usually takes very little effort, it is almost always a good idea.

• I'm sure @ArnoldNeumaier has more to say on this topic. – Christian Clason May 25 '13 at 17:26
• The condition number of the matrices I am using (unscaled variables) is ~1.25. Does this seem reasonable? This is calculated using the 2-norm method (docs.scipy.org/doc/numpy/reference/generated/…). – boyfarrell May 26 '13 at 2:20
• @boy, if there's nothing else multiplying the 1.25, then that's an okay condition number, especially if you recall that unitary matrices have $\kappa_2=1$. – J. M. May 26 '13 at 5:00
• @boyfarrell: I routinely work with condition numbers as large as 10^7, with acceptable results. I wouldn't accept condition numbers much higher than 10^9, however. – jvriesem Aug 14 '15 at 15:07

PDEs of which the solutions have sharp boundaries pose problems that go beyond being able to represent the solution in floating point. This is especially true when solutions have a certain physical meaning, e.g., a density (that per se cannnot be smaller than 0). Consider, for example, $$-\varepsilon \Delta u + u = 0 \text{ on } \Omega,\\ u = 1 \text{ on } \partial\Omega.$$ The exact solution is pointwise positive, but when discretizing with Finite Elements, the may no longer be true. A similar question came up recently.

That said, there is no scaling of variables or domains that removes this difficulty.

One sensible reason for scaling out physical quantities is to have the actual equation as simple as possible, i.e., as few parameters as possible involved (that you'd have to play around with). For example, instead of looking for the solution $u_{\alpha}$ of $$-\alpha^2\Delta u = f_\alpha \text{ on } \alpha\Omega$$ (with some boundary conditions) for many different $\alpha$, you could just go ahead and find the solution $u_1$ of $$-\Delta u = f \text{ on } \Omega.$$ You deduce immediately that $u_{\alpha}(\mathbf{x}):=u_1(\mathbf{x}/\alpha)$ for any $\alpha$. This shows that scaling out the parameter $\alpha$ doesn't really influence the behavior of the solution at all -- a fact that is quite obvious here, but may be quite difficult to see in other cases. You can construct quite complex examples with more parameters in the same way, e.g., Navier-Stokes and its nondimensionalization into the Reynolds-number-formulation.

• And the remaining parameters must then be essential for determining the qualitative behavior of the solution - this is why the Reynolds number is so important in fluid dynamics. This process is called Nondimensionalization. – Christian Clason May 25 '13 at 17:37
• of course, finding such parameter equivalences is essentially the problem of finding symmetry groups of the PDE, a problem which is hard in general – lurscher May 26 '13 at 5:04

Dealing with floating point numbers can be trick with regards to subtraction of small numbers from larger numbers, as well as with many other aspects. I would recommend reading John D. Cooks blog posts on them, such as

as well as Oracle's

What Every Computer Scientist Should Know About Floating-Point Arithematic

Also certain numerical algorithms for minimization, or maximization require normalization for numerical stability.