LouisDesdoigts / dLux

Differentiable optical models as parameterised neural networks in Jax using Zodiax
https://louisdesdoigts.github.io/dLux/
BSD 3-Clause "New" or "Revised" License
43 stars 6 forks source link

Fresnel Tests. #203

Closed Jordan-Dennis closed 1 year ago

Jordan-Dennis commented 1 year ago

Hi all, I will be referring back to this issue and making comments as I construct tests for the fresnel code as part of merging it into main. Regards Jordan.

Jordan-Dennis commented 1 year ago

So the Fresnel code is really out of date in terms of not using @property tags. This will need to be implemented. Additionally, I will attempt side effects in the setters (perhaps this is a bad idea), it should compile out.

Jordan-Dennis commented 1 year ago

So, we need to add the name parameter to the GaussianWavefront constructor. We also need to add error checking into the constructor.

Jordan-Dennis commented 1 year ago

I think that I have a nice way of doing fixtures. I ihave a default_gaussian_wavefront_parameters fixture and a default_gaussian_wavefront fixture. The scope of each will be module, which is safe because equinox objects are immutable. It makes sense for the default_gaussian_wavefront_parameters to be a dictionary, but a list indexed by constants is always tempting. Sigh, I'll do a dictionary.

Jordan-Dennis commented 1 year ago

Hmmmm... discovering new things as I go. Firstly, the syntax I suggested in my previous comment is not too bad, but it is very long given my choice of name (default_gaussian_wavefront_params). I think this name needs to be shortened, perhaps, (def_gauss_wave_params). Moreover, dictionaries were not a good choice in data structure since they are not immutable. This means that a tuple will probably be the best data structure to use.

Jordan-Dennis commented 1 year ago

I could use additional fixtures for non_scalar_param. I use this in two places so I think it is probably justified.

Jordan-Dennis commented 1 year ago

I need to make sure that pixel_scale is correctly accessed, since it can change. This is a problem with the way equinox deals with inheritance and immutability. Basically in dLux::Wavefront pixel_scale is a parameterso GaussianWavefront must also have it as a parameter. By accessing it through a getter method previously I was able to work around this, however, I do not want to re-implement pixel_positions in GaussianWavefront. This is perhaps silly of me, but I am not sure. Actually, I believe it is necessary.

Jordan-Dennis commented 1 year ago

The pixel scale various depending on the units of the wavefront (see #204). As a result a great many functions have an implicit conditional in them. I do not think this branch needs testing as each method is assumed to be an independent unit. However, we have the time honoured question, how does one verify correctness in a way that isn't simply checking one equals one. Some simple questions are:

Is the shape correct? Are the primary features correct? (e.g. Is most of the power in the center, is it an inncreasing function of ...)

These are the ways that I think we can actually create useful scientific tests. Of course within science there is the unanswered question on the value of unit tests. For now it is a requirement of the package, so I might as well do them well. Regards Jordan.

Jordan-Dennis commented 1 year ago

Perhaps an interesting case study is the humble pixel_coordinates function. There are four constructive ways you can test the correctness without explicitly repeating the calculation.

def test_pixel_coordinates_is_linear() -> None:
    pass

def test_pixel_coordinates_is_increasing() -> None:
    pass 

def test_pixel_coordinates_x_is_broadcast() -> None:
    pass

def test_pixel_coordinates_y_is_broadcast() -> None:
    pass

My logic is the following, if you are repeating the original calculation to make sure that you are correct then you already new you were correct. This still has value for integration testing, but that process can be fully automated. These tests are more useful in that they abstractly capture the behaviour of the result without resorting to the specifics of the implementation.

LouisDesdoigts commented 1 year ago

I would say the best way to test the Fresnel is that it matches the other propagators in the focal plane after n intermediary propagation

Jordan-Dennis commented 1 year ago

Definitely for integration tests :+1:.

Jordan-Dennis commented 1 year ago

Gradients are another great way of testing function's behaviour. For example, in quadratic_phase we can test that in log space the second difference is approximately linear, and that the first difference changes sign through the centre of the array. I'm not one-hundred percent sure that this will work, but I think it is an interesting concept.

LouisDesdoigts commented 1 year ago

Yeah we definitely have a lot of room for generating much more detailed and analytic tests. I think for the time being however just testing that stuff runs and doesn't give nans is as detailed as we need to go.

Jordan-Dennis commented 1 year ago

So I've been thinking about how to test radially symmetric functions. There is obviously the special cases along the axis, but how does one broadcast this around the circle. If we are interested in questions like is the function radially increasing, we can do this by calculating averages in circular slices. I love the way this ends up looking very similar to calculus and surface integrals.

LouisDesdoigts commented 1 year ago

Not sure about doing it with calculus but if you want to test radially increasing values you can always convert the values to radii from whatever the center is and then sort the function output by those radial values and check that the next values if >= than the last

Jordan-Dennis commented 1 year ago

Another way we can do our tests is to look for asymptotic behaviour. It is also worth working in the natural units of the problem, which for the GaussianWavefront happens to be the waist_radius and the rayleigh_distanc.

Jordan-Dennis commented 1 year ago

Not sure about doing it with calculus but if you want to test radially increasing values you can always convert the values to radii from whatever the center is and then sort the function output by those radial values and check that the next values if >= than the last

Is what I was describing. My calculus analogy comes from a rather obscure method of solving the wave equation using the mean value theorem. It works with surface integrals of arbitrary balls around a point having the same value as the point (the proof is very confusing). The idea of taking slices by the symmetry of the problem and looking at averages in that symmetry has commonality.

Jordan-Dennis commented 1 year ago

It should be possible to abstract a lot of this testing information. However, the structure of the abstraction will determine the readability of the testing code. For example,

@test_function_property(is_increasing)
def waist_radius(default_gaussian_wavefront: fixture[object]) -> None:
    return default_gaussian_wavefront.waist_radius

Although this presents its own problems. Regards Jordan.

LouisDesdoigts commented 1 year ago

Given the major re-work in 0.13, this issue should be able to be closed as the tests will need to be built from the ground up.