численное дифференциирование

Если важна точность, то лучше найти такое решение где не требуется дифференцировать.

Увы, тут или очень большая неуниверсальность или дифференциирование.

Впрочем, 5 точек + 16 байтного long double хватает

численно дифференциировать не по 3-м, а по пяти точкам?

а кстати как ? Интерполировать - интерполировал..но чтобы дифференцировать по 3-5 точкам - ни разу не пришлось.

Что имеется в виду ? поиск N-ной производной в точке по K ближайшим измерениям ??

Аппроксимация кривой не параболой, а полиномом 4 порядка.

Коэффициенты для 3 точек: -1/(2*H), 0, 1/(2*H)

Коэффициенты для 5 точек: 1/(12*H), -2/(3*H), 0, 2/(3*H), -1/(12*H)

3-х точечной производной хватает для 3-х оборотов вокруг черной дыры по фотонной сфере. 5-точечной для 9 оборотов (тут еще long double понадобился). Использовал рунге-кутта 4 порядка

можешь апроксимировать любым полиномом, пока есть физ.смысл.

идёшь к научруку(обкладываешься книгами и решаешь сам) и выясняешь может ли быть сходство аппроксимируемой кривой с полиномом N-ной степени в заданной области и с какими ограничениями. Может у вас там тригонометрические кривые должны выходить - я же незнаю вашей области.

Для данного примера должна быть окружность. Что я и наблюдаю первые несколько оборотов. А потом луч либо падает, либо улетает

Но тут решение нестабильно, так что уточнять можно до бесконечности.

окружность (или около того) - тригонометрическая кривая. Её вроде можно описать/приблизить полиномом только перейдя в сферичекие координаты или в параметрические кривые(в чем неуверен - давненько не занимался подобным).

Я в сферических и считаю. Но любую гладкую кривую можно описать полиномом в малой окрестности некоторой точки. На этом численное дифференциирование и основано by design

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

Можешь попробовать взять библиотеку QD и посчитать с 256-битными float-ами. (http://crd-legacy.lbl.gov/~dhbailey/mpdist/)

выбор дискретности при моделировании - та ещё задачка :)

зы Еще раз намекну - тригонометрические кривые не могут интер/экстра-полироваться степенными полиномами, иначе ты слишком быстро накопишь погрешности, поймаешь циклические микро-процессы(которых на самом деле нет) внезапные особые точки и 100500 всякой всячины. Есть исключения, но очевидно не в вашем случае

зыы не настолько помню физику чтоб понять «3-х точечной производной хватает для 3-х оборотов вокруг черной дыры по фотонной сфере. 5-точечной для 9 оборотов (тут еще long double понадобился)» и хоть убейте меня не знаю понятия 3 или 5 точечной производной.

точности double иногда не хватает

а если использовать длинную арифметику?

и хоть убейте меня не знаю понятия 3 или 5 точечной производной.

По трех точкам:

f'(x) = (f(x+h)-f(x-h))/2h, центральная точка оказывается с нулевым коэффициентом

По пяти точкам: f'(x) = (f(x-2h)-8f(x-h)+8f(x+h)-f(x+2h))/12h, центральная тоже с нулевым коэффициентом

зы Еще раз намекну - тригонометрические кривые не могут интер/экстра-полироваться степенными полиномами

Иначе в общем случае никак.

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

а шаг по времени какой?

По моему параметру 10^-2. Правда это не совсем время.

Хм, мне пока трех хватает, с пятью точками на границах будет больше возни, матрицы получаются не с 3-мя, а с 5-ю диагоналями, в итоге медленнее все будет. Ну а чтобы побороть ограничение на шаг по времени (точнее dt/dx ну или dt/dx^2) использую схемы подобные схеме Кранка-Николсона, то есть точность я повышаю просто уменьшением шага по dx и при этом мне не надо сильно уменьшать dt. Как-то так.

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

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

Численные методы используют те, кто не осилил аналитические.

Множество задач не имеют аналитического решения.

Ну так на самом деле и есть. И из-за шумов луч и должен уходить. Но производная по 5 точкам продлила путь луча по сфере в 3 раза

Вопрос тем, кто здесь занимался численными методами, вообще часто приходится численно дифференциировать не по 3-м, а по пяти точкам? И если часто, то почему об этом даже не упоминают нигде? Хотя вывод формулы по силам первокурснику. Но про три точки везде пишут, а про 5 нигде

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

Чем выше порядок точности схемы, тем медленнее она работает. Поэтому, избыточную точность, естественно, никто не использует.

По трех точкам:

f'(x) = (f(x+h)-f(x-h))/2h, центральная точка оказывается с нулевым коэффициентом

Зачем извращениями страдать? Есть же предиктор-корректор и другие хорошие схемы второго порядка. Центральная разность - плохая схема.

И если часто, то почему об этом даже не упоминают нигде?

http://www.sciencedirect.com/science/article/pii/S002199918571042X

Численные методы используют те, кто не осилил аналитические.

Да, да. Реши мне уравнения Навье-Стокса аналитическим методом, осилятор.

Центральная разность даёт нефизическое решение (в смысле, оно в численной схеме есть, но в физической системе его нет). Поэтому эту схему нигде нынче не используют.

