From d3c8f21023cae166f9e590ada34a9fb3ded53116 Mon Sep 17 00:00:00 2001 From: Kyle Gorkowski Date: Wed, 4 Mar 2026 06:26:54 -0700 Subject: [PATCH 1/3] feat(condensation): implement latent heat step and tests Mirror CondensationIsothermal.step for latent-heat strategies, update particle and gas states, and record per-step latent-heat energy with last_latent_heat_energy. Forward dynamic_viscosity to mass_transfer_rate and add overloads for legacy and data inputs. Add latent-heat step tests covering mass conservation, energy tracking, sign conventions, isothermal parity, gas updates, dynamic viscosity pass-through, and MassCondensation compatibility. Closes #1141 ADW-ID: 2260fbb3 --- .../condensation/condensation_strategies.py | 134 +++++- .../tests/condensation_strategies_test.py | 390 +++++++++++++++++- 2 files changed, 510 insertions(+), 14 deletions(-) diff --git a/particula/dynamics/condensation/condensation_strategies.py b/particula/dynamics/condensation/condensation_strategies.py index 387f53164..98bf0775d 100644 --- a/particula/dynamics/condensation/condensation_strategies.py +++ b/particula/dynamics/condensation/condensation_strategies.py @@ -33,6 +33,7 @@ from particula.dynamics.condensation.mass_transfer import ( get_first_order_mass_transport_k, + get_latent_heat_energy_released, get_mass_transfer, get_mass_transfer_rate, get_mass_transfer_rate_latent_heat, @@ -2036,6 +2037,47 @@ def rate( rates = self._apply_skip_partitioning(rates) return rates + def _update_latent_heat_energy( + self, mass_transfer: NDArray[np.float64], temperature: float + ) -> None: + """Store latent heat energy released for the current step.""" + if self._latent_heat_strategy is None: + self.last_latent_heat_energy = 0.0 + return + + latent_heat = self._latent_heat_strategy.latent_heat(temperature) + self.last_latent_heat_energy = float( + np.sum( + get_latent_heat_energy_released( + mass_transfer, + latent_heat, + ) + ) + ) + + # pylint: disable=too-many-positional-arguments, too-many-arguments + @overload # type: ignore[override] + def step( + self, + particle: ParticleRepresentation, + gas_species: GasSpecies, + temperature: float, + pressure: float, + time_step: float, + dynamic_viscosity: Optional[float] = None, + ) -> Tuple[ParticleRepresentation, GasSpecies]: ... + + @overload + def step( + self, + particle: ParticleData, + gas_species: GasData, + temperature: float, + pressure: float, + time_step: float, + dynamic_viscosity: Optional[float] = None, + ) -> Tuple[ParticleData, GasData]: ... + def step( self, particle: ParticleRepresentation | ParticleData, @@ -2047,7 +2089,13 @@ def step( ) -> ( Tuple[ParticleRepresentation, GasSpecies] | Tuple[ParticleData, GasData] ): - """Advance one condensation step (stub). + """Advance one condensation step with latent-heat diagnostics. + + The mass transfer rate is computed, optional skip-partitioning applied, + and both the particle and gas states are updated while respecting + inventory limits. The per-step latent-heat energy release is stored in + ``last_latent_heat_energy`` (positive for condensation, negative for + evaporation) and overwritten on every call. Args: particle: Particle representation to update. @@ -2059,11 +2107,83 @@ def step( Returns: Tuple containing updated particle and gas species objects. - - Raises: - NotImplementedError: Step is not implemented for latent-heat - condensation. """ - raise NotImplementedError( - "CondensationLatentHeat.step is not implemented." + particle_data, particle_is_legacy = _unwrap_particle(particle) + gas_data, gas_is_legacy = _unwrap_gas(gas_species) + _require_matching_types(particle_is_legacy, gas_is_legacy) + _require_single_box(particle_data.n_boxes, "ParticleData") + _require_single_box(gas_data.n_boxes, "GasData") + + mass_rate = self.mass_transfer_rate( + particle=particle, + gas_species=gas_species, + temperature=temperature, + pressure=pressure, + dynamic_viscosity=dynamic_viscosity, + ) + + if isinstance(mass_rate, (int, float)): + mass_rate_array = np.array([mass_rate]) + else: + mass_rate_array = mass_rate + + mass_rate_array = self._apply_skip_partitioning(mass_rate_array) + + gas_mass_array: NDArray[np.float64] = np.atleast_1d( + np.asarray(gas_data.concentration[0], dtype=np.float64) + ) + volume = particle_data.volume[0] + norm_conc = particle_data.concentration[0] / volume + mass_transfer = get_mass_transfer( + mass_rate=mass_rate_array, + time_step=time_step, + gas_mass=gas_mass_array, + particle_mass=particle_data.masses[0], + particle_concentration=norm_conc, ) + + self._update_latent_heat_energy(mass_transfer, temperature) + + species_mass = particle_data.masses[0] + if species_mass.ndim == 1: + species_count = 1 + else: + species_count = species_mass.shape[1] + if mass_transfer.ndim == 1: + mass_transfer = mass_transfer.reshape(-1, species_count) + elif mass_transfer.shape[1] > species_count: + mass_transfer = mass_transfer[:, :species_count] + elif mass_transfer.shape[1] < species_count: + mass_transfer = np.broadcast_to( + mass_transfer, (mass_transfer.shape[0], species_count) + ) + + if particle_is_legacy: + particle.add_mass(added_mass=mass_transfer) # type: ignore[union-attr] + else: + per_particle = np.divide( + mass_transfer, + norm_conc[:, np.newaxis] + if mass_transfer.ndim == 2 + else norm_conc, + out=np.zeros_like(mass_transfer), + where=norm_conc[:, np.newaxis] != 0 + if mass_transfer.ndim == 2 + else norm_conc != 0, + ) + particle_data.masses[0] = np.maximum( + particle_data.masses[0] + per_particle, 0.0 + ) + if self.update_gases: + if gas_is_legacy: + gas_species.add_concentration( # type: ignore[union-attr] + added_concentration=-mass_transfer.sum(axis=0) + ) + else: + gas_data.concentration[0] -= mass_transfer.sum(axis=0) + if particle_is_legacy: + return cast( + Tuple[ParticleRepresentation, GasSpecies], + (particle, gas_species), + ) + return cast(Tuple[ParticleData, GasData], (particle_data, gas_data)) diff --git a/particula/dynamics/condensation/tests/condensation_strategies_test.py b/particula/dynamics/condensation/tests/condensation_strategies_test.py index 0384c5764..33b46990f 100644 --- a/particula/dynamics/condensation/tests/condensation_strategies_test.py +++ b/particula/dynamics/condensation/tests/condensation_strategies_test.py @@ -22,6 +22,7 @@ _unwrap_gas, _unwrap_particle, ) +from particula.dynamics.condensation.mass_transfer import get_mass_transfer from particula.gas.gas_data import GasData, from_species from particula.gas.latent_heat_strategies import ConstantLatentHeat from particula.gas.species import GasSpecies @@ -1914,13 +1915,6 @@ def test_base_class_params_forwarded(self): self.assertEqual(cond.accommodation_coefficient, 0.7) self.assertFalse(cond.update_gases) - def test_stubs_raise_not_implemented(self): - """step() remains stubbed for latent heat strategy.""" - cond = CondensationLatentHeat(molar_mass=0.018) - - with self.assertRaises(NotImplementedError): - cond.step(None, None, 298.15, 101325, time_step=1.0) - def test_mass_transfer_rate_isothermal_parity(self): """Latent heat strategy falls back to isothermal rates.""" latent_cond = CondensationLatentHeat( @@ -2229,3 +2223,385 @@ def test_inherited_methods_work(self): self.assertGreater(mean_free_path, 0.0) self.assertGreater(knudsen_number, 0.0) + + def test_step_mass_conservation_single_species(self): + """Single-species step conserves total mass.""" + gas_species, particle, _, _ = self._build_single_species_fixture() + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + ) + + initial_mass = ( + particle.get_mass_concentration() + + gas_species.get_concentration().sum() + ) + particle_new, gas_new = cond.step( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + final_mass = ( + particle_new.get_mass_concentration() + + gas_new.get_concentration().sum() + ) + + relative_error = abs(final_mass - initial_mass) / initial_mass + self.assertLess(relative_error, 1e-12) + + def test_step_mass_conservation_multi_species(self): + """Multi-species step conserves total mass.""" + particle = copy.deepcopy(self.particle) + gas_species = copy.deepcopy(self.gas_species) + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + ) + + initial_mass = ( + particle.get_mass_concentration() + + gas_species.get_concentration().sum() + ) + particle_new, gas_new = cond.step( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + final_mass = ( + particle_new.get_mass_concentration() + + gas_new.get_concentration().sum() + ) + + relative_error = abs(final_mass - initial_mass) / initial_mass + self.assertLess(relative_error, 1e-12) + + def test_step_mass_conservation_data_only(self): + """Data-only step conserves total mass.""" + particle_data, gas_data = self._make_data_inputs() + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + activity_strategy=self.activity_strategy, + surface_strategy=self.surface_strategy, + vapor_pressure_strategy=self.vapor_pressure_strategy, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + ) + + initial_mass = np.sum( + particle_data.masses[0] + * particle_data.concentration[0][:, np.newaxis] + / particle_data.volume[0] + ) + np.sum(gas_data.concentration[0]) + particle_new, gas_new = cond.step( + particle=particle_data, + gas_species=gas_data, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + final_mass = np.sum( + particle_new.masses[0] + * particle_new.concentration[0][:, np.newaxis] + / particle_new.volume[0] + ) + np.sum(gas_new.concentration[0]) + + relative_error = abs(final_mass - initial_mass) / initial_mass + self.assertLess(relative_error, 1e-12) + + def test_step_energy_tracking_matches_dm_times_l(self): + """Latent heat energy matches dm × L for a step.""" + particle = copy.deepcopy(self.particle) + gas_species = copy.deepcopy(self.gas_species) + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + ) + + particle_data, _ = _unwrap_particle(particle) + gas_data, _ = _unwrap_gas(gas_species) + mass_rate = cond.mass_transfer_rate( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + ) + mass_rate_array = ( + np.array([mass_rate]) + if isinstance(mass_rate, (int, float)) + else mass_rate + ) + mass_rate_array = cond._apply_skip_partitioning(mass_rate_array) + norm_conc = particle_data.concentration[0] / particle_data.volume[0] + mass_transfer = get_mass_transfer( + mass_rate=mass_rate_array, + time_step=1.0, + gas_mass=np.asarray(gas_data.concentration[0], dtype=np.float64), + particle_mass=particle_data.masses[0], + particle_concentration=norm_conc, + ) + + cond.step( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + latent_heat = cond._latent_heat_strategy.latent_heat(self.temperature) + + np.testing.assert_allclose( + cond.last_latent_heat_energy, + np.sum(mass_transfer * latent_heat), + rtol=1e-14, + ) + + def test_step_energy_sign_convention(self): + """Condensation is positive energy; evaporation is negative.""" + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + ) + gas_species, particle, _, _ = self._build_single_species_fixture() + + cond.step( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + self.assertGreater(cond.last_latent_heat_energy, 0.0) + + vp_water = par.gas.VaporPressureFactory().get_strategy("water_buck") + mm_water = 18.015e-3 + gas_sub = ( + par.gas.GasSpeciesBuilder() + .set_name("H2O") + .set_molar_mass(np.array([mm_water]), "kg/mol") + .set_vapor_pressure_strategy([vp_water]) + .set_concentration(np.array([1e-30]), "kg/m^3") + .set_partitioning(True) + .build() + ) + particle_evap = copy.deepcopy(particle) + cond.step( + particle=particle_evap, + gas_species=gas_sub, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + self.assertLess(cond.last_latent_heat_energy, 0.0) + + def test_step_isothermal_parity(self): + """Latent heat fallback matches isothermal step outputs.""" + particle_latent = copy.deepcopy(self.particle) + gas_latent = copy.deepcopy(self.gas_species) + particle_iso = copy.deepcopy(self.particle) + gas_iso = copy.deepcopy(self.gas_species) + latent_cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + ) + iso_cond = CondensationIsothermal( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + ) + + particle_out, gas_out = latent_cond.step( + particle=particle_latent, + gas_species=gas_latent, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + particle_iso_out, gas_iso_out = iso_cond.step( + particle=particle_iso, + gas_species=gas_iso, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + + np.testing.assert_allclose( + particle_out.get_mass(), particle_iso_out.get_mass(), rtol=1e-15 + ) + np.testing.assert_allclose( + gas_out.get_concentration(), + gas_iso_out.get_concentration(), + rtol=1e-15, + ) + self.assertEqual(latent_cond.last_latent_heat_energy, 0.0) + + def test_step_returns_correct_types_legacy(self): + """Legacy inputs return legacy types.""" + cond = CondensationLatentHeat(molar_mass=self.molar_mass) + + particle_new, gas_new = cond.step( + particle=copy.deepcopy(self.particle), + gas_species=copy.deepcopy(self.gas_species), + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + + self.assertIsInstance(particle_new, ParticleRepresentation) + self.assertIsInstance(gas_new, GasSpecies) + + def test_step_returns_correct_types_data(self): + """Data-only inputs return data types.""" + particle_data, gas_data = self._make_data_inputs() + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + activity_strategy=self.activity_strategy, + surface_strategy=self.surface_strategy, + vapor_pressure_strategy=self.vapor_pressure_strategy, + ) + + particle_new, gas_new = cond.step( + particle=particle_data, + gas_species=gas_data, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + + self.assertIsInstance(particle_new, ParticleData) + self.assertIsInstance(gas_new, GasData) + + def test_step_updates_gas_when_enabled(self): + """Gas concentrations decrease by mass transfer when enabled.""" + particle = copy.deepcopy(self.particle) + gas_species = copy.deepcopy(self.gas_species) + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + update_gases=True, + ) + + particle_data, _ = _unwrap_particle(particle) + gas_data, _ = _unwrap_gas(gas_species) + mass_rate = cond.mass_transfer_rate( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + ) + mass_rate_array = ( + np.array([mass_rate]) + if isinstance(mass_rate, (int, float)) + else mass_rate + ) + mass_rate_array = cond._apply_skip_partitioning(mass_rate_array) + norm_conc = particle_data.concentration[0] / particle_data.volume[0] + expected_mass_transfer = get_mass_transfer( + mass_rate=mass_rate_array, + time_step=1.0, + gas_mass=np.asarray(gas_data.concentration[0], dtype=np.float64), + particle_mass=particle_data.masses[0], + particle_concentration=norm_conc, + ) + expected_gas = ( + gas_species.get_concentration() - expected_mass_transfer.sum(axis=0) + ) + + _, gas_new = cond.step( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + + np.testing.assert_allclose( + gas_new.get_concentration(), expected_gas, rtol=1e-15 + ) + + def test_step_dynamic_viscosity_pass_through(self): + """Dynamic viscosity override is forwarded to mass_transfer_rate.""" + particle = copy.deepcopy(self.particle) + gas_species = copy.deepcopy(self.gas_species) + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + ) + override = 1.23e-5 + + with patch.object( + cond, "mass_transfer_rate", wraps=cond.mass_transfer_rate + ) as spy: + cond.step( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + dynamic_viscosity=override, + ) + + self.assertEqual(spy.call_args.kwargs["dynamic_viscosity"], override) + + def test_step_mass_condensation_runnable_compatibility(self): + """MassCondensation executes with latent heat strategies.""" + particle = copy.deepcopy(self.particle) + gas_species = copy.deepcopy(self.gas_species) + atmosphere = par.gas.Atmosphere( + temperature=self.temperature, + total_pressure=self.pressure, + partitioning_species=gas_species, + gas_only_species=[], + ) + aerosol = par.Aerosol(atmosphere=atmosphere, particles=particle) + + cond_latent = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + ) + condensation = par.dynamics.MassCondensation( + condensation_strategy=cond_latent + ) + aerosol_out = condensation.execute(aerosol, time_step=1.0) + + self.assertIsInstance(aerosol_out.particles, ParticleRepresentation) + self.assertIsInstance( + aerosol_out.atmosphere.partitioning_species, GasSpecies + ) + + aerosol_iso = par.Aerosol( + atmosphere=copy.deepcopy(atmosphere), + particles=copy.deepcopy(self.particle), + ) + cond_iso = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + ) + condensation_iso = par.dynamics.MassCondensation( + condensation_strategy=cond_iso + ) + aerosol_iso_out = condensation_iso.execute(aerosol_iso, time_step=1.0) + + self.assertIsInstance(aerosol_iso_out.particles, ParticleRepresentation) + self.assertIsInstance( + aerosol_iso_out.atmosphere.partitioning_species, GasSpecies + ) From d3b8fd035e9063e46c7dd9db85cc1b43b48c3ea4 Mon Sep 17 00:00:00 2001 From: Kyle Gorkowski Date: Wed, 4 Mar 2026 07:00:10 -0700 Subject: [PATCH 2/3] docs(condensation): document latent heat step updates Update feature docs, plans, and README to reflect the implemented CondensationLatentHeat step, energy diagnostics, and viscosity override notes. Expand the strategy system guide with latent heat details and example usage. Closes #1141 ADW-ID: 2260fbb3 --- adw-docs/dev-plans/README.md | 2 +- .../epics/E5-non-isothermal-condensation.md | 3 +- ...E5-F3-condensation-latent-heat-strategy.md | 5 ++- adw-docs/dev-plans/features/index.md | 2 +- docs/Features/condensation_strategy_system.md | 45 ++++++++++++++++++- docs/index.md | 4 +- .../condensation/condensation_strategies.py | 16 +++++-- readme.md | 3 +- 8 files changed, 67 insertions(+), 13 deletions(-) diff --git a/adw-docs/dev-plans/README.md b/adw-docs/dev-plans/README.md index 969118530..77a2ea268 100644 --- a/adw-docs/dev-plans/README.md +++ b/adw-docs/dev-plans/README.md @@ -77,7 +77,7 @@ and rollout. - Scope: Thermal resistance factor and non-isothermal mass transfer rate pure functions with energy tracking. - [E5-F3: CondensationLatentHeat Strategy Class][e5-f3] — Status: In Progress - (P1, #1139) + (P3, #1141) - Scope: New condensation strategy with latent heat correction and energy diagnostics. - [E5-F4: Builder, Factory, and Exports][e5-f4] — Status: Planning diff --git a/adw-docs/dev-plans/epics/E5-non-isothermal-condensation.md b/adw-docs/dev-plans/epics/E5-non-isothermal-condensation.md index 71dc38072..7b3bd769b 100644 --- a/adw-docs/dev-plans/epics/E5-non-isothermal-condensation.md +++ b/adw-docs/dev-plans/epics/E5-non-isothermal-condensation.md @@ -284,7 +284,7 @@ L -> 0 as T -> T_c. Used in engineering thermodynamics and EOS-based models. relative), reduced rate when L > 0 (thermal resistance always slows condensation), array shapes for single and multi-species -- [ ] **E5-F3-P3**: Implement `step()` for particle-resolved with energy +- [x] **E5-F3-P3**: Implement `step()` for particle-resolved with energy tracking and tests - `step()` follows `CondensationIsothermal.step()` flow (lines 922-1040) with two additions: @@ -623,3 +623,4 @@ class CondensationLatentHeat(CondensationStrategy): | 2026-03-02 | Initial epic creation | ADW | | 2026-03-02 | Split E5-F1-P3 into P3 (builders) + P4 (factory+exports); split E5-F3-P3 into P3 (particle-resolved step) + P4 (discrete+continuous) + P5 (data-only parity); added missing details: function signatures, file references, thermal conductivity source, vapor_pressure_surface parameter, test tolerances, literature targets | ADW | | 2026-03-04 | Noted E5-F3-P1 issue #1139 and logging expectations | ADW | +| 2026-03-04 | Marked E5-F3-P3 complete for issue #1141 | ADW | diff --git a/adw-docs/dev-plans/features/E5-F3-condensation-latent-heat-strategy.md b/adw-docs/dev-plans/features/E5-F3-condensation-latent-heat-strategy.md index 97bf96004..c98bae8d3 100644 --- a/adw-docs/dev-plans/features/E5-F3-condensation-latent-heat-strategy.md +++ b/adw-docs/dev-plans/features/E5-F3-condensation-latent-heat-strategy.md @@ -146,9 +146,9 @@ additions: - Tests: numerical parity with `CondensationIsothermal` when L=0 (< 1e-15 relative), reduced rate when L > 0, array shapes for single/multi-species -- [ ] **E5-F3-P3**: Implement `step()` for particle-resolved with energy +- [x] **E5-F3-P3**: Implement `step()` for particle-resolved with energy tracking and tests - - Issue: TBD | Size: M (~100 LOC) | Status: Not Started + - Issue: #1141 | Size: M (~100 LOC) | Status: Complete - `step()` follows `CondensationIsothermal.step()` flow (lines 922-1040) with two additions: 1. Uses `mass_transfer_rate` from P2 (includes thermal correction) @@ -241,3 +241,4 @@ additions: |------|--------|--------| | 2026-03-02 | Initial feature document created from E5 epic | ADW | | 2026-03-04 | Marked P1 in progress for issue #1139 | ADW | +| 2026-03-04 | Completed P3 step() implementation and tests (issue #1141) | ADW | diff --git a/adw-docs/dev-plans/features/index.md b/adw-docs/dev-plans/features/index.md index 2151d86ce..88d2749bb 100644 --- a/adw-docs/dev-plans/features/index.md +++ b/adw-docs/dev-plans/features/index.md @@ -28,7 +28,7 @@ work, typically ~100 LOC per phase, that deliver user-facing functionality. |----|------|--------|----------|--------| | E5-F1 | [Latent Heat Strategy Pattern](E5-F1-latent-heat-strategy.md) | In Progress | P1 | 4 | | E5-F2 | [Non-Isothermal Mass Transfer Functions](E5-F2-non-isothermal-mass-transfer.md) | Planning | P1 | 3 | -| E5-F3 | [CondensationLatentHeat Strategy Class](E5-F3-condensation-latent-heat-strategy.md) | In Progress | P1 | 5 | +| E5-F3 | [CondensationLatentHeat Strategy Class](E5-F3-condensation-latent-heat-strategy.md) | In Progress (P3) | P1 | 5 | | E5-F4 | [Builder, Factory, and Exports](E5-F4-builder-factory-exports.md) | Planning | P1 | 2 | | E5-F5 | [Validation and Integration Tests](E5-F5-validation-integration-tests.md) | Planning | P1 | 2 | | E5-F6 | [Documentation and Examples](E5-F6-documentation-examples.md) | Planning | P2 | 3 | diff --git a/docs/Features/condensation_strategy_system.md b/docs/Features/condensation_strategy_system.md index c2e24b88e..3f5067d50 100644 --- a/docs/Features/condensation_strategy_system.md +++ b/docs/Features/condensation_strategy_system.md @@ -12,12 +12,20 @@ into runnable pipelines. You can choose simultaneous isothermal updates or staggered two-pass Gauss-Seidel sweeps with theta-controlled first-pass fractions, batch partitioning, and gas-field updates for stability. +CondensationLatentHeat mirrors the isothermal workflow for particle-resolved +runs while adding latent-heat-aware mass transfer and per-step energy +diagnostics. The step now supports a `dynamic_viscosity` override and records +`last_latent_heat_energy` for parity with the isothermal API plus energy +tracking. + This feature is built around user-facing APIs exposed via `particula.dynamics`: - `CondensationStrategy` – abstract base defining `rate` / `step`. - `CondensationIsothermal` – simultaneous isothermal mass transfer. - `CondensationIsothermalStaggered` – two-pass staggered (Gauss-Seidel) update with theta modes (`"half"`, `"random"`, `"batch"`) and batching. +- `CondensationLatentHeat` – latent-heat-corrected rate with per-step energy + diagnostics. - `CondensationIsothermalBuilder`, `CondensationIsothermalStaggeredBuilder` – fluent builders with validation and unit handling. - `CondensationFactory` – factory selecting a condensation strategy by name. @@ -38,6 +46,8 @@ This feature is built around user-facing APIs exposed via `particula.dynamics`: fields are depleted. - **Pipeline-ready**: Use `MassCondensation` with `sub_steps` for tight coupling to other runnables in a single pipeline. +- **Latent heat diagnostics**: Track per-step energy release for + particle-resolved runs with `CondensationLatentHeat`. ## Who It's For @@ -65,6 +75,7 @@ par.dynamics.CondensationStrategy # Concrete implementations par.dynamics.CondensationIsothermal par.dynamics.CondensationIsothermalStaggered +par.dynamics.CondensationLatentHeat # Builders and factory par.dynamics.CondensationIsothermalBuilder @@ -131,6 +142,36 @@ particle, gas = iso.step( ) ``` +### CondensationLatentHeat (energy diagnostics) + +`CondensationLatentHeat` mirrors the isothermal step but applies a latent-heat +correction when a latent heat strategy (or scalar fallback) is provided. It +records `last_latent_heat_energy` each step (positive for condensation, +negative for evaporation) and accepts a `dynamic_viscosity` override for +particle-resolved workflows. + +```python +latent = par.dynamics.CondensationLatentHeat( + molar_mass=0.018, + diffusion_coefficient=2e-5, + accommodation_coefficient=1.0, + latent_heat=2.4e6, # J/kg fallback +) +particle, gas = latent.step( + particle=particle, + gas_species=gas, + temperature=298.15, + pressure=101325.0, + time_step=1.0, + dynamic_viscosity=1.8e-5, +) +energy_released = latent.last_latent_heat_energy +``` + +When no latent heat strategy is configured (or a nonpositive scalar is +provided), the step follows the isothermal path and reports +`last_latent_heat_energy = 0.0`. + ### CondensationIsothermalStaggered (two-pass Gauss-Seidel) `CondensationIsothermalStaggered` splits each timestep into two passes. Theta @@ -476,6 +517,8 @@ aerosol, gas = selective.step( | `accommodation_coefficient` | Mass accommodation coefficient (unitless). | `1.0` | | `update_gases` | Whether to deplete gas concentrations during step. | `True` | | `skip_partitioning_indices` | Species indices to exclude from partitioning. | `None` | +| `latent_heat_strategy` | Optional latent heat strategy. | `None` | +| `latent_heat` | Scalar fallback latent heat [J/kg]. | `0.0` | | `theta_mode` | Staggered theta selection: `"half"`, `"random"`, `"batch"`. | `"half"` (staggered) | | `num_batches` | Gauss-Seidel batch count (clipped to particle count). | `1` | | `shuffle_each_step` | Shuffle particle order each step (staggered). | `True` | @@ -500,7 +543,7 @@ aerosol, gas = selective.step( - Staggered solver is Gauss-Seidel only; other solvers are not exposed. - Factory supports `"isothermal"` and `"isothermal_staggered"` only. -- No latent-heat or temperature feedback; condensation is isothermal. +- No temperature feedback; latent heat is diagnostic only. - Minimum-radius clamp (1e-10 m) enforces continuum validity; sub-continuum physics is out of scope. diff --git a/docs/index.md b/docs/index.md index abee15f70..2fea125dd 100644 --- a/docs/index.md +++ b/docs/index.md @@ -31,8 +31,8 @@ Whether you’re a researcher, educator, or industry expert, Particula is design - **Supporting non-isothermal condensation** with thermal resistance, latent-heat mass transfer rate utilities, latent-heat energy release bookkeeping, and the `CondensationLatentHeat` strategy with - latent-heat-corrected `mass_transfer_rate()`/`rate()` (step integration - pending). + latent-heat-corrected `mass_transfer_rate()`/`rate()` and an implemented + `step()` that tracks per-step latent heat energy. - **Interrogating your experimental data** to validate and expand your impact. - **Fostering open-source collaboration** to share ideas and build on each other’s work. diff --git a/particula/dynamics/condensation/condensation_strategies.py b/particula/dynamics/condensation/condensation_strategies.py index 98bf0775d..c19b44ffc 100644 --- a/particula/dynamics/condensation/condensation_strategies.py +++ b/particula/dynamics/condensation/condensation_strategies.py @@ -1762,7 +1762,9 @@ class CondensationLatentHeat(CondensationStrategy): Attributes: latent_heat_strategy_input: Strategy input provided at initialization. latent_heat_input: Raw latent heat value provided at initialization. - last_latent_heat_energy: Diagnostic latent heat energy tracker. + last_latent_heat_energy: Net latent heat energy released in the most + recent step [J]. Positive values indicate condensation and + negative values indicate evaporation. Overwritten each call. """ # pylint: disable=R0913, R0917 @@ -2040,7 +2042,12 @@ def rate( def _update_latent_heat_energy( self, mass_transfer: NDArray[np.float64], temperature: float ) -> None: - """Store latent heat energy released for the current step.""" + """Store latent heat energy released for the current step. + + Args: + mass_transfer: Per-particle mass transfer in kg for this step. + temperature: System temperature in Kelvin. + """ if self._latent_heat_strategy is None: self.last_latent_heat_energy = 0.0 return @@ -2093,7 +2100,7 @@ def step( The mass transfer rate is computed, optional skip-partitioning applied, and both the particle and gas states are updated while respecting - inventory limits. The per-step latent-heat energy release is stored in + inventory limits. The per-step latent heat release is stored in ``last_latent_heat_energy`` (positive for condensation, negative for evaporation) and overwritten on every call. @@ -2103,7 +2110,8 @@ def step( temperature: System temperature in Kelvin. pressure: System pressure in Pascals. time_step: Integration timestep in seconds. - dynamic_viscosity: Optional dynamic viscosity override. + dynamic_viscosity: Optional dynamic viscosity override used in the + mass-transfer calculation. Returns: Tuple containing updated particle and gas species objects. diff --git a/readme.md b/readme.md index dacb07973..e08978f18 100644 --- a/readme.md +++ b/readme.md @@ -99,7 +99,8 @@ particula/ `particula.dynamics.get_latent_heat_energy_released` - **Condensation Strategies** — `CondensationIsothermal` plus `CondensationLatentHeat` with latent-heat-corrected - `mass_transfer_rate()`/`rate()` (step integration pending) + `mass_transfer_rate()`/`rate()` and `step()` energy tracking via + `last_latent_heat_energy`, with optional `dynamic_viscosity` override - **Latent Heat Factories** — Build constant, linear, and power-law latent heat strategies via `particula.gas.LatentHeatFactory` with unit-aware builders and gas-phase exports for upcoming non-isothermal workflows From db7b2526e116203d4d3266523e3b194d3e85d716 Mon Sep 17 00:00:00 2001 From: Kyle Gorkowski Date: Wed, 4 Mar 2026 18:50:34 -0700 Subject: [PATCH 3/3] fix(validate): address validation gaps for #1148 Successfully fixed: - Refactored CondensationLatentHeat.step into helpers to reduce complexity and keep behavior intact. - Sanitized non-finite pressure_delta values to mirror isothermal behavior. - Validated volume/concentration inputs for norm_conc and raise clear errors. - Normalized mass_transfer shape before latent heat energy accounting. - Clamped gas concentration updates to nonnegative values. - Documented latent heat energy units and expanded latent heat step tests. Remaining issues: none. Closes #1148 ADW-ID: 31535132 --- adw-docs/dev-plans/README.md | 2 +- .../epics/E5-non-isothermal-condensation.md | 3 +- ...E5-F3-condensation-latent-heat-strategy.md | 5 +- adw-docs/dev-plans/features/index.md | 2 +- docs/Features/condensation_strategy_system.md | 5 +- .../condensation/condensation_strategies.py | 149 ++++++--- .../tests/condensation_strategies_test.py | 303 +++++++++++++++++- 7 files changed, 410 insertions(+), 59 deletions(-) diff --git a/adw-docs/dev-plans/README.md b/adw-docs/dev-plans/README.md index 77a2ea268..969118530 100644 --- a/adw-docs/dev-plans/README.md +++ b/adw-docs/dev-plans/README.md @@ -77,7 +77,7 @@ and rollout. - Scope: Thermal resistance factor and non-isothermal mass transfer rate pure functions with energy tracking. - [E5-F3: CondensationLatentHeat Strategy Class][e5-f3] — Status: In Progress - (P3, #1141) + (P1, #1139) - Scope: New condensation strategy with latent heat correction and energy diagnostics. - [E5-F4: Builder, Factory, and Exports][e5-f4] — Status: Planning diff --git a/adw-docs/dev-plans/epics/E5-non-isothermal-condensation.md b/adw-docs/dev-plans/epics/E5-non-isothermal-condensation.md index 7b3bd769b..71dc38072 100644 --- a/adw-docs/dev-plans/epics/E5-non-isothermal-condensation.md +++ b/adw-docs/dev-plans/epics/E5-non-isothermal-condensation.md @@ -284,7 +284,7 @@ L -> 0 as T -> T_c. Used in engineering thermodynamics and EOS-based models. relative), reduced rate when L > 0 (thermal resistance always slows condensation), array shapes for single and multi-species -- [x] **E5-F3-P3**: Implement `step()` for particle-resolved with energy +- [ ] **E5-F3-P3**: Implement `step()` for particle-resolved with energy tracking and tests - `step()` follows `CondensationIsothermal.step()` flow (lines 922-1040) with two additions: @@ -623,4 +623,3 @@ class CondensationLatentHeat(CondensationStrategy): | 2026-03-02 | Initial epic creation | ADW | | 2026-03-02 | Split E5-F1-P3 into P3 (builders) + P4 (factory+exports); split E5-F3-P3 into P3 (particle-resolved step) + P4 (discrete+continuous) + P5 (data-only parity); added missing details: function signatures, file references, thermal conductivity source, vapor_pressure_surface parameter, test tolerances, literature targets | ADW | | 2026-03-04 | Noted E5-F3-P1 issue #1139 and logging expectations | ADW | -| 2026-03-04 | Marked E5-F3-P3 complete for issue #1141 | ADW | diff --git a/adw-docs/dev-plans/features/E5-F3-condensation-latent-heat-strategy.md b/adw-docs/dev-plans/features/E5-F3-condensation-latent-heat-strategy.md index c98bae8d3..97bf96004 100644 --- a/adw-docs/dev-plans/features/E5-F3-condensation-latent-heat-strategy.md +++ b/adw-docs/dev-plans/features/E5-F3-condensation-latent-heat-strategy.md @@ -146,9 +146,9 @@ additions: - Tests: numerical parity with `CondensationIsothermal` when L=0 (< 1e-15 relative), reduced rate when L > 0, array shapes for single/multi-species -- [x] **E5-F3-P3**: Implement `step()` for particle-resolved with energy +- [ ] **E5-F3-P3**: Implement `step()` for particle-resolved with energy tracking and tests - - Issue: #1141 | Size: M (~100 LOC) | Status: Complete + - Issue: TBD | Size: M (~100 LOC) | Status: Not Started - `step()` follows `CondensationIsothermal.step()` flow (lines 922-1040) with two additions: 1. Uses `mass_transfer_rate` from P2 (includes thermal correction) @@ -241,4 +241,3 @@ additions: |------|--------|--------| | 2026-03-02 | Initial feature document created from E5 epic | ADW | | 2026-03-04 | Marked P1 in progress for issue #1139 | ADW | -| 2026-03-04 | Completed P3 step() implementation and tests (issue #1141) | ADW | diff --git a/adw-docs/dev-plans/features/index.md b/adw-docs/dev-plans/features/index.md index 88d2749bb..2151d86ce 100644 --- a/adw-docs/dev-plans/features/index.md +++ b/adw-docs/dev-plans/features/index.md @@ -28,7 +28,7 @@ work, typically ~100 LOC per phase, that deliver user-facing functionality. |----|------|--------|----------|--------| | E5-F1 | [Latent Heat Strategy Pattern](E5-F1-latent-heat-strategy.md) | In Progress | P1 | 4 | | E5-F2 | [Non-Isothermal Mass Transfer Functions](E5-F2-non-isothermal-mass-transfer.md) | Planning | P1 | 3 | -| E5-F3 | [CondensationLatentHeat Strategy Class](E5-F3-condensation-latent-heat-strategy.md) | In Progress (P3) | P1 | 5 | +| E5-F3 | [CondensationLatentHeat Strategy Class](E5-F3-condensation-latent-heat-strategy.md) | In Progress | P1 | 5 | | E5-F4 | [Builder, Factory, and Exports](E5-F4-builder-factory-exports.md) | Planning | P1 | 2 | | E5-F5 | [Validation and Integration Tests](E5-F5-validation-integration-tests.md) | Planning | P1 | 2 | | E5-F6 | [Documentation and Examples](E5-F6-documentation-examples.md) | Planning | P2 | 3 | diff --git a/docs/Features/condensation_strategy_system.md b/docs/Features/condensation_strategy_system.md index 3f5067d50..347556acb 100644 --- a/docs/Features/condensation_strategy_system.md +++ b/docs/Features/condensation_strategy_system.md @@ -165,13 +165,16 @@ particle, gas = latent.step( time_step=1.0, dynamic_viscosity=1.8e-5, ) -energy_released = latent.last_latent_heat_energy +energy_released = latent.last_latent_heat_energy # total energy [J] ``` When no latent heat strategy is configured (or a nonpositive scalar is provided), the step follows the isothermal path and reports `last_latent_heat_energy = 0.0`. +`last_latent_heat_energy` records the total latent heat released per step +(sum of dm × L), not an energy density. + ### CondensationIsothermalStaggered (two-pass Gauss-Seidel) `CondensationIsothermalStaggered` splits each timestep into two passes. Theta diff --git a/particula/dynamics/condensation/condensation_strategies.py b/particula/dynamics/condensation/condensation_strategies.py index c19b44ffc..243a9ccb4 100644 --- a/particula/dynamics/condensation/condensation_strategies.py +++ b/particula/dynamics/condensation/condensation_strategies.py @@ -1762,7 +1762,7 @@ class CondensationLatentHeat(CondensationStrategy): Attributes: latent_heat_strategy_input: Strategy input provided at initialization. latent_heat_input: Raw latent heat value provided at initialization. - last_latent_heat_energy: Net latent heat energy released in the most + last_latent_heat_energy: Total latent heat energy released in the most recent step [J]. Positive values indicate condensation and negative values indicate evaporation. Overwritten each call. """ @@ -1955,11 +1955,9 @@ def mass_transfer_rate( pressure_delta = self.calculate_pressure_delta( particle, gas_species, temperature, radius_with_fill ) - if not np.all(np.isfinite(pressure_delta)): - raise ValueError( - "Non-finite pressure_delta computed for latent-heat " - "condensation." - ) + pressure_delta = np.nan_to_num( + pressure_delta, posinf=0.0, neginf=0.0, nan=0.0 + ) if self._latent_heat_strategy is None: return get_mass_transfer_rate( @@ -2062,6 +2060,98 @@ def _update_latent_heat_energy( ) ) + def _calculate_norm_conc( + self, volume: float, concentration: NDArray[np.float64] + ) -> NDArray[np.float64]: + """Validate inputs and compute normalized concentration. + + Args: + volume: Particle volume for the single box [m^3]. + concentration: Particle concentration array. + + Returns: + Normalized concentration (concentration / volume). + + Raises: + ValueError: If volume or concentration inputs are invalid. + """ + if not np.isfinite(volume) or volume <= 0.0: + raise ValueError("volume must be finite and positive.") + if not np.all(np.isfinite(concentration)) or np.any( + concentration < 0.0 + ): + raise ValueError( + "concentration must be finite and nonnegative for norm_conc." + ) + return concentration / volume + + def _normalize_mass_transfer_shape( + self, + mass_transfer: NDArray[np.float64], + species_mass: NDArray[np.float64], + ) -> NDArray[np.float64]: + """Normalize mass transfer shape to match species count. + + Args: + mass_transfer: Per-particle mass transfer array. + species_mass: Particle mass array used to infer species count. + + Returns: + Mass transfer array with consistent species dimension. + """ + species_count = 1 if species_mass.ndim == 1 else species_mass.shape[1] + if mass_transfer.ndim == 1: + return mass_transfer.reshape(-1, species_count) + if mass_transfer.shape[1] > species_count: + return mass_transfer[:, :species_count] + if mass_transfer.shape[1] < species_count: + return np.broadcast_to( + mass_transfer, (mass_transfer.shape[0], species_count) + ) + return mass_transfer + + def _apply_mass_transfer_to_particles( + self, + particle_data: ParticleData, + mass_transfer: NDArray[np.float64], + norm_conc: NDArray[np.float64], + ) -> None: + """Apply mass transfer to particle data in place.""" + if mass_transfer.ndim == 2: + nonzero_mask = norm_conc != 0 + denom = norm_conc[:, np.newaxis] + where_mask = nonzero_mask[:, np.newaxis] + else: + denom = norm_conc + where_mask = norm_conc != 0 + per_particle = np.divide( + mass_transfer, + denom, + out=np.zeros_like(mass_transfer), + where=where_mask, + ) + np.add( + particle_data.masses[0], + per_particle, + out=particle_data.masses[0], + ) + np.maximum(particle_data.masses[0], 0.0, out=particle_data.masses[0]) + + def _apply_mass_transfer_to_gas( + self, gas_data: GasData, mass_transfer: NDArray[np.float64] + ) -> None: + """Apply mass transfer to gas concentrations in place.""" + np.subtract( + gas_data.concentration[0], + mass_transfer.sum(axis=0), + out=gas_data.concentration[0], + ) + np.maximum( + gas_data.concentration[0], + 0.0, + out=gas_data.concentration[0], + ) + # pylint: disable=too-many-positional-arguments, too-many-arguments @overload # type: ignore[override] def step( @@ -2101,8 +2191,8 @@ def step( The mass transfer rate is computed, optional skip-partitioning applied, and both the particle and gas states are updated while respecting inventory limits. The per-step latent heat release is stored in - ``last_latent_heat_energy`` (positive for condensation, negative for - evaporation) and overwritten on every call. + ``last_latent_heat_energy`` [J] (positive for condensation, negative + for evaporation) and overwritten on every call. Args: particle: Particle representation to update. @@ -2130,10 +2220,7 @@ def step( dynamic_viscosity=dynamic_viscosity, ) - if isinstance(mass_rate, (int, float)): - mass_rate_array = np.array([mass_rate]) - else: - mass_rate_array = mass_rate + mass_rate_array = np.atleast_1d(np.asarray(mass_rate, dtype=np.float64)) mass_rate_array = self._apply_skip_partitioning(mass_rate_array) @@ -2141,7 +2228,8 @@ def step( np.asarray(gas_data.concentration[0], dtype=np.float64) ) volume = particle_data.volume[0] - norm_conc = particle_data.concentration[0] / volume + concentration = particle_data.concentration[0] + norm_conc = self._calculate_norm_conc(volume, concentration) mass_transfer = get_mass_transfer( mass_rate=mass_rate_array, time_step=time_step, @@ -2150,37 +2238,20 @@ def step( particle_concentration=norm_conc, ) - self._update_latent_heat_energy(mass_transfer, temperature) + mass_transfer = self._normalize_mass_transfer_shape( + mass_transfer, + particle_data.masses[0], + ) - species_mass = particle_data.masses[0] - if species_mass.ndim == 1: - species_count = 1 - else: - species_count = species_mass.shape[1] - if mass_transfer.ndim == 1: - mass_transfer = mass_transfer.reshape(-1, species_count) - elif mass_transfer.shape[1] > species_count: - mass_transfer = mass_transfer[:, :species_count] - elif mass_transfer.shape[1] < species_count: - mass_transfer = np.broadcast_to( - mass_transfer, (mass_transfer.shape[0], species_count) - ) + self._update_latent_heat_energy(mass_transfer, temperature) if particle_is_legacy: particle.add_mass(added_mass=mass_transfer) # type: ignore[union-attr] else: - per_particle = np.divide( + self._apply_mass_transfer_to_particles( + particle_data, mass_transfer, - norm_conc[:, np.newaxis] - if mass_transfer.ndim == 2 - else norm_conc, - out=np.zeros_like(mass_transfer), - where=norm_conc[:, np.newaxis] != 0 - if mass_transfer.ndim == 2 - else norm_conc != 0, - ) - particle_data.masses[0] = np.maximum( - particle_data.masses[0] + per_particle, 0.0 + norm_conc, ) if self.update_gases: if gas_is_legacy: @@ -2188,7 +2259,7 @@ def step( added_concentration=-mass_transfer.sum(axis=0) ) else: - gas_data.concentration[0] -= mass_transfer.sum(axis=0) + self._apply_mass_transfer_to_gas(gas_data, mass_transfer) if particle_is_legacy: return cast( Tuple[ParticleRepresentation, GasSpecies], diff --git a/particula/dynamics/condensation/tests/condensation_strategies_test.py b/particula/dynamics/condensation/tests/condensation_strategies_test.py index 33b46990f..dff15bd2f 100644 --- a/particula/dynamics/condensation/tests/condensation_strategies_test.py +++ b/particula/dynamics/condensation/tests/condensation_strategies_test.py @@ -2126,7 +2126,7 @@ def test_rate_skip_partitioning_out_of_bounds_raises(self): ) def test_nonfinite_pressure_delta_sanitized(self): - """Non-finite pressure deltas raise an error.""" + """Non-finite pressure deltas are sanitized to zero.""" gas_species, particle, _, _ = self._build_single_species_fixture( n_particles=3 ) @@ -2137,13 +2137,50 @@ def test_nonfinite_pressure_delta_sanitized(self): "calculate_pressure_delta", return_value=np.array([np.nan, np.inf, -np.inf]), ): - with self.assertRaises(ValueError): - cond.mass_transfer_rate( - particle=particle, - gas_species=gas_species, - temperature=self.temperature, - pressure=self.pressure, - ) + rates = cond.mass_transfer_rate( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + ) + self.assertTrue(np.all(np.isfinite(rates))) + self.assertTrue(np.allclose(rates, 0.0)) + + def test_nonfinite_pressure_delta_matches_isothermal(self): + """Latent heat sanitization mirrors isothermal behavior.""" + gas_species, particle, _, _ = self._build_single_species_fixture( + n_particles=3 + ) + latent_cond = CondensationLatentHeat(molar_mass=self.molar_mass) + iso_cond = CondensationIsothermal(molar_mass=self.molar_mass) + nonfinite = np.array([np.nan, np.inf, -np.inf]) + + with ( + patch.object( + CondensationLatentHeat, + "calculate_pressure_delta", + return_value=nonfinite, + ), + patch.object( + CondensationIsothermal, + "calculate_pressure_delta", + return_value=nonfinite, + ), + ): + latent_rates = latent_cond.mass_transfer_rate( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + ) + iso_rates = iso_cond.mass_transfer_rate( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + ) + + np.testing.assert_allclose(latent_rates, iso_rates, rtol=1e-14) def test_nonfinite_latent_heat_raises(self): """Non-finite latent heat raises an error.""" @@ -2315,6 +2352,40 @@ def test_step_mass_conservation_data_only(self): relative_error = abs(final_mass - initial_mass) / initial_mass self.assertLess(relative_error, 1e-12) + def test_step_rejects_nonpositive_volume(self): + """Step rejects nonpositive volumes for norm_conc.""" + particle_data, gas_data = self._make_data_inputs() + particle_data.volume[0] = 0.0 + + with self.assertRaisesRegex(ValueError, "volume must be finite"): + self.data_strategy.step( + particle=particle_data, + gas_species=gas_data, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + + def test_step_rejects_invalid_concentration_for_norm_conc(self): + """Step rejects non-finite or negative concentrations for norm_conc.""" + particle_data, gas_data = self._make_data_inputs() + + for bad_value in (-1.0, np.nan): + with self.subTest(bad_value=bad_value): + particle_copy = copy.deepcopy(particle_data) + particle_copy.concentration[0][0] = bad_value + + with self.assertRaisesRegex( + ValueError, "concentration must be finite" + ): + self.data_strategy.step( + particle=particle_copy, + gas_species=gas_data, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + def test_step_energy_tracking_matches_dm_times_l(self): """Latent heat energy matches dm × L for a step.""" particle = copy.deepcopy(self.particle) @@ -2334,11 +2405,58 @@ def test_step_energy_tracking_matches_dm_times_l(self): temperature=self.temperature, pressure=self.pressure, ) - mass_rate_array = ( - np.array([mass_rate]) - if isinstance(mass_rate, (int, float)) - else mass_rate + if isinstance(mass_rate, (int, float)): + mass_rate_array = np.array([mass_rate]) + else: + mass_rate_array = mass_rate + mass_rate_array = cond._apply_skip_partitioning(mass_rate_array) + norm_conc = particle_data.concentration[0] / particle_data.volume[0] + mass_transfer = get_mass_transfer( + mass_rate=mass_rate_array, + time_step=1.0, + gas_mass=np.asarray(gas_data.concentration[0], dtype=np.float64), + particle_mass=particle_data.masses[0], + particle_concentration=norm_conc, + ) + + cond.step( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + latent_heat = cond._latent_heat_strategy.latent_heat(self.temperature) + + np.testing.assert_allclose( + cond.last_latent_heat_energy, + np.sum(mass_transfer * latent_heat), + rtol=1e-14, + ) + + def test_step_energy_tracking_update_gases_false(self): + """Energy tracking works when gas updates are disabled.""" + gas_species, particle, _, _ = self._build_single_species_fixture() + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + update_gases=False, + ) + + particle_data, _ = _unwrap_particle(particle) + gas_data, _ = _unwrap_gas(gas_species) + mass_rate = cond.mass_transfer_rate( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, ) + if isinstance(mass_rate, (int, float)): + mass_rate_array = np.array([mass_rate]) + else: + mass_rate_array = mass_rate mass_rate_array = cond._apply_skip_partitioning(mass_rate_array) norm_conc = particle_data.concentration[0] / particle_data.volume[0] mass_transfer = get_mass_transfer( @@ -2348,6 +2466,17 @@ def test_step_energy_tracking_matches_dm_times_l(self): particle_mass=particle_data.masses[0], particle_concentration=norm_conc, ) + species_mass = particle_data.masses[0] + species_count = 1 if species_mass.ndim == 1 else species_mass.shape[1] + if mass_transfer.ndim == 1: + mass_transfer = mass_transfer.reshape(-1, species_count) + elif mass_transfer.shape[1] > species_count: + mass_transfer = mass_transfer[:, :species_count] + elif mass_transfer.shape[1] < species_count: + mass_transfer = np.broadcast_to( + mass_transfer, (mass_transfer.shape[0], species_count) + ) + initial_conc = gas_species.get_concentration().copy() cond.step( particle=particle, @@ -2356,8 +2485,128 @@ def test_step_energy_tracking_matches_dm_times_l(self): pressure=self.pressure, time_step=1.0, ) + + np.testing.assert_allclose( + gas_species.get_concentration(), initial_conc, rtol=1e-15 + ) + latent_heat = cond._latent_heat_strategy.latent_heat(self.temperature) + np.testing.assert_allclose( + cond.last_latent_heat_energy, + np.sum(mass_transfer * latent_heat), + rtol=1e-14, + ) + + def test_step_energy_tracking_data_inputs(self): + """Latent heat energy tracks on data container inputs.""" + particle_data, gas_data = self._make_data_inputs() + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + activity_strategy=self.activity_strategy, + surface_strategy=self.surface_strategy, + vapor_pressure_strategy=self.vapor_pressure_strategy, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + ) + + mass_rate = cond.mass_transfer_rate( + particle=particle_data, + gas_species=gas_data, + temperature=self.temperature, + pressure=self.pressure, + ) + if isinstance(mass_rate, (int, float)): + mass_rate_array = np.array([mass_rate]) + else: + mass_rate_array = mass_rate + mass_rate_array = cond._apply_skip_partitioning(mass_rate_array) + norm_conc = particle_data.concentration[0] / particle_data.volume[0] + mass_transfer = get_mass_transfer( + mass_rate=mass_rate_array, + time_step=1.0, + gas_mass=np.asarray(gas_data.concentration[0], dtype=np.float64), + particle_mass=particle_data.masses[0], + particle_concentration=norm_conc, + ) + species_mass = particle_data.masses[0] + species_count = 1 if species_mass.ndim == 1 else species_mass.shape[1] + if mass_transfer.ndim == 1: + mass_transfer = mass_transfer.reshape(-1, species_count) + elif mass_transfer.shape[1] > species_count: + mass_transfer = mass_transfer[:, :species_count] + elif mass_transfer.shape[1] < species_count: + mass_transfer = np.broadcast_to( + mass_transfer, (mass_transfer.shape[0], species_count) + ) + + cond.step( + particle=particle_data, + gas_species=gas_data, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + latent_heat = cond._latent_heat_strategy.latent_heat(self.temperature) + np.testing.assert_allclose( + cond.last_latent_heat_energy, + np.sum(mass_transfer * latent_heat), + rtol=1e-14, + ) + + def test_step_energy_tracking_skip_partitioning(self): + """Skip-partitioning affects latent heat bookkeeping.""" + particle = copy.deepcopy(self.particle) + gas_species = copy.deepcopy(self.gas_species) + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + skip_partitioning_indices=[1], + ) + particle_data, _ = _unwrap_particle(particle) + gas_data, _ = _unwrap_gas(gas_species) + mass_rate = cond.mass_transfer_rate( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + ) + if isinstance(mass_rate, (int, float)): + mass_rate_array = np.array([mass_rate]) + else: + mass_rate_array = mass_rate + mass_rate_array = cond._apply_skip_partitioning(mass_rate_array) + norm_conc = particle_data.concentration[0] / particle_data.volume[0] + mass_transfer = get_mass_transfer( + mass_rate=mass_rate_array, + time_step=1.0, + gas_mass=np.asarray(gas_data.concentration[0], dtype=np.float64), + particle_mass=particle_data.masses[0], + particle_concentration=norm_conc, + ) + species_mass = particle_data.masses[0] + species_count = 1 if species_mass.ndim == 1 else species_mass.shape[1] + if mass_transfer.ndim == 1: + mass_transfer = mass_transfer.reshape(-1, species_count) + elif mass_transfer.shape[1] > species_count: + mass_transfer = mass_transfer[:, :species_count] + elif mass_transfer.shape[1] < species_count: + mass_transfer = np.broadcast_to( + mass_transfer, (mass_transfer.shape[0], species_count) + ) + + cond.step( + particle=particle, + gas_species=gas_species, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + + latent_heat = cond._latent_heat_strategy.latent_heat(self.temperature) np.testing.assert_allclose( cond.last_latent_heat_energy, np.sum(mass_transfer * latent_heat), @@ -2534,6 +2783,36 @@ def test_step_updates_gas_when_enabled(self): gas_new.get_concentration(), expected_gas, rtol=1e-15 ) + def test_step_clamps_gas_concentration_nonnegative(self): + """Gas concentrations are clamped to nonnegative values.""" + particle_data, gas_data = self._make_data_inputs() + gas_data.concentration[0] = np.array([1e-12, 2e-12]) + cond = CondensationLatentHeat( + molar_mass=self.molar_mass, + diffusion_coefficient=self.diffusion_coefficient, + accommodation_coefficient=self.accommodation_coefficient, + activity_strategy=self.activity_strategy, + surface_strategy=self.surface_strategy, + vapor_pressure_strategy=self.vapor_pressure_strategy, + latent_heat_strategy=ConstantLatentHeat(latent_heat_ref=2.26e6), + ) + mass_transfer = np.full_like(particle_data.masses[0], 1.0) + + with patch( + "particula.dynamics.condensation.condensation_strategies." + "get_mass_transfer", + return_value=mass_transfer, + ): + _, gas_new = cond.step( + particle=particle_data, + gas_species=gas_data, + temperature=self.temperature, + pressure=self.pressure, + time_step=1.0, + ) + + np.testing.assert_allclose(gas_new.concentration[0], 0.0, rtol=1e-15) + def test_step_dynamic_viscosity_pass_through(self): """Dynamic viscosity override is forwarded to mass_transfer_rate.""" particle = copy.deepcopy(self.particle)