Radiation forces blocked by particles

[AstroEloy]

I think it would be very interesting if the radiation forces implemented in REBOUNDx would take into account other particles (at least the active bodies) and be blocked according to their radii and geometry with the source. In other words, not allowing the radiation to penetrate the section of certain particles.

[dtamayo]

Hi Eloy,

Thanks for the comment. This is a bit different from what REBOUNDx currently does (just calculating accelerations on individual particles), and presumably would require some specialized ray-tracing methods (?). It's probably outside the scope of what we are likely to implement in REBOUNDx, but I'm happy to answer questions if you wanted to try and implement this for yourself!

[AstroEloy]

One possible approach might incorporate a shadowing effect. This effect would cause the radiation pressure to drop to zero as particles move 'behind' an active particle within the system. This could be achieved by introducing a cylindrical shape with a radius matching that of the shadowing object (the active particle that blocks the radiation), as the particles are spherical in REBOUND. The cylinder's axis would intersect with the source particle, and it would be assumed that this cylinder extends infinitely from the source, in a direction opposite to that of the source.

Radiation forces should only be canceled when the particle is farther away from the source than the shadowing particle. Then it should be checked if the distance from the particle to the axis of the cylinder is less than or equal to the radius of the cylinder (which is the radius of the shadowing particle).

Something like that in Python:

def is_in_shadow(particle_pos, source_pos, shadowing_pos, shadowing_radius):

    height = np.linalg.norm(shadowing_pos - source_pos)

    axis = (shadowing_pos - source_pos) / height
    particle_to_base = particle_pos - source_pos

    distance_to_axis = np.linalg.norm(np.cross(particle_to_base, axis))

    if distance_to_axis <= shadowing_radius:
        position_along_axis = np.dot(particle_to_base, axis)

        if height <= position_along_axis:
            return True
        else:
            return False
    else:
        return False

[dtamayo]

That seems like a reasonable approach depending on the application, but things also get complicated right? Does partial shadowing / the penumbra matter? You also potentially run into issues with the integration timestep to resolve these shadow crossings. If you or someone else wants to do a lot of work with this and would be willing to write it in C and test it, we could add a flag to the radiation_forces effect that turns this check on and off. I'll close this for now, but feel free to ask questions if you are interested in taking this on!