Drop GSL dependency

[krystophny]

Addresses #134.

KAMEL is MIT-licensed but linked GSL (GPL-2.0+) via cmake/FetchGSL.cmake. This PR removes GSL entirely: every call site now uses a Fortran intrinsic, the libneo math kit, or a small in-tree replacement. No shipped target links GSL.

Depends on the libneo math kit: itpplasma/libneo#290 must merge first, since KIM now uses neo_bessel_i from LIBNEO::math. Draft until that lands; see also the golden-record note below.

Changes

KIM (Fortran)

• gsl_sf_erf bindings replaced by the Fortran 2008 erf intrinsic (integrands.f90, integrands_adaptive.f90, numerics_utils.f90).
• gsl_sf_bessel_In replaced by bessel_in from libneo neo_bessel_i (FLR2_asymptotics.f90); KIM_lib links LIBNEO::math.
• gsl_mod.f90 deleted.

KiLCA / QL-Balance (C++), new modules under KiLCA/math/numeric/:

• adaptive_quad: adaptive Gauss-Kronrod quadrature with QUADPACK QAG semantics (GK15/21/31). Nodes, weights, and error tests are the public-domain Netlib QUADPACK data (dqk15/21/31, dqage), with attribution in the source. Replaces gsl_integration_qag in four_transf.cpp, hyper1F1.cpp, and the calc_I_array_drift test driver.
• quadpack_qagi: infinite-range integration as a faithful translation of the public-domain QUADPACK dqagie, dqk15i, dqelg, dqpsrt, including the epsilon-algorithm extrapolation. Replaces gsl_integration_qagi/qagiu in vel_integral.cpp. Matches gsl_integration_qagi to 2.6e-16 worst-case over 216 parameter sets spanning the resonant vel_integral regime (omega_E through zero, small collisionality).
• ode_rk8pd: embedded Runge-Kutta 8(7) with the RK8(7)13M tableau of Prince and Dormand (J. Comp. Appl. Math. 7, 1981) and the gsl_odeiv step-control semantics: scale eps_abs + eps_rel*|y|, safety 0.9, growth limits 0.2 to 5, truncated final-step suggestion, acceptance when a rejected step cannot make progress. Replaces gsl_odeiv rk8pd/control/evolve in calc_back.cpp, incompressible.cpp, compressible_flow.cpp. Single steps reproduce GSL bitwise in 98 percent of probes; the rest differ by 1 to 2 ulp from compiler-specific FMA contraction.
• brent_root: Brent-Dekker bracketing solver (Brent 1973). Replaces the gsl_root_fsolver brent loop in calc_mode.cpp.
• deriv_central: adaptive 5-point central difference. Replaces gsl_deriv_central.
• sort_index_doubles (std::sort over an index array, in core/shared.cpp) replaces gsl_heapsort_index in adaptive_grid.cpp, calc_cond.cpp, sysmat_profs.cpp.
• The eigenmode determinant zero search (calc_eigmode.cpp) uses a plain 2D Newton iteration with the existing finite-difference Jacobian instead of gsl_multiroot_fdfsolver_hybridsj. Stopping criteria unchanged (same residual and step-size tests).
• Dead GSL-only code deleted: the Levin-acceleration variant in hyper1F1.cpp (marked "do not use" upstream), the never-called calc_circle_integral_gkq, and the unused analytic Jacobians of the explicit rk8pd stepper.
• All replacements keep state on the stack or in caller-owned workspaces; no globals, safe under OpenMP.

Build system

• cmake/FetchGSL.cmake deleted, include(FetchGSL) removed from cmake/Dependencies.cmake.
• gsl/gslcblas removed from KiLCA, KIM, and QL-Balance link lists, include paths, setup scripts, and READMEs.

Verification

All on macOS arm64, Homebrew clang 21 / gfortran, GSL 2.8 as test oracle only.

Before (main): GSL fetched and linked. After (this branch): full build succeeds, no GSL symbol in any shipped binary.

$ for b in KIM.x KiLCA_Normal_... ql-balance.x fouriermodes.x; do otool -L install/bin/$b | grep -ci gsl; done
0 0 0 0

ctest: 17/18 pass. The one failure (test_rhs_balance, gamma_a_e expected -0.0 vs 2.78e-14*1e13) fails identically on unmodified main built with the same toolchain; it is fp-association noise in a pure-Fortran flux test untouched by this PR.

main:    17/18 passed, test_rhs_balance FAILED (gamma_a_e mismatch 0.2775...)
branch:  17/18 passed, test_rhs_balance FAILED (gamma_a_e mismatch 0.2775...)

Oracle tests against GSL at the real call-site parameters (gate 1e-12 rel):

OK qag GK31 real/imag part          rel=0.0
OK qag GK21 tight tol               rel=1.3e-15
OK qagi vel_integral ind=1..3       rel<=1.6e-16
OK qagiu smooth                     rel=3.1e-16
OK brent root                       rel=0.0 (matches GSL and exact root)
OK deriv_central (5 points)         rel=0.0
OK rk8pd eps=1e-16 / 1e-8           rel<=6.7e-15, accepted steps equal (76/76, 15/15)
qagi stress: 216 cases, 0 mismatches, worst rel diff 2.6e-16
worst rel diff erf intrinsic vs gsl_sf_erf:    1.7e-16
worst rel diff bessel_in vs gsl_sf_bessel_In:  9.2e-16

Golden record (ql-balance.x vs golden.h5 from main+GSL on the same machine; on this platform the run produces LinearProfiles 0-1 and KinProfiles 1000 on main and branch alike, so 33 of 114 quantities are comparable):

init_params (7), T_EM_phi_e/i, KinProfiles/1000 (4): PASS
dql*, Br* at steps 0,1 (18): FAIL at rtol=1e-8, max rel 1.6e-6 (dql), 4.4e-5 (Br_Im)
main vs main rerun (same binary):                    0 value failures

The 18 failures are a one-time transition artifact, bisected to calc_back.cpp: the background equilibrium ODE integrates at eps_abs = eps_rel = 1e-16, where the step controller operates at the roundoff floor and accepts or rejects steps on the last ulp of the error estimate. The resulting u profile differs from the GSL build by at most 2.9e-15 relative (measured by instrumenting both builds), and the resonant layer of the linear mode solve amplifies that to the observed 1e-6-level diffs in Br and dql. At this tolerance even two GSL builds with different compilers produce different step sequences, so no independent implementation can hold rtol 1e-8 on this path; once this PR is on main, the golden reference regenerates GSL-free and the test is consistent again. All other golden quantities, and the full pipeline up to the linear profiles, match the GSL build to better than 1e-8.

[krystophny]

Superseded by #139, which drops GSL by routing all of KAMEL (KIM Fortran, KiLCA/QL-Balance C/C++) onto the fortnum numerical core directly rather than reimplementing the routines in-tree and pulling Bessel from libneo. #139 builds GSL-free with 18/18 ctests passing and zero gsl_ symbols in every shipped binary. Leaving this open for a human to decide on closing.

[krystophny]

Superseded by #139, which drops GSL by depending on fortnum (lazy-fortran/fortnum) directly instead of routing through libneo math-kit. #139 builds green with zero gsl symbols in every binary and ctest 18/18. Closing in favor of #139.