Previous Page Next Page Contents

numeric::fsolve -- search for a numerical root of a system of equations

Introduction

numeric::fsolve(eqs, ..) returns a numerical approximation of a solution of the system of equations eqs.

Call(s)

numeric::fsolve(eq, x  <, Options>)
numeric::fsolve(eq, x = a  <, Options>)
numeric::fsolve(eq, x = a..b  <, Options>)
numeric::fsolve(eqs, [x1, x2, ..]  <, Options>)
numeric::fsolve(eqs, {x1, x2, ..}  <, Options>)
numeric::fsolve(eqs, [x1 = a1, x2 = a2, ..]  <, Options>)
numeric::fsolve(eqs, {x1 = a1, x2 = a2, ..}  <, Options>)
numeric::fsolve(eqs, [x1 = a1..b1, x2 = a2..b2, ..]  <, Options>)
numeric::fsolve(eqs, {x1 = a1..b1, x2 = a2..b2, ..}  <, Options>)

Parameters

eq - an arithmetical expression or an equation in one indeterminate x. An expression eq is interpreted as the equation eq = 0.
x - an identifier or an indexed identifier to be solved for.
a - real or complex numerical starting value for the internal search. Typically, a crude approximation of the solution.
a..b - a range of numerical values defining a search interval for the numerical root.
eqs - a list or a set of expressions or equations in several indeterminates x1, x2, ... Expressions are interpreted as homogeneous equations.
x1, x2, .. - identifiers or indexed identifiers to be solved for.
a1, a2, .. - real or complex numerical starting values for the internal search. Typically, crude approximations of solution.
a1..b1, a2..b2, .. - ranges of numerical values defining search intervals for the numerical root.

Options

RestrictedSearch - makes numeric::fsolve return only numerical roots in the user-defined search range x = a..b and [x1 = a1..b1, x2 = a2..b2, ..], respectively. This is the default search strategy, if a search range with real range parameters is specified for at least one of the unknowns.
UnrestrictedSearch - allows numeric::fsolve to find and return solutions outside the specified search range. With this option, the search range is only used to choose random starting points for the internal numerical search.
MultiSolutions - makes numeric::fsolve return all solutions found in the internal search
Random - With this option, several calls to numeric::fsolve with the same input parameters may produce different roots.

Returns

A single numerical root is returned as a list of equations [x = value] or [x1 = value1, x2 = value2, ..], respectively. FAIL is returned, if no solution is found. With the option MultiSolutions, sequences of solutions may be returned.

Side Effects

The function is sensitive to the environment variable DIGITS, which determines the numerical working precision.

Related Functions

linsolve, numeric::linsolve, numeric::realroot, numeric::realroots, numeric::polyroots, numeric::polysysroots, numeric::solve, polylib::realroots, solve

Details

Option: RestrictedSearch

Option: UnrestrictedSearch

Option: MultiSolutions

Option: Random

Example 1

We compute roots of the sine function:

>> numeric::fsolve(sin(x) = 0, x)
                            [x = -226.1946711]

With the option Random, several calls may result in different roots:

>> numeric::fsolve(sin(x), x, Random)
                             [x = 97.38937226]
>> numeric::fsolve(sin(x), x, Random)
                             [x = 53.40707511]

Particular solutions can be chosen by an appropriate starting point close to the wanted solution, or by a search interval:

>> numeric::fsolve(sin(x), x = 3), 
   numeric::fsolve(sin(x), x = -4..-3)
                   [x = 3.141592653], [x = -3.141592654]

The solutions found by numeric::fsolve can be used in subs and assign to substitute or assign the indeterminates:

>> eqs := [x^2 = sin(y), y^2 = cos(x)]:
>> solution := numeric::fsolve(eqs, [x, y])
                   [x = -0.8517004887, y = 0.8116062152]
>> eval(subs(eqs, solution)) 
        [0.7253937224 = 0.7253937224, 0.6587046485 = 0.6587046485]
>> assign(solution): x, y
                        -0.8517004887, 0.8116062152
>> delete eqs, solution, x, y:

Example 2

We demonstrate the use of search ranges. The following system has solutions with positive and negative x. The solution with x > 0 is obtained with the search interval x = 0..infinity:

>> numeric::fsolve([x^2 = exp(x*y), x^2 = y^2], 
                   [x = 0..infinity, y])
                   [x = 0.7530891649, y = -0.753089165]
