Difference between revisions of "NumLib Documentation"

From Free Pascal wiki
Jump to navigationJump to search
(Initial version)
 
Line 17: Line 17:
 
* <code>f</code> is the function for which the root is to be determined. It must be a function of one floating point argument (type <code>ArbFloat</code>. The type of the function, <code>rfunc1r</code>, is declared in unit <code>typ</code>.
 
* <code>f</code> is the function for which the root is to be determined. It must be a function of one floating point argument (type <code>ArbFloat</code>. The type of the function, <code>rfunc1r</code>, is declared in unit <code>typ</code>.
 
*<code>a</code> and <code>b</code> are the endpoints of the test interval. The root must be located between these two values, i. e. and the function values ''f(a)'' and ''f(b)'' must have different signs.
 
*<code>a</code> and <code>b</code> are the endpoints of the test interval. The root must be located between these two values, i. e. and the function values ''f(a)'' and ''f(b)'' must have different signs.
*<code>ae</code> and <code>re</code> determine the absolute and relative precision, respectively, with which the root will be determined. <code>re</code> is relative to the maxiumum of ''abs(a)'' and ''abs(b)''. Note that precision and speed are conflicting issues. Highest accuracy is achieved if <code>ae</code> is given as <code>MachEps</code> (see unit <code>typ</code>). Both parameters must not be negative.
+
*<code>ae</code> and <code>re</code> determine the absolute and relative precision, respectively, with which the root will be determined. <code>re</code> is relative to the maximum of ''abs(a)'' and ''abs(b)''. Note that precision and speed are conflicting issues. Highest accuracy is achieved if <code>ae</code> is given as <code>MachEps</code> (see unit <code>typ</code>). Both parameters must not be negative.
 
*<code>x</code> returns the value of the found root.
 
*<code>x</code> returns the value of the found root.
 
*<code>term</code> returns whether the process has been successful:
 
*<code>term</code> returns whether the process has been successful:

Revision as of 20:31, 15 April 2017

UNDER CONSTRUCTION

Introduction

Data types and declarations (unit "typ")

Finding the roots of a function (unit "roo")

The roots are the x values at which a function f(x) is zero.

Bisection

In the bisection method two x values a and b are guessed which enclose the root and at which the function has opposite signs. The center point of the interval is determined, and the search contines with that subinterval where the function values of the end points have oppositive signs again.

In NumLib, this approach is supported by the procedure roof1r:

procedure roof1r(f: rfunc1r; a, b, ae, re: ArbFloat; var x: ArbFloat; var term: ArbInt);
  • f is the function for which the root is to be determined. It must be a function of one floating point argument (type ArbFloat. The type of the function, rfunc1r, is declared in unit typ.
  • a and b are the endpoints of the test interval. The root must be located between these two values, i. e. and the function values f(a) and f(b) must have different signs.
  • ae and re determine the absolute and relative precision, respectively, with which the root will be determined. re is relative to the maximum of abs(a) and abs(b). Note that precision and speed are conflicting issues. Highest accuracy is achieved if ae is given as MachEps (see unit typ). Both parameters must not be negative.
  • x returns the value of the found root.
  • term returns whether the process has been successful:
    • 1 - successful termination, a zero point has been found with an absolute accuracy ae, or a relative precision re
    • 2 - the required accuracy of the root could not be reached; However, the value of x is called the "best achievable" approach
    • 3 - error in input parameters: ae < 0 or re < 0, or f(a)*f(b) > 0

Example

The following program determines the square root of 2. This is the number at which the function f(x) = x^2 - 2 is zero. Since f(1) = 1^2 - 2 = -^< 0 and f(2) = 2^2 - 2 = 2 > 0 we can assue a and b to be 1 and 2, respectively.

program bisection_demo;

uses
  typ, roo;

function f(x: ArbFloat): ArbFloat;
begin
  Result := x*x - 2;
end;

var
  x: ArbFloat = 0.0;
  term : ArbInt;

begin
  roof1r(@f, 1.0, 2.0, 1e-9, 0, x, term);
  WriteLn('Bisection result ', x);
  WriteLn('sqrt(2)          ', sqrt(2.0));
end.

Special functions (unit "spe")

Solving linear systems of equations (unit "sle")