INVERSE FUNCTION PAPER

ITERATION ROUTINE

by Jerry Hughes
March 18, 1997

Overview:

This paper describes an iteration process which is useful for determining the domain value of a given range value. An inverse function routine where,

y = f(x) and its inverse x = f(y).
The process is a slight modification of the FOR-NEXT or DO-LOOP originally used to iterate, with the desired degree of accuracy, toward solutions to such problems.
The process described has 3 advantages over the original method and the Newtonian approximation method:

1. More efficient, less than 11 iterations for each significant digit required plus some number of iterations which is a function of the initial guess and increment value.

2. Does not require the derviative of the function or use of numerical differenation routines.

3. Is not sensitive to the closeness of the initial guess or largness of initial increment value.

    This process, with variations in "testing for closeness," is also applicable for determining: local maximum and minimum points, inflection points, intersection of two functions.
    Although this process is described for a function in the form of a simple equation, it is also applicable for iteration of numeric models in the form of a group of equations or other logical group of computer code.
    One restriction is that the overall results of the model, equations or code should be reasonably smooth and continuous.

Original process:

In the original form of the FOR-NEXT or DO-LOOP iteration process the principle components consisted of,

1. Define initial values: initial guess and increment value (xinc),

2. Step through the loop and test for closeness,

3. Drop out of the loop when the test passes.

The increment value is the desired degree of accuracy,usually 1 x 10^n where "n" is the position of the least significant digit required. Examples:

1 x 10^2 yields dddd00

1 x 10^-6 yields ddd.dddddd

The efficiency of this process is dependant upon the initial starting value or guess and the increment value. With an initial value of "0" and incremental value of "1" the first example above requires thousands of iterations, the second 100's of millions.
The test for closeness usually consists of evaluating the difference between either two consecutive domain values or the given range value and currently guessed range value.

Modified process:

In the modified iteration process the components consist of,

1. Define initial values: domain limits (x minimum and x maximum), increment and tolerance values, initialize solution counter, and set iteration time limit,

2. Step through the loop and test for closeness,

3. If close then test for degree of accuracy,

4. If within degree of accuracy then reset increment value to its original magnitude and continue loop,

5. If not then a) refine increment value, b) back up one guess including refinement, c) continue loop,

6. If maximum domain limit is reached then drop out of loop.

     The changes in step 1, as compared to the original, allows for the occurance of more than one solution as well as none. Specifing these domain limits also defines the initial guess (x minimum) and a loop drop-out test value (x maximum).
     Also, seperating the dual purpose of the increment value, incrementation and tolerance test, allows the refinement of the iteration guesses.
     Since the modified process now allows for multiple solutions within a specified domain, a counter is used to note how many were found. If an array is used to save the solutions the counter can be tested against array limits. If a file is used it documents the number of solutions in the file.
     In order to account for the possiblity of never passing the accuracy test, example 0.9999999 does not equal 1.000000, a time limit is employed. This consists of obtaining the clock value immediately before entering the loop, adding some number of seconds to it and then during the iteration loop compare this "time limit" with the clock.
     The efficiency of this process is improved by changing the test for closeness and refining the increment value, steps 2 and 5. These two changes allow for initially large increment values and thus fewer iterations.
     First, the test for closeness is simplified to one of testing to see if the given range value is between the currently guessed range value and the next adjecant value.
     If true then the process needs to back up one guess plus a little and then use more refined guesses. The extra little is to cause the process to repeat the previous guess so that there are no potential gaps in the iteration.
     Second, the refinement of guessing is accomplished by simply decrementing the increment value by an order of magnitude. Thus,
xinc = xinc * .1
     After this first occurance of refinement, the iteration process becomes one of determining subsequent significant digits till the least required is obtained.
     Specifing the accuracy is the same as that for the increment value in step 1 of the original process. Except a more appropriate term would be "tolerance value" along with a different variable name, xtol.
     Third, the test for accuracy can be the same as those for "closeness" in the original process. There is also the option of simply testing the current increment value to see if it is less than or equal to the tolerance value. Thus,
IF xinc <= xtol THEN quit

Variations:

     The variations, as stated before, consist of changes in testing for closeness.
     For local maximum or minimum points in a function the test for closeness is to sense a change in slope. This is accomplished by numerical differenation of the current domain value and the next incremental value then comparing the sign (slope) with the previously determined sign.
     The components of the process are,

1. Define initial values and initialize switches and counters, and set time limit,

2. Determine the initial slope,

3. Begin the loop,

4. Determine slope for current domain value,

5. Test for sign change in slope (test for closeness),

6. If changed then back up 2 points and decrement the increment value,

7. Test for accuracy,

