2014 dxdy logo

Научный форум dxdy

Математика, Физика, Computer Science, Machine Learning, LaTeX, Механика и Техника, Химия,
Биология и Медицина, Экономика и Финансовая Математика, Гуманитарные науки




Начать новую тему Ответить на тему На страницу 1, 2  След.
 
 Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 11:50 


10/11/08
35
Одесса, ОНУ, ИМЭМ
Возникла задача расчётов интегралов такого типа:
$\int\limits_{0}^{\pi}\ln\left\lvert T_n(t)\right\rvert dt$
Где $T_n(t)$ тригонометрический полином большой степени (и с большими коэффициентами) с некоторым числом нулей на отрезке, часто включая его границы. Примерное их расположение известно.

Проблема заключается в том, что при $n$ порядка тысяч NIntegrate начинает считать очень долго, и не всегда даёт ответ.
Пробую подбирать параметры примерно таким образом,
PrecisionGoal -> 25, WorkingPrecision -> 200, Method -> \
{"GlobalAdaptive", "SingularityHandler" -> "DoubleExponential"}
иногда добавляя Maxrecursions, но это не всегда помогает, и интегралы часто считаются по полчаса и больше.

Подскажите, какие параметры помогут подсчитать интегралы максимально эффективно (достаточно в принципе 10 гарантированно точных знаков).

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 11:58 
Аватара пользователя


11/06/12
10390
стихия.вздох.мюсли
Для начала непонятно следующее: зачем выставлять WorkingPrecision -> 200, если вам не нужна такая точность?

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 12:20 


10/11/08
35
Одесса, ОНУ, ИМЭМ
Мне нужна итоговая точность в 10 знаков. Значения подынтегральной функции сильно колеблются. Разница в значении коэффициентов как раз все эти 200 знаков и могут накрывать спокойно. Проблема в сингулярностях, где всё сокращается. Там нужна точность очень большая, максимальные коэффициенты полинома ведут себя как $C_{2n}^n$. Хочется как-то указать, на какие места стоит обратить внимание и каким образом.

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 12:30 
Аватара пользователя


11/06/12
10390
стихия.вздох.мюсли
Почитайте вот это. Авось извлечёте что-то полезное. Я сейчас ещё попробую в справке порыться. Кажись, видел там когда-то что-то этакое.

UPD. В справке нашёл только это.

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 12:35 


10/11/08
35
Одесса, ОНУ, ИМЭМ
Aritaborian в сообщении #980463 писал(а):
Почитайте вот это. Авось извлечёте что-то полезное. Я сейчас ещё попробую в справке порыться. Кажись, видел там когда-то что-то этакое.

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

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 12:42 
Аватара пользователя


11/06/12
10390
стихия.вздох.мюсли
Xaliuss, на всякий случай посоветую вам также обратиться за помощью в VK-сообщество. Там тусуются отличные специалисты.

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 12:53 


10/11/08
35
Одесса, ОНУ, ИМЭМ
Aritaborian, написал. Параметров к NIntegrate очень много, так просто по справке понять какие нужны не выходит.

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 15:28 
Заслуженный участник


25/02/11
1786
Насколько много там нулей и реально ли найти их все с большой точностью? Если не очень много и да, можно попробовать разбить отрезока на интервалы между соседними нулями. Есть еще такой прием: сделать замену, чтобы особенность исчезла. Например, $\int_0^1 \ln x\,dx$ после замены $x=y^2$ становится $4\int_0^1y\ln y\,dy$.

Второй способ: воспользоваться тем, что под логарифмом тригонометрический многочлен. Записать через мнимые экспоненты и разложить на множители. Типа так:
Код:
p[y_] = T[x] /. {Cos[a$_. x] -> (y^a$ + y^-a$)/2} // Together
pn[y_] = Numerator@p[y]
pd[y_] = Denominator@p[y]
NSolve[pn[y] == 0, y, WorkingPrecision -> 20]

