diff --git a/exm/ball.ml b/exm/ball.ml index c7bb991..542f883 100644 --- a/exm/ball.ml +++ b/exm/ball.ml @@ -12,7 +12,7 @@ let of_array = Bigarray.Array1.of_array Bigarray.Float64 Bigarray.c_layout type state = { zin : Zls.zarray; - l_x : Zls.carray; (* [h';h;x';x] *) + lx : Zls.carray; (* [h';h;x';x] *) i : bool } let g = -9.81 @@ -30,16 +30,16 @@ let bouncing_ball () = begin yd.{0} <- 0.0; yd.{1} <- 0.0; yd.{2} <- 0.0; yd.{3} <- 0.0 end in let fzer _ (y: Zls.carray) (zout: Zls.carray) = zout.{0} <- -. y.{1} in let fout _ (y: Zls.carray): Zls.carray = of_array [| y.{1} |] in - let step ({ zin; l_x; _ } as s) _ : Zls.carray * state = - let l_x = if zin.{0} = 1l then - of_array [| -. 0.8 *. l_x.{0}; 0.0; l_x.{2}; l_x.{3} |] else l_x in - of_array [| s.l_x.{1} |], { zin = zfalse; l_x; i = false } in - let cget s = s.l_x in - let cset s l_x = { s with l_x } in + let step ({ zin; lx; _ } as s) _ : Zls.carray * state = + let lx = if zin.{0} = 1l then + of_array [| -. 0.8 *. lx.{0}; 0.0; lx.{2}; lx.{3} |] else lx in + of_array [| s.lx.{1} |], { zin = zfalse; lx; i = false } in + let cget s = s.lx in + let cset s lx = { s with lx } in let zset s zin = { s with zin } in let yd = Zls.cmake 4 in let zout = Zls.cmake 1 in - let state = { zin = zfalse; l_x = of_array [|y'0;y0;x'0;x0|]; i = true } in + let state = { zin = zfalse; lx = of_array [|y'0;y0;x'0;x0|]; i = true } in let reset _ _ = state in let jump _ = true in let zsize = 1 in diff --git a/exm/vdp.ml b/exm/vdp.ml new file mode 100644 index 0000000..2175bd4 --- /dev/null +++ b/exm/vdp.ml @@ -0,0 +1,59 @@ + +open Hsim.Types +open Solvers + +let of_array a = Bigarray.Array1.of_array Bigarray.Float64 Bigarray.c_layout a + +type state = { lx : Zls.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 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 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 jump _ = true + +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_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 5a95924..14a8235 100644 --- a/src/bin/main.ml +++ b/src/bin/main.ml @@ -9,9 +9,15 @@ let stop = ref 30.0 let greedy = ref false 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 set_model s = + if !hasmodel then raise (Arg.Bad "Too many arguments.") else + model := s; hasmodel := true let opts = [ "-sample", Arg.Int (gt0i sample), "n \tSample count (default=10)"; @@ -22,16 +28,19 @@ let opts = [ "-steps", Arg.Int (gt0i steps), "n \tSplit into [n] steps (default=1)"; ] -let errmsg = "Usage: " ^ Sys.executable_name ^ " [OPTIONS]\nOptions are:" +let errmsg = "Usage: " ^ Sys.executable_name ^ " [OPTIONS] MODEL\nOptions are:" -let () = try Arg.parse (Arg.align opts) (fun _ -> ()) errmsg with _ -> exit 2 +let () = try Arg.parse (Arg.align opts) set_model errmsg with _ -> exit 2 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 -let m = Ball.bouncing_ball () +let m = match !model with + | "ball" -> Ball.bouncing_ball () + | "vdp" -> Vdp.van_der_pol () + | _ -> raise (Arg.Bad (Format.sprintf "Unknown model: %s." !model)) let state = if !inplace then (module State.InPlaceSimState : State.SimState) else (module State.FunctionalSimState : State.SimState)