(**************************************************************************)
(*                                                                        *)
(*                                Zelus                                   *)
(*               A synchronous language for hybrid systems                *)
(*                       http://zelus.di.ens.fr                           *)
(*                                                                        *)
(*                    Marc Pouzet and Timothy Bourke                      *)
(*                                                                        *)
(*  Copyright 2012 - 2019. All rights reserved.                           *)
(*                                                                        *)
(*  This file is distributed under the terms of the CeCILL-C licence      *)
(*                                                                        *)
(*  Zelus is developed in the INRIA PARKAS team.                          *)
(*                                                                        *)
(**************************************************************************)

(* Type declarations and values that must be linked with *)
(* the generated code *)
type 'a continuous = { mutable pos: 'a; mutable der: 'a }

type ('a, 'b) zerocrossing = { mutable zin: 'a; mutable zout: 'b }

type 'a signal = 'a * bool
type zero = bool

(* a synchronous stream function with type 'a -D-> 'b *)
(* is represented by an OCaml value of type ('a, 'b) node *)
type ('s, 'a, 'b) node_rec = {
  alloc : unit -> 's; (* allocate the state *)
  step : 's -> 'a -> 'b; (* compute a step *)
  reset : 's -> unit; (* reset/inialize the state *)
}

type ('a, 'b) node =
    Node: ('s, 'a, 'b) node_rec -> ('a, 'b) node

(* the same with a method copy *)
type ('a, 'b) cnode =
    Cnode:
      { alloc : unit -> 's; (* allocate the state *)
        copy : 's -> 's -> unit; (* copy the source into the destination *)
        step : 's -> 'a -> 'b; (* compute a step *)
        reset : 's -> unit; (* reset/inialize the state *)
      } -> ('a, 'b) cnode

open Bigarray

type time    = float
type cvec    = (float, float64_elt, c_layout) Array1.t
type dvec    = (float, float64_elt, c_layout) Array1.t
type zinvec  = (int32, int32_elt,   c_layout) Array1.t
type zoutvec = (float, float64_elt, c_layout) Array1.t

(* The interface with the ODE solver *)
type cstate =
  { mutable dvec : dvec; (* the vector of derivatives *)
    mutable cvec : cvec; (* the vector of positions *)
    mutable zinvec : zinvec; (* the vector of boolean; true when the
                                solver has detected a zero-crossing *)
    mutable zoutvec : zoutvec; (* the corresponding vector that define
                                  zero-crossings *)
    mutable cindex : int; (* the position in the vector of positions *)
    mutable zindex : int; (* the position in the vector of zero-crossings *)
    mutable cend : int; (* the end of the vector of positions *)
    mutable zend : int; (* the end of the zero-crossing vector *)
    mutable cmax : int; (* the maximum size of the vector of positions *)
    mutable zmax : int; (* the maximum number of zero-crossings *)
    mutable horizon : float; (* the next horizon *)
    mutable major : bool; (* integration iff [major = false] *)
  }

(* The interface with the ODE solver (new Zélus version). *)
type cstate_new =
  { mutable dvec    : dvec;          (* Derivative vector.                    *)
    mutable cvec    : cvec;          (* Position vector.                      *)
    mutable zinvec  : zinvec;        (* Zero-crossing result vector.          *)
    mutable zoutvec : zoutvec;       (* Zero-crossing value vector.           *)
    mutable cindex  : int;           (* Position in position vector.          *)
    mutable zindex  : int;           (* Position in zero-crossing vector.     *)
    mutable cend    : int;           (* End of position vector.               *)
    mutable zend    : int;           (* End of zero-crossing vector.          *)
    mutable cmax    : int;           (* Maximum size of position vector.      *)
    mutable zmax    : int;           (* Maximum size of zero-crossing vector. *)
    mutable horizon : float;         (* Next horizon.                         *)
    mutable major   : bool;          (* Step mode: major <-> discrete.        *)
    mutable time    : float;         (* Simulation time.                      *)
  }

(* A hybrid node is a node that is parameterised by a continuous state *)
(* all instances points to this global parameter and read/write on it *)
type ('a, 'b) hnode = cstate -> (time * 'a, 'b) node

type 'b hsimu =
    Hsim:
      { alloc : unit -> 's;
        (* allocate the initial state *)
        maxsize : 's -> int * int;
        (* returns the max length of the *)
        (* cvector and zvector *)
        csize : 's -> int;
        (* returns the current length of the continuous state vector *)
        zsize : 's -> int;
        (* returns the current length of the zero-crossing vector *)
        step : 's -> cvec -> dvec -> zinvec -> time -> 'b;
        (* computes a step *)
        derivative : 's -> cvec -> dvec -> zinvec -> zoutvec -> time -> unit;
        (* computes the derivative *)
        crossings : 's -> cvec -> zinvec -> zoutvec -> time -> unit;
        (* computes the zero-crossings *)
        reset : 's -> unit;
        (* resets the state *)
        horizon : 's -> time;
        (* gives the next time horizon *)
      } -> 'b hsimu

(* a function with type 'a -C-> 'b, when given to a solver, is *)
(* represented by an OCaml value of type ('a, 'b) hsnode *)
type ('a, 'b) hsnode =
    Hnode:
      { state : 's;
        (* the discrete state *)
        zsize : int;
        (* the maximum size of the zero-crossing vector *)
        csize : int;
        (* the maximum size of the continuous state vector (positions) *)
        derivative : 's -> 'a -> time -> cvec -> dvec -> unit;
        (* computes the derivative *)
        crossing : 's -> 'a -> time -> cvec -> zoutvec -> unit;
        (* computes the derivative *)
        output : 's -> 'a -> cvec -> 'b;
        (* computes the zero-crossings *)
        setroots : 's -> 'a -> cvec -> zinvec -> unit;
        (* returns the zero-crossings *)
        majorstep : 's -> time -> cvec -> 'a -> 'b;
        (* computes a step *)
        reset : 's -> unit;
        (* resets the state *)
        horizon : 's -> time;
        (* gives the next time horizon *)
      } -> ('a, 'b) hsnode

 (* An idea suggested by Adrien Guatto, 26/04/2021 *)
 (* provide a means to the type for input/outputs of nodes *)
 (* express them with GADT to ensure type safety *)
 (* type ('a, 'b) node =
  |  Fun : { step : 'a -> 'b;
             typ_arg: 'a typ;
             typ_return: 'b typ
           } -> ('a, 'b) node
  | Node :
      { state : 's; step : 's -> 'a -> 'b * 's;
        typ_arg: 'a typ;
        typ_state : 's typ;
        typ_return: 'b typ } -> ('a, 'b) node

and 'a typ =
  | Tunit : unit typ
  | Tarrow : 'a typ * 'b typ -> ('a * 'b) typ
  | Tint : int -> int typ
  | Ttuple : 'a typlist -> 'a typ
  | Tnode : 'a typ * 'b typ -> ('a,'b) node typ

and 'a typlist =
  | Tnil : unit typlist
  | Tpair : 'a typ * 'b typlist -> ('a * 'b) typlist

Q1: do it for records? sum types ? How?
Q2: provide a "type_of" function for every introduced type?
*)