PARALLEL.RU

Дискуссионный клуб по параллельным вычислениям
Текущее время: 22 окт 18 19:18

Часовой пояс: UTC + 4 часа [ Летнее время ]




Начать новую тему Ответить на тему  [ Сообщений: 9 ] 
Автор Сообщение
СообщениеДобавлено: 7 июл 09 22:36 
Не в сети

Зарегистрирован: 27 мар 08 16:39
Сообщения: 11
Есть банальная задача расчета суммы ряда с заданной точностью. В последовательном цикле это реализуется очень просто - если следующий член ряда дает вклад в сумму меньше заданой точности, то прекращаем считать. Как сделать тоже самое, но в параллельном цикле? Должно ж быть решение, очень уж часто такие задачи в природе встречаются...


Вернуться к началу
 Профиль  
 
СообщениеДобавлено: 8 июл 09 15:33 
Не в сети

Зарегистрирован: 13 сен 08 18:39
Сообщения: 74
Откуда: Москва
Задача, кстати говоря, очень не банальная в общем случае, и алгоритм "если следующий член ряда меньше tol, то стоп" работает только для сходящихся знакоперменных монотонно убывающих рядов (рядов Лейбница), а вот насчет знакопостоянных рядов, скорее всего, нужно применять другие методы оценки остатка ряда. Или речь идет о каких-то конкретных видах рядов? Это всё в качестве урпражнения по программированию или надо решать какие-то прикладные задачи?

_________________
Дмитрий О. Коломиец.
IBM // МГУ, физфак, каф. математики.


Вернуться к началу
 Профиль  
 
СообщениеДобавлено: 8 июл 09 23:37 
Не в сети

Зарегистрирован: 27 мар 08 16:39
Сообщения: 11
Задача у меня очень прикладная, очень хочется ее наконец сделать масштабируемой... Ряд вроде бы у меня сходящийся достаточно монотонно, но знак не меняется, в фортране последовательный цикл я делаю так:

COMMON m
.....
.....
PRECIS = 0.0001
M_max = 100

err =100.0
sum =0.0
S_old=0.0

m=0

do while ((err>PRECIS).AND.(m<M_max))

CALL F(S)
sum = sum + S

err=dabs(S)/dabs(S_old+S)
S_old=S_old+S
m=m+1

enddo

Т.е. выходим из цикла, если достигли точности или лимита итераций.
Процедура F зависит от m (используется COMMON, т.к в жизни я считаю сумму ряда интегралов, которые используют IMSL процедуры, а они не воспринимают ф-ции с несколькими переменным).
Имно относительная ошибка err=dabs(S)/dabs(S_old+S) - как раз хорошо должна подходить для оценки точности.


Вернуться к началу
 Профиль  
 
СообщениеДобавлено: 9 июл 09 17:53 
Не в сети

Зарегистрирован: 13 сен 08 18:39
Сообщения: 74
Откуда: Москва
А, собственно, S как считается? из каких параметов?

Кстати, к вышесказанному, такой алгоритм не позволяет гарантировать каую-либо наперед заданную точность результата, т.е. если сказать, что res = sum +/- PRECIS, то это будет неверно в общем случае.

_________________
Дмитрий О. Коломиец.
IBM // МГУ, физфак, каф. математики.


Вернуться к началу
 Профиль  
 
СообщениеДобавлено: 10 июл 09 4:54 
Не в сети

Зарегистрирован: 27 мар 08 16:39
Сообщения: 11
Вот примерный кусок кода как все в натуре выглядит:
Код:
COMMON a,Di,m
!------------------
      Nker = OMP_GET_NUM_PROCS( )
      call OMP_SET_NUM_THREADS(Nker)
      iter = 10

sum = 0.0d0
err1 = 99999.0

!$OMP PARALLEL DO ORDERED DEFAULT(shared) PRIVATE(m,Y1,sinsin) REDUCTION(+:sum) SCHEDULE (dynamic , iter) IF (Nker .GT. 1)
       do m=1,M_max,1
!$OMP ORDERED
         sinsin = dsin(0.5d0*m*Di/a)/(0.5d0*m*Di/a)
         if (err1>errrel) then
                               CALL DQDAG (funY  ,0.d0,1.0d0,errabs,errrel,1,Y1,errest)
                               Y1=sinsin*Y1
                               sum  = sum  + Y1                             
                               err1=dabs(Y1/sum)
         endif             
!$OMP END ORDERED
       enddo
!$OMP END PARALLEL DO

!------------------
...........
...........

real*8 function funY (x)
implicit none

COMMON a,Di,m
real*8 a,Di
integer m

real*8,  INTENT(IN) ::  x

real*8 .., .., ...

     funY = (... ... ..)
END


DQDAG - ф-ция из IMSL.
Собственно проблема с ним одна - он ерунду считает, когда цикл параллельный, когда последовательный - все ОК.

Так, а чем плохо для монотонного знакопостоянного ряда оценка ошибки в виде err1=dabs(Y1/sum)? т.е. отшение последнего члена суммы к сумме ряда? Для знакопеременного ряда - это не очень годится, т.к. надо смотреть на разницу между эл-ми ряда, которые по сравнению рядом по отдельности могут быть вполне значительными, а для постоянного, имно очень даже годится...


