[antoinepado]

@farhi @willend @ebknudsen

Two main fixes of the bending magnet are proposed. The first concerns the flux and the second the sampling probability densities.

Flux fix

I first started working on the flux fix 2.5 years ago. For archival purposes, work can be found here.
The work presented in this pull request continues this previous work mainly in regards to the sampling fix.

The current McXtrace component uses the spectral flux per $$mrad^2$$ solid angle and per $$0.1%$$ photon energy bandwidth. It is based on equation (28) from B.D. Patterson, Am. J. Phys. 79, 1046 (2011)

Remark : The equation from the paper has a typo and shows $$10^{33}$$ instead of $$10^{13}$$.

This is the spectral flux density in the horizontal plane. This approximation holds when rays from the bending magnet are close to 0 vertically.

Rather than the above, the full angular distribution of synchrotron radiation is used in the fixed component. See equation (1) from Kim in the X-Ray Data Booklet, Section 2.1 Characteristics of Synchrotron Radiation.

Notice that equation (28) is simply is the full angular distribution of synchrotron radiation equation (1) with the vertical angle $$\psi=0$$.

I show below one way I justify this change.
An example of a real source can be found in a paper from Brunelle, Journal of Synchrotron Radiation, 2023 showing the new superbend of SOLEIL's ROCK beamline. The paper describes a slit of aperture $$13mm \times 2.4mm$$ placed 8.7m away from the source.

A simulation of ROCK's bending magnet is done with the synchrotron source software SPECTRA and with McXtrace. This comparison is done in order to validate McXtrace's bending magnet.
The following graph shows the flux in units of photons/sec/0.1%bw.
A simulation with and without the slit is run for three different cases:

• old bending magnet
• fixed bending magnet
• spectra's bending magnet

Flux (photons/sec/0.1%bandwidth), normalized (to see the shape at first):

The graph below shows the same as the above with the flux in units of photons/sec , normalized (to see the shape at first).

In the graphs above we can see SPECTRA and McXtrace correspond for the new fixed bending magnet. The old bending magnet does not correspond to SPECTRA. Additionally, with or without the slit the flux does not change at all for the old bending magnet (the yellow and purple graphs are the same).

Below is the graph showing SPECTRA and the new fixed bending magnet in absolute units instead of normalized:

They correspond in absolute terms also.

Sampling probability densities fix

The fixed component also properly takes into account the sampling probability densities for both the energy and angles (vertical and horizontal).

The correct sampling densities are used whether the sampling is done uniformly or in a truncated gaussian way.

Here is a graph showing the flux (photons/sec/0.1%bandwidth) in absolute value of the new fixed bending magnet for the three different sampling methods (no focus, focus gaussian, focus uniform).

With slit:

Without slit:

The three different sampling methods show the same flux which is to be expected.

McXtrace instruments and simulation data, SPECTRA's configuration files and data, and the plotting script can be found in this uploaded folder in this pull request.
Please let me know if you agree with my analysis. If I have misunderstood something or made mistakes please let me know.
I can go into more detail if needed.
Thank you.
tests_bending_magnet_fixed.zip

[ebknudsen]

Excellent work if I may say so @antoinepado . This looks really compelling.

Haven't gone through the math so can't comment on that, but two practical things.

1. The tests failing and being off by a large amount is not surprising as you definitely change the absolute value of whatever intensity is emitted. Meaning the tests should be updated. Especially if they can be verified to match with SPECTRA.
2. If possible, it'd be useful if we we could avoid all these whitespace edits in a commit. It just makes it more difficult to see what changes were actually made. If style changes are to be made, it's best to keep them in separate commits, which just say formatting or something. This also makes the case for creating a specific mccode clang-format file that could be used automatically.

more to come...

[farhi]

Hi @antoinepado

Very nice work !! especially the comparison with Spectra.

The failing tests still show very large changes