>> numeric::fsolve([x^2 = exp(x*y), x^2 = y^2], 
                   [x = -infinity..0, y])
                   [x = -0.753089165, y = 0.7530891649]

Example 3

Multiple roots can only be computed with a restricted precision:

>> numeric::fsolve(expand((x - 1/3)^5), x = 0.3)
                            [x = 0.3333929968]

Example 4

The following system of equations is degenerate and has a 1-parameter family of solutions. Each call to numeric::fsolve picks out one random solution:

>> numeric::fsolve([x^2 - y^2, x^2 - y^2], [x, y], Random) $ i=1..3
      [x = -140.1698476, y = 140.1698476],
      
         [x = 34.70258251, y = 34.70258251],
      
         [x = -29.16650501, y = 29.16650501]

The equation may also be specified as an underdetermined system:

>> numeric::fsolve([x^2 - y^2], [x, y])
                    [x = -140.1698476, y = 140.1698476]

Example 5

The following equation has no real solution. Consequently, the numerical search with real starting values fails:

>> numeric::fsolve(sin(x) + cos(x)^2 = 3, x)
                                   FAIL

With a complex starting value, a solution is found:

>> numeric::fsolve(sin(x) + cos(x)^2 = 3, x = I)
                    [x = 0.2972513613 + 1.128383965 I]

Also complex search ranges may be specified. In the following, the internal starting point is a random value on the line from 2 + I to 3 + I:

>> numeric::fsolve(sin(x) + cos(x)^2 = 3, x = 2 + I..3 + I)
                     [x = 2.844341292 + 1.128383965 I]

Example 6

Starting values and search intervals can be mixed:

>> numeric::fsolve([x^2 + y^2 = 1, y^2 + z^2 = 1, x^2 + z^2 = 1],
                   [x = 1, y = 0..10, z])
          [x = 0.7071067812, y = 0.7071067812, z = 0.7071067812]

Example 7

With UnrestrictedSearch, search intervals are only used for choosing starting values for the internal Newton search. The numerical iteration may drift towards a solution outside the search range:

>> eqs := [x*sin(10*x) = y^3, y^2 = exp(-2*x/3)]: 
>> numeric::fsolve(eqs, [x = 0..1, y = -1..0], 
                   UnrestrictedSearch)
                    [x = 1.232766202, y = -0.663038602]

With the default strategy RestrictedSearch, only solutions inside the search range are accepted:

>> numeric::fsolve(eqs, [x = 0..1, y = -1..0])
                   [x = 0.9816416007, y = -0.7209295436]

In the last search, also the previous solution outside the search range was found. With the option MultiSolutions, numeric::fsolve returns a sequence of all solutions that were found in the internal search:

>> numeric::fsolve(eqs, [x = 0..1, y = -1..0], MultiSolutions)
      [x = 0.9816416007, y = -0.7209295436],
      
         [x = 1.232766202, y = -0.663038602]
>> delete eqs:

Example 8

Usually, most of the time is spent internally searching for some (crude) approximations of the root. If high precision roots are required, it is recommended to compute first approximations with moderate values of DIGITS and use them as starting values for a refined search:

>> eq := exp(-x) = x:
>> DIGITS := 10: firstApprox := numeric::fsolve(eq, x)
                            [x = 0.5671432904]

This output is suitable as input defining a starting value for x:

>> DIGITS := 1000: numeric::fsolve(eq, firstApprox)
      [x = 0.5671432904097838729999686622103555497538...]
>> delete eq, firstApprox, DIGITS:

Example 9

Specifying starting values for the indeterminates launches a single Newton iteration. This may fail, if the starting values are not sufficiently close to the solution:

>> eq := [x*y = x + y - 4, x/y = x - y + 4]:
>> numeric::fsolve(eq, [x = 1, y = 1])
                                   FAIL

If a search range is specified for at least one of the unknowns, then several Newton iterations with random starting values in the search range are used, until a solution is found or until numeric::fsolve gives up:

>> numeric::fsolve(eq, [x = 1, y = 0..10])
                      [x = 4.026449604e-14, y = 4.0]
>> delete eq:

Background

Changes




Do you have questions or comments?


Copyright © SciFace Software GmbH & Co. KG 2000