KineticForces - BUG - nutype="zero" crashes with BoundsError when the resonant integrand is non-finite and no pole was subtracted

Aperta
#428 0 commenti 0 reazioni 0 assegnatari Vedi su GitHub

I maintainer di solito rispondono entro 3 giorni

Nessuno ha ancora preso questa issue.

Valutazione

Difficoltà
4/5
Tempo stimato
3-5 giorni
Idoneità per principianti
48/100
Tipo di issue
Bug
Chiarezza
Abbastanza chiara
Stato di attività
Tranquilla
Stack tecnologico
julia
Ambito
backend

Direzione di ricerca

Read src/KineticForces/EnergyIntegration.jl, especially _resonant_integrand and the npole == 0 integration branch. Trace how non-finite values are handled when no poles were found, and clarify the intended physical value for that case before changing the fallback. Then add a regression case using a DIII-D-like equilibrium with nutype = "zero"; verify the run completes while the existing solovev_kinetic_nuzero case still passes.

Scritto dal modello di indicizzazione a partire dal testo della issue.

Descrizione

Symptom

A DIII-D-like kinetic-calculated run with [KineticForces] nutype = "zero" (collisionless) at nominal rotation dies during kinetic matrix formation:

ERROR: TaskFailedException
  nested task error: BoundsError: attempt to access 0-element Vector{ComplexF64} at index [1]
   [3] ... src/KineticForces/EnergyIntegration.jl:227 [inlined]
   [4] ... src/KineticForces/EnergyIntegration.jl:268
   [5] QuadGK evalrule ...

Reproduced on both a plain run and a variant; the same deck with nutype = "harmonic" runs fine, and nutype = "zero" plus toroidal_rotation_factor = 0.2 also runs fine — so it needs the collisionless operator at nominal rotation.

Mechanism

_resonant_integrand (EnergyIntegration.jl:~215-234) has a non-finite-value fallback that substitutes the nearest pole's analytic regular part:

if !isfinite(val)
    k = 1
    for j in 2:npole
        abs(x - real(x_poles[j])) < abs(x - real(x_poles[k])) && (k = j)
    end
    val = _real_pole_regular_part(real(x_poles[k]), ...)   # x_poles[1] on an empty vector

The fallback assumes at least one pole was found. But the caller builds the same closure for the npole == 0 branch (:~268) and integrates it:

integrand = x -> _resonant_integrand(x, p, residues, x_poles, npole, leff, wb, n, wd)
if npole == 0
    val, _ = quadgk(integrand, 0.0, X_ENERGY_MAX; ...)

With nutype = "zero" there is no collisional broadening, so Ω can round to exactly 0 at a quadrature node even when no pole was registered in [0, X_ENERGY_MAX] — val becomes non-finite, the fallback runs with npole = 0, and x_poles[1] throws.

Suggested fix

Guard the fallback on npole > 0, and decide what the physical value is when the integrand is non-finite with no subtracted pole (returning 0 for that node, or excluding it, is a modelling choice worth stating in the annotation rather than implying). A regression case exercising nutype = "zero" on a DIII-D-like equilibrium would cover this — the existing solovev_kinetic_nuzero case passes, so Solovev alone does not reach the branch.

Not introduced by the issue #376 kinetic stack: EnergyIntegration.jl has no commits on that branch relative to develop.

Lingua principale
Julia
Stelle
11
Fork
7
Merge medio
3g 22h
PR unite (30g)
70

Preparare l'ambiente

  • Nessun Dockerfile né file Docker Compose
  • Ha un modello di pull request
  • Nessuna guida per i contributori

Come iniziare

  1. Leggi tutta la issue e poi la guida ai contributi del progetto.
  2. Commenta sulla issue per dire che te ne occupi tu — evita che due persone facciano lo stesso lavoro.
  3. Fai un fork del repository e lavora su un branch.
  4. Apri una pull request che faccia riferimento al numero della issue.

Altre issue di OpenFUSIONToolkit/GPEC

Tutte le issue di OpenFUSIONToolkit/GPEC

Issue simili

Altre issue su Julia

Ricevi le nuove issue nella tua casella

Un breve riepilogo di issue GitHub adatte ai principianti.