Running tests / getting status...
Czerny_Turner         :   Display OK (0s)
Czerny_Turner         :   3.03    [val: 0.0739757 / 1.10881e-06 = 6671630 %] <--- BIG DISCREPANCY??
Czerny_Turner_2       :  58.45    [scan: 21/21 points within 20%, worst 101 %, 21/21 points run]
SOLEIL_DIFFABS        :   Display OK (0s)
SOLEIL_DIFFABS        :   NO TEST
SOLEIL_DISCO          :   Display OK (1s)
SOLEIL_DISCO          :   1.29    [val: 1.44686e+16 / 208292000000.0 = 6946306 %] <--- BIG DISCREPANCY??
SOLEIL_DISCO_2        :  21.19    [scan: 19/21 points within 20%, worst 0 %, 21/21 points run] <--- BIG DISCREPANCY??
SOLEIL_MARS           :   Display OK (0s)
SOLEIL_MARS           :  15.18    [val: 130261000.0 / 2560000000.0 = 5 %] <--- BIG DISCREPANCY??
SOLEIL_ROCK           :   Display OK (1s)
SOLEIL_ROCK           :   2.66    [val: 140982000000.0 / 2000000000.0 = 7049 %] <--- BIG DISCREPANCY??
SOLEIL_SAMBA          :   Display OK (0s)
SOLEIL_SAMBA          :   1.15    [val: 52266000000.0 / 12396800000.0 = 422 %] <--- BIG DISCREPANCY??
SOLEIL_SIRIUS         :   Display OK (0s)
SOLEIL_SIRIUS         :   1.30    [val: 7399170000.0 / 21452200000.0 = 34 %] <--- BIG DISCREPANCY??
Test_grating_reflect  :   Display OK (0s)
Test_grating_reflect  :   1.01    [val: 2289840000.0 / 3793780.0 = 60358 %] <--- BIG DISCREPANCY??
Test_BM               :   Display OK (0s)
Test_BM               :   0.97    [val: 920888000000000.0 / 29313800000000.0 = 3141 %] <--- BIG DISCREPANCY??
Test_BM_2             :   0.97    [val: 1288060000000000.0 / 1285490000000000.0 = 100 %]

The flux can change dramatically compared to the current reference values, and the BL. This is expected, however I would recommend a careful check of these, esp. ROCK, Test_BM_1, DISCO_1, to understand for each the origin of the change.

[willend]

@antoinepado agreed, really nice work. :-)

And for @ebknudsen:

This also makes the case for creating a specific mccode clang-format file that could be used automatically.
This sort of already exists in ./devel/bin/mccode-clangformat

[ebknudsen]

@antoinepado agreed, really nice work. :-)

And for @ebknudsen:

This also makes the case for creating a specific mccode clang-format file that could be used automatically.
This sort of already exists in ./devel/bin/mccode-clangformat
@willend Ah - goody! Apologies - I had not noticed that.

[farhi]

I just pushed some fixes to this PR. Mostly the 'lambda' branch is not initialized properly. I removed it to avoid undefined 'qe'.

[antoinepado]

@farhi @willend @ebknudsen Thank you.
I will look at the failing tests and update them.
I'll make sure to keep the same formatting in the future.
Thank you for the fixes @farhi .

[farhi]

Just put back lambda mode, and add checks for E0-dE>0.

[ebknudsen]

No biggie certainly and not even sure it did before - but a note in the header about lambda-mode is perhaps warranted wrt that sampling is uniform in the corresponding energy-interval, not in lambda as such.

[willend]

Thanks @antoinepado, really nice work. I've rerun every affected %Example with this PR's Bending_magnet.comp (mccode-dev, macOS arm64, CPU only). The new references are the mean of 5 seeds at 1e7 rays, pushed in 1c7352f and refined in 18107d4 (8 instr files plus their READMEs). They agree with the CI numbers above within statistics.

Updated %Example values

Instrument	Old	New (5-seed mean, 1e7)	seed spread
Czerny_Turner w_monitor_I	1.10881e-06	0.0721	2.2%
Test_grating_reflect psd_monitor_I	3.79378e+06	3.53e+09	3.7%
SOLEIL_MARS MARS_detector_diffraction_I	2.56e+09	1.341e+08	0.2%
SOLEIL_ROCK ROCK_fluo_monitor_I	2e+09	1.61e+11	18%
SOLEIL_SAMBA SAMBA_det_fluo_I	1.23968e+10	5.69e+10	6.9%
SOLEIL_SIRIUS SIRIUS_def_loc_I	2.14522e+10	6.69e+09	12%
SOLEIL_DISCO DISCO_second_HFP_I	2.08292e+11	1.445e+16	<0.1%
Test_BM SOURCE=0 emon_I	2.93138e+13	9.212e+14	<0.1%
Test_BM SOURCE=1 emon_I	1.28549e+15	unchanged