А потом считать для каждого корня $c_i$ слагаемое $\int_0^\pi\ln|e^{ix}-c_i|\,dx$.

Если многочлен только от косинусов, т.е. разложится по степеням $e^{2ix}$, некоторые интегралы, похоже, будут давать ноль: $\int_0^\pi\ln|e^{2ix}-a|\,dx=\frac12\int_0^{2\pi}\ln|e^{ix}-a|\,dx=0$ при $|a|\le1$. Именно те, что имеют особенность.

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 16:08 


10/11/08
35
Одесса, ОНУ, ИМЭМ
Нулей не очень много, считает их долго. В одном важном случае только в 0 и $\pi$. Пробовал разбивать в случае одного промежуточного нуля, кардинально ситуацию не улучшило. Замену так просто сделать в моём случае нельзя.

Второй приём это по сути формула Иенсена. Найти все корни многочлена степени 1000+ (с БОЛЬШИМИ коэффициентами) с нормальной точностью - не вариант. Но с разложением на множители идея, может можно будет найти асимптотику около нулей, и в окрестности них делить на неё.

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

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 16:34 
Заслуженный участник


25/02/11
1786
Так многочлен по косинусам? И коэффициенты точно заданы (целочисленные)?
Xaliuss в сообщении #980509 писал(а):
Интересно есть ли метод, который позволит указать напрямую как считать в окрестности нулей, так как поведение около них предсказуемое.

Типа $\ln|x-c_i|$? Можно просто вычесть, а интеграл с логарифмом посчитать отдельно.

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 16:50 


10/11/08
35
Одесса, ОНУ, ИМЭМ
Есть и по косинусам, есть и по синусам. Коэффициенты целочисленные, но толку от этого если они как $C_{2n}^n$ примерно растут.

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

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 16:52 
Аватара пользователя


11/06/12
10390
стихия.вздох.мюсли
Возможно я сейчас глупость скажу, но что, если сначала найти нули, а затем явно их указать, разбивая область на промежутки?

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 17:04 


10/11/08
35
Одесса, ОНУ, ИМЭМ
Aritaborian в сообщении #980516 писал(а):
Возможно я сейчас глупость скажу, но что, если сначала найти нули, а затем явно их указать, разбивая область на промежутки?


Этого мало. Делал так. Хотя вроде помогает, но нужно что-то добавить.

Вот так примерно выглядит график логарифма для одного частного случая

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 17:15 
Заслуженный участник


25/02/11
1786
Xaliuss в сообщении #980509 писал(а):
Замену так просто сделать в моём случае нельзя.

Когда только на концах нули, разбить пополам и в каждой половине замену. Если только из-за особенности все тормозит, то должно помочь. Аналогично можно разбивать отрезок между нулями.

Еще можно попробовать изоображать на комплессной плоскости множества нулей для небольших $n$. Мб это что подскажет. Если будут только четные степени $y=e^{i x}$, то остается вклад только тех нулей, которые вне единичного круга. Вдруг их мало, симметрии какие есть или асимптотически куда-то кучкуются.

 Профиль  
                  
 
 Re: Расчет сингулярных интегралов в Mathematica 10
Сообщение20.02.2015, 21:44 
Аватара пользователя


15/01/06
200
Xaliuss, а можно привести пример кода хотя бы одной функции, которую нужно проинтегрировать и которая долго интегрируется?

 Профиль  
                  
Показать сообщения за:  Поле сортировки  
Начать новую тему Ответить на тему  [ Сообщений: 20 ]  На страницу 1, 2  След.

Модераторы: Karan, Toucan, PAV, maxal, Супермодераторы



Кто сейчас на конференции

Сейчас этот форум просматривают: нет зарегистрированных пользователей


Вы не можете начинать темы
Вы не можете отвечать на сообщения
Вы не можете редактировать свои сообщения
Вы не можете удалять свои сообщения
Вы не можете добавлять вложения

Найти:
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group