Analytic gradient computation Aᵤ

[smartalecH]

Currently, we use a finite difference approximation to compute the gradient of the operator (A) w.r.t the design variables (u). This approach is convenient and flexible but has its limitations. Namely, it's rather hard to choose a sufficient finite-difference step size when β approaches ∞. (see #1854)

The proper implementation should involve the analytic "backpropagation" through the averaging routine, which involves another integral.

Not only would this approach improve the accuracy of the gradients, but it would also speed the gradient computation up by a factor of 2.

[smartalecH]

(cc @oskooi )

We also may be able to circumvent the brute force quadrature altogether (at least when β=∞). The averaging integral has a nice analytic form that can be simplified by parts.

[stevengj]

One thing we might consider is to always do the subpixel averaging as if β = infinity, regardless of the actual value of β.

Rationale:

1. The β = infinity averaging is easy. You just use the analytical formula for the volume of a spherical cap to compute the fraction of the sphere occupied by each material, and then the effective ε formula is simple (even if we do full "Kottke" averaging for the anisotropic case, which we aren't doing yet here).
2. Subpixel averaging only really matters for accuracy when the interface is discontinuous (or approximately so on the scale of a pixel), in which case the β=infinity approximation should be good. When β is small and the interface is smooth, you don't need to do any averaging at all — you get second-order accuracy even just by sampling the center of the pixel. So I'm hopeful that you would still get second-order accuracy from the fictitious "β=infinity" average in this case, though that would have to be checked more carefully.

Given the fill fraction, however, you then have to take the finite β into account when defining the materials that multiply the two fill fractions. You will want the materials that correspond to actually sampling ε(P(u)) at points in the sphere (say the center and the edge)?