zfp & fpzip: Floating Point Compression

zfp and Derivatives

Many analysis tasks require computing spatial derivatives or other derived fields, e.g. to estimate vorticity in a turbulence simulation, or to perform topological Morse segmentation from gradients.  When lossy compression is used to store the original field, such derived field computations may amplify any compression-induced errors, resulting in visual artifacts or numerical inaccuracies that were not readily apparent in the original field.  Spatial derivatives are analogous to edge enhancement in image processing, and tend to reveal artifacts in lossily compressed fields.

The images below visualize the divergence of a velocity field (i.e. the sum of three partial derivatives) estimated using central differencing.  The divergence field itself is not compressed; the visible artifacts are due to using lossy compression to represent the original velocity field.  In spite of zfp using the least storage among all methods, its divergence field exhibits no artifacts and is visually identical to the uncompressed field.

Higher-order derivatives tend to be even more sensitive to compression-induced errors.  Below are comparisons of the Laplacian fxx + fyy + fzz (sum of second partial derivatives) of an analytical function f : [-1, 1]3R, given by

f(x, y, z) = c (x4 + y4 + z4) - x y z - log(1 + x2 + y2 + z2)

where c = 104 / 75.  Consequently, the Laplacian is given by

fxx + fyy + fzz = 2 [6 c r2 - (3 + r2) / (1 + r2)2]

where r2 = x2 + y2 + z2.  The function f was chosen so that the contours of its Laplacian are concentric spheres.

As in the case of the divergence computation above, the Laplacian is estimated using second-order central differences from the reconstructed field f that has undergone lossy compression.  The bottom row below shows the zero level sets for f and its Laplacian, the latter which is a sphere of radius r = 1/2.  The dark bands are highlight lines that emphasize the quality of the surface normal, which here depends on the third partial derivatives of f.  As is evident, zfp provides the highest quality in the Laplacian of all compressors while using the least amount of storage.  Indeed, using less than 6 bits/value, zfp achieves higher quality than 32-bit IEEE single precision floating point.

Derivatives estimated using finite differences are subject to two types of error: truncation error, which has to do with the finite support of derivative stencils (i.e. truncation of the local Taylor expansion), and roundoff error, which is due to finite precision.  When estimating derivatives on a uniform grid, the truncation error decreases with the grid spacing.  The roundoff error, on the other hand, increases with decreasing grid spacing.  This is because differences between nearby grid samples become smaller in relation to the function values themselves, and cannot be resolved well using finite precision arithmetic.  For a given precision and fine enough grid, consecutive values become identical, falsely suggesting that most derivatives are zero.

The plot above shows the RMS error in partial derivatives estimated using sixth-order central differences.  The 2D function being differentiated is sinusoidal and has analytic derivatives.  The kinks in the IEEE floating-point curves are points where truncation and roundoff errors are balanced.  Where this occurs depends on the numerical precision available.

Paradoxically, derivatives computed from zfp’s representation are often more accurate than the same derivatives computed from IEEE floating point, even when zfp uses fewer bits of storage.  The zfp error curves decrease monotonically with finer grid spacing, and eventually match the IEEE double-precision curve.  The reason for this is that zfp is far less sensitive to roundoff error.  As the function is sampled more and more finely, it also locally behaves more “smoothly,” which makes it easier to compress.  zfp exploits that when nearby samples share many leading bits, these redundant bits need not all be stored.  This frees up space for storing trailing mantissa bits instead, which improve accuracy and reduce roundoff error.  zfp uses 62 mantissa bits internally compared to 24 respectively 53 for IEEE single and double precision, and can therefore produce more accurate derivatives (and function values) than double precision.  Because zfp’s decompressor returns IEEE double-precision values from which the derivatives are estimated, it is limited by the accuracy afforded by double precision.  However, more accurate derivatives could in principle be computed directly from zfp’s compressed representation.