Documentation
¶
Index ¶
Examples ¶
Constants ¶
This section is empty.
Variables ¶
This section is empty.
Functions ¶
Types ¶
type DerivativeApprox ¶
type DerivativeApprox struct {
// Function of from which the jacobian should be calculated
F VecFunc
// Stepsize used to calculate the inner product between the jacobian and a vector v
// J.dot(v) = (F(x + Eps*v) - F(x))/Eps
Eps float64
// X is the position at which the jacobian should be calculated
X []float64
// contains filtered or unexported fields
}
DerivativeApprox is a type that is used to approximate the dot product between a search direction and the jacobian of a function
func NewCentral ¶
func NewCentral(F VecFunc, Eps float64) DerivativeApprox
NewCentral returns a derivative approximator of using the central differences. Jv = (F(x + eps*v) - F(x - eps*v)/(2*eps), where J is the jacobian
func NewEightPoint ¶
func NewEightPoint(F VecFunc, Eps float64) DerivativeApprox
NewEightPoint returns a derivative approximator of using the central eight point differences. Jv = (-F(x-3*eps*v) + 9*F(x-2*eps*v) - 45*F(x-eps*v)
+45*F(x-eps*v) - 9*F(x+2*eps*v) + F(x+3*v))/(60*eps)
where J is the jacobian
func NewFourPoint ¶
func NewFourPoint(F VecFunc, Eps float64) DerivativeApprox
NewFourPoint returns a derivative approximator of using the central four point differences. Jv = (F(x-2*eps*v) - 8*F(x-eps*v) + 8*F(x+eps*v) - F(x+2*v))/(12*eps) where J is the jacobian
func NewSixPoint ¶
func NewSixPoint(F VecFunc, Eps float64) DerivativeApprox
NewSixPoint returns a derivative approximator of using the central six point differences. Jv = (-F(x-3*eps*v) + 9*F(x-2*eps*v) - 45*F(x-eps*v)
+45*F(x-eps*v) - 9*F(x+2*eps*v) + F(x+3*v))/(60*eps)
where J is the jacobian
type NewtonKrylov ¶
type NewtonKrylov struct {
// Tolerance for convergence. The method terminates when InfNorm(F(x)) + InfNorm(dx)
// is less than this value. Here, dx represents the amound the position changed on the
// current Newton step
Tol float64
// Stepsize used in finite difference approximation of the Jacobian.
StepSize float64
// Maximum number of outer iterations (e.g. the number of Newton steps). If not given
// a default value of 10000 will be used
Maxiter int
// Number of points used to approximate the jacobian. If not given (or set to zero)
// a four point stencil will be used
Stencil int
// InnerMethod is the linear solver used in the "inner" iterations of the NewtonKrylov
// method. If not given, GMRES with the default parameters in Gonum will be used
InnerMethod linsolve.Method
// Settings passed to linsolve.Iterative.
// See https://godoc.org/github.com/gonum/exp/linsolve#Iterative
// for further details
InnerSettings *linsolve.Settings
}
NewtonKrylov solve a non-linear system of equation using the conjugate gradient Jacobian Free Newton Krylov method. It is Newton based method, but the Jacobian is never explicitly constructed, since we only needs its action onto a vector (e.g. Jv, where J is the jacobian matrix and v is a vector). This, makes the method memory efficient. The jacobian matrix along a particular direction is however estimates using finithe differences. Jv = (F(x + eps*v) - F(x))/eps where eps is the Stepsize. F is the non-linear (vector-valued) function of which one seeks the folusion F(x) = 0. The method will try at most Maxiter Newton steps before it terminates. In case, the maximum number of iterations is reached the Converged attribute of the Result struct will be false.
Example ¶
package main
import (
"fmt"
"math"
"github.com/davidkleiven/gononlin/nonlin"
)
func main() {
// This example shows how one can use NewtonKrylov to solve the
// system of equations
// (x-1)^2*(x - y) = 0
// (x-2)^3*cos(2*x/y) = 0
problem := nonlin.Problem{
F: func(out, x []float64) {
out[0] = math.Pow(x[0]-1.0, 2.0) * (x[0] - x[1])
out[1] = math.Pow(x[1]-2.0, 3.0) * math.Cos(2.0*x[0]/x[1])
},
}
solver := nonlin.NewtonKrylov{
// Maximum number of Newton iterations
Maxiter: 1000,
// Stepsize used to appriximate jacobian with finite differences
StepSize: 1e-2,
// Tolerance for the solution
Tol: 1e-7,
}
x0 := []float64{0.0, 3.0}
res := solver.Solve(problem, x0)
fmt.Printf("Root: (x, y) = (%.2f, %.2f)\n", res.X[0], res.X[1])
fmt.Printf("Function value: (%.2f, %.2f)\n", res.F[0], res.F[1])
}
Output: Root: (x, y) = (1.00, 2.00) Function value: (-0.00, 0.00)
func (*NewtonKrylov) Solve ¶
func (nk *NewtonKrylov) Solve(p Problem, x []float64) Result
Solve solves the non-linear system of equations. The method terminates when the inifinity norm of F + the inifinity norm of dx is less than the tolerance. dx is the change in x between two sucessive iterations.
type Problem ¶
type Problem struct {
// F is a function where we seek the solution of F(x) = 0. The value of the function
// should be given to out. out will always be the same length as x. Note that F must
// not modify x.
F func(out, x []float64)
}
Problem to be solved
type Result ¶
type Result struct {
// X is the solution vector at termination
X []float64
// F is the value the function at termination
F []float64
// MaxF is the infinity norm of F
MaxF float64
// Converges is True if the convergence criteria are met, otherwise it is False
Converged bool
}
Result is a result type returned by the solvers.
Source Files