Чем лучше посчитать четырехэтажный неопределенный интеграл?

[AD]

Учебная задачка по численным методам доехала до решения краевой задачи для не очень уж глубокомысленного уравнения $x''''(t)=f(t)$ с краевыми условиями $x(1)=x'(0)=x''(1)=x'''(0)=0$.

Ясно, что решение выписывается в виде

$$\def\d{\,\mathrm{d}} x(t)=\int_t^1\int_0^s\int_r^1\int_0^q f(p)\d p\d q\d r\d s$$

Функция $f:[0,1]\to\mathbb{R}$ бесконечно гладкая, знакопостоянная и монотонная на $[0,1]$, и вычисляется мгновенно.

Ну, собственно, вот и интересуюсь, как нормальные люди считают такие интегралы? Пришло в голову несколько вариантов - и интересно советы послушать. То есть пробовать все методы подряд не хочется, хотя, возможно, было бы полезно.

Мои мысли в хронологическом порядке:

Первая: так же, как и предыдущую задачку из этого практикума - Рунге-Куттом, а потом пристреливаться.

Вторая: Просто прогнать четыре раза какую-нибудь квадратурную формулу.

Третья: Считать этот интеграл кратным, и шарахнуть Монте-Карлой.

Четвертая: Сочинить что-нибудь сеточное типа
$x(t_{i+2})-4x(t_{i+1})+6x(t_i)-4x(t_{i-1})+x(t_{i-2})=f(x_i)$, $i=3,\ldots,N-2$
(ну узлы, скажем, равномерные, и в этом случае я тут просто какой-то множитель забыл типа $N^{-4}$)
$x(t_N)=0$, $x(t_{N-2})-2x(t_{N-1})+x(t_N)=0$,
$x(t_2)-x(t_1)=0$, $x(t_4)-3x(t_3)+3x(t_2)-x(t_1)=0$.

Какой из методов предпочтительнее? Хотелось бы что-нибудь не очень заумное, но чтобы и точность была не совсем уж бестолковая.

[Gafield]

Ну, конктретно такой интеграл имело бы смысл преобразовать к виду

$$ x(t)=\int_0^1K(t,p)f(p)\, dp, $$

где $K(t,p)$ имеет (скорее всего, но я не проверял) простой вид.
Например, если бы все условия были бы с одной стороны, то

$$ \int_0^t\int_0^s\int_0^r\int_0^q f(p)d pd qd rd s=\frac16\int_0^t(t-p)^3f(p)\, dp. $$

[AD]

Ух ты как клёво

Вот что вышло:

$$x(t)=(1-t)\int_0^1\frac{1-r^2}2 f(r)\,dr-\frac{(1-t)^3}6\int_0^1f(p)\,dp+\int_t^1\frac{(s-t)^3}6f(s)\,ds$$

Ну я еще попроверяю, конечно. Прежде всего, при помощи сравнения с каким-нибудь глупым методом, который уже написан

То есть, фактически, надо просто посчитать несколько интегралов вида $\int z^k f(z)\,dz$. И все они разные - одни от нуля, другие от $t$. Но все равно приятно, что погрешность не накапливается, да?

А такие интегралы уже как-нибудь тупо считать, да? Или для них вообще уже всё известно, наверное, есть готовые формулы хорошие?

[TOTAL]

$h=1/N$
$x'''_{i+1}=x'''_{i}+hf(x_{i+1/2})$
$x''_{i-1}=x''_{i}-h(x'''_{i-1}+x'''_{i})/2$
$x'_{i+1}=x'_{i}+h(x''_{i+1}+x''_{i})/2$
$x_{i-1}=x_{i}-h(x'_{i-1}+x'_{i})/2$

[AD]

TOTAL, ну это ж как-то тупо-прямолинейно очень ... Это я и [почти] сам написал При $N=100000$ погрешность чуть ли не в третьем знаке после запятой, а время становится ощутимое. Предыдущую задачу до точности $10^{-7}$ доводил почти мгновенно, а она была сложнее (в смысле уравнение кривее) ...

[TOTAL]

AD писал(а):

При $N=100000$ погрешность чуть ли не в третьем знаке после запятой, а время становится ощутимое.

Ради спортивного интереса и для проверки правильности программы. Какого порядка точности получается метод? Во сколько раз уменьшается погрешность при увеличении $N$ в два раза? (Должна в 4 раза.) Не надо такие большие $N$.

И еще. Пусть этот метод кажется хуже другого из другой задачи. А если сравнивать с другим методом для этой же задачи?

[AD]

TOTAL Н-даа ... Это я, собственно, брал не средние значения, а левые или правые. Теперь при $N=5$ уже неплохо

Добавлено спустя 12 минут 51 секунду:

Сейчас еще какой-нибудь метод напишу что-ли ...

Добавлено спустя 1 минуту 9 секунд:

Метод действительно второго порядка , однако в первой задаче нам велели не ниже четвертого

Добавлено спустя 1 минуту 15 секунд:

А правильно ли я понимаю, что погрешность на предыдущих шагах будет $\overline{\overline{o}}$ от погрешности на последнем шаге, да?

[TOTAL]

AD писал(а):

А правильно ли я понимаю, что погрешность на предыдущих шагах будет $\overline{\overline{o}}$ от погрешности на последнем шаге, да?

На первом шаге погрешность $O(h^2)$. Она так и остается до самого конца.
На первом шаге несложно и более высокий порядок сделать, но вот сохранить его на следующих шагах уже тяжелее. Лично я (меня не заставляют считать с таким-то порядкам) считаю, что методы второго порядка в большинстве случаев самые хорошие. Они несравненно лучше методов первого порядка, а при переходе со второго на третий улучшение уже не столь ощутимое, хотя усложнение алгоритма, наоборот, существенное.

Программируйте альтернативный более крутой метод и устраивайте между ними соревнование.

[AD]

Так, у меня некоторое непонимание. Вот у нас формула Симпсона

$$\int_{x_{2i-2}}^{x_{2i}}f(p)\,dp\sim \frac{x_{2i}-x_{2i-2}}6\bigl(f(x_{2i-2})+4f(x_{2i-1})+f(x_{2i})\bigr)$$

обещает четвертый порядок. Пишу:

// данные берутся из f1, кладутся в f2
void integrate_lr(double *f1, double *f2, long n) {
   long i,k;
   f2[0]=0.;
   for (i=1; i<=n/2; ++i) {
      f2[i]=f2[i-1]+(f1[2*i-2]+4.*f1[2*i-1]+f1[2*i])/(double)(3*(n-1));
   }
}

Запускаю - порядок явно как будто первый. Че-то я не так пишу или так вообще нельзя?

Добавлено спустя 35 минут 38 секунд:

Ну, ладно, надеюсь, дальше сам разберусь, поотлаживаю еще.

Добавлено спустя 49 минут 9 секунд:

Короче, запрогал методом Gafieldа в форме меня, то есть по формуле

Цитата:

$$x(t)=(1-t)\int_0^1\frac{1-r^2}2 f(r)\,dr-\frac{(1-t)^3}6\int_0^1f(p)\,dp+\int_t^1\frac{(s-t)^3}6f(s)\,ds$$

, считая четыре однократных интеграла $\int z^kf(z)\,dz$, $k=0,1,2,3$, Рунге-Куттом 6-го порядка на шаге (то есть пятого в целом). Ответ сошелся .

То есть всем спасибо, буду просветляться дальше.