exoplanet-dev / jaxoplanet

Astronomical time series analysis with JAX
https://jax.exoplanet.codes
MIT License
32 stars 11 forks source link

Intent of the `central` kwarg in `System.add_body` #190

Closed lgrcia closed 1 month ago

lgrcia commented 2 months ago

@dfm, what is the intent of the central argument in the System.add_body method?

Here is the method:

def add_body(
    self,
    body: Body | None = None,
    central: Central | None = None,
    **kwargs: Any,
) -> "System":
    """Add a body to the system and return a new system

    Args:
        body (Body | None, optional): body to add. Defaults to None.
        central (Central | None, optional): TODO. Defaults to None.

    Returns:
        System: :py:class:`~jaxoplanet.orbits.keplerian.System` with the added body
    """
    body_: Body | OrbitalBody | None = body
    if body_ is None:
        body_ = Body(**kwargs)
    if central is not None:
        body_ = OrbitalBody(central, body_)
    return System(central=self.central, bodies=self.bodies + (body_,))
dgegen commented 2 months ago

On that note: I found it a bit confusing sometime ago that a stack of bodies possibly orbiting different central objects are allowed within a class whose instances have a unique central attribute. Concretely, initialisations like the following:

from jaxoplanet.orbits.keplerian import Body, Central, System, OrbitalBody

central = Central(mass=1.1, radius=0.9)
body = Body(period=10)
orbital_body = OrbitalBody(central, body)

system = System(bodies=(body, orbital_body))
print([b.central == system.central for b in system.bodies])  # [True, False]

Unless there is a specific reason for this, it might make sense to enforce consistency such that all bodies orbit the central attribute of the instance:

class System(eqx.Module):
    """A Keplerian orbital system"""

    central: Central
    _body_stack: ObjectStack[OrbitalBody]

    def __init__(
        self,
        central: Central | None = None,
        *,
        bodies: Iterable[Body | OrbitalBody] = (),
    ):
        self.central = Central() if central is None else central
        self._body_stack = ObjectStack(
            *(
                b if isinstance(b, OrbitalBody) else OrbitalBody(self.central, b)
                for b in bodies
            )
        )
        self.validate_bodies()  # New

    def validate_bodies(self):
        if not all(body.central == self.central for body in self.bodies):
            raise ValueError("All bodies must orbit the provided central body.")

This would then of course also make the central argument in add_body superfluous.

lgrcia commented 2 months ago

@dgegen, actually users are not supposed to use OrbitalBody, it's only a computational tool. Maybe we can turn it into a private class. I'm working on documenting that more properly but I think this is a separate issue.

dfm commented 1 month ago

The main situation where you might want to have a system where each body has a different "effective" central is when you want to use the photoeccentric effect. In that case, you use approximate circular orbits for each body in the system, but each body has mutually inconsistent period/semi-major axis or, equivalently, a different central density.

In early versions of exoplanet I didn't allow users to have different centrals for each body, but it ended up being a practically useful feature, so I wanted it to be possible with jaxoplanet from the start!

lgrcia commented 1 month ago

I think that's great! Thanks for clarifying that.