From d398989ece505d00070a97c79d191b8a5416adc6 Mon Sep 17 00:00:00 2001 From: Henri Saudubray Date: Tue, 29 Apr 2025 15:33:57 +0200 Subject: [PATCH] feat (exm): sincos example --- exm/ball.ml | 4 +++ exm/sincos.ml | 44 ++++++++++++++++++++++++++ exm/vdp.ml | 68 ++++++++++++++++------------------------- src/bin/main.ml | 20 +++++++----- src/lib/common/debug.ml | 5 +-- src/lib/hsim/sim.ml | 4 --- 6 files changed, 87 insertions(+), 58 deletions(-) create mode 100644 exm/sincos.ml diff --git a/exm/ball.ml b/exm/ball.ml index 542f883..97da4d4 100644 --- a/exm/ball.ml +++ b/exm/ball.ml @@ -52,3 +52,7 @@ let bouncing_ball () = horizon = (fun _ -> max_float); cset; cget; zset; reset; jump; zsize } +let errmsg = "Too many arguments for the model (needed: 0)" +let bouncing_ball = function + | [] -> bouncing_ball () + | _ -> raise (Invalid_argument errmsg) diff --git a/exm/sincos.ml b/exm/sincos.ml new file mode 100644 index 0000000..be1d6c6 --- /dev/null +++ b/exm/sincos.ml @@ -0,0 +1,44 @@ + +open Hsim.Types +open Solvers.Zls + +let of_array a = Bigarray.Array1.of_array Bigarray.Float64 Bigarray.c_layout a + +type state = { si : bool; sx : carray } + +let yd = cmake 3 +let zout = cmake 1 +let zsize = 1 + +let fzer _ _ _ = zout +let fout _ _ y = of_array [| y.{0}; y.{1}; y.{2} |] +let cget s = s.sx +let cset s lx = { s with sx = lx } +let zset s _ = s +let jump _ = true +let horizon _ = max_float + +let sinus_cosinus theta0 omega = + let sin0 = Float.sin theta0 in + let cos0 = Float.cos theta0 in + let fder _ _ y = + yd.{0} <- omega *. y.{1}; yd.{1} <- -.omega *. y.{0}; yd.{2} <- 1.0; yd in + let step { si; sx } _ = + let sx = if si then of_array [| sin0; cos0; 0.0 |] else sx in + of_array [| sx.{0}; sx.{1}; sx.{2} |], { sx; si = false } in + let state = { sx = of_array [| sin0; cos0; 0.0 |]; si = true } in + let reset _ _ = state in + HNode { + state; fder; fzer; fout; step; horizon; cset; cget; zset; zsize; reset; jump + } + +let errmsg_invalid = "Invalid arguments to model (needed: 2 floats)" +let errmsg_few = "Too few arguments to model (needed: 2 floats)" +let errmsg_many = "Too many arguments to model (needed: 2 floats)" +let sinus_cosinus = function + | [t0; om] -> + let t0, om = try float_of_string t0, float_of_string om + with Failure _ -> raise (Invalid_argument errmsg_invalid) in + sinus_cosinus t0 om + | [] | [_] -> raise (Invalid_argument errmsg_few) + | _ -> raise (Invalid_argument errmsg_many) diff --git a/exm/vdp.ml b/exm/vdp.ml index 2175bd4..568e9c5 100644 --- a/exm/vdp.ml +++ b/exm/vdp.ml @@ -1,59 +1,43 @@ open Hsim.Types -open Solvers +open Solvers.Zls -let of_array a = Bigarray.Array1.of_array Bigarray.Float64 Bigarray.c_layout a +let of_array a : carray = + Bigarray.Array1.of_array Bigarray.Float64 Bigarray.c_layout a -type state = { lx : Zls.carray; i : bool } +type state = { lx : carray; i : bool } let mu = 5.0 let x0 = 1.0 let y0 = 1.0 -let fder _ y yd = - yd.{0} <- y.{1}; yd.{1} <- (mu *. (1.0 -. (y.{0} *. y.{0})) *. y.{1}) -. y.{0} -let fzero _ _ _ = () -let fout _ y = of_array [| y.{0}; y.{1} |] +let zsize = 1 +let yd = cmake 2 +let zout = cmake 1 + +let fder _ _ y = + yd.{0} <- y.{1}; + yd.{1} <- (mu *. (1.0 -. (y.{0} *. y.{0})) *. y.{1}) -. y.{0}; + yd +let fzer _ _ _ = zout +let fout _ _ y = of_array [| y.{0}; y.{1} |] let step { i; lx } _ = let lx = if i then of_array [| x0; y0 |] else lx in of_array [| lx.{0}; lx.{1} |], { lx; i = false } - let cget s = s.lx -let cset s lx = { s with lx } +let cset s lx = { s with lx } let zset s _ = s -let yd = Zls.cmake 2 -let zout = Zls.cmake 1 -let zsize = 1 -let s_init = { lx = of_array [| x0; y0 |]; i = true } -let reset _ _ = s_init +let state = { lx = of_array [| x0; y0 |]; i = true } +let reset _ _ = state let jump _ = true +let horizon _ = max_float -let van_der_pol () = - HNode { state = s_init; - fder = (fun _ _ y -> fder 0.0 y yd; yd); - fzer = (fun _ _ y -> fzero 0.0 y zout; zout); - fout = (fun s _ y -> fout s y); - step; - horizon = (fun _ -> max_float); - cset; cget; zset; zsize; reset; jump } +let van_der_pol () = HNode { + state; fder; fzer; fout; step; reset; horizon; jump; cset; cget; zset; zsize +} + +let errmsg = "Too many arguments for the model (needed: 0)" +let van_der_pol = function + | [] -> van_der_pol () + | _ -> raise (Invalid_argument errmsg) -(* let van_der_pol_prop_record () = *) -(* let s_init = *) -(* { lx = of_array [| x0; y0 |]; i = true } in *) -(* { name = "van_der_pol_prop"; *) -(* s = s_init; *) -(* fder = (fun s a y -> fder 0.0 y yd; yd); *) -(* fzero = (fun s a y -> fzero 0.0 y zout; zout); *) -(* fout = (fun s { lx } y -> 1.0 /. (abs_float (s.lx.{0} -. lx.{0}))); *) -(* fstep = (fun s { lx } -> *) -(* let v, s = fstep s () in *) -(* s.lx.{0} -. lx.{0}, s); *) -(* horizon = (fun s -> max_float); *) -(* cset; cget; zset; csize; zsize; reset; jump } *) -(* let van_der_pol_with_assert () = *) -(* Fun_hybrid *) -(* { body = van_der_pol_record (); *) -(* assertions = *) -(* [Fun_hybrid { body = *) -(* van_der_pol_prop_record (); *) -(* assertions = [] }] } *) diff --git a/src/bin/main.ml b/src/bin/main.ml index 14a8235..8c8f920 100644 --- a/src/bin/main.ml +++ b/src/bin/main.ml @@ -11,13 +11,14 @@ let inplace = ref false let steps = ref 1 let model = ref "" -let gt0i v i = v := if i <= 0 then 1 else i -let gt0f v f = v := if f <= 0.0 then 1.0 else f -let hasmodel = ref false +let gt0i v i = v := if i <= 0 then 1 else i +let gt0f v f = v := if f <= 0.0 then 1.0 else f +let hasmodel = ref false +let modelargs = ref [] let set_model s = - if !hasmodel then raise (Arg.Bad "Too many arguments.") else - model := s; hasmodel := true + if !hasmodel then modelargs := s :: !modelargs + else model := s; hasmodel := true let opts = [ "-sample", Arg.Int (gt0i sample), "n \tSample count (default=10)"; @@ -37,10 +38,13 @@ let output = Output.print !sample let c = StatefulRK45.(if !inplace then InPlace.csolve else Functional.csolve) let z = StatefulZ.(Functional.zsolve) let s = Solver.solver_c c z +open Format let m = match !model with - | "ball" -> Ball.bouncing_ball () - | "vdp" -> Vdp.van_der_pol () - | _ -> raise (Arg.Bad (Format.sprintf "Unknown model: %s." !model)) + | "ball" -> Ball.bouncing_ball + | "vdp" -> Vdp.van_der_pol + | "sincos" -> Sincos.sinus_cosinus + | _ -> eprintf "Unknown model: %s\n" !model; exit 2 +let m = try m !modelargs with Invalid_argument s -> eprintf "%s\n" s; exit 2 let state = if !inplace then (module State.InPlaceSimState : State.SimState) else (module State.FunctionalSimState : State.SimState) diff --git a/src/lib/common/debug.ml b/src/lib/common/debug.ml index fdd37e7..e34b3cf 100644 --- a/src/lib/common/debug.ml +++ b/src/lib/common/debug.ml @@ -1,7 +1,4 @@ let debug = ref false -let fmt () = if !debug then Format.std_formatter - else Format.make_formatter (fun _ _ _ -> ()) (fun () -> ()) - -let print = Format.fprintf (fmt ()) "%s" +let print s = if !debug then Format.printf "%s\n" s else () diff --git a/src/lib/hsim/sim.ml b/src/lib/hsim/sim.ml index ba03a34..b0d7599 100644 --- a/src/lib/hsim/sim.ml +++ b/src/lib/hsim/sim.ml @@ -82,8 +82,6 @@ module LazySim (S : SimState) = (** Run the model on multiple inputs. *) let run_on_n model solver inputs use = ignore @@ List.fold_left (fun (DNode sim) i -> - Common.Debug.print - (Format.sprintf "New input: %.10e\t%.10e\n" i.start i.length); let state = match sim.step sim.state (Some i) with | None, s -> s | _ -> assert false in let rec loop (DNode s) = @@ -185,8 +183,6 @@ module GreedySim (S : SimState) = (** Run the model on multiple inputs. *) let run_on_n model solver inputs use = let o, _ = List.fold_left (fun (acc, DNode sim) i -> - Common.Debug.print - (Format.sprintf "new input: %.10e\t%.10e\n" i.start i.length); let o, state = sim.step sim.state i in o::acc, DNode { sim with state }) ([], run model solver) inputs in List.iter use (List.concat (List.rev o))