Skip to main content

alice_physics/
phase_change.rs

1//! Phase Change Simulation Modifier
2//!
3//! Temperature-driven phase transitions between solid, liquid, and gas:
4//! - **Melting** (solid → liquid): Surface softens, material flows downward
5//! - **Vaporization** (liquid → gas): Material expands, then dissipates
6//! - **Solidification** (liquid → solid): Material grows rigid, surface freezes
7//! - **Condensation** (gas → liquid): Gas contracts back to liquid
8//!
9//! # Physics Model
10//!
11//! Each cell tracks: temperature, phase state, and a latent-heat buffer.
12//! Transitions use the enthalpy method (Voller & Cross 1981; Carslaw &
13// LIMITATION(COV-THERM-069): with the heat capacity normalised to 1 the latent heats are in kelvin-equivalent units
14//! Jaeger ch. XI for the Stefan problem): with the heat capacity normalised
15//! to 1 the latent heats are in kelvin-equivalent units, and each `update`
16//! moves temperature above a transition point into the buffer (or releases
17//! it back on cooling) until the buffer reaches the latent heat. The cell
18//! therefore sits **on the plateau** `T = T_m` while it melts and
19//! `T + latent` is conserved exactly by the transition step; the melt time
20//! is set by how fast diffusion / applied heat feed the cell, not by the
21//! step size. Phase transitions modify the SDF: liquid flows downward
22//! (mass-conserving transfer), gas expands outward.
23//!
24//! Author: Moroya Sakamoto
25
26use crate::coupled_field::{CoupledField, CoupledFieldError, CoupledScalar};
27use crate::math::Fix128;
28use crate::sim_field::ScalarField3D;
29use crate::sim_modifier::PhysicsModifier;
30use crate::sim_modifier::{observe_max, observe_sum, StateReader, StateWriter};
31use crate::world_participant::{
32    ObservationSink, Participant, ParticipantFault, ParticipantKind, StateError, SubstepCtx,
33};
34
35#[cfg(not(feature = "std"))]
36use alloc::vec::Vec;
37
38// ============================================================================
39// Phase State
40// ============================================================================
41
42/// Phase of matter
43#[derive(Clone, Copy, Debug, PartialEq, Eq)]
44#[repr(u8)]
45pub enum Phase {
46    /// Rigid — no SDF modification
47    Solid = 0,
48    /// Flows downward under gravity
49    Liquid = 1,
50    /// Expands outward, dissipates over time
51    Gas = 2,
52}
53
54// ============================================================================
55// Configuration
56// ============================================================================
57
58/// Phase change modifier configuration
59#[derive(Clone, Copy, Debug, PartialEq)]
60pub struct PhaseChangeConfig {
61    /// Temperature at which solid melts to liquid
62    pub melt_temperature: f32,
63    /// Temperature at which liquid vaporizes to gas
64    pub boil_temperature: f32,
65    /// Latent heat of fusion in kelvin-equivalent units (`L_f / c_p`):
66    /// the temperature rise a cell must supply, and hold at
67    /// `melt_temperature`, before it turns liquid.
68    pub latent_heat_fusion: f32,
69    /// Latent heat of vaporization in kelvin-equivalent units (`L_v / c_p`).
70    pub latent_heat_vaporization: f32,
71    /// Thermal diffusion rate
72    pub diffusion_rate: f32,
73    /// Ambient temperature
74    pub ambient_temperature: f32,
75    /// Cooling rate toward ambient
76    pub cooling_rate: f32,
77    /// Liquid flow speed (gravity-driven downward): fraction of a liquid
78    /// cell's offset moved to the cell below per second, and (× 0.1) the
79    /// rate at which liquid softens the surface (offset per second).
80    pub liquid_flow_speed: f32,
81    /// Gas expansion rate (SDF offset per second)
82    pub gas_expansion_rate: f32,
83    /// Gas dissipation rate (material removal per second)
84    pub gas_dissipation_rate: f32,
85    /// Maximum SDF offset from phase changes
86    pub max_offset: f32,
87}
88
89impl Default for PhaseChangeConfig {
90    fn default() -> Self {
91        Self {
92            melt_temperature: 200.0,
93            boil_temperature: 500.0,
94            latent_heat_fusion: 50.0,
95            latent_heat_vaporization: 100.0,
96            diffusion_rate: 0.5,
97            ambient_temperature: 20.0,
98            cooling_rate: 0.1,
99            liquid_flow_speed: 1.0,
100            gas_expansion_rate: 0.5,
101            gas_dissipation_rate: 0.3,
102            max_offset: 3.0,
103        }
104    }
105}
106
107// ============================================================================
108// Phase Change Modifier
109// ============================================================================
110
111/// Phase-change SDF modifier with solid/liquid/gas transitions
112pub struct PhaseChangeModifier {
113    /// Configuration
114    pub config: PhaseChangeConfig,
115    /// Temperature field
116    pub temperature: ScalarField3D,
117    /// Phase state per cell (stored as f32: 0=Solid, 1=Liquid, 2=Gas)
118    pub phase: ScalarField3D,
119    /// Latent heat buffer (energy absorbed/released during phase transition)
120    pub latent_heat: ScalarField3D,
121    /// Accumulated SDF offset (positive = material removed)
122    pub sdf_offset: ScalarField3D,
123    /// Whether modifier is enabled
124    pub enabled: bool,
125}
126
127impl PhaseChangeModifier {
128    /// Create a new phase change modifier
129    #[must_use]
130    pub fn new(
131        config: PhaseChangeConfig,
132        resolution: usize,
133        min: (f32, f32, f32),
134        max: (f32, f32, f32),
135    ) -> Self {
136        Self {
137            config,
138            temperature: ScalarField3D::new_filled(
139                resolution,
140                resolution,
141                resolution,
142                min,
143                max,
144                config.ambient_temperature,
145            ),
146            phase: ScalarField3D::new(resolution, resolution, resolution, min, max),
147            latent_heat: ScalarField3D::new(resolution, resolution, resolution, min, max),
148            sdf_offset: ScalarField3D::new(resolution, resolution, resolution, min, max),
149            enabled: true,
150        }
151    }
152
153    /// Apply heat at a position
154    pub fn apply_heat_at(&mut self, x: f32, y: f32, z: f32, amount: f32, radius: f32) {
155        self.temperature.splat(x, y, z, amount, radius);
156    }
157
158    /// Get temperature at a point
159    #[must_use]
160    pub fn temperature_at(&self, x: f32, y: f32, z: f32) -> f32 {
161        self.temperature.sample(x, y, z)
162    }
163
164    /// Get phase at a point
165    #[must_use]
166    pub fn phase_at(&self, x: f32, y: f32, z: f32) -> Phase {
167        let v = self.phase.sample(x, y, z);
168        if v < 0.5 {
169            Phase::Solid
170        } else if v < 1.5 {
171            Phase::Liquid
172        } else {
173            Phase::Gas
174        }
175    }
176
177    /// Process phase transitions (enthalpy method).
178    ///
179    /// The latent buffer holds the latent enthalpy stored in the cell, in
180    /// kelvin-equivalent units: `0..L_f` while a solid melts,
181    /// `L_f..L_f + L_v` while a liquid boils, `L_f + L_v` for gas. Heat above
182    /// a transition temperature is moved from `temperature` into the buffer
183    /// (melting / boiling) and back out of it on cooling (freezing /
184    /// condensing), so `temperature + latent_heat` is conserved cell by cell
185    /// and the cell sits on the plateau until the buffer fills or empties.
186    /// Independent of `dt`: the step only redistributes enthalpy already
187    /// present in the cell.
188    fn process_transitions(&mut self) {
189        let melt_t = self.config.melt_temperature;
190        let boil_t = self.config.boil_temperature;
191        let lh_fusion = crate::det_math::max(self.config.latent_heat_fusion, 0.0);
192        let lh_vaporization = crate::det_math::max(self.config.latent_heat_vaporization, 0.0);
193        let lh_total = lh_fusion + lh_vaporization;
194
195        let n = self.temperature.cell_count();
196        for i in 0..n {
197            let mut temp = self.temperature.data[i];
198            let mut phase_val = self.phase.data[i];
199            let mut lh = self.latent_heat.data[i];
200
201            // Solid → Liquid: superheat above T_m fills the fusion buffer
202            if phase_val < 0.5 && temp > melt_t {
203                let absorbed =
204                    crate::det_math::max(crate::det_math::min(temp - melt_t, lh_fusion - lh), 0.0);
205                temp -= absorbed;
206                lh += absorbed;
207                if lh >= lh_fusion {
208                    phase_val = 1.0;
209                    lh = lh_fusion;
210                }
211            }
212
213            // Liquid → Gas: superheat above T_b fills the vaporization buffer
214            if (0.5..1.5).contains(&phase_val) && temp > boil_t {
215                let absorbed =
216                    crate::det_math::max(crate::det_math::min(temp - boil_t, lh_total - lh), 0.0);
217                temp -= absorbed;
218                lh += absorbed;
219                if lh >= lh_total {
220                    phase_val = 2.0;
221                    lh = lh_total;
222                }
223            }
224
225            // Gas → Liquid: subcooling below T_b drains the vaporization buffer
226            if phase_val >= 1.5 && temp < boil_t && lh > lh_fusion {
227                let released = crate::det_math::min(boil_t - temp, lh - lh_fusion);
228                temp += released;
229                lh -= released;
230                if lh <= lh_fusion {
231                    phase_val = 1.0;
232                    lh = lh_fusion;
233                }
234            }
235
236            // Liquid: a partly-boiled cell that cools first gives back the
237            // vaporization progress, then subcooling below T_m drains the
238            // fusion buffer until the cell freezes
239            if (0.5..1.5).contains(&phase_val) {
240                if temp < boil_t && lh > lh_fusion {
241                    let released = crate::det_math::min(boil_t - temp, lh - lh_fusion);
242                    temp += released;
243                    lh -= released;
244                }
245                if temp < melt_t && lh > 0.0 {
246                    let released = crate::det_math::min(melt_t - temp, lh);
247                    temp += released;
248                    lh -= released;
249                }
250                if temp < melt_t && lh <= 0.0 {
251                    phase_val = 0.0;
252                    lh = 0.0;
253                }
254            }
255
256            // Solid: a partly-melted cell that cools gives its progress back
257            if phase_val < 0.5 && temp < melt_t && lh > 0.0 {
258                let released = crate::det_math::min(melt_t - temp, lh);
259                temp += released;
260                lh -= released;
261            }
262
263            self.temperature.data[i] = temp;
264            self.phase.data[i] = phase_val;
265            self.latent_heat.data[i] = lh;
266        }
267    }
268
269    /// Accumulate SDF offsets based on phase state
270    fn accumulate_offsets(&mut self, dt: f32) {
271        let liquid_flow = self.config.liquid_flow_speed;
272        let gas_expand = self.config.gas_expansion_rate;
273        let gas_dissipate = self.config.gas_dissipation_rate;
274        let max_offset = self.config.max_offset;
275
276        let n = self.sdf_offset.cell_count();
277        for i in 0..n {
278            let phase_val = self.phase.data[i];
279
280            if phase_val >= 1.5 {
281                // Gas: expand and dissipate (increase SDF = material removed)
282                self.sdf_offset.data[i] = crate::det_math::min(
283                    (gas_expand + gas_dissipate).mul_add(dt, self.sdf_offset.data[i]),
284                    max_offset,
285                );
286            } else if phase_val >= 0.5 {
287                // Liquid: slight material softening
288                self.sdf_offset.data[i] = crate::det_math::min(
289                    (liquid_flow * 0.1).mul_add(dt, self.sdf_offset.data[i]),
290                    max_offset,
291                );
292            }
293            // Solid: no offset change
294        }
295
296        // Liquid gravity flow: move a fraction `liquid_flow · dt` (≤ 1) of
297        // each liquid cell's offset to the cell below. The amount removed
298        // from the source equals the amount the destination accepts (its
299        // `max_offset` head-room), so the column total is conserved.
300        if liquid_flow > 0.0 {
301            let fraction = crate::det_math::min(liquid_flow * dt, 1.0);
302            let nx = self.sdf_offset.nx;
303            let ny = self.sdf_offset.ny;
304            let nz = self.sdf_offset.nz;
305            for iz in 0..nz {
306                for ix in 0..nx {
307                    for iy in 1..ny {
308                        let idx = self.phase.index(ix, iy, iz);
309                        let below_idx = self.phase.index(ix, iy - 1, iz);
310                        if !(0.5..1.5).contains(&self.phase.data[idx]) {
311                            continue;
312                        }
313                        let source = self.sdf_offset.data[idx];
314                        if source <= 0.0 {
315                            continue;
316                        }
317                        let room =
318                            crate::det_math::max(max_offset - self.sdf_offset.data[below_idx], 0.0);
319                        let transfer = crate::det_math::min(source * fraction, room);
320                        if transfer > 0.0 {
321                            self.sdf_offset.data[below_idx] += transfer;
322                            self.sdf_offset.data[idx] = source - transfer;
323                        }
324                    }
325                }
326            }
327        }
328    }
329}
330
331impl PhysicsModifier for PhaseChangeModifier {
332    #[inline]
333    fn modify_distance(&self, x: f32, y: f32, z: f32, original_dist: f32) -> f32 {
334        if !self.enabled {
335            return original_dist;
336        }
337
338        let offset = self.sdf_offset.sample(x, y, z);
339        original_dist + offset
340    }
341
342    fn update(&mut self, dt: f32) {
343        if !self.enabled {
344            return;
345        }
346
347        // 1. Diffuse temperature
348        self.temperature.diffuse(dt, self.config.diffusion_rate);
349
350        // 2. Cool toward ambient
351        self.temperature.decay_toward(
352            self.config.ambient_temperature,
353            self.config.cooling_rate,
354            dt,
355        );
356
357        // 3. Process phase transitions (enthalpy redistribution, dt-free)
358        self.process_transitions();
359
360        // 4. Accumulate SDF offsets
361        self.accumulate_offsets(dt);
362    }
363
364    fn name(&self) -> &'static str {
365        "phase_change"
366    }
367
368    fn is_active(&self) -> bool {
369        self.enabled
370    }
371}
372
373// ============================================================================
374// Coupling channel
375// ============================================================================
376
377/// The temperature field is shared through a deterministic `Fix128` channel.
378///
379/// The counterpart of the impl on [`crate::thermal::ThermalModifier`]: the two
380/// modifiers hold independent copies of the same physical quantity, and
381/// [`crate::coupled_field::reconcile_mean`] is what makes a point have one
382/// temperature instead of two.
383impl CoupledScalar for PhaseChangeModifier {
384    fn coupled_name(&self) -> &'static str {
385        "temperature"
386    }
387
388    fn coupled_channel(&self) -> Result<CoupledField, CoupledFieldError> {
389        CoupledField::try_matching(&self.temperature)
390    }
391
392    fn publish(&self, out: &mut CoupledField) -> Result<(), CoupledFieldError> {
393        out.copy_from_f32(&self.temperature)
394    }
395
396    fn adopt(&mut self, src: &CoupledField) -> Result<(), CoupledFieldError> {
397        src.write_to_f32(&mut self.temperature)
398    }
399
400    fn check_adopt(&self, src: &CoupledField) -> Result<(), CoupledFieldError> {
401        src.require_same_grid_f32(&self.temperature)
402    }
403}
404
405// ============================================================================
406// World participant
407// ============================================================================
408
409impl PhaseChangeModifier {
410    /// Snapshot tag of this type as a world participant: `"PHAS"`, the four
411    /// ASCII bytes read big endian. Never changes.
412    pub const PARTICIPANT_KIND: ParticipantKind =
413        ParticipantKind::new(u32::from_be_bytes(*b"PHAS"));
414
415    fn decode_state(bytes: &[u8]) -> Result<Self, StateError> {
416        let mut r = StateReader::new(bytes)?;
417        let config = PhaseChangeConfig {
418            melt_temperature: r.f32()?,
419            boil_temperature: r.f32()?,
420            latent_heat_fusion: r.f32()?,
421            latent_heat_vaporization: r.f32()?,
422            diffusion_rate: r.f32()?,
423            ambient_temperature: r.f32()?,
424            cooling_rate: r.f32()?,
425            liquid_flow_speed: r.f32()?,
426            gas_expansion_rate: r.f32()?,
427            gas_dissipation_rate: r.f32()?,
428            max_offset: r.f32()?,
429        };
430        let enabled = r.bool()?;
431        let temperature = r.field()?;
432        let phase = r.field()?;
433        let latent_heat = r.field()?;
434        let sdf_offset = r.field()?;
435        r.finish()?;
436        Ok(Self {
437            config,
438            temperature,
439            phase,
440            latent_heat,
441            sdf_offset,
442            enabled,
443        })
444    }
445}
446
447/// One `update(h)` per world substep (`h` converted with
448/// [`Fix128::to_f32`]); see the module documentation of
449/// [`crate::sim_modifier`] for the payload and what is not coupled yet.
450///
451/// Observations: channel 0 the highest cell temperature (absent for a grid
452/// without cells), channel 1 the number of cells that are not solid
453/// (stored phase ≥ 0.5), channel 2 the total SDF offset.
454impl Participant for PhaseChangeModifier {
455    fn kind(&self) -> ParticipantKind {
456        Self::PARTICIPANT_KIND
457    }
458
459    fn substep(&mut self, _ctx: &mut SubstepCtx<'_>, h: Fix128) -> Result<(), ParticipantFault> {
460        PhysicsModifier::update(self, h.to_f32());
461        Ok(())
462    }
463
464    fn observe(&self, out: &mut ObservationSink) {
465        observe_max(out, 0, &self.temperature);
466        let fluid = self.phase.data.iter().filter(|&&p| p >= 0.5).count();
467        out.push(1, Fix128::from_int(fluid as i64));
468        observe_sum(out, 2, &self.sdf_offset);
469    }
470
471    fn write_state(&self, out: &mut Vec<u8>) {
472        let mut w = StateWriter::new(out);
473        let c = &self.config;
474        for v in [
475            c.melt_temperature,
476            c.boil_temperature,
477            c.latent_heat_fusion,
478            c.latent_heat_vaporization,
479            c.diffusion_rate,
480            c.ambient_temperature,
481            c.cooling_rate,
482            c.liquid_flow_speed,
483            c.gas_expansion_rate,
484            c.gas_dissipation_rate,
485            c.max_offset,
486        ] {
487            w.f32(v);
488        }
489        w.bool(self.enabled);
490        w.field(&self.temperature);
491        w.field(&self.phase);
492        w.field(&self.latent_heat);
493        w.field(&self.sdf_offset);
494    }
495
496    fn check_state(&self, bytes: &[u8]) -> Result<(), StateError> {
497        Self::decode_state(bytes).map(drop)
498    }
499
500    fn read_state(&mut self, bytes: &[u8]) {
501        match Self::decode_state(bytes) {
502            Ok(state) => *self = state,
503            Err(e) => panic!("read_state called with a payload check_state refuses: {e:?}"),
504        }
505    }
506}
507
508// ============================================================================
509// Tests
510// ============================================================================
511
512#[cfg(test)]
513mod tests {
514    use super::*;
515
516    #[test]
517    fn test_phase_change_no_heat() {
518        let config = PhaseChangeConfig::default();
519        let modifier = PhaseChangeModifier::new(config, 8, (-2.0, -2.0, -2.0), (2.0, 2.0, 2.0));
520
521        // No heat: everything is solid, no offset
522        let d = modifier.modify_distance(0.0, 0.0, 0.0, -1.0);
523        assert!(
524            (d - (-1.0)).abs() < 0.01,
525            "No heat should not change SDF, got {d}"
526        );
527        assert_eq!(modifier.phase_at(0.0, 0.0, 0.0), Phase::Solid);
528    }
529
530    #[test]
531    fn test_melting_transition() {
532        let config = PhaseChangeConfig {
533            melt_temperature: 100.0,
534            latent_heat_fusion: 10.0,
535            ambient_temperature: 20.0,
536            cooling_rate: 0.0,   // No cooling
537            diffusion_rate: 0.0, // No diffusion (prevent re-solidification)
538            ..Default::default()
539        };
540        let mut modifier = PhaseChangeModifier::new(config, 8, (-2.0, -2.0, -2.0), (2.0, 2.0, 2.0));
541
542        // Apply heat well above melt temperature (large radius for low-res grid)
543        modifier.apply_heat_at(0.0, 0.0, 0.0, 2000.0, 1.5);
544
545        // Simulate enough steps for transition
546        for _ in 0..100 {
547            modifier.update(0.016);
548        }
549
550        let phase = modifier.phase_at(0.0, 0.0, 0.0);
551        assert!(
552            phase == Phase::Liquid || phase == Phase::Gas,
553            "High heat should melt material, got {phase:?}"
554        );
555    }
556
557    #[test]
558    fn test_vaporization_transition() {
559        let config = PhaseChangeConfig {
560            melt_temperature: 50.0,
561            boil_temperature: 100.0,
562            latent_heat_fusion: 1.0,
563            latent_heat_vaporization: 1.0,
564            ambient_temperature: 20.0,
565            cooling_rate: 0.0,
566            diffusion_rate: 0.0,
567            ..Default::default()
568        };
569        let mut modifier = PhaseChangeModifier::new(config, 8, (-2.0, -2.0, -2.0), (2.0, 2.0, 2.0));
570
571        // Apply extreme heat (large radius for low-res grid)
572        modifier.apply_heat_at(0.0, 0.0, 0.0, 5000.0, 1.5);
573
574        for _ in 0..200 {
575            modifier.update(0.016);
576        }
577
578        let phase = modifier.phase_at(0.0, 0.0, 0.0);
579        assert_eq!(
580            phase,
581            Phase::Gas,
582            "Extreme heat should vaporize, got {phase:?}"
583        );
584
585        // Gas should have increased SDF offset
586        let offset = modifier.sdf_offset.sample(0.0, 0.0, 0.0);
587        assert!(offset > 0.0, "Gas should add SDF offset, got {offset}");
588    }
589
590    #[test]
591    fn test_solidification() {
592        let config = PhaseChangeConfig {
593            melt_temperature: 50.0,
594            latent_heat_fusion: 1.0,
595            ambient_temperature: 20.0,
596            cooling_rate: 0.0,
597            diffusion_rate: 0.0,
598            ..Default::default()
599        };
600        let mut modifier = PhaseChangeModifier::new(config, 8, (-2.0, -2.0, -2.0), (2.0, 2.0, 2.0));
601
602        // Heat to melt (large radius for low-res grid)
603        modifier.apply_heat_at(0.0, 0.0, 0.0, 2000.0, 1.5);
604        for _ in 0..50 {
605            modifier.update(0.016);
606        }
607
608        // Should be liquid by now
609        let phase_hot = modifier.phase_at(0.0, 0.0, 0.0);
610        assert!(
611            phase_hot == Phase::Liquid || phase_hot == Phase::Gas,
612            "Should be liquid/gas after heating, got {phase_hot:?}"
613        );
614
615        // Enable fast cooling and cool down
616        modifier.config.cooling_rate = 5.0;
617        for _ in 0..500 {
618            modifier.update(0.016);
619        }
620
621        let phase_cold = modifier.phase_at(0.0, 0.0, 0.0);
622        assert_eq!(
623            phase_cold,
624            Phase::Solid,
625            "Should solidify after cooling, got {phase_cold:?}"
626        );
627    }
628
629    fn still_config() -> PhaseChangeConfig {
630        PhaseChangeConfig {
631            melt_temperature: 200.0,
632            boil_temperature: 500.0,
633            latent_heat_fusion: 50.0,
634            latent_heat_vaporization: 100.0,
635            diffusion_rate: 0.0,
636            ambient_temperature: 20.0,
637            cooling_rate: 0.0,
638            liquid_flow_speed: 0.0,
639            gas_expansion_rate: 0.5,
640            gas_dissipation_rate: 0.25,
641            max_offset: 3.0,
642        }
643    }
644
645    /// oracle: a gas cell holding the full buffer `L_f + L_v = 150` at 400
646    /// (100 below boiling) condenses: it releases `min(100, 150 − 50) = 100`,
647    /// so it ends liquid at 500 with buffer `L_f = 50`.
648    #[test]
649    fn gas_condenses_to_liquid_closed_form() {
650        let mut m = PhaseChangeModifier::new(still_config(), 1, (0.0, 0.0, 0.0), (1.0, 1.0, 1.0));
651        m.phase.data[0] = 2.0;
652        m.latent_heat.data[0] = 150.0;
653        m.temperature.data[0] = 400.0;
654        m.process_transitions();
655        assert_eq!(m.phase.data[0], 1.0);
656        assert_eq!(m.temperature.data[0], 500.0);
657        assert_eq!(m.latent_heat.data[0], 50.0);
658    }
659
660    /// oracle: a gas cell grows its offset by `(expansion + dissipation)·dt =
661    /// 0.75·(1/2) = 3/8` per call and stops at `max_offset`; a solid cell
662    /// gets nothing.
663    #[test]
664    fn gas_offset_grows_and_is_capped() {
665        let mut m = PhaseChangeModifier::new(still_config(), 2, (0.0, 0.0, 0.0), (1.0, 1.0, 1.0));
666        m.phase.data[0] = 2.0;
667        m.accumulate_offsets(0.5);
668        assert_eq!(m.sdf_offset.data[0], 0.375);
669        assert_eq!(m.sdf_offset.data[1], 0.0);
670        m.accumulate_offsets(100.0);
671        assert_eq!(m.sdf_offset.data[0], 3.0);
672    }
673
674    /// oracle: with every cell liquid, flow speed 1 and `dt = 1/2`, each cell
675    /// first gains `0.1·1·(1/2) = 1/20`; then a fraction `1·(1/2)` of the
676    /// upper cell flows down: lower `1/20 + 1/40 = 3/40`, upper `1/40`. In a
677    /// column whose lower cell is at `max_offset` there is no room, and an
678    /// upper cell with a negative offset has nothing to give.
679    #[test]
680    fn liquid_flows_down_with_room_and_source() {
681        let mut m = PhaseChangeModifier::new(
682            PhaseChangeConfig {
683                liquid_flow_speed: 1.0,
684                ..still_config()
685            },
686            2,
687            (0.0, 0.0, 0.0),
688            (1.0, 1.0, 1.0),
689        );
690        m.phase.data.fill(1.0);
691        let full_below = m.sdf_offset.index(1, 0, 0);
692        let full_above = m.sdf_offset.index(1, 1, 0);
693        m.sdf_offset.data[full_below] = 3.0;
694        let neg_below = m.sdf_offset.index(0, 0, 1);
695        let neg_above = m.sdf_offset.index(0, 1, 1);
696        m.sdf_offset.data[neg_above] = -1.0;
697        m.accumulate_offsets(0.5);
698        let below = m.sdf_offset.index(0, 0, 0);
699        let above = m.sdf_offset.index(0, 1, 0);
700        assert!((m.sdf_offset.data[below] - 0.075).abs() < 1e-6);
701        assert!((m.sdf_offset.data[above] - 0.025).abs() < 1e-6);
702        assert_eq!(m.sdf_offset.data[full_below], 3.0);
703        assert!((m.sdf_offset.data[full_above] - 0.05).abs() < 1e-6);
704        assert!((m.sdf_offset.data[neg_above] + 0.95).abs() < 1e-6);
705        assert!((m.sdf_offset.data[neg_below] - 0.05).abs() < 1e-6);
706    }
707
708    /// oracle: a disabled modifier leaves the distance and its fields as
709    /// they are.
710    #[test]
711    fn disabled_phase_change_is_inert() {
712        let mut m = PhaseChangeModifier::new(still_config(), 1, (0.0, 0.0, 0.0), (1.0, 1.0, 1.0));
713        m.sdf_offset.data[0] = 0.5;
714        m.temperature.data[0] = 900.0;
715        m.enabled = false;
716        assert!(!m.is_active());
717        assert_eq!(m.modify_distance(0.5, 0.5, 0.5, 1.0), 1.0);
718        m.update(1.0);
719        assert_eq!(m.temperature.data[0], 900.0);
720        assert_eq!(m.phase.data[0], 0.0);
721    }
722
723    fn sample_modifier() -> PhaseChangeModifier {
724        let mut m = PhaseChangeModifier::new(
725            PhaseChangeConfig {
726                melt_temperature: 150.0,
727                boil_temperature: 450.0,
728                latent_heat_fusion: 40.0,
729                latent_heat_vaporization: 80.0,
730                diffusion_rate: 0.125,
731                ambient_temperature: 22.5,
732                cooling_rate: 0.25,
733                liquid_flow_speed: 0.5,
734                gas_expansion_rate: 0.375,
735                gas_dissipation_rate: 0.625,
736                max_offset: 2.5,
737            },
738            2,
739            (-1.0, -2.0, -3.0),
740            (1.0, 2.0, 3.0),
741        );
742        for (i, v) in m.temperature.data.iter_mut().enumerate() {
743            *v = 100.0 + i as f32;
744        }
745        m.phase.data[1] = 1.0;
746        m.phase.data[2] = 2.0;
747        m.latent_heat.data[1] = 40.0;
748        m.sdf_offset.data[2] = 0.75;
749        m.enabled = false;
750        m
751    }
752
753    /// oracle: the participant kind is `"PHAS"` read big endian; a payload
754    /// written by `write_state` is accepted and read back into a different
755    /// modifier, which then has the original config, flag and the four
756    /// fields and writes the same bytes. Observations: channel 0 the hottest
757    /// cell (107), channel 1 the liquid or gas cells (2), channel 2 the total
758    /// offset (3/4).
759    #[test]
760    fn participant_state_round_trip_and_observations() {
761        assert_eq!(
762            PhaseChangeModifier::PARTICIPANT_KIND,
763            ParticipantKind::new(0x5048_4153)
764        );
765        let src = sample_modifier();
766        assert_eq!(
767            Participant::kind(&src),
768            PhaseChangeModifier::PARTICIPANT_KIND
769        );
770        let mut bytes = Vec::new();
771        src.write_state(&mut bytes);
772        assert_eq!(src.check_state(&bytes), Ok(()));
773        let mut dst = PhaseChangeModifier::new(
774            PhaseChangeConfig::default(),
775            1,
776            (0.0, 0.0, 0.0),
777            (1.0, 1.0, 1.0),
778        );
779        dst.read_state(&bytes);
780        assert_eq!(dst.config, src.config);
781        assert!(!dst.enabled);
782        assert_eq!(dst.temperature.data, src.temperature.data);
783        assert_eq!(dst.phase.data, src.phase.data);
784        assert_eq!(dst.latent_heat.data, src.latent_heat.data);
785        assert_eq!(dst.sdf_offset.data, src.sdf_offset.data);
786        let mut again = Vec::new();
787        dst.write_state(&mut again);
788        assert_eq!(again, bytes);
789
790        let mut sink = ObservationSink::new();
791        dst.observe(&mut sink);
792        assert_eq!(
793            sink.values(),
794            &[
795                (0, Fix128::from_int(107)),
796                (1, Fix128::from_int(2)),
797                (2, Fix128::from_ratio(3, 4))
798            ]
799        );
800    }
801
802    /// oracle: the enabled flag sits after the version (4) and eleven `f32`s
803    /// (44), at offset 48; a value other than 0/1 is refused, as is a
804    /// payload one byte short or long.
805    #[test]
806    fn participant_state_refusals() {
807        let src = sample_modifier();
808        let mut bytes = Vec::new();
809        src.write_state(&mut bytes);
810        assert_eq!(bytes[48], 0);
811        let mut bad = bytes.clone();
812        bad[48] = 5;
813        assert_eq!(src.check_state(&bad), Err(StateError::InvalidValue));
814        assert!(matches!(
815            src.check_state(&bytes[..bytes.len() - 1]),
816            Err(StateError::Length { .. })
817        ));
818        let mut long = bytes;
819        long.push(0);
820        assert!(matches!(
821            src.check_state(&long),
822            Err(StateError::Length { .. })
823        ));
824    }
825
826    #[test]
827    #[should_panic(expected = "read_state called with a payload check_state refuses")]
828    fn participant_read_state_panics_on_a_refused_payload() {
829        let mut m = sample_modifier();
830        m.read_state(&[1, 0, 0]);
831    }
832
833    /// oracle: in a world of one substep `dt = 1/4`, a one-cell solid at 220
834    /// (melt 200, fusion buffer 50, no diffusion or cooling) absorbs
835    /// `min(20, 50) = 20` into its buffer and stays solid at 200: channel 0
836    /// reads 200, channel 1 (fluid cells) 0, channel 2 (offset) 0.
837    #[test]
838    fn participant_substep_in_a_world() {
839        let mut m = PhaseChangeModifier::new(still_config(), 1, (0.0, 0.0, 0.0), (1.0, 1.0, 1.0));
840        m.temperature.data[0] = 220.0;
841        let mut world = crate::solver::PhysicsWorld::new(crate::solver::SolverConfig {
842            substeps: 1,
843            ..Default::default()
844        });
845        world.add_participant(Box::new(m)).expect("register");
846        world.step(Fix128::from_ratio(1, 4));
847        let Some(crate::world_participant::Observed::Exact(sink)) = world.observe_participant(0)
848        else {
849            panic!("observation");
850        };
851        assert_eq!(
852            sink.values(),
853            &[
854                (0, Fix128::from_int(200)),
855                (1, Fix128::ZERO),
856                (2, Fix128::ZERO)
857            ]
858        );
859    }
860}