Это свойство называется «немонотонностью» (пара разных определений) или «положительностью». Обладают ей все линейные схемы выше первого порядка точности.

Это какие именно, где про них почитать? И чем они лучше расчета по 5 точкам?

Если что, вот мой код расчета:

tensor2 dGi(Lvector p, int n)
{
	int L = p.dim();
	tensor2 res(L);
	double coef[5] = {1/(12*H), -2/(3*H), 0, 2/(3*H), -1/(12*H)};
	Lvector  h(L);
	int i;
	for (i = 0; i < L; i++)
	{
		h[i] = (i==n)*H;
	}

	for (i = 0; i < 5; i++)
	{
                if (i == 2)
                   continue;
		Lvector dp = h * (i-2);
		tensor2 metr = Metric(p + dp);
		res += metr*(coef[i]);
	}


	return res;
}

По моему параметру 10^-2. Правда это не совсем время.

уравнение какое?

http://upload.wikimedia.org/math/d/f/9/df9964e9250e597ed0f1f23f3d1ddb21.png

только там t - не время

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

Чем выше порядок точности схемы, тем медленнее она работает. Поэтому, избыточную точность, естественно, никто не использует.

ЕМНИП есть теорема, которая указывает, что при порядке выше четвертого численный алгоритм будет обязательно нестабильным. Если интересно, могу найти.

Центральная разность даёт нефизическое решение (в смысле, оно в численной схеме есть, но в физической системе его нет).

фигня какая-то.

Поэтому эту схему нигде нынче не используют.

может быть и не используют, но уж точно не поэтому.

А что это за Г?

символ Кристофеля

http://upload.wikimedia.org/math/5/1/3/513b1509cc346cb0522b374c4916b676.png

а почему не взять какой-нибудь из спектральных методов вместо конечных разностей? у вас же на окружности все происходит? периодические ГУ => сам бог велел зафигачить производную через фурье. при переходе от трех-точечной схемы к 5-точечной вы получите всего +2 порядка точности...

найдите. мне интересно, потому что мне так сдается, что это неправда.

Это какие именно, где про них почитать?

Второго порядка схемы? Да их полно. Какую именно использовать, зависит от конкретной задачи.

И чем они лучше расчета по 5 точкам?

Тем что они 2-го порядка точности. Т.е. меньше вычислений требуют.

Не в лотерею, а в карты, но теорема Годунова всё же.

может быть и не используют, но уж точно не поэтому.

В газодинамике именно поэтому. Центральную разность (leap frog) использовали довольно долго, но решение приходилось фильтровать. Сейчас есть более хорошие схемы.

Если интересно, могу найти.

Буду благодарен.

Не в лотерею, а в карты

что? о_0

вот полез в вики. там написано:

The theorem states that:

Linear numerical schemes for solving partial differential equations (PDE's), having the property of not generating new extrema (monotone scheme), can be at most first-order accurate.

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

у вас же на окружности все происходит?

Это очень частный случай, выбранный для проверки точности вычислений.

В общем случае все происходит в 4-мерном кривом пространстве-времени

что? о_0

Классика же

— Правда ли, что шахматист Петросян выиграл в лотерею тысячу рублей?

— Правда, только не шахматист Петросян, а футболист «Арарата» Акопян, и не тысячу, а десять тысяч, и не рублей, а долларов, и не в лотерею, а в карты, и не выиграл, а проиграл.

Соотносится это очевидно по причинам идиотичности утверждения дикого, которому про схемы Рабинович напел.

кстати, отвечая на вопрос — почему про 5-точечные схемы нигде не пишут. наверное, правильный ответ примерно такой: 3 точки достаточно, чтобы изложить концепцию конечных разностей, показать, в чем состоит аппроксимация и пояснить, откуда берется ошибка и насколько она велика. 5-точечная схема концептуально ничем не отличается, уравнения в конечных разностях выписываются точно также, ошибка считается точно также (и — как вы и отметили — вывод под силу первокуру). так зачем тогда про эту схему писать?

Пример можешь привести? Гугл мне упорно пытается подсунуть вычисление производных второго порядка.

найдите. мне интересно, потому что мне так сдается, что это неправда.

действительно - неправда.

там про ограничение Dahlquist'а. Про ограничение на порядок приближения производной в зависимости от числа «шагов» алгоритма.

В газодинамике именно поэтому. Центральную разность (leap frog) использовали довольно долго, но решение приходилось фильтровать.

а дифур какой был? Навье стокса? Так наверное потому фильтровать приходилось, что харакеристики дифура лежали за передлами численных характеристик метода.

понятно. тем не менее если для вас точность настолько важна, все-таки обратите внимание на спектральные методы. они сформулированы и для не-периодических ГУ... (хотя мне по-прежнему кажется, что по угловой переменной они у вас должны быть периодическими :))

adams-bashforth-moulton?

тока сча дошло :)

Гугл мне упорно пытается подсунуть вычисление производных второго порядка.

Так это оно и есть. Целое семейство схем предиктор-корректор есть, всякие upwind схемы. Каждая хороша для своей задачи.