"""Exact-volume response, independent refinement, and retained fit failures. An empirical hull of angular/finite-difference estimates is numerical evidence, not a formal certificate. Fit-versus-exact disagreement is recorded separately and prevents a categorical reader label near the approximate threshold. """ import numpy as np from scipy.optimize import brentq from independent_volume import independent_log_derivative, independent_radius from orbit_reference import constant_gamma_closed from roche_reference import eggleton, radius_log_derivative, volume_response def orbital_response(q, retained, gamma): return 2*(retained*q-1)+(1-retained)*(2*gamma+1)*q/(1+q) def response_estimates(q): estimates = { "implicit_x_axis_128": volume_response(q, 128)[1], "implicit_x_axis_256": volume_response(q, 256)[1], "independent_z_axis_128_logstep_003": independent_log_derivative(q, 128, .003), "independent_z_axis_256_logstep_003": independent_log_derivative(q, 256, .003), "independent_z_axis_256_logstep_0015": independent_log_derivative(q, 256, .0015), } values = list(estimates.values()) fitted = float(radius_log_derivative(q)) return {"q": q, "exact_volume_log_response_estimates": estimates, "empirical_exact_response_hull": [min(values), max(values)], "fitted_radius_log_response": fitted, "fit_minus_implicit_response": fitted-estimates["implicit_x_axis_256"], "fit_relative_derivative_difference": fitted/estimates["implicit_x_axis_256"]-1, "numerical_hull_width": max(values)-min(values)} def feedback_label(record, stellar_response, retained=1., gamma=2.): q = record["q"] base = orbital_response(q, retained, gamma) multiplier = 1+retained*q exact_low, exact_high = [base+multiplier*value for value in record["empirical_exact_response_hull"]] fitted = base+multiplier*record["fitted_radius_log_response"] reader_low, reader_high = min(exact_low, fitted), max(exact_high, fitted) numerical_label = ("local-detachment" if stellar_response > exact_high else "growing-overfill" if stellar_response < exact_low else "numerically-unresolved") label = ("near-threshold-sensitive-to-radius-approximation" if reader_low <= stellar_response <= reader_high else numerical_label) return {"stellar_response": stellar_response, "retained_fraction": retained, "gamma_J": gamma, "exact_lobe_response_hull": [exact_low, exact_high], "fitted_lobe_response": fitted, "exact_numerical_feedback": numerical_label, "reader_label": label, "signed_margin_from_fit_and_numerical_hull": stellar_response-reader_high if stellar_response > reader_high else stellar_response-reader_low if stellar_response < reader_low else 0.} def finite_removal_case(q, response, removed): md, ma = q-removed, 1+removed new_q = md/ma distance = float(constant_gamma_closed(q, 1, 1, removed, 1, 2)) radius_ratio = (md/q)**response actual = radius_ratio*independent_radius(q, 384)/(distance*independent_radius(new_q, 384))-1 fitted = radius_ratio*eggleton(q)/(distance*eggleton(new_q))-1 labels = feedback_label(response_estimates(q), response) assert actual > 0 and fitted < 0 assert labels["reader_label"] == "near-threshold-sensitive-to-radius-approximation" return {"q": q, "retained": 1, "stellar_response": response, "removed_mass": removed, "independent384_overfill_relative_change": float(actual), "fitted_overfill_relative_change": float(fitted), "feedback": labels, "old_fit_classification_failure_retained": True} def run_stability_reference(): # These are explicit sampled q values, not a certified continuous error bound. sampled = sorted(set([float(f"{q:.12g}") for q in np.geomspace(.05, 20, 13)] + [.1, .3, .4, .6323329047973826, 1., 1.5, 3., 10.])) records = [response_estimates(q) for q in sampled] assert max(row["numerical_hull_width"] for row in records) < 2e-7 implicit_root = brentq(lambda q: orbital_response(q, 1, 2)+(1+q)*volume_response(q, 256)[1]+1/3, .55, .75, xtol=2e-12) independent_root = brentq(lambda q: orbital_response(q, 1, 2)+(1+q)*independent_log_derivative(q, 256, .002)+1/3, .55, .75, xtol=2e-9) fitted_root = brentq(lambda q: orbital_response(q, 1, 2)+(1+q)*radius_log_derivative(q)+1/3, .55, .75, xtol=2e-12) assert abs(implicit_root-independent_root) < 2e-8 failures = [finite_removal_case(1., .460, .0001), finite_removal_case(.6323329047973826, -1/3, .6323329047973826*.0001)] presets = [] for name, q, adiabatic, thermal in [("lighter-donor-detaches", .4, -1/3, 0.), ("heavier-convective-donor", 1.5, -1/3, 0.), ("adiabatic-detachment-thermal-instability", 1.5, 4., .8)]: estimates = response_estimates(q) ad, th = feedback_label(estimates, adiabatic), feedback_label(estimates, thermal) assert "near-threshold" not in ad["reader_label"] and "near-threshold" not in th["reader_label"] presets.append({"preset": name, "q": q, "adiabatic": ad, "thermal": th}) return {"status": "PASS exact-volume response with retained approximation failures", "sampled_mass_ratio_responses": records, "largest_sampled_numerical_hull_width": max(row["numerical_hull_width"] for row in records), "largest_sampled_relative_fit_derivative_bias": max(abs(row["fit_relative_derivative_difference"]) for row in records), "retained_sign_flip_counterexamples": failures, "conditional_core_free_polytrope_qcrit": {"zeta_ad": -1/3, "implicit_volume_256": implicit_root, "independent_volume_256_finite_difference": independent_root, "eggleton_fit": fitted_root, "not_universal_stellar_threshold": True}, "proposed_three_transfer_presets": presets, "ce_preset_stability_label": "Not a local transfer criterion; separate alpha endpoint", "label_rule": "Use exact-volume response estimates. If the stellar response lies within the empirical numerical/fit threshold hull, refuse a categorical reader label and expose the sensitivity. No arbitrary error inflation or formal certificate.", "continuous_control_boundary": "Sampled comparisons do not certify interpolation over all q; compute/refine each selected q or restrict labels to verified presets."}