diff --git a/garnets/garnets.py b/garnets/garnets.py index 1209454..9a6792a 100644 --- a/garnets/garnets.py +++ b/garnets/garnets.py @@ -205,7 +205,10 @@ def coalesce_planetesimals(disk, planets, canidate, do_moons): and canidate.mass > .0001 * earth_mass \ and planet.mass_of_moons < planet.mass * .05 \ and planet.mass > canidate.mass: - # TODO: Remove planet.mass > canidate.mass distinction, just switch the canidate and planet! + # TODO: Refactor moon capture logic. If `canidate.mass > planet.mass`, + # swap `candidate` and `planet` to ensure the more massive body is treated + # as the primary object capturing the smaller one as a moon. + # This avoids the `planet.mass > canidate.mass` check. planet.add_moon( convert_planetesimal_to_protomoon( canidate, planet)) @@ -219,9 +222,11 @@ def coalesce_planetesimals(disk, planets, canidate, do_moons): finished = True break else: - # TODO: Reasons. + # TODO: Enhance logging for failed moon capture. Specify which condition(s) + # (e.g., candidate mass too high/low, planet's existing moon mass too high, + # candidate more massive than planet) prevented capture. logging.info( - "Did not capture potential moon at %s. Collision imminent." % quantity_repr( + "Did not capture potential moon at %s. Conditions not met or collision imminent." % quantity_repr( planet.orbit.a, au) ) @@ -234,7 +239,11 @@ def coalesce_planetesimals(disk, planets, canidate, do_moons): planet.gas_mass = planet.gas_mass + canidate.gas_mass # + new_gas finished = True - # Accrete MORE DUST! TODO: Refactor to this, mark the planet as needing to re-accrete. + # TODO: Consider refactoring the post-collision accretion step. + # Instead of immediate re-accretion via `disk.accrete_dust(planet)`, + # explore marking the planet as 'dirty' or 'needs_re_accretion' to + # handle this in a subsequent simulation phase. This could allow for + # more complex interactions or batch processing of accretion events. disk.accrete_dust(planet) logging.info("Conglomerate is now %s at %s." % ( @@ -243,9 +252,11 @@ def coalesce_planetesimals(disk, planets, canidate, do_moons): )) if not finished: - # TODO: Extra info. - logging.info("New Protoplanet at %s." % - quantity_repr(canidate.orbit.a, au)) + # TODO: Enhance logging for new protoplanet formation. Include + # additional information such as its mass and potentially key + # compositional details (e.g., dust/gas ratio if available). + logging.info("New Protoplanet formed at %s with mass %s." % + (quantity_repr(canidate.orbit.a, au), mass_repr(canidate.mass))) planets.append(convert_planetesimal_to_protoplanet(canidate)) @@ -259,7 +270,11 @@ def calculate_gases(star, planet, planet_id): for i in range(len(GASES)): - # TODO(woursler): WTH is this? + # TODO: Clarify the derivation and meaning of the variable `yp` + # and the formula used for its calculation. This formula appears to adjust + # boiling points based on pressure. Document the source of the constants + # (e.g., -5050.5, 373., 0.001) and explain its physical basis + # (e.g., Clausius-Clapeyron approximation for water, adapted for other gases). yp = GASES[i].boil / \ (373. * ((log((pressure) + 0.001) / -5050.5) + (1.0 / 373.))) @@ -274,7 +289,11 @@ def calculate_gases(star, planet, planet_id): fract = 1.0 pres2 = 1.0 - # TODO(woursler): This needs to go in the chemtable somehow. + # TODO: Refactor gas reactivity logic. Move the hardcoded rules + # for specific gases (Ar, He, O, O2, CO2) and their dependence on + # stellar age and pressure into `chemtable.py`, possibly by adding + # reactivity parameters or functions to the `Gas` class or a + # related data structure. This would centralize chemical properties. if GASES[i].symbol == "Ar": react = .15 * star.age / 4e9 / year @@ -383,7 +402,10 @@ def roche_limit(planet, moon): def hill_sphere(planet, star): - # TODO(woursler): Confirm this + # TODO: The formula used is a common approximation for the Hill sphere + # radius ( R_H = a * (m / (3*M))^(1/3) ). Confirm if this level of precision + # is adequate for the simulation or if more exact formulations are needed for + # specific scenarios (e.g., highly eccentric orbits, larger moon/planet mass ratios). return planet.orbit.a * pow( (planet.mass / (3.0 * star.mass)), (1.0 / 3.0)) @@ -431,12 +453,19 @@ def generate_planet(protoplanet, planet.radius = volume_radius(planet.mass, planet.density) planet.surf_accel = acceleration(planet.mass, planet.radius) - planet.surf_grav = planet.surf_accel # TODO(woursler): Remove ambiguity. + # TODO: Clarify the distinction between `surf_accel` (calculated + # from M, R) and `surf_grav`. If they are intended to be identical + # (gravitational acceleration at the surface), consider using only one + # variable (e.g., `surface_gravity`). If `surf_grav` might later include + # other effects (e.g., centrifugal due to rotation), document this. + planet.surf_grav = planet.surf_accel planet.molec_weight = min_molec_weight(planet) if ((planet.mass > 1 * earth_mass) - # TODO(woursler): Understand where this comes from. + # TODO: Document the source or rationale for the 0.05 (5%) + # gas mass fraction threshold used in classifying gas giant types. + # Cite relevant planetary science literature or model assumptions. and ((planet.gas_mass / planet.mass) > 0.05) and (planet.molec_weight <= lookup_gas('He').weight)): @@ -555,6 +584,11 @@ def generate_planet(protoplanet, ), atm, )''' + # TODO: Review the surface pressure calculation for rocky planets. + # Currently, it's set to a fixed `1 * atm` before `iterate_surface_temp`. + # Evaluate if the commented-out detailed `pressure(...)` calculation + # (based on `volatile_gas_inventory`) should be re-enabled and integrated + # into the iteration loop for a self-consistent atmosphere model. if planet.surf_pressure == 0 * atm: planet.boil_point = 0 * K @@ -576,8 +610,12 @@ def generate_planet(protoplanet, # Next we assign a type to the planet. - # TODO(woursler): Figure out if the pressure unit is atm or millibar. - if (planet.surf_pressure < 1 * millibar): + # TODO: Critical: Ensure consistent units for `planet.surf_pressure`. + # It's initialized in `atm` (e.g. `1 * atm` or `INCREDIBLY_LARGE_NUMBER * atm`) + # but compared with `millibar` values in type classification (e.g. `< 1 * millibar`). + # Convert values to a consistent unit (e.g., millibar or atm) before + # comparison to ensure correct planet type assignment. + if (planet.surf_pressure < 1 * millibar): # Assuming surf_pressure is converted to mbar for this block if (not is_moon) and planet.mass < ASTEROID_MASS_LIMIT: planet.type = PlanetType.ASTERIODS else: @@ -615,7 +653,10 @@ def generate_planet(protoplanet, elif (planet.surf_temp < FREEZING_POINT_OF_WATER): planet.type = PlanetType.ICE else: - # TODO(woursler): Consider throwing an error here. + # TODO: Implement error handling for unclassified planets. + # If planet type remains `UNKNOWN` after classification logic, + # raise an error or log a critical warning to identify and address + # gaps in the classification rules. planet.type = PlanetType.UNKNOWN '''if (flag_verbose & 0x0001) fprintf (stderr, "%12s\tp=%4.2Lf\tm=%4.2Lf\tg=%4.2Lf\tt=%+.1Lf\t%s\t Unknown %s\n", @@ -656,9 +697,12 @@ def generate_planet(protoplanet, moon.orbit = Orbit(a=moon_a, e=moon_e) else: - # TODO(woursler): This seems like an error? - # The moon has no stable space to orbit. - moon.orbit = Orbit(a=0, e=0) + # TODO: Address unstable moon orbits. If + # `roche_limit_r * 1.5 >= hill_sphere_r / 2`, setting orbit to a=0, e=0 + # is not a physical solution. Consider alternative outcomes: moon is + # destroyed (becomes rings), ejected, or crashes into the planet. + # Log this event clearly. + moon.orbit = Orbit(a=0*km, e=0) # Representing no stable orbit found planet.moons.append(moon) @@ -713,7 +757,10 @@ def generate_planet(protoplanet, print(system) - # Output SVG TODO(woursler): Configure as flags? + # Output SVG + # TODO: Make SVG output configurable. Add command-line flags + # (e.g., using `argparse`) to control whether the SVG visualization is + # generated and to specify the output path. max_x = 1500 max_y = 120 diff --git a/garnets/stellar_system.py b/garnets/stellar_system.py index 368a832..569d4cd 100644 --- a/garnets/stellar_system.py +++ b/garnets/stellar_system.py @@ -41,7 +41,10 @@ def mass(self): @property # Approximates the luminosity of the star. - # TODO: express only as ratio? + # TODO: Clarify if this property should return a pre-calculated ratio or if the + # current calculation (which results in a de facto ratio to Solar luminosity) + # is the intended design. If the latter, consider renaming to + # `calculated_luminosity_ratio` or adding a comment explaining units. # Source: http://en.wikipedia.org/wiki/Mass%E2%80%93luminosity_relation def luminosity_ratio(self): if (self.mass_ratio < .43): @@ -55,12 +58,18 @@ def luminosity_ratio(self): return 3200 * self.mass_ratio @property - # Source: StarGen, TODO Verify against current data. + # Source: StarGen + # TODO: Verify the formula `200 * (self.mass_ratio**(1/3)) * au` for stellar + # dust limit against current astrophysical models and data, as its source is StarGen. def stellar_dust_limit(self): return 200 * (self.mass_ratio**(1/3)) * au @property - def r_ecosphere(self): # Source: StarGen, TODO Name? Value? Possible habitable zone? + # Source: StarGen + # TODO: Verify the term 'r_ecosphere' and its calculation + # `sqrt(self.luminosity_ratio) * au` against modern definitions of the habitable + # zone center. Confirm if this simplified model is appropriate for the simulation's scope. + def r_ecosphere(self): return sqrt(self.luminosity_ratio) * au @property @@ -72,7 +81,11 @@ def max_r_ecosphere(self): return sqrt(self.luminosity_ratio / 0.48) * au @property - def life(self): # Source: StarGen, TODO Name? Value? + # Source: StarGen + # TODO: Verify the formula `10**10 * (self.mass_ratio / self.luminosity_ratio)` + # for stellar lifetime against current astrophysical models. Ensure the + # implicit unit (years) is clear or made explicit by returning a Quantity. + def life(self): return 10**10 * (self.mass_ratio / self.luminosity_ratio) @property @@ -112,7 +125,11 @@ def apoapsis(self): return (1 + self.e) * self.a -STAR_MASS = 1 * solar_mass # TODO(woursler): Use the actual star mass +# TODO: Critical: Replace global `STAR_MASS` constant with the actual mass of +# the star associated with the planetoid (e.g., `self.disk.star.mass` for +# Planetesimals or `self.star.mass` for Protoplanets) in the `reduced_mass` +# calculation for accuracy. This constant is a placeholder. +STAR_MASS = 1 * solar_mass @attrs @@ -130,19 +147,24 @@ def reduced_mass(self): # To understand what this is all about... # http://spiff.rit.edu/classes/phys440/lectures/reduced/reduced.html # But some sort of 3 body case, see dole. - # TODO: Understand better? - # TODO(woursler): Actually use the mass of the star? + # TODO: Investigate and document the origin and physical basis of the + # `reduced_mass` formula `(self.mass / (STAR_MASS + self.mass))**0.25`, + # particularly the 0.25 exponent. Cite Dole or other relevant sources. + # This is also linked to the TODO regarding the global `STAR_MASS`. return (self.mass / (STAR_MASS + self.mass))**0.25 @property def inner_effect_limit(self): - # TODO(woursler): Should extend further in, include dust destabilized - # by indirect gravitational effects. + # TODO: Enhance `inner_effect_limit` model. Consider extending the limit + # further inward to account for dust destabilization by indirect gravitational + # effects, beyond the current geometric calculation. return self.orbit.a * (1.0 - self.orbit.e) / (1.0 + DISK_ECCENTRICITY) @property def outer_effect_limit(self): - # TODO(woursler): See note on inner_effect_limit + # TODO: Enhance `outer_effect_limit` model, similar to `inner_effect_limit`. + # Consider effects beyond the current geometric calculation (e.g., indirect + # gravitational perturbations) that might extend the zone of influence. return self.orbit.a * (1.0 + self.orbit.e) / (1.0 - DISK_ECCENTRICITY) @@ -154,11 +176,14 @@ class Planetesimal(Planetoid): def critical_mass(self): perihelion_dist = self.orbit.a * (1 - self.orbit.e) temp = perihelion_dist * sqrt(self.disk.star.luminosity_ratio) - # TODO(woursler): Understand the basis for this. + # TODO: Document the origin and theoretical basis of the `critical_mass` + # formula `B * (dimensionless_with_units(temp, au)**-0.75) * solar_mass`. + # Cite the relevant planetary formation theory or source (e.g., Dole, Hayashi, etc.) + # for this specific calculation involving constant B. return B * (dimensionless_with_units(temp, au)**-0.75) * solar_mass -# TODO(woursler): Migrate to xatu core... +# TODO: Refactor `mass_repr` by migrating its functionality to `xatu.core`. # quantity_repr(mass, {lunar_mass, earth_mass, jupiter_mass, solar_mass}) # or even quantity_repr(mass, celestial_mass_units) def mass_repr(mass) -> str: @@ -204,9 +229,20 @@ class Protomoon(Planetoid): protoplanet = attr() -# TODO(woursler): Go over these with a fine tooth comb. Many are not relevant, or only relevant during initialization. -# Many should be properties. -# Many should be initialized differently. +# TODO: Major refactor for `Planet` class attributes: +# 1. Review each attribute currently defaulted to zero or an empty factory +# (e.g., `moon_a`, `core_radius`, `density`, `orb_period`, temperatures, etc.). +# 2. Identify attributes that should be calculated as properties from other +# more fundamental attributes (e.g., `density` from mass and radius). +# 3. Determine which attributes are true state variables requiring storage +# versus derived values that can be computed on-demand. +# 4. For attributes primarily used during the generation/initialization phase, +# assess if they need to be stored on the `Planet` instance long-term or +# if they can be transient. +# 5. Ensure appropriate default values (e.g., `None` or unit-consistent zeros +# like `0*km` or `0*K`) for attributes that are not immediately known or set +# at instantiation, rather than using plain scalar `0` where it might be +# misleading for quantity types. @attrs(repr=False) class Planet(): # Orbital details. @@ -219,7 +255,6 @@ class Planet(): moons = attr(factory=list) - # ZEROES start here -- TODO(woursler): A bunch of these should be other Zero-like types. gas_giant = attr(default=False) # TRUE if the planet is a gas giant # semi-major axis of lunar orbit moon_a = attr(factory=lambda: 0*au, repr=quantity_formatter(au)) @@ -268,8 +303,6 @@ class Planet(): atmosphere = attr(default=None) type = attr(default=PlanetType.UNKNOWN) # Type code - # ZEROES end here - def __repr__(self): if self.atmosphere is None: atmosphere_string = "No Atmosphere"