Improve (interpretation of) per-component diagnostics?

[maxentile]

Following up on #922 (review)

In principle, per-component decomposition can be misleading. Here is a contrived example:

import numpy as np
np.random.seed(0)

print(
"""
Problem:
    Interpreting per-component delta_us can be misleading,
    since per-component delta_us can be arbitrarily large
    even in a case where u(x, lam) == 0 for all x, lam.

Construction:
    u(x, lam) = u_a(x, lam) + u_b(x, lam)
    where
    u_a(x, lam) is an arbitrary function
    u_b(x, lam) = - u_a(x, lam)
"""
)

def u_a(x, lam):
    return 1_000_000 * (np.sin(1_000_000 * (x + lam)) + lam)

def u_b(x, lam):
    return - u_a(x, lam)

def u(x, lam):
    """u(x, lam) == 0 for all x and all lam"""
    return u_a(x, lam) + u_b(x, lam)

lam_i = 0.0
lam_j = 1.0

x = np.linspace(-1, 1, 1_000_000)

delta_us = u(x, lam_j) - u(x, lam_i)
print(f'overall delta_us: mean = {np.mean(delta_us)}, stddev = {np.std(delta_us)}')

delta_us_a = u_a(x, lam_j) - u_a(x, lam_i)
print(f'component a delta_us: mean = {np.mean(delta_us_a)}, stddev = {np.std(delta_us_a)}')

delta_us_b = u_b(x, lam_j) - u_b(x, lam_i)
print(f'component b delta_us: mean = {np.mean(delta_us_b)}, stddev = {np.std(delta_us_b)}')

But in practice, this kind of diagnostic is still informative.

To interpret these plots more confidently, perhaps it is sufficient to assert some further assumptions about the input (unsure what these assumptions are, maybe just that the energy components are "independent enough"?), or perhaps these plots can be augmented in some way.