Вероятно, вам стоит задать этот вопрос на стековой бирже Physics or Mathematics.
Я думаю, что это очень сложный вопрос до такой степени, что великие математики изучали стабильность Солнечной системы. С проблемой двух тел все «легко», но проблема трех тел печально известна своим хаотичным поведением (Пуанкаре внимательно изучил ее и в процессе заложил основы качественной теории динамических систем). Если я не ошибаюсь (не стесняйтесь проверить это онлайн), нестабильность орбитальной динамики большого количества тел (большое означает три и более) - это состояние, вероятность возникновения которого очень высока. Между тем, наличие стабильных конфигураций имеет очень низкую вероятность.
Теперь, для так называемых интегрируемых систем («точно решаемых»), таких как n копий несвязанных моделей солнечной / звездной системы, небольшие возмущения с большей вероятностью дадут стабильную динамику из-за теорема Колмогорова-Арнольда-Мозера. Так что я могу сказать, что для вас более вероятно, что вы столкнетесь со стабильностью, если вы сначала настроите тела в своей модели как сравнительно небольшие источники гравитации, вращающиеся вокруг одного значительно более крупного гравитационного источника. Каждое тело имеет одну доминирующую силу от большого источника и множество гораздо меньших возмущений от остальных тел (или усредненных источников вашего алгоритма Barnes-Hut). Если вы рассмотрите только доминирующую силу и отключите возмущения, у вас будет солнечная система с n разделенными системами двух тел (каждое тело следует эллиптическому движению вокруг общего центра тяжести). Если включить возмущения, эта динамика изменится, но имеет тенденцию очень медленно отклоняться от невозмущенной, и, скорее всего, будет стабильной. Итак, начните с упорядоченной динамики и начните немного изменять массы тела, их положения и скорости. Вы можете проследить, как меняется динамика, когда вы меняете параметры и начальные условия.
Еще одна вещь: всегда полезно поместить инерциальную систему координат, относительно которой представлены положения и скорости тел, в центр масс группы тел. Это более или менее гарантировано, когда сумма начальных импульсов равна нулевому вектору. Эта установка дает центр масс системы всегда фиксирован в некоторой точке пространства, поэтому простой перенос переместит его в начало системы координат.