home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 9 / 09.iso / l / l042 / 1.ddi / CHAP2.ARC / ROOTSEQU.PAS < prev    next >
Encoding:
Pascal/Delphi Source File  |  1987-12-30  |  21.4 KB  |  347 lines
  1. unit RootsEqu;
  2.  
  3. {----------------------------------------------------------------------------}
  4. {-                                                                          -}
  5. {-     Turbo Pascal Numerical Methods Toolbox                               -}
  6. {-     Copyright (c) 1986, 87 by Borland International, Inc.                -}
  7. {-                                                                          -}
  8. {-  This unit provides procedures for finding the roots                     -}
  9. {-  of a single equation in one variable.                                   -}
  10. {-                                                                          -}
  11. {----------------------------------------------------------------------------}
  12.  
  13. {$I Float.inc} { Determines the setting of the $N compiler directive }
  14.  
  15. interface
  16.  
  17. {$IFOPT N+}
  18. type
  19.   Float = Double; { 8 byte real, requires 8087 math chip }
  20.  
  21. const
  22.   TNNearlyZero = 1E-015;
  23. {$ELSE}
  24. type
  25.   Float = real;   { 6 byte real, no math chip required }
  26.  
  27. const
  28.   TNNearlyZero = 1E-07;
  29. {$ENDIF}
  30.  
  31.   TNArraySize  = 30;   { maximum size of vectors }
  32.  
  33. type
  34.   TNvector     = array[0..TNArraySize] of Float;
  35.   TNIntVector  = array[0..TNArraySize] of integer;
  36.   TNcomplex    = record
  37.                    Re, Im : Float;
  38.                  end;
  39.   TNCompVector = array[0..TNArraySize] of TNcomplex;
  40.  
  41. procedure Bisect(LeftEnd  : Float;
  42.                  RightEnd : Float;
  43.                  Tol      : Float;
  44.                  MaxIter  : integer;
  45.              var Answer   : Float;
  46.              var fAnswer  : Float;
  47.              var Iter     : integer;
  48.              var Error    : byte;
  49.                  FuncPtr  : Pointer);
  50.  
  51.  
  52. {----------------------------------------------------------------------------}
  53. {-                                                                          -}
  54. {- Input:  LeftEnd, RightEnd, Tol, MaxIter                                  -}
  55. {- Output: Answer, fAnswer, Iter, Error                                     -}
  56. {-                                                                          -}
  57. {-          Purpose: This unit provides a procedure for finding a root      -}
  58. {-                   of a user specified function, for a user specified     -}
  59. {-                   interval, [a,b], where f(a) and f(b) are of opposite   -}
  60. {-                   signs.  The algorithm successively bisects the         -}
  61. {-                   interval, closing in on the root.  The user must       -}
  62. {-                   supply the desired tolerance to which the root should  -}
  63. {-                   be found.                                              -}
  64. {-                                                                          -}
  65. {- Global Variables: LeftEnd  : real    left endpoint                       -}
  66. {-                   RightEnd : real    right endpoint                      -}
  67. {-                   Tol      : real    tolerance of error                  -}
  68. {-                   MaxIter  : real    maximum number of iterations        -}
  69. {-                   Answer   : real    root of the function TNTargetF      -}
  70. {-                   fAnswer  : real    TNTargetF(Answer)                   -}
  71. {-                                      (should be close to 0)              -}
  72. {-                   Iter     :integer  number of iterations                -}
  73. {-                   Error    : byte    flags if something went wrong       -}
  74. {-                                                                          -}
  75. {-           Errors: 0: No error                                            -}
  76. {-                   1: maximum number of iterations exceeded               -}
  77. {-                   2: f(a) and f(b) are not of opposite signs             -}
  78. {-                   3: Tol <= 0                                            -}
  79. {-                   4: MaxIter < 0                                         -}
  80. {-                                                                          -}
  81. {----------------------------------------------------------------------------}
  82.  
  83. procedure Newton_Raphson(Guess    : Float;
  84.                          Tol      : Float;
  85.                          MaxIter  : integer;
  86.                      var Root     : Float;
  87.                      var Value    : Float;
  88.                      var Deriv    : Float;
  89.                      var Iter     : integer;
  90.                      var Error    : byte;
  91.                          FuncPtr1 : Pointer;
  92.                          FuncPtr2 : Pointer);
  93.  
  94. {----------------------------------------------------------------------------}
  95. {-                                                                          -}
  96. {-             Input: Guess, Tol, MaxIter                                   -}
  97. {-            Output: Root, Value, Deriv, Iter, Error                       -}
  98. {-                                                                          -}
  99. {-           Purpose: This unit provides a procedure for finding a single   -}
  100. {-                    real root of a user specified function with a known   -}
  101. {-                    continuous first derivative, given a user             -}
  102. {-                    specified initial guess.  The procedure implements    -}
  103. {-                    Newton-Raphson's algorithm for finding a single       -}
  104. {-                    zero of a function.                                   -}
  105. {-                    The user must specify the desired tolerance           -}
  106. {-                    in the answer.                                        -}
  107. {-                                                                          -}
  108. {-  Global Variables: Guess   : real;    user's estimate of root            -}
  109. {-                    Tol     : real;    tolerance in answer                -}
  110. {-                    MaxIter : integer; maximum number of iterations       -}
  111. {-                    Root    : real;    real part of calculated roots      -}
  112. {-                    Value   : real;    value of the polynomial at root    -}
  113. {-                    Deriv   : real;    value of the derivative at root    -}
  114. {-                    Iter    : real;    number of iterations it took       -}
  115. {-                                       to find root                       -}
  116. {-                    Error   : byte;    flags if something went wrong      -}
  117. {-                                                                          -}
  118. {-            Errors: 1: Iter >= MaxIter                                    -}
  119. {-                    2: The slope was zero at some point                   -}
  120. {-                    3: Tol <= 0                                           -}
  121. {-                    4: MaxIter < 0                                        -}
  122. {-                                                                          -}
  123. {----------------------------------------------------------------------------}
  124.  
  125. procedure Secant(Guess1  : Float;
  126.                  Guess2  : Float;
  127.                  Tol     : Float;
  128.                  MaxIter : integer;
  129.              var Root    : Float;
  130.              var Value   : Float;
  131.              var Iter    : integer;
  132.              var Error   : byte;
  133.                  FuncPtr : Pointer);
  134. {----------------------------------------------------------------------------}
  135. {-                                                                          -}
  136. {-             Input: Guess1, Guess2, Tol, MaxIter                          -}
  137. {-            Output: Root, Value, Iter, Error                              -}
  138. {-                                                                          -}
  139. {-           Purpose: This unit provides a procedure for finding a single   -}
  140. {-                    real root of a user specified function, given a       -}
  141. {-                    specified initial guess.  The procedure implements    -}
  142. {-                    the secant method for finding a single                -}
  143. {-                    root of a function.                                   -}
  144. {-                    The user must specify the desired tolerance           -}
  145. {-                    in the answer.                                        -}
  146. {-                                                                          -}
  147. {-  Global Variables: Guess1  : real;    initial approx #1                  -}
  148. {-                    Guess2  : real;    initial approx #2                  -}
  149. {-                    Tol     : real;    tolerance in answer                -}
  150. {-                    MaxIter : integer; maximum number of iterations       -}
  151. {-                    Root    : real;    real part of calculated roots      -}
  152. {-                    Value   : real;    value of the polynomial at root    -}
  153. {-                    Iter    : real;    number of iterations it took       -}
  154. {-                                       to find root                       -}
  155. {-                    Error   : byte;    flags if something went wrong      -}
  156. {-                                                                          -}
  157. {-            Errors: 1: Iter >= MaxIter                                    -}
  158. {-                    2: The slope was zero at some point                   -}
  159. {-                    3: Tol <= 0                                           -}
  160. {-                    4: MaxIter < 0                                        -}
  161. {-                                                                          -}
  162. {----------------------------------------------------------------------------}
  163.  
  164. procedure Newt_Horn_Defl(InitDegree : integer;
  165.                          InitPoly   : TNvector;
  166.                          Guess      : Float;
  167.                          Tol        : Float;
  168.                          MaxIter    : integer;
  169.                      var Degree     : integer;
  170.                      var NumRoots   : integer;
  171.                      var Poly       : TNvector;
  172.                      var Root       : TNvector;
  173.                      var Imag       : TNvector;
  174.                      var Value      : TNvector;
  175.                      var Deriv      : TNvector;
  176.                      var Iter       : TNIntVector;
  177.                      var Error      : byte);
  178.  
  179. {----------------------------------------------------------------------------}
  180. {-                                                                          -}
  181. {-             Input: InitDegree, InitPoly, Guess, Tol, MaxIter             -}
  182. {-            Output: Degree, NumRoots, Poly, Root, Imag, Value, Deriv      -}
  183. {-                    Iter, Error                                           -}
  184. {-                                                                          -}
  185. {-           Purpose: This unit provides a procedure for finding several    -}
  186. {-                    roots of a user specified polynomial given a user     -}
  187. {-                    specified initial guess.  The procedure implements    -}
  188. {-                    Newton-Horner's algorithm for finding a single        -}
  189. {-                    root of a polynomial and deflation techniques for     -}
  190. {-                    reducing a polynomial given one of its roots.         -}
  191. {-                    Should the polynomial contain no more than two        -}
  192. {-                    complex roots, they may also be determined.           -}
  193. {-                    The user must specify the desired tolerance in the    -}
  194. {-                    answer and the maximum number of iterations.          -}
  195. {-                                                                          -}
  196. {- Pre-Defined Types: TNvector    = array[0..TNArraySize] of real;          -}
  197. {-                    TNIntVector = array[0..TNArraySize] of integer;       -}
  198. {-                                                                          -}
  199. {-  Global Variables: InitDegree : integer;  degree of user's polynomial    -}
  200. {-                    InitPoly   : TNvector; coefficients of user's         -}
  201. {-                                           polynomial where InitPoly[n]   -}
  202. {-                                           is the coefficient of X^n      -}
  203. {-                    Guess      : real;     user's estimate of root        -}
  204. {-                    Tol        : real;     tolerance in answer            -}
  205. {-                    Degree     : real;     degree of reduced polynomial   -}
  206. {-                                           left when procedure is done    -}
  207. {-                                           (>0 if all the roots were      -}
  208. {-                                           not Found)                     -}
  209. {-                    Poly       : TNvector; coefficients of reduced poly   -}
  210. {-                    NumRoots   : integer;  number of roots calculated     -}
  211. {-                    Root       : TNvector; real parts of calculated roots -}
  212. {-                    Imag       : TNvector; imaginary parts of roots (non- -}
  213. {-                                           zero for no more than 2)       -}
  214. {-                    Value      : TNvector; value of the polynomial at     -}
  215. {-                                           each root                      -}
  216. {-                    Deriv      : TNvector; value of the derivative at     -}
  217. {-                                           each root                      -}
  218. {-                    Iter       : TNIntVector; number of iterations it     -}
  219. {-                                              took to find each root      -}
  220. {-                    Error      : byte;     flags if something went wrong  -}
  221. {-                                                                          -}
  222. {-            Errors: 0: No error                                           -}
  223. {-                    1: Iter >= MaxIter                                    -}
  224. {-                    2: The slope was zero at some point                   -}
  225. {-                    3: Degree <= 0                                        -}
  226. {-                    4: Tol <= 0                                           -}
  227. {-                    5: MaxIter < 0                                        -}
  228. {-                                                                          -}
  229. {----------------------------------------------------------------------------}
  230.  
  231. procedure Muller(Guess   : TNcomplex;
  232.                  Tol     : Float;
  233.                  MaxIter : integer;
  234.              var Answer  : TNcomplex;
  235.              var yAnswer : TNcomplex;
  236.              var Iter    : integer;
  237.              var Error   : byte;
  238.                  FuncPtr : Pointer);
  239.  
  240. {----------------------------------------------------------------------------}
  241. {-                                                                          -}
  242. {-               Input: Guess, Tol, MaxIter                                 -}
  243. {-              Output: Answer, yAnswer, Iter, Error                        -}
  244. {-                                                                          -}
  245. {-             Purpose: This program uses Muller's method to find a root    -}
  246. {-                      of a user defined function Y=TNTargetF given an     -}
  247. {-                      initial approximation.  The root may be complex.    -}
  248. {-                                                                          -}
  249. {-                                                                          -}
  250. {-   User-Defined                                                           -}
  251. {-          Procedures: TNTargetF(X : TNcomplex; VAR Y : TNcomplex);        -}
  252. {-                                                                          -}
  253. {-  User-Defined Types: TNcomplex = record                                  -}
  254. {-                                    Re, Im : real;                        -}
  255. {-                                  end;                                    -}
  256. {-                                                                          -}
  257. {-    Global Variables: Guess   : real;            initial guess            -}
  258. {-                      Tol     : real;            tolerance in the         -}
  259. {-                                                 answer                   -}
  260. {-                      MaxIter : integer;         maximum number of        -}
  261. {-                                                 iterations               -}
  262. {-                      Answer  : TNcomplex;       a root of the            -}
  263. {-                                                 polynomial               -}
  264. {-                      yAnswer : TNcomplex;       value of the             -}
  265. {-                                                 polynomial at the        -}
  266. {-                                                 root (close to zero)     -}
  267. {-                      Iter    : integer;         number of iterations     -}
  268. {-                                                 it took to find root     -}
  269. {-                      Error   : byte;            flags an error           -}
  270. {-                                                                          -}
  271. {-              Errors: 0: No errors                                        -}
  272. {-                      1: Iter > MaxIter                                   -}
  273. {-                      2: parabola could not                               -}
  274. {-                         be formed                                        -}
  275. {-                      3: Tol <= 0                                         -}
  276. {-                      4: MaxIter < 0                                      -}
  277. {-                                                                          -}
  278. {----------------------------------------------------------------------------}
  279.  
  280. procedure Laguerre(var Degree    : integer;
  281.                    var Poly      : TNCompVector;
  282.                        InitGuess : TNcomplex;
  283.                        Tol       : Float;
  284.                        MaxIter   : integer;
  285.                    var NumRoots  : integer;
  286.                    var Roots     : TNCompVector;
  287.                    var yRoots    : TNCompVector;
  288.                    var Iter      : TNIntVector;
  289.                    var Error     : byte);
  290.  
  291. {----------------------------------------------------------------------------}
  292. {-                                                                          -}
  293. {-            Input: Degree, Poly, InitGuess, Tol, MaxIter                  -}
  294. {-           Output: Degree, Poly, NumRoots, Roots, yRoots, Iter, Error     -}
  295. {-                                                                          -}
  296. {-          Purpose: This unit provides a procedure for finding all the     -}
  297. {-                   roots (real and complex) to a polynomial.              -}
  298. {-                   Laguerre's method with deflation is used.              -}
  299. {-                   The user must input the coefficients of the quadratic  -}
  300. {-                   and the tolerance in the answers generated.            -}
  301. {-                                                                          -}
  302. {-  Pre-defined Types: TNcomplex    = record                                -}
  303. {-                                      Re, Im : real;                      -}
  304. {-                                    end;                                  -}
  305. {-                                                                          -}
  306. {-                     TNIntVector  = array[0..TNArraySize] of integer;     -}
  307. {-                     TNCompVector = array[0..TNArraySize] of TNcomplex;   -}
  308. {-                                                                          -}
  309. {- Global Variables: Degree    : integer;      degree of deflated           -}
  310. {-                                             polynomial                   -}
  311. {-                   Poly      : TNCompVector; coefficients of deflated     -}
  312. {-                                             polynomial where Poly[n] is  -}
  313. {-                                             the coefficient of X^n       -}
  314. {-                   InitGuess : TNcomplex;    initial guess to a root      -}
  315. {-                                             (may be very crude)          -}
  316. {-                   Tol       : real;         tolerance in the answer      -}
  317. {-                   MaxIter   : integer;      number of iterations         -}
  318. {-                   NumRoots  : integer;      number of roots calculated   -}
  319. {-                   Roots     : TNCompVector; the roots calculated         -}
  320. {-                   yRoots    : TNCompVector; the value of the function    -}
  321. {-                                             at the calculated roots      -}
  322. {-                   Iter      : TNIntVector;  number iteration it took to  -}
  323. {-                                             find each root               -}
  324. {-                   Error     : byte;         flags an error               -}
  325. {-                                                                          -}
  326. {-           Errors: 0: No error                                            -}
  327. {-                   1: Iter > MaxIter                                      -}
  328. {-                   2: Degree <= 0                                         -}
  329. {-                   3: Tol <= 0                                            -}
  330. {-                   4: MaxIter < 0                                         -}
  331. {-                                                                          -}
  332. {----------------------------------------------------------------------------}
  333.  
  334. implementation
  335.  
  336. {$F+}
  337. {$L RootsEqu.OBJ} { Link in external routines }
  338.  
  339. function UserFunction(X : Float; ProcAddr : Pointer) : Float; external;
  340.  
  341. procedure UserProcedure(X : TNcomplex; var Y : TNcomplex; ProcAddr : Pointer); external;
  342. {$F-}
  343.  
  344. {$I RootsEqu.inc}  { Include procedure code }
  345.  
  346. end. { RootsEqu }
  347.