Оптимизация интерполяции в Mathematica

В рамках моей работы мне часто приходится визуализировать сложные трехмерные плотности. Один набор программ, с которым я работаю, выводит радиальную составляющую плотности в виде набора из 781 точки на логарифмической сетке, ri = (Rmax/Rstep)^((i-1)/(pts-1), умноженной на сферическую гармонику. Для систем с низкой симметрией число сферических гармоник может быть достаточно большим для обеспечения точности, например. одна система требует 49 гармоник, соответствующих lmax = 6. Таким образом, чтобы использовать эти данные в Mathematica, у меня будет сумма до 49 интерполированных функций, каждая из которых умножается на другую сферическую гармонику. При использовании v.6 и построении интерполированных радиальных функций с использованием Interpolation и установки r = Sqrt(x^2 + y^2 + z^2) я останавливал ContourPlot3D спустя более часа, когда ничего не отображалось. Это включало уменьшение как InterpolationOrder, так и MaxRecursion до 1.

Появилось несколько альтернатив:

1. Оцените функцию плотности на фиксированной сетке и используйте вместо нее ListContourPlot.
2. Или линейно сплайнируйте радиальную функцию и используйте Piecewise, чтобы сшить их вместе. (Это представилось само собой, так как я мог использовать упрощение, чтобы уменьшить сложность результирующей функции.)

В итоге я использовал обе, так как InterpolatingFunction дает заметную задержку в оценке, а при оценке до 49 интерполируемых функций любая задержка может стать заметной. Кроме того, ContourPlot3D был быстрее со сплайном, но это не дало мне желаемого ускорения.

Я честно признаюсь, что я не пробовал Interpolation на v.7, и я не пробовал это на моем обновленном оборудовании (G4 против Intel Core i5). Тем не менее, я ищу альтернативы моей текущей схеме; желательно тот, где я могу использовать ContourPlot3D напрямую. Я мог бы попробовать другую форму сплайна, например B-spline, и, возможно, объединить что с UnitBox вместо использования Piecewise.

Редактировать: Просто чтобы уточнить, моя текущая реализация включает в себя создание сплайна первого порядка для каждой радиальной части, умножение каждого из них на соответствующую сферическую гармонику, суммирование и Simplifyобработка уравнений для каждого радиального интервала, а затем использование Piecewise связать их в одну функцию. Итак, моя реализация является полуаналитической в ​​том смысле, что сферические гармоники точны, а числовой является только радиальная часть. Это одна из причин, по которой я хотел бы использовать ContourPlot3D, чтобы воспользоваться полуаналитическим характером данных. Следует отметить, что радиальная сетка достаточно мелкая, чтобы создать хорошее представление радиальной части и ее можно плавно интерполировать. Хотя это дало мне значительное ускорение, когда я писал код, он все еще был медленным для оборудования, которое я использовал в то время.

Итак, вместо использования ContourPlot3D я бы сначала сгенерировал функцию, как указано выше, а затем оценил бы ее на декартовой сетке 80³. Именно данные этого шага я использовал в ListContourPlot3D. Так как это не адаптивная сетка, местами это было слишком конечно, и мне не хватало фич.