Kirill Pankratov (neznaika_nalune) wrote,
Kirill Pankratov
neznaika_nalune

Category:

Немного прикладной математики

Довольно регулярно - в среднем 1-2 раза в месяц - я до сих пор получаю е-мейлы от пользователей моих математических программ написанных на Матлабе, на досуге, около 15 лет назад. У некоторых из них было десятки тысяч (или больше) пользователей, они разошлись по миру, кем-то переписывались и видоизменялись. Но большинство из них уже не работают в их изначальной версии - из-за несовместимости с нынешней версией Матлаба. Основные вопросы которые мне задают, связаны с этим.

Но интересно другое - то что нередко мои программы написанные лет 15 назад, остаются state-of-the-art - всё равно лучше по многим параметрам, чем почти всё что написано на эту тему на Матлабе с тех пор. На днях получил е-мейл насчёт программы для численного интегрирования, которую я запостил в матлабовскую ньюсгруппу ещё в 1995, потом кто-то в 1999 залил на сайт Mathworks:

Dear Dr. Kirill Pankratov,
I am a PhD student at MIT .

I was looking at your quads.m and quadr.m files online from
http://www.mathworks.com/matlabcentral/linkexchange/links/41-quadrmquadsm

I get the following error message when I try to run then, I am not sure why they are not working.

[Q,xo,fo,s] = quads('1/sqrt(x)/(x-2)/(x-4)^2',0,5);
??? Subscript indices must either be real positive integers or logicals.
Error in ==> lininsrt at 75
xo(mask,:) = x;
Error in ==> quads at 160
[x_o,mask] = lininsrt(x_o,F);

Is it because you wrote these functions in 1999 and now I am using a newer Matlab version?
I really appreciate your help,
Thank you very much,
G*******


Чуть-чуть покопался, изменил что надо - работает. Послал просившему, тот оказался вполне удовлетворён.
Заодно посмотрел как она ведёт себя по сравнению с тем что в Матлабе имеется сейчас. Выяснилось, что, хотя за эти годы программы численного интегрирования улучшились, они недотягивают до того что я написал 15 лет назад.
Мои программы были особенно полезны для "трудных" функций - в том числе с высокочастотными колебаниями, очень широким динамическим диапазоном и сингулярностями.

Как правило мои программа одновременно быстрее, точнее и надёжнее чем собственная матлабовская
Например, моя функция quadr запросто может проинтегрировать 10 тысяч синусоид (точный ответ - 2.0):

>> quadr('sin(x)',0,10001*pi,1e-6)
Recursion level limit reached. Singularity likely.

ans =
2.0001

В то время как родную матлабовскую quad заносит на виражах в бесконечность даже на 1/25 этого диапазона:
>> quad('sin(x)',0,401*pi,1e-6)
Warning: Maximum function count exceeded; singularity likely.
> In quad at 92

ans =
-105.5199

Моя функция может даже автоматически найти, проанализировать и, если возможно, проинтегрировать сингулярности. Например, она может определить является ли сингулярность интегрируемой, интегрируемой в смысле главного значения (типа 1/х) или неинтегрируемой. Ничего подобного матлабовские родные программы не делают даже близко.

Из описания прилагаемого к моим программам:

>> [Q,xo,fo,s] = quads('1/x',-1,2);
Singularity at x = -0.000000, integrable in a principal value sense

Check the accuracy:

>> Q/log(2)-1

ans =
-1.1324e-06

(LOG(2) is a value of integral from 1 to 2. Integral from -1
to 1 is zero in a principal value sense).
...
Let's try a function with several singularities of different
types:

>> [Q,xo,fo,s] = quads('1/sqrt(x)/(x-2)/(x-4)^2',0,5);
Singularity at x = 0.000008, integrable
Singularity at x = 2.000000, integrable in a principal value sense
Singularity at x = 4.000000, non-integrable

Now try a function which is not strictly polynomial or
rational:

>> [Q,xo,fo,s] = quads('x*tan(x)',0,pi*4);
Singularity at x = 1.570796, integrable in a principal value sense
Singularity at x = 4.712389, integrable in a principal value sense
Singularity at x = 7.853982, integrable in a principal value sense
Singularity at x = 10.995574, integrable in a principal value sense

One can easily check that singularities are at the points
pi/2, 3*pi/2, 5*pi/2, 7*pi/2:

>> s(1,:)/pi

ans =
0.5000 1.5000 2.5000 3.5000


По-моему, круто! :)
Subscribe

  • Пермь больного человека

    Заметка в "Ленте" Добыча нефти в США обернулась экологической катастрофой по мотивам статьи в "Гардиан" How Texas’s zombie oil wells are creating…

  • "Going postal" - Кучно пошло

    Очередные пострелушки - теперь в Индианаполисе, в офисе Fedex, 8 убитых.…

  • Апрельская зима по расписанию

    16 апреля 2021 - более 10 см снега. Повторяется прошлогодняя история почти день в день, тогда было 18 апреля. И примерно так почти каждый год в…

  • Post a new comment

    Error

    default userpic

    Your reply will be screened

    Your IP address will be recorded 

    When you submit the form an invisible reCAPTCHA check will be performed.
    You must follow the Privacy Policy and Google Terms of use.
  • 21 comments

  • Пермь больного человека

    Заметка в "Ленте" Добыча нефти в США обернулась экологической катастрофой по мотивам статьи в "Гардиан" How Texas’s zombie oil wells are creating…

  • "Going postal" - Кучно пошло

    Очередные пострелушки - теперь в Индианаполисе, в офисе Fedex, 8 убитых.…

  • Апрельская зима по расписанию

    16 апреля 2021 - более 10 см снега. Повторяется прошлогодняя история почти день в день, тогда было 18 апреля. И примерно так почти каждый год в…