Вернуться к началу
 Профиль  
 
СообщениеДобавлено: 10 июл 09 16:44 
Не в сети

Зарегистрирован: 13 сен 08 18:39
Сообщения: 74
Откуда: Москва
Pir0texnik писал(а):
Так, а чем плохо для монотонного знакопостоянного ряда оценка ошибки в виде err1=dabs(Y1/sum)? т.е. отшение последнего члена суммы к сумме ряда? Для знакопеременного ряда - это не очень годится, т.к. надо смотреть на разницу между эл-ми ряда, которые по сравнению рядом по отдельности могут быть вполне значительными, а для постоянного, имно очень даже годится...


Ну это банальный матан... ошибкой в данном случае является сумма остатка ряда (то есть сумма всех членов, начиная с некоторого). Так вот, есть теорема Лейбница, в которой показано, что сумма остатка знакопеременного убывающего (|a(n)| > |a(n+1)|) ряда меньше модуля первого члена остатка. Вот для такого случая ваш алгоритм подходит.

Теперь посмотрим, что получится, если таким образом посчитать сумму монотонно убывающего положительного ряда.
Набросал простенький пример на мэпле:
Код:
> restart;
# Просим считать с точностью 20 знаков после запятой
> Digits:= 20;
# находим аналитически сумму ряда 1/k^2 для k от 1 до бесконечности.
> exact := sum(1/k^2, k=1..infinity);

                                        2
                                      Pi
                             exact := ---
                                       6

# Сумма равна Pi^2/6
# Что в десятичном виде с 20 знаками после запятой выглядит вот так:
> evalf(exact);

                        1.6449340668482264365

# Выберем желаемую точность для вашего алгоритма
> tol:=0.001;

                             tol := 0.001

# Пусть на пером шаге "относительная ошибка" равна 1.0
> err:=1.0;

                              err := 1.0

# Частичную сумму ряда для начала инициализируем нулем
> total:=0;

                              total := 0

# Далее по алгоритму, останавливаемся, когда последний член суммы
# отнесенный к сумме менше tol
> for i from 1 by 1 while err > tol do
>     total := total + evalf(1/i^2):
>     err := evalf((1/i^2))/total:
> end do:
# получаем результат
> total;

                        1.6057234035910008572

# сотрим, насколько он далеко от точного
> realError := evalf(exact - total);

                  realError := 0.0392106632572255793

> видим что, реальная точность в 39 раз хуже tol


Если задать tol := 0.0001, то реальная точность будет уже в 125 раз хуже tol.

Надеюсь, понятно, что таким алгоритмом считать нельзя?

_________________
Дмитрий О. Коломиец.
IBM // МГУ, физфак, каф. математики.


Вернуться к началу
 Профиль  
 
СообщениеДобавлено: 10 июл 09 19:26 
Не в сети

Зарегистрирован: 27 мар 08 16:39
Сообщения: 11
Да, я понимаю, что у меня будет ухудшаться отношение между абсолютной и относительной ошибкой, но во-1, мне неизвестно(будь оно известно зачем тогда считать ряд :) ) абсолютное значение суммы ряда, поэтому оно мне бесполезно, во-2, мне нужно как-то ограничить суму ряда, я же не могу считать до бесконечности, с этой ошибкой у меня есть хоть какой-то критерий точности, который в принципе я могу как-то соотнести потом с абс. ошибкой набрав побольше информации о ряде, во-3 с увеличением относительной точности кол-во членов ряда растет => абсолютная ошибка как ни крути - уменьшается. Или есть более адекватный способ оценки абс. точности?


Вернуться к началу
 Профиль  
 
СообщениеДобавлено: 13 июл 09 12:28 
Не в сети

Зарегистрирован: 13 сен 08 18:39
Сообщения: 74
Откуда: Москва
Не хочу никого расстраивать, но даже для "хорошего" ряда, увеличение количества членов частичной суммы может приводить к ухудшению абсолюбной точности, очень легко могу соорудить пример )

_________________
Дмитрий О. Коломиец.
IBM // МГУ, физфак, каф. математики.


Вернуться к началу
 Профиль  
 
СообщениеДобавлено: 13 июл 09 22:02 
Не в сети

Зарегистрирован: 27 мар 08 16:39
Сообщения: 11
Что ж делать тогда? :)
Понятно, что ряд может, например, периодически возрастать и убывать - надо смотреть, когда его обрывать, но все равно чем больше - тем лучше.
Как бы там ни было, сумма моего ряда получается с удовлетворительной точностью. Даже с такой оценкой. Вот только параллельно не считаеся...


Вернуться к началу
 Профиль  
 
Показать сообщения за:  Поле сортировки  
Начать новую тему Ответить на тему  [ Сообщений: 9 ] 

Часовой пояс: UTC + 4 часа [ Летнее время ]


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

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


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

Найти:
Перейти:  
cron
Создано на основе phpBB® Forum Software © phpBB Group
Русская поддержка phpBB