diff --git a/src/FSharp/FindRoots.fs b/src/FSharp/FindRoots.fs index fb6304d3..30e4c30e 100644 --- a/src/FSharp/FindRoots.fs +++ b/src/FSharp/FindRoots.fs @@ -66,6 +66,10 @@ module FindRoots = | true, root -> Some root | false, _ -> None + let broyden maxIterations accuracy guess (f:double[]->double[]) = + match Broyden.TryFindRoot(tobcl f, guess, accuracy, maxIterations) with + | true, root -> Some root + | false, _ -> None // simple usage diff --git a/src/FSharpUnitTests/FindRootsTests.fs b/src/FSharpUnitTests/FindRootsTests.fs index 393eb79d..a74e5895 100644 --- a/src/FSharpUnitTests/FindRootsTests.fs +++ b/src/FSharpUnitTests/FindRootsTests.fs @@ -9,6 +9,13 @@ module FindRootsTests = let f x = (x - 3.0)*(x - 4.0) let df x = 2.0*x - 7.0 + let g (xa:float[]) = + let x = xa.[0]; + let T = xa.[1]; + let k = 0.12 * Math.Exp(12581.0 * (T - 298.0) / (298.0 * T)); + let fx = 120.0 * x - 75.0 * k * (1.0 - x); + let fT = -x * (873.0 - T) + 11.0 * (T - 300.0); + [|fx; fT|] [] let ``Bisection should find both roots of (x - 3) * (x - 4)``() = @@ -35,9 +42,19 @@ module FindRootsTests = (f, df) ||> FindRoots.newtonRaphsonRobust 100 20 1e-14 -5.0 3.5 |> should equal (Some 3.0) (f, df) ||> FindRoots.newtonRaphsonRobust 100 20 1e-14 3.2 5.0 |> should equal (Some 4.0) + [] + let ``Bryoden should find both roots of (x - 3) * (x - 4) and Twoeq2``() = + f |> (fun g (x:float[]) -> [|g x.[0]|]) |> FindRoots.broyden 100 1e-14 [|1.0;|] |> should equal (Some [|3.0|]) + f |> (fun g (x:float[]) -> [|g x.[0]|]) |> FindRoots.broyden 100 1e-14 [|9.0;|] |> should equal (Some [|4.0|]) + g |> FindRoots.broyden 100 1e-6 [|1.0; 400.0;|] |> (fun x -> + match x with + | Some x -> should (equalWithin 1e-1) [|0.9638680512795; 346.16369814640|] + | None -> failwith "The element in array is not equal.") |> ignore + [] let ``Simple method should find both roots of (x - 3) * (x - 4)``() = f |> FindRoots.ofFunction -5.0 3.5 |> Option.get |> should (equalWithin 1e-8) 3.0 f |> FindRoots.ofFunction 3.2 5.0 |> Option.get |> should (equalWithin 1e-8) 4.0 (f, df) ||> FindRoots.ofFunctionDerivative -5.0 3.5 |> Option.get |> should (equalWithin 1e-8) 3.0 (f, df) ||> FindRoots.ofFunctionDerivative 3.2 5.0 |> Option.get |> should (equalWithin 1e-8) 4.0 +