Test_BM SOURCE=1 is the Wiggler, so it's unaffected. The %Scan references are unchanged too: I reran the DISCO scan locally and all 21 points match within ~1%. SOLEIL_DIFFABS has no %Example.

For ROCK, SAMBA and SIRIUS, the scatter between seeds is several times larger than the error each run reports. At the 1e6 rays the tests use, ROCK, Test_grating_reflect and SIRIUS scatter by 16–25%, so they may still occasionally fail the 20% tolerance.

Review comments

1. E0 vs lambda0: the header says E0 overrides lambda0, but the code overwrites E0/dE whenever lambda0/dlambda are set. Either guard it with if (!E0 && lambda0 > 0 ...) or fix the doc.
2. dE=0 now aborts. Since the flux is integrated over the band, that's physically reasonable, but please say so in %D. Monochromatic use was allowed before.
3. Output unit: please state in %D that p is now photons/s into the sampled solid angle and energy band. Before, it was effectively per mrad² per 0.1% bw. That change affects how every downstream monitor should be read.
4. No-focus mode samples a hard-coded ±2.5 mrad in both planes, so the total (e.g. Test_BM SOURCE=0) is the flux into a 5 mrad horizontal fan. Please document this. A Gaussian with σ = 2.5 mrad cut at ±1σ is nearly flat (0.61 at the edge), so uniform sampling would be just as efficient and simpler:

   theta = randpm1()*angle; p *= 2*angle;

   That lets you drop cnstc, qtheta/qpsi and both do/while loops. The same holds for gauss_t with s1x = focus_xw/2.
5. Leftovers: oneoverk, F1, l, the commented alpha/MELE/s1x/s1y lines, the redundant focus_xw!=0 &&, and the stale comments on s1y and qe ("no lambda handling here").
6. Formatting: the whitespace changes make the diff harder to read (git diff -w shows the actual change), but that's fine for this PR. For future changes, running devel/bin/mccode-clangformat or putting formatting in its own commit would help reviewers.

I checked the flux expression against Kim's eq. (1), including the angular Jacobian dist/(dist²+dx²), the 2dE/E bandwidth factor and the truncated-Gaussian inverse densities, and it all looks right to me.

Not tested: GPU, Linux, Windows.

[farhi]

amazing !
Thanks.

[willend]

Update: we've gone ahead and addressed review points 1–5 above in fa469fa, so no action is needed on those, @antoinepado.

1. E0 vs lambda0: E0 now overrides lambda0 as documented; lambda mode is only used when E0 is unset. The error message now spells out what is required.
2. dE>0: documented in %D.
3. Output unit: %D now says p is photons/s into the sampled solid angle and the band [E0-dE, E0+dE]. It also cites Kim's eq. (1), with Patterson's eq. (28) as its ψ=0 limit.
4. No-focus mode: θ and ψ are now sampled uniformly in ±2.5 mrad (p *= 2*angle) instead of with the 1σ-truncated Gaussian, and the ±2.5 mrad / 5 mrad fan is documented. The expected values are the same, but the RNG is used differently, so seeded results move within statistics.
5. Leftovers: removed (oneoverk, F1, l, cnstc, qtheta/qpsi, the commented-out code, the redundant check, and the stale comments).

Separately, 8645e0e raises ncount for two noisy tests. The SOLEIL_DISCO %Scan now uses -n1e6: at the default of 1e5 per point, two points got no rays on w_monitor, so Wavelength came out as 0. Test_grating_reflect now uses -n 1e7, because it scatters by about 25% between seeds at 1e6.

A CI-style local check (seed 1000, 2 MPI ranks, macOS arm64) after these changes:

• all %Example tests are within 14% of the references
• both %Scan tests are within 1.2% at all 21 points
• Test_BM SOURCE=1 and MARS are unchanged

Not tested: GPU, Linux, Windows (CI will cover the last two).

[antoinepado]

@willend Thank you very much.