From 2b0197dda7f18a1bd0d36da325b43ebd43e46b61 Mon Sep 17 00:00:00 2001 From: Tristan Gingold Date: Wed, 20 Jul 2022 19:44:39 +0200 Subject: grt: add analog_solver (work in progress) --- src/grt/grt-analog_solver.adb | 84 +++++++++++++++++++++++++++++++++++++++++++ src/grt/grt-analog_solver.ads | 69 +++++++++++++++++++++++++++++++++++ src/grt/grt-processes.adb | 43 +++++++++++++++++----- src/grt/grt-processes.ads | 10 ++++++ 4 files changed, 197 insertions(+), 9 deletions(-) create mode 100644 src/grt/grt-analog_solver.adb create mode 100644 src/grt/grt-analog_solver.ads diff --git a/src/grt/grt-analog_solver.adb b/src/grt/grt-analog_solver.adb new file mode 100644 index 000000000..c4d2c3b39 --- /dev/null +++ b/src/grt/grt-analog_solver.adb @@ -0,0 +1,84 @@ +-- GHDL Run Time (GRT) - analog solver for sundials. +-- Copyright (C) 2022 Tristan Gingold +-- +-- This program is free software: you can redistribute it and/or modify +-- it under the terms of the GNU General Public License as published by +-- the Free Software Foundation, either version 2 of the License, or +-- (at your option) any later version. +-- +-- This program is distributed in the hope that it will be useful, +-- but WITHOUT ANY WARRANTY; without even the implied warranty of +-- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +-- GNU General Public License for more details. +-- +-- You should have received a copy of the GNU General Public License +-- along with this program. If not, see . +-- +-- As a special exception, if other files instantiate generics from this +-- unit, or you link this unit with other files to produce an executable, +-- this unit does not by itself cause the resulting executable to be +-- covered by the GNU General Public License. This exception does not +-- however invalidate any other reasons why the executable file might be +-- covered by the GNU Public License. + +with Grt.Errors; use Grt.Errors; + +package body Grt.Analog_Solver is + function Sundials_Init (Sz : Ghdl_I32) return Ghdl_I32; + pragma Import (C, Sundials_Init, "grt_sundials_init"); + + function Sundials_Start return Ghdl_I32; + pragma Import (C, Sundials_Start, "grt_sundials_start"); + + procedure Sundials_Solve (T0 : Ghdl_F64; Tn : in out Ghdl_F64; + Res : out Integer); + pragma Import (C, Sundials_Solve, "grt_sundials_solve"); + + procedure Init (Size : Ghdl_I32) + is + Res : Ghdl_I32; + begin + Res := Sundials_Init (Size); + if Res < 0 then + Internal_Error ("sundials initialization failure"); + end if; + end Init; + + function Sundials_Get_Yy_Vec return F64_C_Arr_Ptr; + pragma Import (C, Sundials_Get_Yy_Vec, "grt_sundials_get_yy_vec"); + + function Get_Init_Val_Ptr return F64_C_Arr_Ptr is + begin + return Sundials_Get_Yy_Vec; + end Get_Init_Val_Ptr; + + function Sundials_Get_Yp_Vec return F64_C_Arr_Ptr; + pragma Import (C, Sundials_Get_Yp_Vec, "grt_sundials_get_yp_vec"); + + function Get_Init_Der_Ptr return F64_C_Arr_Ptr is + begin + return Sundials_Get_Yp_Vec; + end Get_Init_Der_Ptr; + + procedure Start is + begin + if Sundials_Start /= 0 then + Internal_Error ("sundials start"); + end if; + end Start; + + procedure Set_Root_Size (Size : Ghdl_I32) is + begin + Internal_Error ("sundials set_root_size"); + end Set_Root_Size; + + procedure Solve (T : Ghdl_F64; Tn : in out Ghdl_F64; Res : out Integer) is + begin + Sundials_Solve (T, Tn, Res); + end Solve; + + procedure Finish is + begin + Internal_Error ("sundials finish"); + end Finish; +end Grt.Analog_Solver; diff --git a/src/grt/grt-analog_solver.ads b/src/grt/grt-analog_solver.ads new file mode 100644 index 000000000..84490f16f --- /dev/null +++ b/src/grt/grt-analog_solver.ads @@ -0,0 +1,69 @@ +-- GHDL Run Time (GRT) - analog solver. +-- Copyright (C) 2022 Tristan Gingold +-- +-- This program is free software: you can redistribute it and/or modify +-- it under the terms of the GNU General Public License as published by +-- the Free Software Foundation, either version 2 of the License, or +-- (at your option) any later version. +-- +-- This program is distributed in the hope that it will be useful, +-- but WITHOUT ANY WARRANTY; without even the implied warranty of +-- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +-- GNU General Public License for more details. +-- +-- You should have received a copy of the GNU General Public License +-- along with this program. If not, see . +-- +-- As a special exception, if other files instantiate generics from this +-- unit, or you link this unit with other files to produce an executable, +-- this unit does not by itself cause the resulting executable to be +-- covered by the GNU General Public License. This exception does not +-- however invalidate any other reasons why the executable file might be +-- covered by the GNU Public License. + +with Grt.Types; use Grt.Types; + +package Grt.Analog_Solver is + type F64_Array is array (Natural range <>) of Ghdl_F64; + subtype F64_Fat_Array is F64_Array (Natural); + + -- A pointer to an F64 array. + type F64_C_Arr_Ptr is access all F64_Fat_Array; + pragma Convention (C, F64_C_Arr_Ptr); + + -- Initialize the analog solver, SIZE is the number of scalar quantities. + procedure Init (Size : Ghdl_I32); + + -- Return the address of the initial values vector and derivative vector. + function Get_Init_Val_Ptr return F64_C_Arr_Ptr; + function Get_Init_Der_Ptr return F64_C_Arr_Ptr; + + -- Finish initialization. + procedure Start; + + procedure Residues (T : Ghdl_F64; + Y : F64_C_Arr_Ptr; + Yp : F64_C_Arr_Ptr; + Res : F64_C_Arr_Ptr); + pragma Import (C, Residues, "grt__analog_solver__residues"); + + procedure Set_Root_Size (Size : Ghdl_I32); + + procedure Roots (T : Ghdl_F64; + Y : F64_C_Arr_Ptr; + Yp: F64_C_Arr_Ptr; + Res : F64_C_Arr_Ptr); + pragma Import (C, Roots, "grt__analog_solver__roots"); + + procedure Set_Values (Y : F64_C_Arr_Ptr; + Yp: F64_C_Arr_Ptr); + pragma Import (C, Set_Values, "grt__analog_solver__set_values"); + + -- Return status: + -- 0: Ok, Tn reached + -- 1: Stopped due to zero crossing + -- 2: failure + procedure Solve (T : Ghdl_F64; Tn : in out Ghdl_F64; Res : out Integer); + + procedure Finish; +end Grt.Analog_Solver; diff --git a/src/grt/grt-processes.adb b/src/grt/grt-processes.adb index 9eba3b677..6133343d8 100644 --- a/src/grt/grt-processes.adb +++ b/src/grt/grt-processes.adb @@ -22,8 +22,6 @@ -- covered by the GNU Public License. with Grt.Table; with Ada.Unchecked_Deallocation; -with System.Storage_Elements; -- Work around GNAT bug. -pragma Unreferenced (System.Storage_Elements); with Grt.Disp; with Grt.Astdio; with Grt.Astdio.Vhdl; use Grt.Astdio.Vhdl; @@ -39,6 +37,7 @@ with Grt.Stats; with Grt.Threads; use Grt.Threads; pragma Elaborate_All (Grt.Table); with Grt.Stdio; +with Grt.Analog_Solver; package body Grt.Processes is Last_Time : constant Std_Time := Std_Time'Last; @@ -880,12 +879,20 @@ package body Grt.Processes is -- - The time of the next simulation cycle (which in this case is the -- first simulation cycle), Tn, is calculated according to the rules -- of step f of the simulation cycle, below. - Next_Time := Compute_Next_Time; - if Next_Time /= 0 then - if Has_Callbacks (Hooks.Cb_Last_Known_Delta) then - Call_Callbacks (Hooks.Cb_Last_Known_Delta); - Flush_Active_Chain; - Next_Time := Compute_Next_Time; + + -- LRM 1076.1-2007 + -- - The time of the next simulation cycle (which in this case is the + -- first simulation cycle), Tn, is set to 0.0 + if Flag_AMS then + Next_Time := 0; + else + Next_Time := Compute_Next_Time; + if Next_Time /= 0 then + if Has_Callbacks (Hooks.Cb_Last_Known_Delta) then + Call_Callbacks (Hooks.Cb_Last_Known_Delta); + Flush_Active_Chain; + Next_Time := Compute_Next_Time; + end if; end if; end if; @@ -896,12 +903,26 @@ package body Grt.Processes is -- Launch a simulation cycle. function Simulation_Cycle return Integer is + use Grt.Options; Tn : Std_Time; + Tn_AMS : Ghdl_F64; Status : Integer; begin -- LRM08 14.7.5.3 Simulation cycle (ex LRM93 12.6.4) -- A simulation cycle consists of the following steps: -- + + -- LRM 1076.1-2007 12.6.4 Simulation cycle + -- a) The analog solver is executed + if Flag_AMS and Next_Time > Current_Time then + Current_Time_AMS := Ghdl_F64 (Current_Time) * Time_Phys_To_Real; + Tn_AMS := Ghdl_F64 (Next_Time) * Time_Phys_To_Real; + Grt.Analog_Solver.Solve (Current_Time_AMS, Tn_AMS, Status); + if Status /= 0 then + Internal_Error ("simulation_cycle - analog_solver"); + end if; + end if; + -- a) The current time, Tc is set equal to Tn. Simulation is complete -- when Tn = TIME'HIGH and there are no active drivers or process -- resumptions at Tn. @@ -1013,7 +1034,11 @@ package body Grt.Processes is if Options.Flag_Stats then Stats.Start_Next_Time; end if; - Tn := Compute_Next_Time; + if Flag_AMS and Break_Flag then + Tn := Current_Time; + else + Tn := Compute_Next_Time; + end if; -- h) If the next simulation cycle will be a delta cycle, the remainder -- of the step is skipped. Otherwise the following actions occur diff --git a/src/grt/grt-processes.ads b/src/grt/grt-processes.ads index 80df93641..da820a8b0 100644 --- a/src/grt/grt-processes.ads +++ b/src/grt/grt-processes.ads @@ -70,6 +70,16 @@ package Grt.Processes is -- Number of non-delta cycles. Nbr_Cycles : Ghdl_I64; + -- If True, do AMS simulation + Flag_AMS : Boolean := False; + + -- The break flag of AMS simulation + -- + -- LRM 1076.1-2017 14.7.5.2 Initialization + -- At the beginning of initialization, [...] and the break flag is assumed + -- to be cleared. + Break_Flag : Boolean := False; + type Process_Type is private; -- type Process_Acc is access all Process_Type; -- cgit v1.2.3