@@ -208,7 +208,140 @@ class _Eval:
208208
209209
210210@needs_eos
211- @pytest .mark .smoke
211+ def test_entropy_staggered_and_temperature_staggered_accessors_for_extended_state ():
212+ """The state accessors ``entropy_staggered`` and
213+ ``temperature_staggered`` have three structurally distinct
214+ branches: gradient (covered by
215+ ``test_entropy_solver_gradient_smoke``), extended-state
216+ (energy_balance / bower2018, lines 1750-1751, 1761-1762), and
217+ plain quasi_steady (lines 1752, 1764). The latter two need
218+ explicit accessor calls.
219+
220+ Discriminator: the accessors must return the entropy / temperature
221+ block in a shape that downstream PROTEUS code (atmosphere
222+ coupling) can index. A regression that lost the extended-state
223+ slice would surface as a length-N+1 array where N+1 includes
224+ the trailing dSdr_cmb / T_core component; the slice must drop it.
225+ """
226+ from aragog .eos .entropy import EntropyEOS
227+ from aragog .solver .entropy_solver import EntropySolver
228+
229+ eos = EntropyEOS (EOS_DIR )
230+
231+ # ── Energy-balance: extended state (length N+1) ─────────────
232+ parameters = _build_params (inner_bc = 2 , core_bc = 'energy_balance' )
233+ solver = EntropySolver (parameters , entropy_eos = eos )
234+ solver .initialize ()
235+ solver .set_initial_entropy (3050.0 )
236+ solver .solve ()
237+ n_stag = solver ._n_stag
238+ s_block = solver .entropy_staggered
239+ assert s_block .shape [0 ] == n_stag , (
240+ f'extended-state entropy_staggered shape { s_block .shape } ; '
241+ f'expected first axis = n_stag = { n_stag } '
242+ )
243+ t_block = solver .temperature_staggered
244+ assert t_block .shape [0 ] == n_stag
245+
246+ # ── Quasi-steady: plain state (length N) ────────────────────
247+ parameters = _build_params (inner_bc = 2 , core_bc = 'quasi_steady' )
248+ solver = EntropySolver (parameters , entropy_eos = eos )
249+ solver .initialize ()
250+ solver .set_initial_entropy (3050.0 )
251+ solver .solve ()
252+ n_stag = solver ._n_stag
253+ s_block = solver .entropy_staggered
254+ assert s_block .shape [0 ] == n_stag
255+ t_block = solver .temperature_staggered
256+ assert t_block .shape [0 ] == n_stag
257+
258+
259+ @pytest .mark .unit
260+ def test_temperature_staggered_const_properties_branch ():
261+ """Lines 1768-1769: when ``entropy_eos is None`` (const_properties
262+ path), ``temperature_staggered`` must return the analytic
263+ ``T_ref * exp((S - S_ref) / Cp)`` form rather than calling EOS.
264+
265+ Discriminator: directly compare against the const_properties
266+ analytic formula. A regression that always called the EOS
267+ branch would surface as an AttributeError on ``None.temperature``.
268+ """
269+ from unittest .mock import MagicMock
270+
271+ import numpy as np
272+
273+ from aragog .solver .entropy_solver import EntropySolver
274+
275+ solver = EntropySolver .__new__ (EntropySolver )
276+ solver ._core_bc = 'quasi_steady'
277+ solver ._n_stag = 5
278+ solver .entropy_eos = None # const_properties
279+ pm = MagicMock ()
280+ pm .const_T_ref = 3000.0
281+ pm .const_S_ref = 3000.0
282+ pm .const_Cp = 1000.0
283+ solver .parameters = MagicMock ()
284+ solver .parameters .phase_mixed = pm
285+ fake_sol = MagicMock ()
286+ fake_sol .y = np .full ((5 , 1 ), 3050.0 )
287+ solver ._solution = fake_sol
288+
289+ T = solver .temperature_staggered
290+ expected = 3000.0 * np .exp ((3050.0 - 3000.0 ) / 1000.0 )
291+ np .testing .assert_allclose (T , expected , rtol = 1e-12 )
292+
293+
294+ @needs_eos
295+ def test_solve_with_fully_solid_initial_state_uses_relaxed_atol (caplog ):
296+ """Lines 2166-2167: when ``phi0 < 0.01`` (essentially fully
297+ solid IC) the solver relaxes ``atol_scale = 10.0`` and removes
298+ the per-step ``max_step`` cap, since there's no rheological
299+ transition to chase.
300+
301+ Discriminator: the diagnostic log line emitted by ``solve()``
302+ must include ``atol_scale=10.0x``. A regression that lost the
303+ fully-solid branch would still log ``1.0x``, locking the
304+ integrator at unnecessarily tight step bounds.
305+ """
306+ import logging
307+
308+ from aragog .eos .entropy import EntropyEOS
309+ from aragog .solver .entropy_solver import EntropySolver
310+
311+ eos = EntropyEOS (EOS_DIR )
312+ parameters = _build_params (inner_bc = 2 , core_bc = 'quasi_steady' )
313+ solver = EntropySolver (parameters , entropy_eos = eos )
314+ solver .initialize ()
315+ # Pick S well below the solidus at every node so phi0 ~ 0.
316+ P_basic = solver ._P_basic_flat
317+ S_sol_basic = eos .solidus_entropy (P_basic ).ravel ()
318+ S_init = float (S_sol_basic .min ()) - 200.0
319+ solver .set_initial_entropy (S_init )
320+
321+ with caplog .at_level (logging .INFO , logger = 'fwl.aragog.solver.entropy_solver' ):
322+ solver .solve ()
323+ # The relaxed atol_scale value must surface in the
324+ # nondimensionalisation log line.
325+ log_text = '\n ' .join (r .message for r in caplog .records )
326+ assert 'atol_scale=10.0x' in log_text , (
327+ f'expected atol_scale=10.0x in solve() log; full log:\n { log_text } '
328+ )
329+
330+ # Discriminator: get_state() with Phi_global < 0.01 must hit
331+ # the fully-solid rheological-front branch (entropy_solver.py:2903)
332+ # which sets ``rf = r_basic[-1]``. ``RF_depth`` then reduces to
333+ # near zero. A regression that mis-typed the threshold would
334+ # land in the argmin branch and produce a non-zero RF_depth.
335+ out = solver .get_state ()
336+ assert float (out .Phi_global ) < 0.05 , (
337+ f'Phi_global = { float (out .Phi_global ):.4f} not in fully-solid regime'
338+ )
339+ assert float (out .RF_depth ) < 0.01 , (
340+ f'RF_depth = { float (out .RF_depth ):.4f} ; expected near zero in fully-solid regime'
341+ )
342+
343+
344+ @needs_eos
212345def test_phi_global_falls_back_to_mean_when_mass_total_zero ():
213346 """When ``mass_total_for_phi <= 0`` the helper falls back to
214347 ``np.mean(phi_stag)`` (line 2885). This is a pathological edge
0 commit comments