home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 9 / 09.iso / l / l042 / 1.ddi / CHAP5.ARC / ROMBERG.PAS < prev    next >
Encoding:
Pascal/Delphi Source File  |  1987-12-30  |  7.4 KB  |  215 lines
  1. program Romberg_Prog;
  2.  
  3. {----------------------------------------------------------------------------}
  4. {-                                                                          -}
  5. {-     Turbo Pascal Numerical Methods Toolbox                               -}
  6. {-     Copyright (c) 1986, 87 by Borland International, Inc.                -}
  7. {-                                                                          -}
  8. {-           Purpose: This program demonstrates integration with the        -}
  9. {-                    Romberg algorithm.                                    -}
  10. {-                                                                          -}
  11. {-           Unit   : Integrat    procedure Romberg                         -}
  12. {-                                                                          -}
  13. {----------------------------------------------------------------------------}
  14.  
  15. {$I-}      { Disable I/O checking }
  16. {$R+}      { Enable range checking }
  17.  
  18. uses
  19.   Integrat, Dos, Crt, Common;
  20.  
  21. var
  22.   LowerLimit, UpperLimit : Float;   { Limits of integration }
  23.   Tolerance : Float;                { Tolerance in the answer }
  24.   MaxIter : integer;                { Maximum number of iterations }
  25.   Integral : Float;                 { Value of the integral }
  26.   Iter : integer;                   { Number of iterations to find answer }
  27.   Error : byte;                     { Flags if something went wrong }
  28.  
  29. {$F+}
  30. function TNTargetF(X : Float) : Float;
  31.  
  32. {-----------------------------------------------------}
  33. {-         This is the function to integrate         -}
  34. {-----------------------------------------------------}
  35.  
  36. begin
  37.   TNTargetF := Exp(3 * X) + Sqr(X) / 3;
  38. end; { function TNTargetF }
  39. {$F-}
  40.  
  41. procedure Initialize(var LowerLimit : Float;
  42.                      var UpperLimit : Float;
  43.                      var Integral   : Float;
  44.                      var Tolerance  : Float;
  45.                      var MaxIter    : integer;
  46.                      var Iter       : integer;
  47.                      var Error      : byte);
  48.  
  49. {------------------------------------------------------------------}
  50. {- Output: LowerLimit, UpperLimit, Integral,                      -}
  51. {-         Tolerance, MaxIter, Iter, Error                        -}
  52. {-                                                                -}
  53. {- This procedure initializes the above variables to zero         -}
  54. {------------------------------------------------------------------}
  55.  
  56. begin
  57.   LowerLimit := 0;
  58.   UpperLimit := 0;
  59.   Integral := 0;
  60.   Tolerance := 0;
  61.   MaxIter := 0;
  62.   Iter := 0;
  63.   Error := 0;
  64. end; { procedure Initialize }
  65.  
  66. procedure GetData(var LowerLimit : Float;
  67.                   var UpperLimit : Float;
  68.                   var Tolerance  : Float;
  69.                   var MaxIter    : integer);
  70.  
  71. {------------------------------------------------------------}
  72. {- Output: LowerLimit, UpperLimit, Tolerance, MaxIter       -}
  73. {-                                                          -}
  74. {- This procedure assigns values to the above variables     -}
  75. {- from keyboard input                                      -}
  76. {------------------------------------------------------------}
  77.  
  78. procedure GetLimits(var LowerLimit : Float;
  79.                     var UpperLimit : Float);
  80.  
  81. {------------------------------------------------------------}
  82. {- Output: LowerLimit, UpperLimit                           -}
  83. {-                                                          -}
  84. {- This procedure assigns values to the limits of           -}
  85. {- integration from keyboard input                          -}
  86. {------------------------------------------------------------}
  87.  
  88. begin
  89.   repeat
  90.     repeat
  91.       Write('Lower limit of integration? ');
  92.       Readln(LowerLimit);
  93.       IOCheck;
  94.     until not IOerr;
  95.     Writeln;
  96.     repeat
  97.       Write('Upper limit of integration? ');
  98.       Readln(UpperLimit);
  99.       IOCheck;
  100.     until not IOerr;
  101.     if LowerLimit = UpperLimit then
  102.     begin
  103.       Writeln;
  104.       Writeln('       The limits of integration must be different.');
  105.       Writeln;
  106.     end;
  107.   until LowerLimit <> UpperLimit;
  108. end; { procedure GetLimits }
  109.  
  110. procedure GetTolerance(var Tolerance : Float);
  111.  
  112. {--------------------------------------------------}
  113. {- Output: Tolerance                              -}
  114. {-                                                -}
  115. {- This procedure reads in the accepted Tolerance -}
  116. {- from the keyboard.                             -}
  117. {--------------------------------------------------}
  118.  
  119. begin
  120.   Writeln;
  121.   repeat
  122.     Tolerance := 1E-8;
  123.     Write('Tolerance (> 0) ');
  124.     ReadFloat(Tolerance);
  125.     IOCheck;
  126.     if Tolerance <= 0 then
  127.     begin
  128.       IOerr := true;
  129.       Tolerance := 1E-8;
  130.     end;
  131.   until not IOerr;
  132. end; { procedure GetTolerance }
  133.  
  134. procedure GetMaxIter(var MaxIter : integer);
  135.  
  136. {--------------------------------------------------}
  137. {- Output: MaxIter                                -}
  138. {-                                                -}
  139. {- This procedure reads in the accepted MaxIter   -}
  140. {- from the keyboard.                             -}
  141. {--------------------------------------------------}
  142.  
  143. begin
  144.   Writeln;
  145.   repeat
  146.     MaxIter := 100;
  147.     Write('Maximum number of iterations: (> 0) ');
  148.     ReadInt(MaxIter);
  149.     IOCheck;
  150.     if MaxIter <= 0 then
  151.     begin
  152.       IOerr := true;
  153.       MaxIter := 1000;
  154.     end;
  155.   until not IOerr;
  156. end; { procedure GetMaxIter }
  157.  
  158. begin { procedure GetData }
  159.   GetLimits(LowerLimit, UpperLimit);
  160.   GetTolerance(Tolerance);
  161.   GetMaxIter(MaxIter);
  162.   GetOutputFile(OutFile);
  163. end; { procedure GetData }
  164.  
  165. procedure Results(LowerLimit : Float;
  166.                   UpperLimit : Float;
  167.                   Tolerance  : Float;
  168.                   MaxIter    : integer;
  169.                   Integral   : Float;
  170.                   Iter       : integer;
  171.                   Error      : byte);
  172.  
  173. {------------------------------------------------------------}
  174. {- This procedure outputs the results to the device OutFile -}
  175. {------------------------------------------------------------}
  176.  
  177. begin
  178.   Writeln(OutFile);
  179.   Writeln(OutFile);
  180.   Writeln(OutFile, 'Lower Limit:' : 30, LowerLimit : 25);
  181.   Writeln(OutFile, 'Upper Limit:' : 30, UpperLimit : 25);
  182.   Writeln(OutFile, 'Tolerance:' : 30, Tolerance : 25);
  183.   Writeln(OutFile, 'Maximum number of iterations:' : 30, MaxIter : 5);
  184.   Writeln(OutFile, 'Number of iterations:' : 30, Iter : 5);
  185.   Writeln(OutFile);
  186.   if Error = 3 then
  187.     DisplayWarning;
  188.   if Error in [1, 2] then
  189.     DisplayError;
  190.  
  191.   case Error of
  192.     0 : Writeln(OutFile, 'Integral:' : 25, Integral);
  193.  
  194.     1 : Writeln(OutFile, 'The tolerance must be greater than zero.');
  195.  
  196.     2 : Writeln(OutFile,
  197.                 'The maximum number of iterations must be greater than zero.');
  198.     3 : begin
  199.           Writeln(OutFile, 'Tolerance not reached in ', Iter, ' iterations.');
  200.           Writeln(OutFile, 'The last iterate of the integral is:', Integral);
  201.         end;
  202.   end; { case }
  203. end; { procedure Results }
  204.  
  205. begin { program Romberg }
  206.   ClrScr;
  207.   Initialize(LowerLimit, UpperLimit, Integral, Tolerance,
  208.              MaxIter, Iter, Error);
  209.   GetData(LowerLimit, UpperLimit, Tolerance, MaxIter);
  210.   Romberg(LowerLimit, UpperLimit, Tolerance, MaxIter,
  211.           Integral, Iter, Error, @TNTargetF);
  212.   Results(LowerLimit, UpperLimit, Tolerance, MaxIter, Integral, Iter, Error);
  213.   Close(OutFile);
  214. end. { program Romberg }
  215.