На главную Наши проекты:
Журнал   ·   Discuz!ML   ·   Wiki   ·   DRKB   ·   Помощь проекту
ПРАВИЛА FAQ Помощь Участники Календарь Избранное RSS
msm.ru
! Друзья, соблюдайте, пожалуйста, правила форума и данного раздела:
Данный раздел не предназначен для вопросов и обсуждений, он содержит FAQ-заготовки для разных языков программирования. Любой желающий может разместить здесь свою статью. Вопросы же задавайте в тематических разделах!
• Если ваша статья может быть перенесена в FAQ соответствующего раздела, при условии, что она будет оформлена в соответствии с Требованиями к оформлению статей.
• Чтобы остальным было проще понять, указывайте в описании темы (подзаголовке) название языка в [квадратных скобках]!
Модераторы: Модераторы
  
> Численное интегрирование , [Pascal] Решаем интеграл разными способами
    subj in attach (CP-866)

    Че скажите?

    09.05.05: Исправил ошибку в методе Симпсона...
    Сообщение отредактировано: Jin X -

    Прикреплённый файлПрикреплённый файлLAB3_.zip (1.51 Кбайт, скачиваний: 331)
      Теория - это, конечно, хорошо. Но всё же я рекомендую составлять поэтапное решение поставленных задач, опираясь на источники теории, дабы научить народ переводить математические формулы в программный код.
      Важна "история" развития программы: с чего начинать (типы, процедуры, циклы), методы и пути решения, и т.д. ;)
        ok, буду че-то думать ;)
          Долго думал куда это запостить: то ли в "Маразм прогамееров", то ли в топик, где обсуждались самые "корявые исходники" :lol:
          Решил запостить сюда, ведь это все работает(пусть не так оптимально как хотелось бы, но зато мало строк ;) )....
          Тем более
          Цитата Romtek @
          методы и пути решения, и т.д. ;)

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

          З.Ы. че скажите?

          ExpandedWrap disabled
            {
              Writeln('Labinskiy Nikolay aka e-moe (c) 2005');
              Writeln('See you later on forum.sources.ru ;)');
            }
            {$N+,E+,G+}
            const
              a = -1.1;
              b = 0.9;
             
            type
              Tfunc = function(x: extended):extended;
             
            function func(x: extended): extended; far;
            begin
              func:=0.9*cos(-1.77*x)+2.65*exp(-1.57*x)-11.48*x*x*x*x+1.29*x-10.49;
            end;
             
            function Int1(f: TFunc; a,b: extended; n: byte): extended;
            { Правило трапеций }
            var
              res: double;   { времееная переменная, будет хранить результат }
              {если extended, то надо обнулить отдельно, а не в цикле :( TP7.1}
              i: byte;       { счетцик цикла }
            begin
              { на 0 делить нельзя! }
              inc(n,byte(n=0));
              for i:=0 to n do
                res:=res*byte(i<>0)+(byte((i>0)and(i<n))+1)*f(a+i*(b-a)/n)*(b-a)/n*0.5;
                  { интеграл - площадь фигуры, а следовательно не может
                    быть отрицательной, поэтому используем абсолютное значение
                    функции (или по-русски модуль) }
              Int1:=res;
            end;
             
            function Int2(f: TFunc; a,b: extended; n,m: byte): extended;
            begin
              { n<>m && n*m <> 0 }
              inc(n,byte(n=0));
              inc(m,byte(m=0)*(byte(n=1)+1)-byte((m=255) and (n=255))*2+byte(m=n));
              { Для подготовки к вычислению надо в 2 раза больше строк чем
                на само вычисление... Странно, не так ли? :) }
              Int2:=Int1(f,a,b,n)+(Int1(f,a,b,n)-Int1(f,a,b,m))/((n*n)/(m*m)-1)
            end;
             
            function Int3(f: TFunc; a,b: extended; n: byte): extended;
            var
              res: double; {если extended, то надо обнулить отдельно, а не в цикле :( TP7.1}
              i: byte;
            begin
              inc(n,byte(n and 1)+byte(n=0)*2-byte(n=255)*2);  { (n mod 2=0) && n<>0 }
              for i:=0 to n do
                res:=res*byte(n<>0)+(1 shl (byte((i>0)and(i<n))+byte(odd(i)and(i>0)and(i<n))))*f(a+i*(b-a)/n)*(b-a)/n/3;
              Int3:=res;
            end;
             
            function Int4(f: TFunc; a,b: extended; n: byte): extended;
            const
              AiMi : array[2..6,1..12] of extended =
                     ((0.5773502692,1,-0.5773502692,1,0,0,0,0,0,0,0,0),
                     (0.7745966692,0.5555555556,-0.7745966692,0.5555555556,
                       0,0.8888888889,0,0,0,0,0,0),
                     (0.8611363116,0.3478548451,-0.8611363116,0.3478548451,
                       0.3399810436,0.6521451549,-0.3399810436,0.6521451549,0,0,0,0),
                     (0.9061798459,0.2369268851,-0.9061798459,0.2369268851,0.5384693101,
                       0.4786286705,-0.5384693101,0.4786286705,0.0,0.5688888889,0,0),
                     (0.9324695142,0.1713244924,-0.9324695142,0.1713244924,
                       0.6612093865,0.3607615730,-0.6612093865,0.3607615730,
                       0.2386191861,0.4679139346,-0.2386191861,0.4679139346));
            var
              i:byte;
              res: double; {если extended, то надо обнулить отдельно, а не в цикле :( TP7.1}
            begin
              { n>=2 && n <=6 }
              inc(n,byte(n<2)*(2-n)-byte(n>6)*(n-6));
              for i:=1 to n do
                res:=res*byte(i<>1)+AiMi[n,i*2]*((b-a)*f((b-a)*AiMi[n,i*2-1]/2+(b+a)/2)/2);
              Int4:=res;
            end;
             
            begin
              Writeln('     Численное интегрирование:',#10);
              Writeln('f(x) = 0.9*cos(-1.77x) + 2.65*e^(-1.57x) - 11.48*x^4 + 1.29x - 10.49');
              Writeln('Int(-1.1, 0.9)[f(x)]',#10);
             
              Writeln('Правило трапеций:');
              writeln('n=2      I=',Int1(func,a,b,2):5:5);
              writeln('n=8      I=',Int1(func,a,b,8):5:5);
              writeln('n=32     I=',Int1(func,a,b,32):5:5);
              writeln('n=64     I=',Int1(func,a,b,64):5:5);
             
              Writeln('Экстраполяционный переход к пределу::');
              writeln('n=2,8    I=',Int2(func,a,b,2,8):5:5);
              writeln('n=8,16   I=',Int2(func,a,b,8,16):5:5);
              writeln('n=16,32  I=',Int2(func,a,b,16,32):5:5);
             
              Writeln('Правило Симпсона:');
              writeln('n=2      I=',Int3(func,a,b,2):5:5);
              writeln('n=8      I=',Int3(func,a,b,8):5:5);
              writeln('n=32     I=',Int3(func,a,b,32):5:5);
             
              Writeln('Метод Гаусса:');
              writeln('n=2      I=',Int4(func,a,b,2):5:5);
              writeln('n=3      I=',Int4(func,a,b,3):5:5);
              writeln('n=4      I=',Int4(func,a,b,4):5:5);
              writeln('n=5      I=',Int4(func,a,b,5):5:5);
              writeln('n=6      I=',Int4(func,a,b,6):5:5);
             
            end.
            З.Ы.
                 не знаю что на меня нашло :), но я решил макисмально
                 сократить число строк во всех функциях, вот что из
                 всего этого вышло... :)
             
            З.З.Ы.
                  параметры функций:
                    f - ссылка на подинтегральную функцию
                    a,b - координаты начальной и конечной точнек интегрирования
                    n[,m] - чило разбиений.
             
                 в функции можно передавать любые n и m: если нужно,
                 то прога сама их подправит так, как нужно для
                 конкретного метода;)
             
            З.З.З.Ы.
                 как это работает? а кто его знает... я уже и
                 сам не пойму :) а если серьезно, то могу помочь
                 разобраться что тут к чему ;)


          Че-то мне кажется, что более короткой реализации этих процедур еще не делали ;) :lol:
          Сообщение отредактировано: e-moe -
            e-moe, а теперь подробнее про каждый метод: :)
            Что, как и почему.
            Не знаю на кого рассчитан сей код, но я его с трудом понимаю.
              Цитата Romtek @
              e-moe, а теперь подробнее про каждый метод: :)

              Угу, я это и хочу расписать ;) :))
              Только вот не знаю, может сделать отдельную тему где писать все такие "штучки", а тут все-таки написать нормальный исходник по интегрированию?

              Цитата Romtek @
              ...но я его с трудом понимаю.

              Эт точно, если не расписать как это работает, то я и сам скоро забуду.... :lol:
                Цитата e-moe @
                может сделать отдельную тему где писать все такие "штучки", а тут все-таки написать нормальный исходник по интегрированию?
                Зачем плодить лишние темы? "Штучки" можно вполне и здесь обсудить (это и есть цель темы).
                Как только появится что-то оформленное и обсуждённое, можно тогда создавать новую тему, в которую "оно" перемещается.
                  понятно ;)
                    Цитата e-moe @
                    интеграл - площадь фигуры, а следовательно не может быть отрицательной, поэтому используем абсолютное значение функции (или по-русски модуль

                    Ошибка! Вообще-то интеграл запросто может быть отрицательным. И интеграл это не площадь. Наоборот, площадь - это интеграл. А интеграл определяется совсем по-другому.

                    Вообще, не хочу остужать энтузиазм или заниматься саморекламой, но загляните по адресу http://alglib.sources.ru/integral/
                      shadeofgray, спсибо за критику, как появится время исправлю ;)
                      0 пользователей читают эту тему (0 гостей и 0 скрытых пользователей)
                      0 пользователей:


                      Рейтинг@Mail.ru
                      [ Script execution time: 0,0309 ]   [ 17 queries used ]   [ Generated: 14.11.24, 02:02 GMT ]