[question] apparent rarely occurring infinite loop at high energies for certain CPR rootfinder parameters with the combined Kr-C-Morse potential

[drobnyjt]

There seems to be a rare event at high energies (>10 keV) for the Kr-C-Morse combined potential that causes an infinite loop. I suspect a value is going to +/-INF and is not being caught by a NaN check, but more investigation is required.

The occurrence rate is somewhere between 1/250000 ions and 1/10000 ions, but any numerical situation that causes an infinite loop should instead panic. The question is: 1) where is the issue occurring? 2) can it be ameliorated with changes to the CPR rootfinder parameters? 3) what changes should be made to the code?

Example input file that exhibits this issue:

    [options]
    name = "krc_morse_23_default_Es_1_5"
    track_recoils = false
    weak_collision_order = 0
    electronic_stopping_mode = "LOW_ENERGY_NONLOCAL"
    mean_free_path_model = "LIQUID"
    interaction_potential = [[{"KRC_MORSE"={D=6.758838e-20, r0=2.7820000000000003e-10, alpha=14197999999.999998, k=70000000000.0, x0=7.5e-11}}]]
    scattering_integral = [["GAUSS_LEGENDRE"]]
    root_finder = [[{"CPR"={n0=2, nmax=100, epsilon=1E-7, complex_threshold=1E-6, truncation_threshold=1E-9, far_from_zero=1E9, interval_limit=1E-12, derivative_free=true}}]]
    num_threads = 4
    num_chunks = 1

    [particle_parameters]
    length_unit = "ANGSTROM"
    energy_unit = "EV"
    mass_unit = "AMU"
    N = [ 10000 ]
    m = [ 1.008 ]
    Z = [ 1 ]
    E = [ 6189.65818891261 ]
    Ec = [ 0.1 ]
    Es = [ 1.5 ]
    interaction_index = [ 0 ]
    pos = [ [ -4.4, 0.0, 0.0,] ]
    dir = [ [ 0.9999999999984769, 1.7453292519934434e-6, 0.0,] ]

    [geometry_input]
    length_unit = "ANGSTROM"
    electronic_stopping_correction_factor = 1.09
    densities = [ 0.0914 ]

    [material_parameters]
    energy_unit = "EV"
    mass_unit = "AMU"
    Eb = [ 0.0 ]
    Es = [ 5.61 ]
    Ec = [ 3.0 ]
    Z = [ 28 ]
    m = [ 58.69 ]
    interaction_index = [ 0 ]
    surface_binding_model = {"PLANAR"={calculation="INDIVIDUAL"}}
    bulk_binding_model = "AVERAGE"
    

[drobnyjt]

Potential cause: if the CPR rootfinder truncation threshold is too high, there may be an infinite loop caused by truncated terms leading to identical Chebyshev polynomials each iteration and the error not decreasing as n increases, leading to infinite subdivision - although it should hit either nmax or the subdivision interval limit first. To test this, reducing those two parameters should be tried to see if it's not actually an infinite loop, just a very, very long loop caused by many large matrix inversions.

[drobnyjt]

Interestingly, I cannot reproduce this on a different machine. More investigation is needed.

[drobnyjt]

On the original machine with this issue, it appears to be limited to high energy, low nmax, and very small interval limits.