8. If within tolerance then save solution, reset increment value to original magnitude, change current slope sign, continue loop,

9. Do loop drop-out tests.

     Step 1, local maximum, minimum switches serve the purpose of assigning solutions to the appropriate array or file or otherwise identifing and tagging solutions. These switches are initially set false (-1) and when a local maximum or minimum is sensed set true (+1).
     The logic of the slope switches is +1 for positive and -1 for negative.
     Step 2, determine initial slope, is done by assuming a positive slope, calculating difference between the range values of the minimum domain and next "x" value, then resetting the slope switch if the difference is negative.
     A maximum point is determined by a change from positive to negative slope, conversely, a minimum point from negative to positive.
     Steps 4 and 5 are repetitions of step 2 except the current domain value and its next adjacent value are used and the comparison is between current and next slope signs.
     Step 6, if the change is from positive to negative set the local maximum switch to "+1". If it is from negative to positive set the local minimum switch to "+1".
     Step 6 is also where the incremental value is refined for better guesses.
     The iteration is backed up 2 steps since the process is very similar to that for finding inflection points. The exception is that the current domain value and the next 2 increments are used to determine rate of change in slope. With the addition of a few statements, the same code can be used for both local max/min and inflection points.
     Step 7 is the test to see if the "xinc" variable is less than the "xtol" value.
     Step 8 does the solution bookkeeping, updates the logic switches and resets the increment variable to its original large value.
     The drop out tests, step 9, now includes a) "x" greater than maximum domain limit and b) exceeding time limit within the loop.

     For inflection points, adding logic switches, counters and tests for closeness and accuracy is all that is required.
     The logic switches are "rate" and two inflection variables, one for a positive to negative change in rate and the other a negative to positive. The rate switch uses "+1" for an increasing slope and "-1" for a decreasing slope.
     The test for closeness is to compare two adjacent dy/dx slopes and the current rate. If the first dy/dx is greater than the second dy/dx, (decreasing rate), and the current rate was increasing then an inflection point is near. If the first dy/dx is less than the second dy/dx, (increasing rate), and the current rate was decreasing then an inflection point is near.
     When either is true, set the apporpriate inflection switch, back up two places, refine the increment value and continue the loop.
     The remainder of the steps are the same as for local max/min points except the addition of doing the bookkeeping and resetting for inflection variables and switches.

Intersection of two curves:

To find an intersection of two curves, first determine if the range value of the first is larger or smaller than that of the second for the given minimum domain value. Then step through the domain and sense when the first becomes smaller or larger than the second. The components of the routine are,

1. Define domain limits, initialize counters and switches, and set time limit,

2. Determine initial range values and set the test switch,

3. Step through the domain and sense when the range values change relative to each other (test for closeness),

4. If close then back up one step plus a little, refine increment value, continue loop,

5. Test for accuracy {NOT (xinc <= xtol)},

6. If true then save solution, reset test switch, and reset increment value to original magnitude,

7. Do drop out of loop tests.

The logic switch, Test, is set to "-1" for the second function being larger than the first and "+1" for the reverse.
Since an intersection point is to be found, the absolute difference between the two range values is compared to the tolerance value. Otherwise, an acceptable domain tolerance may result in an un-acceptable range tolerance.

Cautions:

     Iteration loop. By providing a 'time out' test all possible causes of an endless loop are accounted for. The amount of time to allow for the process is dependent upon the complexity of the function code, increment and tolerance values and domain limits.
     Simple equations with loose values and small limits should require only a few seconds. Numerical models of mechanisms with large number of lines of code, small values and large limits may require minutes.
     Experience with any given application is the only guidance.
     Initial increment values. For non-cyclic functions this value should be some percentage of the domain to be searched. For cyclic functions the percentage should be of a single cycle.
     The increment should be less than the distance between adjacent points of interest and never large enough to enclose complete cycles. To do so may result in unpredictable behavior of the code.
     Tolerance (accuracy) value. This value determines the least significant digit of solutions and is determined by the application of this iteration routine.
     It is limited by the double-precision of the machine.
     Double-precision. Results of functions, either equations or code, should be assigned to double-precision variables. The tolerance variable is double-precision and requires comparsion to like variables. The same rule applies to other tests for closeness and accuracy.

The Code:

     The code, as found in the companion file ITERCODE.HTM, is written for Microsoft's QB45 syntax. As such it requires that language development environment to operate.
     The source code is designed to be run from within the QB45 environment. The equations or numerical model code are treated as DEF FUNC routines to allow for development of additional code as required for the particular application. That is the user is expected to substitute his own code within the FUNCTIONS to be evaluated.
     Care should be given to avoid variable name conflicts.