На главную
ПРАВИЛА FAQ Помощь Участники Календарь Избранное DigiMania RSS
msm.ru
! Правила раздела Ищу исходники!
Раздел создан для обмена исходниками, как своими, так и чужими. Попытки продажи собственных исходников крайне не приветствуются (для этого у нас есть разделы ПОМОЩЬ СТУДЕНТАМ и Разовые заказы и подработка), а попытки продажи чужих исходников будут наказываться.

Обратите внимание:
1. Прежде чем начать новую тему или отправить сообщение, убедитесь, что Вы не нарушаете правил форума!
2. Обязательно воспользуйтесь поиском. Весьма вероятно, что кто-то до вас уже искал интересующие вас исходники.
3. При создании темы обязательно указывайте интересующий вас язык программирования! Указывать следует в таком формате: "[_язык_] Название темы". Например, "[C++] Ищу исходники простого графического редактора"
4. Не забывайте пользоваться тегом [ code=*** ] ...текст программы... [ /code ] для выделения текста программы подсветкой!

Полезные ссылки:
user posted image Поиск по Разделу
Модераторы: B.V.
  
> У кого нибудь есть QUANC8 на Паскале/Delphi (программа числ.интегрирования), из книги Форсайт, Мальком, Моулер "Машинные методы мат.вычислений".
    Программка полезная, в интернете найти не удается, если ни у кого нет, буду сам переводить с фортрана на паскаль (Delphi).

    У меня есть на паскале процедуры SPLINE и SEVAL из этой книги, если кому надо.
      swe
      Всё нужно.
      Правильный обед должен состоять из 5 блюд приготовленных из 33 ингредиентов.
        У меня получился вот такой код для QUANC8:

        ExpandedWrap disabled
          unit Unit1;
           
          interface
           
          uses
            Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,
            Dialogs, StdCtrls;
           
          type
            TForm1 = class(TForm)
              Button1: TButton;
              procedure Button1Click(Sender: TObject);
            end;
           
          var
            Form1: TForm1;
           
          implementation
           
          {$R *.dfm}
           
          type
            TFun = function(X:real):real;
           
          FUNCTION FUN(X:real):REAL;
          begin
            IF (X = 0) then FUN := 1.0
                       else FUN := SIN(X)/X;
          end;
           
          function iPOW(II,IP:integer):integer; (* для вычисления 2^IP *)
          var I,IR : integer;
          begin
            IR := II; for I := 2 to IP do IR := IR * II; result := IR;
          end;
           
          function MAX(X1,X2:real):real; begin if X2>X1 then MAX:=X2 else MAX:=X1 end;
           
          procedure QUANC8(FUN:TFun;A,B,ABSERR,RELERR:real;var RESULT,ERREST:real;
                           var NOFUN:integer;var FLAG:real);
          // ОЦЕНИТЬ ИНТЕГРАЛ ДЛЯ FUN(X) ОТ А ДО В С ЗАДАННОЙ ПОЛЬЗОВАТЕЛЕМ ТОЧНОСТЬЮ.
          // АВТОМАТИЧЕСКАЯ АДАПТИВНАЯ ПРОГРАММА, ОСНОВАННАЯ НА ФОРМУЛЕ НЬЮТОНА-КОТЕСА
          // 8-ГО ПОРЯДКА
          //
          // ВХОДНАЯ ИНФОРМАЦИЯ..
          //
          // FUN   ИМЯ ФУНКЦИИ FUN{X), РЕАЛИЗУЮЩЕЙ ПОДЫНТЕГРАЛЬНУЮ ФУНКЦИЮ
          //
          // А      НИЖНИЙ ПРЕДЕЛ ИНТЕГРИРОВАНИЯ
          // В      ВЕРХНИЙ ПРЕДЕЛ ИНТЕГРИРОВАНИЯ (В МОЖЕТ БЫТЬ МЕНЬШЕ, ЧЕМ А)
          //
          // RELERR ГРАНИЦА ОТНОСИТЕЛЬНОЙ ПОГРЕШНОСТИ (ДОЛЖНА БЫТЬ НЕОТРИЦАТЕЛЬНА)
          // ABSERR ГРАНИЦА АБСОЛЮТНОЙ ПОГРЕШНОСТИ    (ДОЛЖНА БЫТЬ НЕОТРИЦАТЕЛЬНА)
          //
          //
          // ВЫХОДНАЯ ИНФОРМАЦИЯ..
          //
          // RESULT ПРИБЛИЖЕНИЕ К ИНТЕГРАЛУ, УДОВЛЕТВОРЯЮЩЕЕ, МОЖНО НАДЕЯТЬСЯ, МЕНЕЕ
          //        ЖЕСТКОЙ ИЗ ДВУХ ГРАНИЦ ПОГРЕШНОСТИ.
          // ERREST ОЦЕНКА ВЕЛИЧИНЫ ДЕЙСТВИТЕЛЬНОЙ ОШИБКИ.
          // NOFUN  ЧИСЛО ЗНАЧЕНИЙ ФУНКЦИИ, ИСПОЛЬЗОВАННЫХ ПРИ ВЫЧИСЛЕНИИ RESULT.
          // FLAG   ИНДИКАТОР НАДЕЖНОСТИ. ЕСЛИ FLAG РАВЕН НУЛЮ, ТО RESULT, ВЕРОЯТНО,
          //        УДОВЛЕТВОРЯЕТ ЗАДАННОЙ ГРАНИЦЕ ПОГРЕШНОСТИ.
          //        ЕСЛИ FLAG=XXXYYY, ТО ХХХ-ЧИСЛО ИНТЕРВАЛОВ, ДЛЯ КОТОРЫХ НЕ БЫЛО
          //        СХОДИМОСТИ, А O.YYY= ЧАСТЬ ОСНОВНОГО ИНТЕРВАЛА, ОСТАВШАЯСЯ ДЛЯ
          //        ОБРАБОТКИ В ТОТ МОМЕНТ, КОГДА ПРОГРАММА ПРИБЛИЗИЛАСЬ
          //        К ПРЕДЕЛЬНОМУ ЗНАЧЕНИЮ ДЛЯ NOFUN.
          //
          var
            W0, W1, W2, W3, W4, AREA, X0, F0, STONE, STEP, COR11, TEMP : REAL;
            QPREV, QNOW, QDIFF, QLEFT, ESTERR, TOLERR                  : REAL;
            QRIGHT       : array [1..31] of real;
            F, X         : array [1..16] of real;
            FSAVE, XSAVE : array [1..8,1..30] of real;
            LEVMIN, LEVMAX, LEVOUT, NOMAX, NOFIN, LEV, NIM, I, J, jj : INTEGER;
            QDone,Q72 : boolean; (* лог.переменные для завершения циклов *)
          begin
          //
          // *** ЭТАП 1 *** ПРИСВОЕНИЕ НАЧАЛЬНЫХ ЗНАЧЕНИЙ
          // ПЕРЕМЕННЫМ, НЕ ЗАВИСЯЩИМ ОТ ИНТЕРВАЛА. ГЕНЕРИ-
          // РОВАНИЕ КОНСТАНТ.
          //
            LEVMIN :=    1;
            LEVMAX :=   30;
            LEVOUT :=    6;
            NOMAX  := 5000;
            NOFIN  := NOMAX - 8 * (LEVMAX - LEVOUT + iPOW(2,(LEVOUT + 1)));
          //
          // ЕСЛИ NOFUN ДОСТИГАЕТ ЗНАЧЕНИЯ NOFIN, TO ТРЕВОГА
          //
            W0 :=   3956/14175;
            W1 :=  23552/14175;
            W2 :=  -3712/14175;
            W3 :=  41984/14175;
            W4 := -18160/14175;
          //
          // ПРИСВОИТЬ НУЛЕВЫЕ ЗНАЧЕНИЯ ПЕРЕМЕННЫМ СУММАМ.
          //
            FLAG   := 0;
            RESULT := 0;
            COR11  := 0;
            ERREST := 0;
            AREA   := 0;
            NOFUN  := 0;
            if (A = B) then Exit;
          //
          // ***ЭТАП 2*** ПРИСВОЕНИЕ НАЧАЛЬНЫХ ЗНАЧЕНИЙ
          // ПЕРЕМЕННЫМ, ЗАВИСЯЩИМ ОТ ИНТЕРВАЛА, В
          // СООТВЕТСТВИИ С ПЕРВЫМ ИНТЕРВАЛОМ
          //
            LEV   := 0;
            NIM   := 1;
            X0    := A;
            X[16] := B;
            QPREV := 0;
            F0    := FUN(X0);
            STONE := (B - A)/16;
            X[8]  := (X0  +X[16])/2;
            X[4]  := (X0  +X[8] )/2;
            X[12] := (X[8]+X[16])/2;
            X[2]  := (X0  +X[4] )/2;
            X[6]  := (X[4]+X[8] )/2;
            X[10] := (X[8]+X[12])/2;
            X[14] := (X[12]+X[16])/2;
           
            for jj := 1 to 8 do begin J := jj*2; F[J] := FUN(X[J]); end; NOFUN := 1+8;
          //
          // ***ЭТАП 3*** ОСНОВНЫЕ ВЫЧИСЛЕНИЯ
          // ТРЕБУЮТСЯ QPREV, Х0,Х2,Х4,...,Х16, F0,F2,F4,...,F16.
          // ВЫЧИСЛЯЮТСЯ X1,ХЗ,...,Х15, F1,F3,...,F15, QLEFT, QRIGHT, QNOW, QDIFF, AREA.
          //
          //
           
          /////////////////////////////////////////////////     30
            QDone := false;
            repeat
           
              X[1]  := (X0+X[2])/2;
              F[1]  := FUN(X[1]);
          //  for J := 3 to 15 step 2 do begin
              for jj := 1 to 7 do begin
                J := 2*jj+1;
                X[J] := (X[J-1]+X[J+1])/2;
                F[J] := FUN(X[J]);
              end;
              NOFUN := NOFUN + 8;
              STEP  := (X[16] - X0)/16;
           
              QLEFT := (W0*(F0   + F[8])+
                        W1*(F[1] + F[7])+
                        W2*(F[2] + F[6])+
                        W3*(F[3] + F[5])+
                        W4* F[4]         )*STEP;
              QRIGHT[LEV+1] :=
                       (W0*(F[8] + F[16])+
                        W1*(F[9] + F[15])+
                        W2*(F[10]+ F[14])+
                        W3*(F[11]+ F[13])+
                        W4* F[12]         )*STEP;
              QNOW  := QLEFT + QRIGHT[LEV+1];
              QDIFF := QNOW  - QPREV;
              AREA  := AREA  + QDIFF;
          //
          // ***ЭТАП 4***. ПРОВЕРКА СХОДИМОСТИ ДЛЯ ИНТЕРВАЛА
          //
              ESTERR := ABS(QDIFF)/1023;
              TOLERR := MAX(ABSERR, RELERR*ABS(AREA)) * (STEP / STONE);
           
              IF  (LEV < LEVMIN) or
                  (
                  (Not (LEV >= LEVMAX)) and
                  (Not (NOFUN > NOFIN)) and
                  (Not (ESTERR <= TOLERR))
                  )
              then begin          //   GO TO 50
          /////////////////////////////////////////////// 50
          // ***ЭТАП 5*** СХОДИМОСТИ НЕТ
          // УСТАНОВИТЬ СЛЕДУЮЩИЙ ИНТЕРВАЛ.
          //
                NIM := 2*NIM;
                LEV := LEV+1;
          //
          // ЗАПОМНИТЬ ЭЛЕМЕНТЫ,   ОТНОСЯЩИЕСЯ   К   ПРАВОЙ
          // ПОЛОВИНЕ ИНТЕРВАЛА, ДЛЯ БУДУЩЕГО ИСПОЛЬЗОВАНИЯ.
          //
                for I := 1 to 8 do begin
                   FSAVE[I, LEV] := F[I + 8];
                   XSAVE[I, LEV] := X[I + 8];
                end;
          //
          // СОБРАТЬ ЭЛЕМЕНТЫ, ОТНОСЯЩИЕСЯ К ЛЕВОЙ ПОЛОВИНЕ
          // ИНТЕРВАЛА, ДЛЯ НЕМЕДЛЕННОГО ИСПОЛЬЗОВАНИЯ.
          //
                QPREV := QLEFT;
                for I := 1 to 8 do begin
                  J := -I;
                  F[2*J+18] := F[J + 9];
                  X[2*J+18] := X[J + 9];
                end;
                // вернуться к while Not QDone  =  GO TO 30
              end
              else begin    // else 50
                IF (ESTERR <= TOLERR) OR (NOFUN > NOFIN ) OR  (LEV >= LEVMAX) then begin
                //   = GO TO 70
                  IF (LEV >= LEVMAX) then begin
          //////////////////////////////////////////   62
          //
          // ТЕКУЩЕЕ ПРЕДЕЛЬНОЕ ЗНАЧЕНИЕ ГЛУБИНЫ ДЕЛЕНИЯ
          // ПОПОЛАМ РАВНО LEVMAX
          //
                    FLAG := FLAG + 1.0;
                  end  // 62
                  else begin // 62
                    IF (NOFUN > NOFIN) then begin // GO TO 60
          //////////////////////////////////////////   60
          // ***ЭТАП 6*** "ПОЖАРНЫЙ" РАЗДЕЛ
          // ЧИСЛО ЗНАЧЕНИЙ ФУНКЦИИ БЛИЗКО К ТОМУ, ЧТОБЫ
          // ПРЕВЫСИТЬ УСТАНОВЛЕННЫЙ ПРЕДЕЛ.
          //
                      NOFIN  := 2*NOFIN;
                      LEVMAX := LEVOUT;
                      FLAG   := FLAG + (B - X0) /(B - A);
                    end;  (* 62 *)
             //     GO TO 70
                  end; // else 62
           
          //////////////////////////////////////////// 70
          //
          // ***ЭТАП 7*** СХОДИМОСТЬ ДЛЯ ИНТЕРВАЛА ИМЕЕТ МЕСТО
          // ПРИБАВИТЬ ОЧЕРЕДНЫЕ СЛАГАЕМЫЕ  К ПЕРЕМЕННЫМ СУММАМ.
          //
                  RESULT := RESULT + QNOW;
                  ERREST := ERREST + ESTERR;
                  COR11  := COR11  + QDIFF/1023;
          //
          // УСТАНОВИТЬ СЛЕДУЮЩИЙ ИНТЕРВАЛ.
          //
          ////////////////////////////////////////////  72
                  Q72 := false;
                  repeat
                    Q72 := (NIM = 2*(NIM div 2));
                    if Not Q72 then begin
                      NIM := NIM div 2;
                      LEV := LEV-1;
                    end;
                  until Q72;
          ///////////////////////////////////////////// 75
                  NIM := NIM + 1;
                  IF (LEV <= 0) then QDone := TRUE;  //      GO TO 80
          //
          // СОБРАТЬ ЭЛЕМЕНТЫ, НЕОБХОДИМЫЕ ДЛЯ СЛЕДУЮЩЕГО
          // ИНТЕРВАЛА.
          //
                  QPREV := QRIGHT[LEV];
                  X0    := X[16];
                  F0    := F[16];
                  for I := 1 to 8 do begin
                    F[2*I] := FSAVE[I, LEV];
                    X[2*I] := XSAVE[I, LEV];
                  end;
                END; //  70
              end; //  else 50
            until QDone;  (* QDone = 30 *)
          /////////////////////////////////////////// 80
          //
          // ***ЭТАП 8*** ЗАКЛЮЧИТЕЛЬНЫЕ ОПЕРАЦИИ И ВЫХОД
          //
            RESULT := RESULT + COR11;
          //
          // ОБЕСПЕЧИТЬ, ЧТОБЫ ЗНАЧЕНИЕ ПЕРЕМЕННОЙ ERREST
          // БЫЛО НЕ МЕНЬШЕ УРОВНЯ ОКРУГЛЕНИЙ.
          //
            IF (ERREST = 0) then Exit;
          ////////////////////////////////////////// 82
            TEMP := ABS(RESULT) + ERREST;
            while (TEMP = ABS(RESULT)) do begin
              ERREST := 2.0*ERREST;
              TEMP := ABS(RESULT) + ERREST;
            end;
          END;
           
          (*============================================================================*)
           
          function ISt(I:integer):string;
          begin result := SysUtils.IntToStr(I); end;
           
          function FSt(R:real;N,D:integer):string;(*Nпоз.в строке, D цифр после запятой *)
          begin System.Str(R:N:D,Result); end;
           
          function ESt(R:real;D:integer): string;
          begin SysUtils.FmtStr(Result,'%'+ISt(D)+'e',[R]) end;
           
          procedure Say(S:string);
          Const sHd = 'Сообщение';
          begin
            Windows.MessageBox(GetActiveWindow,PChar(S),PChar(sHd),MB_OK);
          end;
           
          procedure Test_QUANC8;
          var
            fFUN : TFun;
            A, B, ABSERR, RELERR, RESULT, ERREST, FLAG : REAL;
            NOFUN : INTEGER;
          begin
            fFUN := FUN;
            A := 0.0;
            B := 2.0;
            RELERR := 1E-10;
            ABSERR := 0.0;
            QUANC8(fFUN, A, B, ABSERR, RELERR, RESULT, ERREST, NOFUN, FLAG);
           
            Say('Result ='+FSt(RESULT,15,10)+' ERREST='+ESt(ERREST,10)+
            ' NOFUN='+ISt(NOFUN)+' FLAG='+ESt(FLAG,5));
            IF (FLAG <> 0) then
              Say('WARNING.. RESULT MAY BE UNRELIABLE. FLAG = '+FSt(FLAG,6,2));
          end;
           
          procedure TForm1.Button1Click(Sender: TObject);
          begin Test_QUANC8; end;
           
          end.


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


        Рейтинг@Mail.ru
        [ Script Execution time: 0,0842 ]   [ 15 queries used ]   [ Generated: 16.07.19, 20:45 GMT ]