using System.Numerics.Tensors; using System.Runtime.CompilerServices; using Content.Shared.Atmos.Prototypes; using Content.Shared.Atmos.Reactions; using Content.Shared.CCVar; using JetBrains.Annotations; namespace Content.Shared.Atmos.EntitySystems; public abstract partial class SharedAtmosphereSystem { /* Partial class for operations involving GasMixtures. Sometimes methods here are abstract because they need different client/server implementations due to sandboxing. */ /// /// Cached array of molar heat capacities of the gases. /// public float[] GasMolarHeatCapacities => _gasMolarHeatCapacities; private float[] _gasMolarHeatCapacities = new float[Atmospherics.AdjustedNumberOfGases]; /// /// Cached array of gas specific mols /// public float[] GasMolarMasses => _gasMolarMasses; private float[] _gasMolarMasses = new float[Atmospherics.AdjustedNumberOfGases]; /// /// Mask used to determine if a gas is flammable or not. /// /// This is used to quickly determine if a contains any flammable gas. /// When determining flammability, the float is multiplied with the mask and then /// added to see if the mixture is flammable, and how many moles are considered flammable. /// This is done instead of a massive if statement of doom everywhere. /// Say Plasma has the bool set to true. /// Atmospherics will place a 1 in the spot where plasma goes in the masking array. /// Whenever we need to determine if a GasMixture contains fuel gases, we multiply the /// gas array by the mask. Fuel gases will keep their value (being multiplied by one) /// whereas non-fuel gases will be multiplied by zero and be zeroed out. /// The resulting array can be HorizontalAdded, with any value above zero indicating fuel gases. /// This works for multiple fuel gases at the same time, so it's a fairly quick way /// to determine if a mixture has the gases we care about. protected readonly float[] GasFuelMask = new float[Atmospherics.AdjustedNumberOfGases]; /// /// Mask used to determine if a gas is an oxidizer or not. /// Used in the same way as . /// Nothing really super special. /// protected readonly float[] GasOxidizerMask = new float[Atmospherics.AdjustedNumberOfGases]; /// /// Mask used to determine both fuel and oxidizer properties of a gas at the same time. /// Primarily used to quickly report the specific moles in a mixture that caused a flammable reaction to occur. /// protected readonly float[] GasOxidiserFuelMask = new float[Atmospherics.TotalNumberOfGases]; public string?[] GasReagents = new string[Atmospherics.TotalNumberOfGases]; protected readonly GasPrototype[] GasPrototypes = new GasPrototype[Atmospherics.TotalNumberOfGases]; public virtual void InitializeGases() { foreach (var gas in Enum.GetValues()) { var idx = (int)gas; // Log an error if the corresponding prototype isn't found if (!ProtoMan.TryIndex(gas.ToString(), out var gasPrototype)) { Log.Error($"Failed to find corresponding {nameof(GasPrototype)} for gas ID {(int)gas} ({gas}) with expected ID \"{gas.ToString()}\". Is your prototype named correctly?"); continue; } GasPrototypes[idx] = gasPrototype; GasReagents[idx] = gasPrototype.Reagent; } for (var i = 0; i < GasPrototypes.Length; i++) { /* As an optimization routine we pre-divide the specific heat by the heat scale here, so we don't have to do it every time we calculate heat capacity. Most usages are going to want the scaled value anyway. If you would like the unscaled specific heat, you'd need to multiply by HeatScale again. TODO ATMOS: please just make this 2 separate arrays instead of invoking multiplication every time. */ _gasMolarHeatCapacities[i] = GasPrototypes[i].MolarHeatCapacity / HeatScale; _gasMolarMasses[i] = GasPrototypes[i].MolarMass; // """Mask""" built here. Used to determine if a gas is fuel/oxidizer or not decently quickly and clearly. GasFuelMask[i] = GasPrototypes[i].IsFuel ? 1 : 0; // Same for oxidizer mask. GasOxidizerMask[i] = GasPrototypes[i].IsOxidizer ? 1 : 0; // OxidiserFuel mask is just fuel and oxidizer combined, because both are required for a reaction to occur. GasOxidiserFuelMask[i] = GasFuelMask[i] * GasOxidizerMask[i]; } } /// /// Gets only the moles that are considered a fuel and an oxidizer in a . /// /// The to get the flammable moles for. /// A buffer to write the flammable moles into. Must be the same length as the number of gases. /// A of moles where only the flammable and oxidizer moles are returned, and the rest are 0. [PublicAPI] public void GetFlammableMoles(GasMixture mixture, float[] buffer) { TensorPrimitives.Multiply(mixture.Moles, GasOxidiserFuelMask, buffer); } /// /// Determines if a is ignitable or not. /// This is a combination of determining if a mixture both has oxidizer and fuel. /// /// The to determine. /// The minimum amount of moles at which a is /// considered ignitable, for both oxidizer and fuel. /// True if the is ignitable, otherwise, false. [PublicAPI] public bool IsMixtureIgnitable(GasMixture mixture, float epsilon = Atmospherics.Epsilon) { return IsMixtureFuel(mixture, epsilon) && IsMixtureOxidizer(mixture, epsilon); } /// /// Determines if a has fuel gases in it or not. /// /// The to determine. /// The minimum amount of moles at which a /// is considered fuel. /// True if the is fuel, otherwise, false. [PublicAPI] public abstract bool IsMixtureFuel(GasMixture mixture, float epsilon = Atmospherics.Epsilon); /// /// Determines if a has oxidizer gases in it or not. /// /// The to determine. /// The minimum amount of moles at which a /// is considered an oxidizer. /// True if the is an oxidizer, otherwise, false. [PublicAPI] public abstract bool IsMixtureOxidizer(GasMixture mixture, float epsilon = Atmospherics.Epsilon); /// /// Calculates the heat capacity for a . /// /// The to calculate the heat capacity for. /// Whether to apply the heat capacity scaling factor. /// This is an extremely important boolean to consider or else you will get heat transfer wrong. /// See for more info. /// The heat capacity of the . [PublicAPI] public float GetHeatCapacity(GasMixture mixture, bool applyScaling) { var scale = GetHeatCapacityCalculation(mixture.Moles, mixture.Immutable); // By default GetHeatCapacityCalculation() has the heat-scale divisor pre-applied. // So if we want the un-scaled heat capacity, we have to multiply by the scale. return applyScaling ? scale : scale * HeatScale; } /// /// Calculates the thermal energy for a . /// /// The to calculate the thermal /// energy of. /// The 's thermal energy in joules. [PublicAPI] public float GetThermalEnergy(GasMixture mixture) { return mixture.Temperature * GetHeatCapacity(mixture); } /// /// Calculates the thermal energy for a gas mixture, /// using a provided cached heat capacity value. /// /// The to calculate the thermal energy of. /// A cached heat capacity value for the gas mixture, /// to avoid redundant heat capacity calculations. /// The 's thermal energy in joules. [PublicAPI] public float GetThermalEnergy(GasMixture mixture, float cachedHeatCapacity) { return mixture.Temperature * cachedHeatCapacity; } /// /// Merges one into another, modifying the receiver. /// /// The to merge into. This will be modified. /// The to merge from. This will not be modified. [PublicAPI] public void Merge(GasMixture receiver, GasMixture giver) { if (receiver.Immutable) return; if (MathF.Abs(receiver.Temperature - giver.Temperature) > Atmospherics.MinimumTemperatureDeltaToConsider) { var receiverHeatCapacity = GetHeatCapacity(receiver); var giverHeatCapacity = GetHeatCapacity(giver); var combinedHeatCapacity = receiverHeatCapacity + giverHeatCapacity; if (combinedHeatCapacity > Atmospherics.MinimumHeatCapacity) { receiver.Temperature = (GetThermalEnergy(giver, giverHeatCapacity) + GetThermalEnergy(receiver, receiverHeatCapacity)) / combinedHeatCapacity; } } TensorPrimitives.Add(receiver.Moles, giver.Moles, receiver.Moles); } /// /// Performs reactions for a given gas mixture on an optional holder. /// /// The to perform reactions on. /// that holds the . /// used by Atmospherics to determine locality for certain reaction effects. /// The of the reactions performed. [PublicAPI] public abstract ReactionResult React(GasMixture mixture, IGasMixtureHolder? holder); /// /// Gets the heat capacity for a . /// /// The to calculate the heat capacity for. /// The heat capacity of the . /// Note that the heat capacity of the mixture may be slightly different from /// "real life" as we intentionally fake a heat capacity for space in /// in order to allow Atmospherics to cool down space. protected float GetHeatCapacity(GasMixture mixture) { return GetHeatCapacityCalculation(mixture.Moles, mixture.Immutable); } /// /// Gets the mass of a given /// /// The in question /// Returns the volume in kilograms. [PublicAPI] public abstract float GetMass(GasMixture mix); /// [PublicAPI] public abstract float GetMass(float[] moles); /// /// Calculates the amount of volume transferred from one gas mixture to another over time based on flow rate. /// /// /// A /// Another /// The area of transfer, in square meters. One tile of movement is about one square meter. /// delta time, or how much time is passing/has passed. /// Discharge coefficient. An abstract modifier for friction and turbulence. /// /// The volume of gas being moved over dt in Litres. /// If the value is positive it's in the direction of mix1->mix2, /// If it's negative it's in the direction of mix2 -> mix1 /// /// I'm assuming C is always 1 because I'm lazy, you can precalculate it and pass it with the area if you really care. [PublicAPI] public double GetFlowVolume(GasMixture mix1, GasMixture mix2, float area, float dt, float c = 1f) { ArgumentOutOfRangeException.ThrowIfNegativeOrZero(dt); return dt * GetFlowRate(mix1, mix2, area, c); } /// [PublicAPI] public double GetFlowVolume(GasMixture mix1, float deltaP, float area, float dt, float c = 1f) { ArgumentOutOfRangeException.ThrowIfNegativeOrZero(dt); return dt * GetFlowRate(mix1, deltaP, area, c); } /// /// Calculates the volumetric flow rate between two gas mixtures. /// /// A /// Another /// The area of transfer, in square meters. One tile of movement is about one square meter. /// Discharge coefficient. An abstract modifier for friction and turbulence. /// /// The volume of gas being moved in Litres / Second. /// If the value is positive it's in the direction of mix1->mix2, /// If it's negative it's in the direction of mix2 -> mix1 /// /// I'm assuming C is always 1 because I'm lazy, you can precalculate it and pass it with the area if you really care. [PublicAPI] public double GetFlowRate(GasMixture mix1, GasMixture mix2, float area, float c = 1f) { /* Q = C × A × √(2 × ΔP / ρ) Q is the volumetric airflow rate C is the discharge coefficient A is the cross-sectional area ΔP is the measured pressure difference ρ is the air density, adjusted for environmental conditions. We can break this up into Q = A × V where V is the velocity of the gas. */ ArgumentOutOfRangeException.ThrowIfNegativeOrZero(area); return area * GetFlowVelocity(mix1, mix2, c); } /// [PublicAPI] public double GetFlowRate(GasMixture mix1, float deltaP, float area, float c = 1f) { /* Q = C × A × √(2 × ΔP / ρ) Q is the volumetric airflow rate C is the discharge coefficient A is the cross-sectional area ΔP is the measured pressure difference ρ is the air density, adjusted for environmental conditions. We can break this up into Q = A × V where V is the velocity of the gas. */ ArgumentOutOfRangeException.ThrowIfNegativeOrZero(area); return area * GetFlowVelocity(mix1, deltaP, c); } /// /// Calculates the flow velocity between two gas mixtures using Q = C × A × √(2 × ΔP / ρ) but without the A (area) /// Useful for determining flow rate, or how fast a gas is moving. /// /// A /// Another /// Discharge coefficient. An abstract modifier for friction and turbulence. /// /// The velocity of gas movement between two mixtures in Meters / Second. /// If the value is positive it's in the direction of mix1->mix2, /// If it's negative it's in the direction of mix2 -> mix1 /// [PublicAPI] public double GetFlowVelocity(GasMixture mix1, GasMixture mix2, float c = 1f) { if (mix1.Pressure > mix2.Pressure) return GetFlowVelocity(mix1, mix1.Pressure - mix2.Pressure, c); return -GetFlowVelocity(mix2, mix2.Pressure - mix1.Pressure, c); } /// /// Calculates the flow velocity between a gas mixture given a pressure differential. /// /// The mixture which is being allowed to flow /// The difference in pressure between this mixture and where it's flowing to /// Discharge coefficient. An abstract modifier for friction and turbulence. /// /// The velocity of the gas leaving our mixture in Meters / Second. /// [PublicAPI] public double GetFlowVelocity(GasMixture mix1, float deltaP, float c = 1f) { /* V = C × √(2 × ΔP / ρ) V is the velocity of our gas C is the discharge coefficient ΔP is the measured pressure difference ρ is the air density, adjusted for environmental conditions. Density is equivalent to Mass / Volume, so we invert that to divide by density. */ ArgumentOutOfRangeException.ThrowIfNegativeOrZero(deltaP); ArgumentOutOfRangeException.ThrowIfNegativeOrZero(c); return c * Math.Sqrt(2 * deltaP * mix1.Volume / GetMass(mix1)); } /// /// Lets a volume of gas flow throw a constrained area into another volume of gas over a period of time. /// /// Gas volume that is discharging some of its gas. /// Gas volume that is receiving the discharge. /// Time that the discharge occurs in seconds, should be as small as possible since it doesn't use calculus /// Area that our gas is traveling through in m^2, the larger the area the bigger the transfer. /// Default of 2m^2 since that's the area of a single face of an atmos tile. [PublicAPI] public void FlowGas(GasMixture mixture, GasMixture? output, float dt, float area) { FlowGas(mixture, output, mixture.Pressure, dt, area); } /// [PublicAPI] public void FlowGas(GasMixture mixture, GasMixture? output, float pressure, float dt, float area) { if (output == null) { FlowGas(mixture, pressure, dt, area); return; } pressure = Math.Min(pressure, mixture.Pressure - output.Pressure); var removed = FlowGas(mixture, pressure, dt, area); if (removed == null) return; Merge(output, removed); } /// /// Lets a volume of gas flow through constrained area at a constrained pressure delta. /// /// Mixture of gas that is currently flowing /// Pressure our gas is able to flow at. /// Time that the discharge occurs in seconds, should be as small as possible since it doesn't use calculus /// Area that our gas is traveling through in m^2, the larger the area the bigger the transfer. /// Default of 2m^2 since that's the area of a single face of an atmos tile. /// [PublicAPI] public GasMixture? FlowGas(GasMixture mixture, float deltaP, float dt, float area = 2f) { if (deltaP <= 0) return null; return ReleaseGasAt(mixture, (float)GetFlowVolume(mixture, deltaP, area, dt), mixture.Pressure); } /// /// Releases some volume of a gas mixture at a specified pressure. /// /// Mixture which is releasing gas. /// Optional Mixture to receive gas /// Volume we are releasing /// Pressure of the released volume. [PublicAPI] public void ReleaseGasAt(GasMixture mixture, GasMixture? output, float volume, float targetPressure) { if (output == null) { ReleaseGasAt(mixture, volume, targetPressure); return; } targetPressure = Math.Min(targetPressure, mixture.Pressure - output.Pressure); if (targetPressure <= 0) return; var molesNeeded = Math.Min(targetPressure * volume / (Atmospherics.R * mixture.Temperature), MolesToEqualizePressure(mixture, output)); var removed = mixture.Remove(molesNeeded); Merge(mixture, removed); } /// [PublicAPI] public GasMixture? ReleaseGasAt(GasMixture mixture, float volume, float targetPressure) { if (targetPressure <= 0) return null; targetPressure = Math.Min(targetPressure, mixture.Pressure); return RemoveVolumeAtPressure(mixture, volume, targetPressure); } /// /// Removes a specified volume of gas from a mixture, at a specific pressure. /// /// mixture of gas /// volume we're attempting to remove /// pressure that volume will be removed at. public GasMixture RemoveVolumeAtPressure(GasMixture mixture, float volume, float pressure) { var molesNeeded = pressure * volume / (Atmospherics.R * mixture.Temperature); return mixture.Remove(molesNeeded); } /// /// Gets the heat capacity for a . /// /// The moles array of the /// Whether this represents space, /// and thus experiences space-specific mechanics (we cheat and make it a bit cooler). /// See . /// The heat capacity of the . [MethodImpl(MethodImplOptions.AggressiveInlining)] protected abstract float GetHeatCapacityCalculation(float[] moles, bool space); /// /// Calculates the moles that must be transferred from /// to to equalize pressure. /// public float MolesToEqualizePressure(GasMixture gasMixture1, GasMixture gasMixture2) { return gasMixture1.TotalMoles * FractionToEqualizePressure(gasMixture1, gasMixture2); } /// /// Calculates the dimensionless fraction of gas required to equalize pressure between two gas mixtures. /// /// The first gas mixture involved in the pressure equalization. /// This mixture should be the one you always expect to be the highest pressure. /// The second gas mixture involved in the pressure equalization. /// A float (from 0 to 1) representing the dimensionless fraction of gas that needs to be transferred from the /// mixture of higher pressure to the mixture of lower pressure. /// /// /// This properly takes into account the effect /// of gas merging from inlet to outlet affecting the temperature /// (and possibly increasing the pressure) in the outlet. /// /// /// The gas is assumed to expand freely, /// so the temperature of the gas with the greater pressure is not changing. /// /// /// /// If you want to calculate the moles required to equalize pressure between an inlet and an outlet, /// multiply the fraction returned by the source moles. /// public float FractionToEqualizePressure(GasMixture gasMixture1, GasMixture gasMixture2) { /* Problem: the gas being merged from the inlet to the outlet could affect the temp. of the gas and cause a pressure rise. We want the pressure to be equalized, so we have to account for this. For clarity, let's assume that gasMixture1 is the inlet and gasMixture2 is the outlet. We require mechanical equilibrium, so \( P_1' = P_2' \) Before the transfer, we have: \( P_1 = \frac{n_1 R T_1}{V_1} \) \( P_2 = \frac{n_2 R T_2}{V_2} \) After removing fraction \( x \) moles from the inlet, we have: \( P_1' = \frac{(1 - x) n_1 R T_1}{V_1} \) The outlet will gain the same \( x n_1 \) moles of gas. So \( n_2' = n_2 + x n_1 \) After mixing, the outlet temperature will be changed. Denote the new mixture temperature as \( T_2' \). Volume is constant. So we have: \( P_2' = \frac{(n_2 + x n_1) R T_2}{V_2} \) The total energy of the incoming inlet to outlet gas at \( T_1 \) plus the existing energy of the outlet gas at \( T_2 \) will be equal to the energy of the new outlet gas at \( T_2' \). This leads to the following derivation: \( x n_1 C_1 T_1 + n_2 C_2 T_2 = (x n_1 C_1 + n_2 C_2) T_2' \) Where \( C_1 \) and \( C_2 \) are the heat capacities of the inlet and outlet gases, respectively. Solving for \( T_2' \) gives us: \( T_2' = \frac{x n_1 C_1 T_1 + n_2 C_2 T_2}{x n_1 C_1 + n_2 C_2} \) Once again, we require mechanical equilibrium (\( P_1' = P_2' \)), so we can substitute \( T_2' \) into the pressure equation: \( \frac{(1 - x) n_1 R T_1}{V_1} = \frac{(n_2 + x n_1) R}{V_2} \cdot \frac{x n_1 C_1 T_1 + n_2 C_2 T_2} {x n_1 C_1 + n_2 C_2} \) Now it's a matter of solving for \( x \). Not going to show the full derivation here, just steps. 1. Cancel common factor \( R \). 2. Multiply both sides by \( x n_1 C_1 + n_2 C_2 \), so that everything becomes a polynomial in terms of \( x \). 3. Expand both sides. 4. Collect like powers of \( x \). 5. After collecting, you should end up with a polynomial of the form: \( (-n_1 C_1 T_1 (1 + \frac{V_2}{V_1})) x^2 + (n_1 T_1 \frac{V_2}{V_1} (C_1 - C_2) - n_2 C_1 T_1 - n_1 C_2 T_2) x + (n_1 T_1 \frac{V_2}{V_1} C_2 - n_2 C_2 T_2) = 0 \) Divide through by \( n_1 C_1 T_1 \) and replace each ratio with a symbol for clarity: \( k_V = \frac{V_2}{V_1} \) \( k_n = \frac{n_2}{n_1} \) \( k_T = \frac{T_2}{T_1} \) \( k_C = \frac{C_2}{C_1} \) */ // Ensure that P_1 > P_2 so the quadratic works out. if (gasMixture1.Pressure < gasMixture2.Pressure) { (gasMixture1, gasMixture2) = (gasMixture2, gasMixture1); } // Establish the dimensionless ratios. var volumeRatio = gasMixture2.Volume / gasMixture1.Volume; var molesRatio = gasMixture2.TotalMoles / gasMixture1.TotalMoles; var temperatureRatio = gasMixture2.Temperature / gasMixture1.Temperature; var heatCapacityRatio = GetHeatCapacity(gasMixture2) / GetHeatCapacity(gasMixture1); // The quadratic equation is solved for the transfer fraction. var quadraticA = 1 + volumeRatio; var quadraticB = molesRatio - volumeRatio + heatCapacityRatio * (temperatureRatio + volumeRatio); var quadraticC = heatCapacityRatio * (molesRatio * temperatureRatio - volumeRatio); return (-quadraticB + MathF.Sqrt(quadraticB * quadraticB - 4 * quadraticA * quadraticC)) / (2 * quadraticA); } /// /// Determines the fraction of gas to be removed and transferred from a source /// to a target to reach a target pressure /// in the target . /// /// The source that gas will be removed from. /// This should always be of higher pressure than the second . /// The target that will increase in pressure /// to the target pressure. /// The target mixture's desired pressure to target. /// A float representing the dimensionless fraction of gas to transfer from the source /// to the target. This may return negative if you have your mixtures swapped. /// Note that this method doesn't take into account the heat capacity of the /// transferred volume causing a pressure rise in the target . [PublicAPI] public static float FractionToMaxPressure(GasMixture mix1, GasMixture mix2, float targetPressure) { var molesToTransfer = MolesToMaxPressure(mix1, mix2, targetPressure); return molesToTransfer / mix1.TotalMoles; } /// /// Determines the number of moles to be removed and transferred from a source /// to a target to reach a target pressure /// in the target . /// /// The source that gas will be removed from. /// This should always be of higher pressure than the second . /// The target that will increase in pressure /// to the target pressure. /// The target mixture's desired pressure to target. /// The difference in moles required to reach the target pressure. /// Note that this method doesn't take into account the heat capacity of the /// transferred volume causing a pressure rise in the target . [PublicAPI] public static float MolesToMaxPressure(GasMixture mix1, GasMixture mix2, float targetPressure) { /* Calculate the moles required to reach the target pressure. The formula is derived from the ideal gas law and the general Richman's law, under the simplification that all the specific heat capacities are equal. Derivation can also be seen at https://github.com/space-wizards/space-station-14/pull/35211/files/a0ae787fe07a4e792570f55b49d9dd8038eb6e4d#r1961183456 TODO ATMOS Make this properly obey the heat capacity change on the target mixture. Derivation is as follows. Assume A is mix1, B is mix2, C is the combined mixture after transfer. We can express the number of moles in C: n_C = n_A + n_B We can then determine the temperature of C: T_C = \frac{T_A n_A c_A + T_B n_B c_B}{n_A c_A + n_B c_B} Where c_A and c_B are the specific heats of mixtures A and B, respectively. We can then express the pressure of C: P_C = \frac{n_C R T_C}{V_C} Using the above equations, we can express P_C as follows: P_C = \frac{(n_A + n_B) R (\frac{T_a n_A + T_B n_B}{n_A + n_B}}{V_C} Which can be reduced to: P_C = \frac{R (T_A n_A + T_B n_B)}{V_C} Solving for n_A gives: n_A = \frac{P_C V_C - R T_B n_B}{R T_A} Using the ideal gas law to substitute: n_A = \frac{P_C V_C - P_B V_B}{R T_A} The output volume doesn't change: V_B = V_C So: n_A = \frac{(P_C - P_B) V_B}{R T_A} */ var delta = targetPressure - mix2.Pressure; var requiredMoles = (delta * mix2.Volume) / (mix1.Temperature * Atmospherics.R); // Return the fraction of moles to transfer. return requiredMoles; } /// /// Determines the number of moles that need to be removed from a to reach a target pressure threshold. /// /// The gas mixture whose moles and properties will be used in the calculation. /// The target pressure threshold to calculate against. /// The difference in moles required to reach the target pressure threshold. /// The temperature of the gas is assumed to be not changing due to a free expansion. public static float MolesToPressureThreshold(GasMixture gasMixture, float targetPressure) { // Kid named PV = nRT. return gasMixture.TotalMoles - targetPressure * gasMixture.Volume / (Atmospherics.R * gasMixture.Temperature); } }