Numerical gradient for (interferenced?) intensity

[MauriceZeuner]

Is there a package or algorithm you could recommend for calculating numerical gradients of a function that takes the form:

f(pos: Vector3<f64>) -> f64

Strictly speaking, I need it to get the intensity gradient (E-field-squared-gradient) for cases of many-beam interference, where an analytical solution is hard to obtain (Although I will think about it some more... for some special cases (like linear polarized Gaussian beams only).

I've seen https://docs.rs/GSL/4.0.0/rgsl/numerical_differentiation/index.html but I don't know how to "install" it.

Obviously, I could just write a 2,3 or 5 -point algorithm myself, and implement it in maths - but I was wondering if there existed an optimized algorithm for that already.

[ElliotB256]

The question is quite general at the moment. Can you describe the specific case with more detail, as there may be another solution?

For instance, if this is for forces due to interfering beams - can you add the gradient produced by each beam separately? if V_i(x) are summed to give a total V(x) = \sum_i V_i(x), then doesnt d V(x)/ dx = \sum_i d V_i(x)/dx? For most beam shapes it should be possible to determine an analytic form of the gradient (or at least one that can be precomputed)

[MauriceZeuner]

Yes I can: The problem is different to what you described as I don't need the gradient of a summed quantity but the gradient of a complex modulus squared of summed complex quantities. I need: (with c being a few constants that are not important right now)

\grad I = \grad (c * (\sum E_i )^2)

...where the ^2 denotes the modulus squared of the argument arg, i.e. arg.re^2+arg.im^2

For the full Gaussian beam it is of, course still possible to get an analytic form - however we have to find an analytic form that can scale with the number of beams and respects their possibly different polarizations. I'd say this is "hard" to obtain, not impossible but I believe that the resulting analytic form could - in the end - be similarly expensive (in terms of performance) and it's probably an absolute nightmare to maintain the code. I can try to spend some time to derive something for... N gaussian beams so we can actually judge that.

[MauriceZeuner]

If it was linear - well, there would be no point in doing this at all - since interference (like in optical lattices) is a non-linear effect.

[ElliotB256]

Is the solution then to work in the domain of the E-fields (like you previously said), and have an intermediate system which converts E->I?

[ElliotB256]

(helps you for the magnitude but not the gradient I suppose)

[MauriceZeuner]

yes, there will be something like a SampleElectricFieldSystem (where the interference happens) and then another ConvertElectricFieldToIntensityGradientSystem which basically takes the complex modulus and multiplies constants and then computes the gradient of that.

At least that is currently my idea - I am open to better suggestions