diff --git a/exm/ball.ml b/exm/ball.ml index 72e804e..d78109d 100644 --- a/exm/ball.ml +++ b/exm/ball.ml @@ -29,7 +29,7 @@ let fder y yd = else begin yd.{0} <- 0.0; yd.{1} <- 0.0; yd.{2} <- 0.0; yd.{3} <- 0.0 end; yd let fzer y zo = zo.{0} <- -. y.{1}; zo -let fout _ _ y = of_array [| y.{1} |] +let fout _ _ y = of_array [| y.{1}; y.{0} |] let jump _ = true let horizon _ = max_float let cget s = s.lx @@ -38,7 +38,7 @@ let zset s zin = { s with zin } let step ({ zin; lx; _ } as s) zfalse = 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 } + of_array [| s.lx.{1}; s.lx.{0} |], { zin=zfalse; lx; i=false } let bouncing_ball () : (_, _, _, carray, carray, carray, zarray, carray) hrec = let yd = cmake csize in diff --git a/exm/dune b/exm/dune index 4eecb1e..aa40c69 100644 --- a/exm/dune +++ b/exm/dune @@ -1,3 +1,8 @@ +(env + (dev + (flags + (:standard -w -a)))) + (library (name examples) (libraries hsim solvers)) diff --git a/exm/zelus/ballz/ballz.ml b/exm/zelus/ballz/ballz.ml deleted file mode 100644 index 08df6f6..0000000 --- a/exm/zelus/ballz/ballz.ml +++ /dev/null @@ -1,131 +0,0 @@ -(* The Zelus compiler, version 2.2-dev - (2025-06-16-15:24) *) -open Common -open Ztypes -open Solvers - -let (+=) r v = r := !r + v - -let g = 9.81 - -let y0 = 50. - -let y'0 = 0. - -type ball = - { mutable major : bool; - mutable h : float; - mutable init : bool; - mutable z : (bool, float) zerocrossing; - mutable y' : float continuous; - mutable y : float continuous } - -let ball cstate = - - let ball_alloc _ = - cstate.cmax <- cstate.cmax + 2; - cstate.zmax <- cstate.zmax + 1; - { major = false; - h = 42.; - init = false; - z = { zin = false; zout = 1. }; - y' = { pos = 42.; der = 0. }; - y = { pos = 42.; der = 0. } - } in - - let ball_step self (_, ()) = - let cidx = cstate.cindex in let cpos = ref cidx in - let zidx = cstate.zindex in let zpos = ref zidx in - cstate.cindex <- cstate.cindex + 2; - cstate.zindex <- cstate.zindex + 1; - self.major <- cstate.major; - if cstate.major then - for i = cidx to 1 do Zls.set cstate.dvec i 0. done - else begin - self.y'.pos <- Zls.get cstate.cvec !cpos; cpos += 1; - self.y.pos <- Zls.get cstate.cvec !cpos; cpos += 1 - end; - let res0, res1 = - let encore = ref false in - if self.init then self.y'.pos <- y'0; - let last_y' = self.y'.pos in - if self.z.zin then begin - encore := true; self.y'.pos <- -0.8 *. last_y' - end; - self.h <- if !encore then 0. else infinity; - if self.init then self.y.pos <- y0; - self.init <- false; - self.z.zout <- -. self.y.pos; - self.y'.der <- -. g; - self.y.der <- self.y'.pos; - self.y.pos, self.y.der - in - cstate.horizon <- min cstate.horizon self.h; - cpos := cidx; - if cstate.major then begin - Zls.set cstate.cvec !cpos self.y'.pos; cpos += 1; - Zls.set cstate.cvec !cpos self.y.pos; cpos += 1; - self.z.zin <- false - end else begin - self.z.zin <- Zls.get_zin cstate.zinvec !zpos; zpos += 1; - zpos := zidx; - Zls.set cstate.zoutvec !zpos self.z.zout; zpos += 1; - Zls.set cstate.dvec !cpos self.y'.der; cpos += 1; - Zls.set cstate.dvec !cpos self.y.der; cpos += 1 - end; - Bigarray.(Array1.of_array Float64 c_layout [| res0; res1 |]) in - - let ball_reset self = self.init <- true in - - Node { alloc = ball_alloc; step = ball_step ; reset = ball_reset } - -type ('f , 'e , 'd , 'c , 'b , 'a) _main = - { mutable i_73 : 'f ; - mutable major_62 : 'e ; - mutable h_72 : 'd ; - mutable i_70 : 'c ; mutable h_68 : 'b ; mutable t_63 : 'a } - -let main (cstate_80:Ztypes.cstate) = - let Node { alloc = i_73_alloc; step = i_73_step ; reset = i_73_reset } = ball - cstate_80 in - let main_alloc _ = - cstate_80.cmax <- (+) cstate_80.cmax 1; - { major_62 = false ; - h_72 = 42. ; - i_70 = (false:bool) ; - h_68 = (42.:float) ; t_63 = { pos = 42.; der = 0. }; - i_73 = i_73_alloc () (* continuous *) } in - let main_step self ((time_61:float) , ()) = - ((let (cindex_81:int) = cstate_80.cindex in - let cpos_83 = ref (cindex_81:int) in - cstate_80.cindex <- (+) cstate_80.cindex 1 ; - self.major_62 <- cstate_80.major ; - (if cstate_80.major then - for i_1 = cindex_81 to 0 do Zls.set cstate_80.dvec i_1 0. done - else ((self.t_63.pos <- Zls.get cstate_80.cvec !cpos_83 ; - cpos_83 := (+) !cpos_83 1))) ; - (let (result_85) = - let h_71 = ref (infinity:float) in - (if self.i_70 then self.h_68 <- (+.) time_61 0.) ; - (let (z_69:bool) = (&&) self.major_62 ((>=) time_61 self.h_68) in - self.h_68 <- (if z_69 then (+.) self.h_68 0.01 else self.h_68) ; - h_71 := min !h_71 self.h_68 ; - self.h_72 <- !h_71 ; - self.i_70 <- false ; - self.t_63.der <- 1. ; - (let (y_64:float) = (i_73_step self.i_73 (time_61 , ())).{0} in - (begin match z_69 with - | true -> Printf.printf "%.10e\t%.10e\n" self.t_63.pos y_64 - | _ -> () end) ; - Bigarray.(Array1.create Float64 c_layout 0))) in - cstate_80.horizon <- min cstate_80.horizon self.h_72 ; - cpos_83 := cindex_81 ; - (if cstate_80.major then - (((Zls.set cstate_80.cvec !cpos_83 self.t_63.pos ; - cpos_83 := (+) !cpos_83 1))) - else (((Zls.set cstate_80.dvec !cpos_83 self.t_63.der ; - cpos_83 := (+) !cpos_83 1)))) ; result_85))) in - let main_reset self = - ((self.i_70 <- true ; self.t_63.pos <- 0. ; i_73_reset self.i_73 ): - unit) in - Node { alloc = main_alloc; step = main_step ; reset = main_reset } diff --git a/exm/zelus/ballz/ballz.ml.bak b/exm/zelus/ballz/ballz.ml.bak deleted file mode 100644 index d248fc0..0000000 --- a/exm/zelus/ballz/ballz.ml.bak +++ /dev/null @@ -1,137 +0,0 @@ -(* The Zelus compiler, version 2.2-dev - (2025-06-16-15:24) *) -open Common -open Ztypes -open Solvers - -let g = 9.81 - -let y0 = 50. - -let y'0 = 0. - -type ('g , 'f , 'e , 'd , 'c , 'b , 'a) _ball = - { mutable major_50 : 'g ; - mutable h_60 : 'f ; - mutable h_58 : 'e ; - mutable i_56 : 'd ; - mutable x_55 : 'c ; mutable y'_52 : 'b ; mutable y_51 : 'a } - -let ball (cstate_74:Ztypes.cstate) = - - let ball_alloc _ = - cstate_74.cmax <- (+) cstate_74.cmax 2 ; - cstate_74.zmax <- (+) cstate_74.zmax 1; - { major_50 = false ; - h_60 = 42. ; - h_58 = (42.:float) ; - i_56 = (false:bool) ; - x_55 = { zin = false; zout = 1. } ; - y'_52 = { pos = 42.; der = 0. } ; y_51 = { pos = 42.; der = 0. } } in - let ball_step self ((_time_49:float) , ((y0_48:float) , (y'0_47:float))) = - Printf.printf "STEP (%d)\n" cstate_74.cindex; - ((let (cindex_75:int) = cstate_74.cindex in - let cpos_77 = ref (cindex_75:int) in - let (zindex_76:int) = cstate_74.zindex in - let zpos_78 = ref (zindex_76:int) in - cstate_74.cindex <- (+) cstate_74.cindex 2 ; - cstate_74.zindex <- (+) cstate_74.zindex 1 ; - self.major_50 <- cstate_74.major ; - (if cstate_74.major then - for i_1 = cindex_75 to 1 do Zls.set cstate_74.dvec i_1 0. done - else ((self.y'_52.pos <- Zls.get cstate_74.cvec !cpos_77 ; - cpos_77 := (+) !cpos_77 1) ; - (self.y_51.pos <- Zls.get cstate_74.cvec !cpos_77 ; - cpos_77 := (+) !cpos_77 1))) ; - (let (result_79:float) = - let h_59 = ref (infinity:float) in - let encore_57 = ref (false:bool) in - (if self.i_56 then self.y'_52.pos <- y'0_47) ; - (let (l_54:float) = self.y'_52.pos in - (begin match self.x_55.zin with - | true -> - encore_57 := true ; - self.y'_52.pos <- ( *. ) (-0.8) l_54 | _ -> () end) - ; - self.h_58 <- (if !encore_57 then 0. else infinity) ; - h_59 := min !h_59 self.h_58 ; - self.h_60 <- !h_59 ; - (if self.i_56 then self.y_51.pos <- y0_48) ; - self.i_56 <- false ; - self.x_55.zout <- (~-.) self.y_51.pos ; - self.y'_52.der <- (~-.) g ; - self.y_51.der <- self.y'_52.pos ; self.y_51.pos) in - cstate_74.horizon <- min cstate_74.horizon self.h_60 ; - cpos_77 := cindex_75 ; - (if cstate_74.major then - (((Printf.printf "idx: %d\n" !cpos_77; - Zls.set cstate_74.cvec !cpos_77 self.y'_52.pos ; - cpos_77 := (+) !cpos_77 1) ; - (Zls.set cstate_74.cvec !cpos_77 self.y_51.pos ; - cpos_77 := (+) !cpos_77 1)) ; ((self.x_55.zin <- false))) - else (((self.x_55.zin <- Zls.get_zin cstate_74.zinvec !zpos_78 ; - zpos_78 := (+) !zpos_78 1)) ; - zpos_78 := zindex_76 ; - ((Zls.set cstate_74.zoutvec !zpos_78 self.x_55.zout ; - zpos_78 := (+) !zpos_78 1)) ; - ((Zls.set cstate_74.dvec !cpos_77 self.y'_52.der ; - cpos_77 := (+) !cpos_77 1) ; - (Zls.set cstate_74.dvec !cpos_77 self.y_51.der ; - cpos_77 := (+) !cpos_77 1)))) ; result_79)):float) in - let ball_reset self = - (self.i_56 <- true:unit) in - Node { alloc = ball_alloc; step = ball_step ; reset = ball_reset } -type ('f , 'e , 'd , 'c , 'b , 'a) _main = - { mutable main_ball : 'f ; - mutable main_major : 'e ; - mutable h_72 : 'd; mutable i_70 : 'c; mutable h : 'b; - mutable t_63 : 'a } - -let main cs = - let Node - { alloc = ball_alloc; - step = ball_step; - reset = ball_reset } = ball cs in - - let main_alloc _ = - cs.cmax <- cs.cmax + 1; - { main_major = false; - h_72 = 42.0; i_70 = false; h = 42.0; - t_63 = { pos = 42.; der = 0. }; - main_ball = ball_alloc () } in - - let main_step self (time, ()) = - let cindex = cs.cindex in - let cpos = ref cindex in - Printf.printf "main:cindex: %d\n" cs.cindex; - cs.cindex <- cs.cindex + 1; - self.main_major <- cs.major; - if cs.major then for i = cindex to 0 do Zls.set cs.dvec i 0. done - else begin self.t_63.pos <- Zls.get cs.cvec !cpos; cpos := !cpos + 1 end; - let result = - if self.i_70 then self.h <- time; - let z = self.main_major && (time >= self.h) in - if z then self.h <- self.h +. 0.01; - self.h_72 <- min infinity self.h; - self.i_70 <- false; - self.t_63.der <- 1.; - let y_64 = ball_step self.main_ball (time, (y0, y'0)) in - if z then begin - print_float self.t_63.pos; - print_string "\t"; - print_float y_64; - print_newline () - end; - Bigarray.(Array1.create Float64 c_layout 0) in - cs.horizon <- min cs.horizon self.h_72; - cpos := cindex; - if cs.major then Zls.set cs.cvec !cpos self.t_63.pos - else Zls.set cs.dvec !cpos self.t_63.der; - result in - - let main_reset self = - self.i_70 <- true; - self.t_63.pos <- 0.; - ball_reset self.main_ball in - - Node { alloc = main_alloc; step = main_step ; reset = main_reset } diff --git a/exm/zelus/ballz/ballz.ml.bak2 b/exm/zelus/ballz/ballz.ml.bak2 deleted file mode 100644 index 05842bb..0000000 --- a/exm/zelus/ballz/ballz.ml.bak2 +++ /dev/null @@ -1,135 +0,0 @@ -(* The Zelus compiler, version 2.2-dev - (2025-06-16-15:24) *) -open Common -open Ztypes -open Solvers - -let g = 9.81 - -let y0 = 50. - -let y'0 = 0. - -type ('g , 'f , 'e , 'd , 'c , 'b , 'a) _ball = - { mutable major_50 : 'g ; - mutable h_60 : 'f ; - mutable h_58 : 'e ; - mutable i_56 : 'd ; - mutable x_55 : 'c ; mutable y'_52 : 'b ; mutable y_51 : 'a } - -let ball (cstate_74:Ztypes.cstate) = - - let ball_alloc _ = - cstate_74.cmax <- (+) cstate_74.cmax 2 ; - cstate_74.zmax <- (+) cstate_74.zmax 1; - { major_50 = false ; - h_60 = 42. ; - h_58 = (42.:float) ; - i_56 = (false:bool) ; - x_55 = { zin = false; zout = 1. } ; - y'_52 = { pos = 42.; der = 0. } ; y_51 = { pos = 42.; der = 0. } } in - - let ball_step self ((_time_49:float) , ()) = - ((let (cindex_75:int) = cstate_74.cindex in - let cpos_77 = ref (cindex_75:int) in - let (zindex_76:int) = cstate_74.zindex in - let zpos_78 = ref (zindex_76:int) in - cstate_74.cindex <- (+) cstate_74.cindex 2 ; - cstate_74.zindex <- (+) cstate_74.zindex 1 ; - self.major_50 <- cstate_74.major ; - (if cstate_74.major then - for i_1 = cindex_75 to 1 do Zls.set cstate_74.dvec i_1 0. done - else ((self.y'_52.pos <- Zls.get cstate_74.cvec !cpos_77 ; - cpos_77 := (+) !cpos_77 1) ; - (self.y_51.pos <- Zls.get cstate_74.cvec !cpos_77 ; - cpos_77 := (+) !cpos_77 1))) ; - (let (result_79:float) = - let h_59 = ref (infinity:float) in - let encore_57 = ref (false:bool) in - (if self.i_56 then self.y'_52.pos <- y'0) ; - (let (l_54:float) = self.y'_52.pos in - begin match self.x_55.zin with - | true -> - encore_57 := true ; - self.y'_52.pos <- ( *. ) (-0.8) l_54 | _ -> () end; - self.h_58 <- (if !encore_57 then 0. else infinity) ; - h_59 := min !h_59 self.h_58 ; - self.h_60 <- !h_59 ; - (if self.i_56 then self.y_51.pos <- y0) ; - self.i_56 <- false ; - self.x_55.zout <- (~-.) self.y_51.pos ; - self.y'_52.der <- (~-.) g ; - self.y_51.der <- self.y'_52.pos ; self.y_51.pos) in - cstate_74.horizon <- min cstate_74.horizon self.h_60 ; - cpos_77 := cindex_75 ; - (if cstate_74.major then - (((Zls.set cstate_74.cvec !cpos_77 self.y'_52.pos ; - cpos_77 := (+) !cpos_77 1) ; - (Zls.set cstate_74.cvec !cpos_77 self.y_51.pos ; - cpos_77 := (+) !cpos_77 1)) ; ((self.x_55.zin <- false))) - else (((self.x_55.zin <- Zls.get_zin cstate_74.zinvec !zpos_78 ; - zpos_78 := (+) !zpos_78 1)) ; - zpos_78 := zindex_76 ; - ((Zls.set cstate_74.zoutvec !zpos_78 self.x_55.zout ; - zpos_78 := (+) !zpos_78 1)) ; - ((Zls.set cstate_74.dvec !cpos_77 self.y'_52.der ; - cpos_77 := (+) !cpos_77 1) ; - (Zls.set cstate_74.dvec !cpos_77 self.y_51.der ; - cpos_77 := (+) !cpos_77 1)))) ; - Bigarray.(Array1.of_array Float64 c_layout [| result_79 |])))) in - - let ball_reset self = - (self.i_56 <- true:unit) in - - Node { alloc = ball_alloc; step = ball_step ; reset = ball_reset } - -type ('f , 'e , 'd , 'c , 'b , 'a) _main = - { mutable i_73 : 'f ; - mutable major_62 : 'e ; - mutable h_72 : 'd ; - mutable i_70 : 'c ; mutable h_68 : 'b ; mutable t_63 : 'a } - -let main (cstate_80:Ztypes.cstate) = - let Node { alloc = i_73_alloc; step = i_73_step ; reset = i_73_reset } = ball - cstate_80 in - let main_alloc _ = - cstate_80.cmax <- (+) cstate_80.cmax 1; - { major_62 = false ; - h_72 = 42. ; - i_70 = (false:bool) ; - h_68 = (42.:float) ; t_63 = { pos = 42.; der = 0. }; - i_73 = i_73_alloc () (* continuous *) } in - let main_step self ((time_61:float) , ()) = - ((let (cindex_81:int) = cstate_80.cindex in - let cpos_83 = ref (cindex_81:int) in - cstate_80.cindex <- (+) cstate_80.cindex 1 ; - self.major_62 <- cstate_80.major ; - (if cstate_80.major then - for i_1 = cindex_81 to 0 do Zls.set cstate_80.dvec i_1 0. done - else ((self.t_63.pos <- Zls.get cstate_80.cvec !cpos_83 ; - cpos_83 := (+) !cpos_83 1))) ; - (let (result_85) = - let h_71 = ref (infinity:float) in - (if self.i_70 then self.h_68 <- (+.) time_61 0.) ; - (let (z_69:bool) = (&&) self.major_62 ((>=) time_61 self.h_68) in - self.h_68 <- (if z_69 then (+.) self.h_68 0.01 else self.h_68) ; - h_71 := min !h_71 self.h_68 ; - self.h_72 <- !h_71 ; - self.i_70 <- false ; - self.t_63.der <- 1. ; - (let (y_64:float) = (i_73_step self.i_73 (time_61 , ())).{0} in - (begin match z_69 with - | true -> Printf.printf "%.10e\t%.10e\n" self.t_63.pos y_64 - | _ -> () end) ; - Bigarray.(Array1.create Float64 c_layout 0))) in - cstate_80.horizon <- min cstate_80.horizon self.h_72 ; - cpos_83 := cindex_81 ; - (if cstate_80.major then - (((Zls.set cstate_80.cvec !cpos_83 self.t_63.pos ; - cpos_83 := (+) !cpos_83 1))) - else (((Zls.set cstate_80.dvec !cpos_83 self.t_63.der ; - cpos_83 := (+) !cpos_83 1)))) ; result_85))) in - let main_reset self = - ((self.i_70 <- true ; self.t_63.pos <- 0. ; i_73_reset self.i_73 ): - unit) in - Node { alloc = main_alloc; step = main_step ; reset = main_reset } diff --git a/exm/zelus/ballz/ballz.zci b/exm/zelus/ballz/ballz.zci deleted file mode 100644 index 08fa125..0000000 Binary files a/exm/zelus/ballz/ballz.zci and /dev/null differ diff --git a/exm/zelus/ballz/ballz.zls b/exm/zelus/ballz/ballz.zls index fbaef0a..5f7e95f 100644 --- a/exm/zelus/ballz/ballz.zls +++ b/exm/zelus/ballz/ballz.zls @@ -2,18 +2,21 @@ let g = 9.81 let y0 = 50.0 let y'0 = 0.0 -let hybrid ball (y0, y'0) = y where +let hybrid ball (y0, y'0) = (y, y', z) where rec der y = y' init y0 and der y' = -. g init y'0 reset z -> -0.8 *. (last y') and z = up(-. y) let hybrid main () = let der t = 1.0 init 0.0 in - let y = ball (y0, y'0) in - let z = period(0.01) in - present z -> ( + let rec der p = 1.0 init -0.01 reset s -> -0.01 + and s = up(p) in + let (y, y', z) = ball (y0, y'0) in + present z | s -> ( print_float t; print_string "\t"; print_float y; + print_string "\t"; + print_float y'; print_newline () ); () diff --git a/exm/zelus/ballz/dune b/exm/zelus/ballz/dune new file mode 100644 index 0000000..81087c1 --- /dev/null +++ b/exm/zelus/ballz/dune @@ -0,0 +1,6 @@ +(rule + (targets ballz.ml ballz.zci) + (deps + (:zl ballz.zls)) + (action + (run zeluc %{zl}))) diff --git a/exm/zelus/sincos/dune b/exm/zelus/sincos/dune new file mode 100644 index 0000000..ef360bb --- /dev/null +++ b/exm/zelus/sincos/dune @@ -0,0 +1,6 @@ +(rule + (targets sincosz.ml sincosz.zci) + (deps + (:zl sincosz.zls)) + (action + (run zeluc %{zl}))) diff --git a/exm/zelus/sincos/sincosz.ml b/exm/zelus/sincos/sincosz.ml deleted file mode 100644 index ebbd09d..0000000 --- a/exm/zelus/sincos/sincosz.ml +++ /dev/null @@ -1,45 +0,0 @@ -(* The Zelus compiler, version 2.2-dev - (2025-06-16-15:24) *) -open Common -open Ztypes -open Solvers - -type ('c , 'b , 'a) _f = - { mutable major_11 : 'c ; mutable sin_13 : 'b ; mutable cos_12 : 'a } - -let f (cstate_14:Ztypes.cstate) = - - let f_alloc _ = - cstate_14.cmax <- (+) cstate_14.cmax 2; - { major_11 = false ; - sin_13 = { pos = 42.; der = 0. } ; cos_12 = { pos = 42.; der = 0. } } in - let f_step self ((_time_10:float) , ()) = - ((let (cindex_15:int) = cstate_14.cindex in - let cpos_17 = ref (cindex_15:int) in - cstate_14.cindex <- (+) cstate_14.cindex 2 ; - self.major_11 <- cstate_14.major ; - (if cstate_14.major then - for i_1 = cindex_15 to 1 do Zls.set cstate_14.dvec i_1 0. done - else ((self.sin_13.pos <- Zls.get cstate_14.cvec !cpos_17 ; - cpos_17 := (+) !cpos_17 1) ; - (self.cos_12.pos <- Zls.get cstate_14.cvec !cpos_17 ; - cpos_17 := (+) !cpos_17 1))) ; - (let (result_19) = - self.cos_12.der <- (~-.) self.sin_13.pos ; - self.sin_13.der <- self.cos_12.pos ; - Bigarray.(Array1.of_array Float64 c_layout - [| self.sin_13.pos; self.cos_12.pos |]) in - cpos_17 := cindex_15 ; - (if cstate_14.major then - (((Zls.set cstate_14.cvec !cpos_17 self.sin_13.pos ; - cpos_17 := (+) !cpos_17 1) ; - (Zls.set cstate_14.cvec !cpos_17 self.cos_12.pos ; - cpos_17 := (+) !cpos_17 1))) - else (((Zls.set cstate_14.dvec !cpos_17 self.sin_13.der ; - cpos_17 := (+) !cpos_17 1) ; - (Zls.set cstate_14.dvec !cpos_17 self.cos_12.der ; - cpos_17 := (+) !cpos_17 1)))) ; result_19))) in - - let f_reset self = - ((self.sin_13.pos <- 0. ; self.cos_12.pos <- 1.):unit) in - Node { alloc = f_alloc; step = f_step ; reset = f_reset } diff --git a/exm/zelus/sincos/sincosz.zci b/exm/zelus/sincos/sincosz.zci deleted file mode 100644 index e531551..0000000 Binary files a/exm/zelus/sincos/sincosz.zci and /dev/null differ diff --git a/exm/ztypes.ml b/exm/ztypes.ml new file mode 100644 index 0000000..567f8ff --- /dev/null +++ b/exm/ztypes.ml @@ -0,0 +1,4 @@ +include Common +include Ztypes +include Solvers + diff --git a/src/bin/lift.ml b/src/bin/lift.ml index 6b4e36d..ca9f336 100644 --- a/src/bin/lift.ml +++ b/src/bin/lift.ml @@ -38,7 +38,6 @@ let lift f = (* the function that compute the zero-crossings *) let fzer { state; _ } input y = - Common.Debug.print "LIFT :: Calling [fzer]"; cstate.major <- false; cstate.zinvec <- no_roots_in; cstate.dvec <- ignore_der; diff --git a/src/bin/main.ml b/src/bin/main.ml index de91d4c..3d6c608 100644 --- a/src/bin/main.ml +++ b/src/bin/main.ml @@ -19,6 +19,7 @@ let maxstep = ref None let mintol = ref None let maxtol = ref None let no_print = ref false +let zelus = 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 @@ -44,6 +45,7 @@ let opts = [ "-mintol", Arg.String (opt mintol), "\tSet minimum solver tolerance"; "-maxtol", Arg.String (opt maxtol), "\tSet maximum solver tolerance"; "-no-print", Arg.Set no_print, "\tDo not print output values"; + "-zelus", Arg.Set zelus, "\tUse the output of the Zélus compiler"; ] let errmsg = "Usage: " ^ Sys.executable_name ^ " [OPTIONS] MODEL\nOptions are:" @@ -53,18 +55,30 @@ let () = try Arg.parse (Arg.align opts) set_model errmsg with _ -> exit 2 let args = List.rev !modelargs let () = ignore lift +let wrap_zelus (HNode m) = + let ret = Bigarray.(Array1.create Float64 c_layout 0) in + let fout s a y = ignore (m.fout s a y); ret in + let step s () = let _, s = m.step s () in ret, s in + HNode { m with fout; step } + let m = - try match !model with - | None -> Format.eprintf "Missing model\n"; exit 2 - | Some "ball" -> Ball.init args - | Some "vdp" -> Vdp.init args - | Some "sincos" -> Sincos.init args - | Some "sqrt" -> Sqrt.init args - | Some "sin1x" -> Sin1x.init args - | Some "sin1xd" -> Sin1x_der.init args - | Some "ballz" -> lift Ballz.ball - | Some "sincosz" -> lift Sincosz.f - | Some s -> Format.eprintf "Unknown model: %s\n" s; exit 2 + try + if !zelus then + match !model with + | None -> Format.eprintf "Missing model\n"; exit 2 + | Some "ballz" -> wrap_zelus (lift Ballz.main) + | Some "sincosz" -> wrap_zelus (lift Sincosz.f) + | Some s -> Format.eprintf "Unknown model: %s\n" s; exit 2 + else + match !model with + | None -> Format.eprintf "Missing model\n"; exit 2 + | Some "ball" -> Ball.init args + | Some "vdp" -> Vdp.init args + | Some "sincos" -> Sincos.init args + | Some "sqrt" -> Sqrt.init args + | Some "sin1x" -> Sin1x.init args + | Some "sin1xd" -> Sin1x_der.init args + | Some s -> Format.eprintf "Unknown model: %s\n" s; exit 2 with Invalid_argument s -> Format.eprintf "%s\n" s; exit 2 let st = if !inplace then (module State.InPlaceSimState : State.SimState) diff --git a/src/lib/hsim/sim.ml b/src/lib/hsim/sim.ml index ac8ed4f..d6dd628 100644 --- a/src/lib/hsim/sim.ml +++ b/src/lib/hsim/sim.ml @@ -7,12 +7,13 @@ module Sim (S : SimState) = struct include S - let step_discrete s step hor fder fzer cget csize zsize jump reset = - Common.Debug.print "SIMU :: DISCRETE :: start"; - let ms, ss = get_mstate s, get_sstate s in + let step_discrete s step hor fder fzer cget zset csize zsize jump reset = + let ms, ss, zin = get_mstate s, get_sstate s, get_zin s in + let ms = match zin with Some z -> zset ms z | None -> ms in let i, now, stop = get_input s, get_now s, get_stop s in let o, ms = step ms (i.u now) in let s = + let s = set_zin None s in let h = hor ms in if h <= 0.0 then set_mstate ms s else if now >= stop then set_idle s @@ -27,30 +28,25 @@ module Sim (S : SimState) = let mode, stop, now = Continuous, i.h, 0.0 in update ms ss (set_running ~mode ~input:i ~stop ~now s) end else set_running ~mode:Continuous s in - Common.Debug.print "SIMU :: DISCRETE :: end"; Utils.dot o, s - let step_continuous s step cset fout zset = - Common.Debug.print "SIMU :: CONTINUOUS :: start"; + let step_continuous s step cset fout hor = let ms, ss = get_mstate s, get_sstate s in let i, now, stop = get_input s, get_now s, get_stop s in - let (h, f, z), ss = step ss stop in + let stop = min stop (hor ms) in + let (h, f, z), ss = step ss (min stop (hor ms)) in let ms = cset ms (f h) in let fy t = f (now +. t) in let fms t = cset ms (fy t) in let fout t = fout ms (i.u (now +. t)) (fy t) in let s, c = match z with | None -> - let s, c = if h >= stop - then set_running ~mode:Discrete ~now:h s, Discontinuous - else set_running ~now: h s, Continuous in - update ms ss s, c - | Some z -> - let s = set_running ~mode:Discrete ~now:h s in - update (zset ms z) ss s, Discontinuous in + if h >= stop + then set_running ~mode:Discrete ~now:h s, Discontinuous + else set_running ~now: h s, Continuous + | Some _ -> set_running ~mode:Discrete ~now:h s, Discontinuous in let h = h -. now in - Common.Debug.print "SIMU :: CONTINUOUS :: end"; - { h; u=fout; c }, s, { h; c; u=fms } + { h; u=fout; c }, update ms ss (set_zin z s), { h; c; u=fms } (** Simulation of a model with any solver. *) let run @@ -59,11 +55,11 @@ module Sim (S : SimState) = : ('p * (('y, 'yder) ivp * ('y, 'zout) zc), 'a, 'b) sim = let state = get_init m.state s.state in let step_discrete st = - let o, s = step_discrete st m.step m.horizon m.fder m.fzer m.cget + let o, s = step_discrete st m.step m.horizon m.fder m.fzer m.cget m.zset m.csize m.zsize m.jump s.reset in Some o, s in let step_continuous st = - let o, s, _ = step_continuous st s.step m.cset m.fout m.zset in + let o, s, _ = step_continuous st s.step m.cset m.fout m.horizon in Some o, s in let step st = function @@ -95,8 +91,8 @@ module Sim (S : SimState) = let step_discrete (st, al) = let m=m.body in let o, st = - step_discrete st m.step m.horizon m.fder m.fzer m.cget m.csize m.zsize - m.jump s.reset in + step_discrete st m.step m.horizon m.fder m.fzer m.cget m.zset m.csize + m.zsize m.jump s.reset in let al = List.map (fun (DNode a) -> let _, state = a.step a.state @@ Some (Utils.dot @@ get_mstate st) in DNode { a with state }) al in @@ -104,7 +100,7 @@ module Sim (S : SimState) = let step_continuous (st, al) = let ({ h; _ } as o), st, u = - step_continuous st s.step m.body.cset m.body.fout m.body.zset in + step_continuous st s.step m.body.cset m.body.fout m.body.horizon in let al = List.map (fun (DNode a) -> (* Step assertions repeatedly until they reach the horizon. *) let rec step s = @@ -140,10 +136,10 @@ module Sim (S : SimState) = = let state = get_init m.state s.state in let step_discrete st = let o, st = step_discrete st m.step m.horizon m.fder m.fzer m.cget - m.csize m.zsize m.jump s.reset in + m.zset m.csize m.zsize m.jump s.reset in Some o, st in let step_continuous st = - let o, st, _ = step_continuous st s.step m.cset m.fout m.zset in + let o, st, _ = step_continuous st s.step m.cset m.fout m.horizon in o, st in let rec step st = function diff --git a/src/lib/hsim/state.ml b/src/lib/hsim/state.ml index 676e12a..29fd8b2 100644 --- a/src/lib/hsim/state.ml +++ b/src/lib/hsim/state.ml @@ -15,59 +15,65 @@ module type SimState = - Idle: waiting for input; - Running: currently integrating; in this case, we have access to the step mode, current input, timestamp and stop time. *) - type ('a, 'ms, 'ss) state + type ('a, 'ms, 'ss, 'zin) state (** Get the model state. *) - val get_mstate : ('a, 'ms, 'ss) state -> 'ms + val get_mstate : ('a, 'ms, 'ss, 'zin) state -> 'ms (** Get the solver state. *) - val get_sstate : ('a, 'ms, 'ss) state -> 'ss + val get_sstate : ('a, 'ms, 'ss, 'zin) state -> 'ss + + (** Get the last zero-crossing value. *) + val get_zin : ('a, 'ms, 'ss, 'zin) state -> 'zin option (** Get the current step mode. ⚠ Should only be called when running (see [is_running]). *) - val get_mode : ('a, 'ms, 'ss) state -> mode + val get_mode : ('a, 'ms, 'ss, 'zin) state -> mode (** Get the current input. ⚠ Should only be called when running (see [is_running]). *) - val get_input : ('a, 'ms, 'ss) state -> 'a value + val get_input : ('a, 'ms, 'ss, 'zin) state -> 'a value (** Get the current timestamp. ⚠ Should only be called when running (see [is_running]). *) - val get_now : ('a, 'ms, 'ss) state -> time + val get_now : ('a, 'ms, 'ss, 'zin) state -> time (** Get the current stop time. ⚠ Should only be called when running (see [is_running]). *) - val get_stop : ('a, 'ms, 'ss) state -> time + val get_stop : ('a, 'ms, 'ss, 'zin) state -> time (** Build an initial state. *) - val get_init : 'ms -> 'ss -> ('a, 'ms, 'ss) state + val get_init : 'ms -> 'ss -> ('a, 'ms, 'ss, 'zin) state (** Is the simulation running or idle ? *) - val is_running : ('a, 'ms, 'ss) state -> bool + val is_running : ('a, 'ms, 'ss, 'zin) state -> bool (** Update the model state. *) - val set_mstate : 'ms -> ('a, 'ms, 'ss) state -> ('a, 'ms, 'ss) state + val set_mstate : 'ms -> ('a, 'ms, 'ss, 'zin) state -> ('a, 'ms, 'ss, 'zin) state (** Update the solver state. *) - val set_sstate : 'ss -> ('a, 'ms, 'ss) state -> ('a, 'ms, 'ss) state + val set_sstate : 'ss -> ('a, 'ms, 'ss, 'zin) state -> ('a, 'ms, 'ss, 'zin) state + + (** Update the zero-crossing value. *) + val set_zin : 'zin option -> ('a, 'ms, 'ss, 'zin) state -> ('a, 'ms, 'ss, 'zin) state (** Update both the solver and model states. *) - val update : 'ms -> 'ss -> ('a, 'ms, 'ss) state -> ('a, 'ms, 'ss) state + val update : 'ms -> 'ss -> ('a, 'ms, 'ss, 'zin) state -> ('a, 'ms, 'ss, 'zin) state (** Update the status to running. *) val set_running : ?mode:mode -> ?input:'a value -> ?now:time -> ?stop:time -> - ('a, 'ms, 'ss) state -> ('a, 'ms, 'ss) state + ('a, 'ms, 'ss, 'zin) state -> ('a, 'ms, 'ss, 'zin) state (** Update the status to idle. *) - val set_idle : ('a, 'ms, 'ss) state -> ('a, 'ms, 'ss) state + val set_idle : ('a, 'ms, 'ss, 'zin) state -> ('a, 'ms, 'ss, 'zin) state end module type SimStateCopy = sig include SimState - val copy : ('a, 'ms, 'ss) state -> ('a, 'ms, 'ss) state + val copy : ('a, 'ms, 'ss, 'zin) state -> ('a, 'ms, 'ss, 'zin) state end module FunctionalSimState : SimState = @@ -88,15 +94,17 @@ module FunctionalSimState : SimState = (** Internal state of the simulation node: model state, solver state and current simulation status. *) - type ('a, 'ms, 'ss) state = + type ('a, 'ms, 'ss, 'zin) state = { status : 'a status; (** Current simulation status. *) mstate : 'ms; (** Model state. *) - sstate : 'ss } (** Solver state. *) + sstate : 'ss; (** Solver state. *) + zin : 'zin option; } (** Last zero-crossing vector *) exception Not_running let get_mstate state = state.mstate let get_sstate state = state.sstate + let get_zin state = state.zin let is_running state = match state.status with Running _ -> true | Idle -> false @@ -120,6 +128,7 @@ module FunctionalSimState : SimState = let set_mstate mstate state = { state with mstate } let set_sstate sstate state = { state with sstate } + let set_zin zin state = { state with zin } let update mstate sstate state = { state with mstate; sstate } @@ -132,7 +141,7 @@ module FunctionalSimState : SimState = let get_stop s = match s.status with Running r -> r.stop | Idle -> raise Not_running - let get_init mstate sstate = { status = Idle; mstate; sstate } + let get_init mstate sstate = { status = Idle; mstate; sstate; zin = None } end module InPlaceSimState : SimState = @@ -146,15 +155,17 @@ module InPlaceSimState : SimState = mutable stop : time; } -> 'a status - type ('a, 'ms, 'ss) state = + type ('a, 'ms, 'ss, 'zin) state = { mutable status : 'a status; mutable mstate : 'ms; - mutable sstate : 'ss } + mutable sstate : 'ss; + mutable zin : 'zin option } exception Not_running let get_mstate state = state.mstate let get_sstate state = state.sstate + let get_zin state = state.zin let is_running state = match state.status with Running _ -> true | Idle -> false @@ -179,6 +190,7 @@ module InPlaceSimState : SimState = let set_mstate mstate state = state.mstate <- mstate; state let set_sstate sstate state = state.sstate <- sstate; state + let set_zin zin state = state.zin <- zin; state let update mstate sstate state = state.mstate <- mstate; state.sstate <- sstate; state @@ -192,6 +204,6 @@ module InPlaceSimState : SimState = let get_stop s = match s.status with Running r -> r.stop | Idle -> raise Not_running - let get_init mstate sstate = { status = Idle; mstate; sstate } + let get_init mstate sstate = { status = Idle; mstate; sstate; zin=None } end diff --git a/src/lib/solvers/illinois.ml b/src/lib/solvers/illinois.ml index 3f5cdb6..7e280ce 100644 --- a/src/lib/solvers/illinois.ml +++ b/src/lib/solvers/illinois.ml @@ -154,7 +154,6 @@ let num_roots { f0; _ } = Zls.length f0 (* f0/t0 take the previous values of f1/t1, f1/t1 are refreshed by g *) let step ({ g; f0 = f0; f1 = f1; t1 = t1; _ } as s) t c = - Common.Debug.print "ZSOL :: Calling [step]"; (* swap f0 and f1; f0 takes the previous value of f1 *) s.f0 <- f1; s.t0 <- t1; diff --git a/src/lib/solvers/statefulRK45.ml b/src/lib/solvers/statefulRK45.ml index 997d0ee..d85658e 100644 --- a/src/lib/solvers/statefulRK45.ml +++ b/src/lib/solvers/statefulRK45.ml @@ -22,8 +22,6 @@ module Functional = { state; vec = init } in let step ({ state ; vec=v } as s) h = - Common.Debug.print "SOLVER STEP"; - Common.Debug.print_entry v; let y_nv = vec v in let h = step state h y_nv in let state = copy state in diff --git a/src/lib/solvers/statefulZ.ml b/src/lib/solvers/statefulZ.ml index f6dc60b..b642a08 100644 --- a/src/lib/solvers/statefulZ.ml +++ b/src/lib/solvers/statefulZ.ml @@ -17,8 +17,6 @@ module Functional = let reset { fzer; init; size } { vec; _ } = let fzer t cvec zout = let zout' = fzer t cvec in blit zout' zout in - Common.Debug.print "ZSolver Reset"; - Common.Debug.print_entry init; { state = initialize size fzer init; vec = if length vec = size then vec else zmake size } in