Численное интегрирование компактной функции на треугольнике

как следует из названия, я пытаюсь вычислить интеграл от компактно поддерживаемой функции (квинтический полином Вендланда) на треугольнике. Обратите внимание, что центр функции находится где-то в трехмерном пространстве. Я интегрирую эту функцию в произвольный, но маленький треугольник ( ). В настоящее время я использую интеграцию, описанную Dunavant, 1985 (p = 19).$area < \frac{(radius/4)^2}{2}$

Однако представляется, что эти квадратурные правила не подходят для задач с компактной поддержкой. Это подтверждается тем фактом, что когда я интегрирую (то есть функцию, которая равна 1 внутри круга радиуса 1) на плоскости, которая дискретизируется с помощью треугольников, мои (нормализованные) результаты находятся между 1,001 и 0,897.$f(r) = [r\leq1]$

Итак, мой вопрос: существует ли специальное квадратурное правило для такого рода проблем? Будет ли лучше работать составное правило интеграции низкого порядка?

К сожалению, эта процедура действительно важна в моем коде, поэтому точность имеет решающее значение. С другой стороны, мне нужно сделать эту интеграцию «пару раз» за один шаг по времени, чтобы вычислительные затраты не были слишком высокими. Распараллеливание не является проблемой, так как я выполню интеграцию сама по себе.

Заранее спасибо за ваши ответы.

РЕДАКТИРОВАТЬ: Quintic многочлен Вендланда дается с и где r_0 - произвольный вектор в \ mathbb {R} ^ 3$W(q) = [q\leq2]\frac{\alpha}{h^3}(1-\frac{q}{2})^4(2q+1)$$\alpha = \frac{21}{16\pi}$$q=\frac{\|r-r_0\|}{h}$$r_0$$\mathbb{R}^3$

РЕДАКТИРОВАТЬ 2: Если $\Delta$ является двумерным треугольником, то я хочу вычислить $\int_\Delta \omega(r) dr$ с $\omega(r) = W(\frac{\|r-r_0\|}{h})$ , Таким образом, $q$ в $W$ никогда не будет меньше 0. Обратите внимание, что интеграл является поверхностным интегралом по 2-D поверхности в $\mathbb{R}^3$

EDIT3: у меня есть аналитическое решение для 1-D (линия) проблемы. Вычисление одного для 2-D (треугольник) также возможно.

Поскольку функция является гладкой в ​​пределах , но не имеет фиксированной степени (то есть в плоскости), я бы предложил использовать простую адаптивную схему, например правило трапеции с методом Ромберга , в обоих измерениях.$q \leq 2$

То есть, если ваш треугольник определен вершинами , и , и у вас есть подпрограмма, которая интегрируется вдоль линии от до , вы можете сделать следующее (в обозначениях Matlab):$x$$y$$z \in \mathbb R^3$romb(f,a,b)fab

int = romb( @(xi) romb( W , xi , y+(z-y)*(xi-x)./(z-x) ) , x , z );

В romb, не используйте фиксированное количество очков, но продолжает расти таблицы , пока разница между двумя последовательными диагоналями не ниже вашего требуемым допуском. Поскольку ваша функция гладкая, это должно быть хорошей оценкой ошибки.

Если части треугольника находятся за пределами области , вы можете попытаться соответствующим образом настроить пределы интегрирования в приведенном выше коде.$W(q)$

Возможно, это не самый эффективный в вычислительном отношении способ решения вашей проблемы, но адаптивность даст вам гораздо большую устойчивость, чем правило с фиксированной степенью.

Для хорошего обзора кубатурных правил см. "Р. Кулс, Энциклопедия кубатурных формул J. Сложность, 19: 445-453, 2003". Использование фиксированного правила может дать вам преимущество, заключающееся в том, что некоторые правила точно интегрируют полиномы (как гауссовская квадратура в одномерном случае).

Cools также является одним из основных авторов CUBPACK , программного пакета для числовой кубатуры.

Правила интегрирования предполагают, что функция локально хорошо аппроксимируется полиномом низкой степени. Ваша проблема не имеет ничего общего с компактной поддержкой. Компактно поддерживаемые радиальные базисные функции являются гладкими на границе опоры, и квадратурные правила вплоть до порядка гладкости могут использоваться без проблем. (Правила более высокого порядка не помогают; поэтому вам, вероятно, не следует использовать правило, которое точно интегрирует полиномы степени 5.)

В вашем случае неточность проистекает из того факта, что предположение о хорошей полиномиальной аппроксимируемости не выполняется в вашем случае для треугольников вблизи , даже если они не содержат .$r_0$$r_0$

$W$ является гладкой как функция от , но является негладкой функцией от с градиентом, который становится бесконечным в пределе . Интеграция по , и составная функция является негладкой функцией .$q$$q$$r$$r\rightarrow r_0$$r$$r$

Если треугольник не содержит , функция является но это не помогает, поскольку старшая производная очень быстро растет вблизи , а ошибка метода высокого порядка пропорциональна производной высокого порядка, следовательно, очень велика !$r_0$$C^inf$$r_0$

Простое решение состоит в том, чтобы разбить каждый треугольник T на число N_T из подтреугольников. Вы можете взять далеко от и близко к . Вы можете в автономном режиме выяснить, насколько большим должен быть для треугольников заданного диаметра и расстояния от чтобы достичь желаемой точности. Кроме того, вы должны использовать только формулы низкого порядка, близкие к .$N_T=1$$r_0$$N_T\gg 1$$r_0$$N_T$$r_0$$r_0$

Поскольку вы интегрируете по треугольнику, но является 3-мерным, треугольник, очевидно, находится в .$r_0$$R^3$

Таким образом, более быстрое лекарство сведет в таблицу интеграл для в зависимости от координат треугольника (нормализуется поворотом его в двумерную плоскость , что одна вершина лежит на оси , и отражением его так, чтобы второй вершина лежит над ней). Эта таблица должна быть достаточно подробной, чтобы сделать линейную или квадратичную интерполяцию достаточно точной. Но вы можете использовать медленный метод, описанный первым, чтобы создать эту таблицу.$r_0=0$$xy$$x$

Другой способ избавиться от этой проблемы - использовать компактную радиальную базисную функцию, которая является полиномом от а не . Это везде гладко и легко интегрируется.$q^2$$q$