summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
Diffstat
-rw-r--r--bin/.gitignore4
-rw-r--r--data/clock.txt23
-rw-r--r--data/column.txt24
-rw-r--r--data/column2.txt22
-rw-r--r--data/corners.txt23
-rw-r--r--data/dripping-pan.txt23
-rw-r--r--data/flat.txt23
-rw-r--r--data/glass-half-empty.txt22
-rw-r--r--data/leidenfrost.txt23
-rw-r--r--data/pour-out.txt19
-rw-r--r--data/tanada.txt19
-rw-r--r--fluid.gpr20
-rw-r--r--license.txt25
-rw-r--r--obj/.gitignore4
-rw-r--r--src/complex_fixed_points.adb529
-rw-r--r--src/complex_fixed_points.ads117
-rw-r--r--src/fluid_simulator.adb197
17 files changed, 1117 insertions, 0 deletions
diff --git a/bin/.gitignore b/bin/.gitignore
new file mode 100644
index 0000000..ea7f887
--- /dev/null
+++ b/bin/.gitignore
@@ -0,0 +1,4 @@
+
+
+*
+!.gitignore
diff --git a/data/clock.txt b/data/clock.txt
new file mode 100644
index 0000000..bfaf5e6
--- /dev/null
+++ b/data/clock.txt
@@ -0,0 +1,23 @@
+ ########################
+ ##xxxxxxxxxxxxxxxxxxxx##
+ ##xxxxxxxxxxxxxxxxxxxx##
+ ##xxxxxxxxxxxxxxxxxxxx##
+ ##xxxxxxxxxxxxxxxxxxxx##
+ ##xxxxxxxxxxxxxxxxxxxx##
+ ###xxxxxxxxxxxxxxxxxx###
+ ###xxxxxxxxxxxxxxxx###
+ ###xxxxxxxxxxxx###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ ## ##
+ ## ##
+ ## ##
+ ## ##
+ ## ##
+ ########################
diff --git a/data/column.txt b/data/column.txt
new file mode 100644
index 0000000..67134a4
--- /dev/null
+++ b/data/column.txt
@@ -0,0 +1,24 @@
+
+
+
+
+
+
+### ............... ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ................. ###
+### ............... ###
+##########################################################################
+ ######################################################################## \ No newline at end of file
diff --git a/data/column2.txt b/data/column2.txt
new file mode 100644
index 0000000..4e871d6
--- /dev/null
+++ b/data/column2.txt
@@ -0,0 +1,22 @@
+
+
+
+
+### .......... ###
+### ........... ###
+### ........... ###
+### ........... ###
+### ........... ###
+### ........... ###
+### ...........
+### ...........
+### ...........
+### ........... ###
+### ........... ###
+### ........... ###
+### ........... ###
+### ........... ###
+### ........... # ###
+### .......... ### ###
+##############################################################################
+ ############################################################################
diff --git a/data/corners.txt b/data/corners.txt
new file mode 100644
index 0000000..45e5dcc
--- /dev/null
+++ b/data/corners.txt
@@ -0,0 +1,23 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+# #
+# #
+# #
+#@ @#
+###############################################################################
diff --git a/data/dripping-pan.txt b/data/dripping-pan.txt
new file mode 100644
index 0000000..2ed702b
--- /dev/null
+++ b/data/dripping-pan.txt
@@ -0,0 +1,23 @@
+
+
+
+
+
+
+
+
+# #
+# #
+# @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ #
+#@ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @ @#
+##################################### #######################################
+
+
+
+
+
+
+
+
+
+
diff --git a/data/flat.txt b/data/flat.txt
new file mode 100644
index 0000000..ad084da
--- /dev/null
+++ b/data/flat.txt
@@ -0,0 +1,23 @@
+
+
+
+
+
+
+
+
+# #
+# #
+# #
+#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#
+###############################################################################
+
+
+
+
+
+
+
+
+
+
diff --git a/data/glass-half-empty.txt b/data/glass-half-empty.txt
new file mode 100644
index 0000000..a329915
--- /dev/null
+++ b/data/glass-half-empty.txt
@@ -0,0 +1,22 @@
+
+
+
+ ###.......................###
+ ###.......................###
+ ###.......................###
+ ###.......................###
+ ###.......................###
+ ###.......................###
+ ###.......................###
+ ###.......................###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ ### ###
+ #############################
+ #############################
diff --git a/data/leidenfrost.txt b/data/leidenfrost.txt
new file mode 100644
index 0000000..c0a365e
--- /dev/null
+++ b/data/leidenfrost.txt
@@ -0,0 +1,23 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+# #
+# #
+# #
+# @ #
+###############################################################################
diff --git a/data/pour-out.txt b/data/pour-out.txt
new file mode 100644
index 0000000..7c9745a
--- /dev/null
+++ b/data/pour-out.txt
@@ -0,0 +1,19 @@
+
+ ### . . . . . . . . . . . ###
+ ###.......................###
+ ### . . . . . . . . . . . ###
+ ###.......................###
+ ### . . . . . . . . . . . ###
+ ###.......................###
+ ### . . . . . . . . . . . ###
+ ###.......................###
+ ### . . . . . . . . . . . ###
+ ###.......................###
+ ### . . . . . . . . . . . ###
+ ###.......................##
+ ### . . . . . . . . . . . ##
+ ###....................... ###
+ ################################# ## ##
+ ############################### ## ##
+ ## ##
+ ###########
diff --git a/data/tanada.txt b/data/tanada.txt
new file mode 100644
index 0000000..62b08f2
--- /dev/null
+++ b/data/tanada.txt
@@ -0,0 +1,19 @@
+
+
+ ##
+ ##. . . . . . . .
+ ##..................
+ ##. . . . . . . . . . ..
+ ##.......................
+ ##. . . . . . . . . . . ..
+ ##........................##
+ ##. . . . . . . . . . . ..##
+ ############################ ##
+ ## ##
+ ############ ##
+ ## ##
+ ############ ##
+ ## ##
+ ############ ##
+ ## ##
+ ############
diff --git a/fluid.gpr b/fluid.gpr
new file mode 100644
index 0000000..1a7a74e
--- /dev/null
+++ b/fluid.gpr
@@ -0,0 +1,20 @@
+
+project Fluid is
+
+ for Languages use ("Ada");
+
+ for Source_Dirs use ("src/**");
+ for Object_Dir use "obj";
+ for Exec_Dir use "bin";
+ for Main use ("fluid_simulator.adb");
+
+ package Builder is
+ for Executable ("fluid_simulator.adb") use "fluid";
+ end Builder;
+
+ package Compiler is
+ for Default_Switches("Ada") use ("-gnaty4aAbcefhiklM100nprt");
+ end Compiler;
+
+end Fluid;
+
diff --git a/license.txt b/license.txt
new file mode 100644
index 0000000..41c1fae
--- /dev/null
+++ b/license.txt
@@ -0,0 +1,25 @@
+
+
+SUNSET LICENSE
+Version 1.0, June 2017
+
+
+1. You may copy, modify, use, sell, or distribute this work, verbatim or
+modified, for any purpose.
+
+2. If you sell or distribute this work, whether verbatim or modified, you must
+include a copy of this license, and you must make the source code available for
+no extra charge.
+
+3. A modified version of this work must be clearly labeled as such.
+
+4. Derivative works must also be licensed under this license or a license of
+equivalent terms. As an exception, linking this work with another, whether
+statically or dynamically, does not impose any license requirements on the
+other work.
+
+5. If a minimum of 15 years have passed since the date of first publishing for
+a part of this work, then that part is placed into the public domain and you
+may do whatever you want with it, regardless of all other clauses.
+
+
diff --git a/obj/.gitignore b/obj/.gitignore
new file mode 100644
index 0000000..ea7f887
--- /dev/null
+++ b/obj/.gitignore
@@ -0,0 +1,4 @@
+
+
+*
+!.gitignore
diff --git a/src/complex_fixed_points.adb b/src/complex_fixed_points.adb
new file mode 100644
index 0000000..a0c7be1
--- /dev/null
+++ b/src/complex_fixed_points.adb
@@ -0,0 +1,529 @@
+
+package body Complex_Fixed_Points is
+
+ -- I have no idea what visibility conflict is requiring these two to be defined
+
+ function "*"
+ (Left, Right : in Real'Base)
+ return Real'Base is
+ begin
+ return Standard."*" (Left, Right);
+ end "*";
+
+ function "/"
+ (Left, Right : in Real'Base)
+ return Real'Base is
+ begin
+ return Standard."/" (Left, Right);
+ end "/";
+
+
+
+ function Re
+ (X : in Complex)
+ return Real'Base is
+ begin
+ return X.Re;
+ end Re;
+
+ function Im
+ (X : in Complex)
+ return Real'Base is
+ begin
+ return X.Im;
+ end Im;
+
+ function Im
+ (X : in Imaginary)
+ return Real'Base is
+ begin
+ return Real'Base (X);
+ end Im;
+
+
+
+ procedure Set_Re
+ (X : in out Complex;
+ Re : in Real'Base) is
+ begin
+ X.Re := Re;
+ end Set_Re;
+
+ procedure Set_Im
+ (X : in out Complex;
+ Im : in Real'Base) is
+ begin
+ X.Im := Im;
+ end Set_Im;
+
+ procedure Set_Im
+ (X : in out Imaginary;
+ Im : in Real'Base) is
+ begin
+ X := Imaginary (Im);
+ end Set_Im;
+
+
+
+ function Cartesian
+ (Re, Im : in Real'Base)
+ return Complex is
+ begin
+ return (Re => Re, Im => Im);
+ end Cartesian;
+
+ function Cartesian
+ (Re : in Real'Base)
+ return Complex is
+ begin
+ return (Re => Re, Im => 0.0);
+ end Cartesian;
+
+ function Cartesian
+ (Im : in Imaginary)
+ return Complex is
+ begin
+ return (Re => 0.0, Im => Real'Base (Im));
+ end Cartesian;
+
+
+
+ function Sqrt
+ (X : in Real'Base)
+ return Real'Base
+ is
+ Old_Result : Real'Base;
+ Result : Real'Base := 1.0;
+ Current_Square : Real'Base := 1.0;
+ Next_L : Real'Base := 3.0;
+ begin
+ -- If we don't check for this an overflow will happen due to taking things to the limit
+ if X = 0.0 then
+ return 0.0;
+ end if;
+
+ -- Find the smallest integer greater than Sqrt (X)
+ while Current_Square < X loop
+ Current_Square := Current_Square + Next_L;
+ Next_L := Next_L + 2.0;
+ Result := Result + 1.0;
+ end loop;
+
+ -- Lucky guess
+ if Current_Square = X then
+ return Result;
+ end if;
+
+ -- Babylonian aka Newton's Method
+ -- To be honest I'm kinda suspicious this may not necessarily terminate
+ loop
+ Old_Result := Result;
+ Result := (Old_Result + X / Old_Result) / 2.0;
+ exit when Real (Result) = Real (Old_Result);
+ end loop;
+ return Result;
+ end Sqrt;
+
+ function Modulus
+ (X : in Complex)
+ return Real'Base is
+ begin
+ return Sqrt (X.Re * X.Re + X.Im * X.Im);
+ end Modulus;
+
+
+
+ function "+"
+ (Right : in Complex)
+ return Complex is
+ begin
+ return Right;
+ end "+";
+
+ function "-"
+ (Right : in Complex)
+ return Complex is
+ begin
+ return (Re => -Right.Re, Im => -Right.Im);
+ end "-";
+
+ function Conjugate
+ (X : in Complex)
+ return Complex is
+ begin
+ return (Re => X.Re, Im => -X.Im);
+ end Conjugate;
+
+
+
+ function "+"
+ (Left, Right : in Complex)
+ return Complex is
+ begin
+ return (Re => Left.Re + Right.Re,
+ Im => Left.Im + Right.Im);
+ end "+";
+
+ function "-"
+ (Left, Right : in Complex)
+ return Complex is
+ begin
+ return (Re => Left.Re - Right.Re,
+ Im => Left.Im - Right.Im);
+ end "-";
+
+ function "*"
+ (Left, Right : in Complex)
+ return Complex is
+ begin
+ return (Re => Left.Re * Right.Re - Left.Im * Right.Im,
+ Im => Left.Re * Right.Im + Left.Im * Right.Re);
+ end "*";
+
+ function "/"
+ (Left, Right : in Complex)
+ return Complex
+ is
+ Denominator : Real'Base := Right.Re * Right.Re + Right.Im * Right.Im;
+ begin
+ return (Re => (Left.Re * Right.Re + Left.Im * Right.Im) / Denominator,
+ Im => (Left.Im * Right.Re - Left.Re * Right.Im) / Denominator);
+ end "/";
+
+
+
+ function "**"
+ (Left : in Complex;
+ Right : in Integer)
+ return Complex
+ is
+ Result : Complex;
+ Exponent : Integer := Right;
+ begin
+ if Exponent >= 1 then
+ Result := Left;
+ while Exponent > 1 loop
+ Result := Result * Left;
+ Exponent := Exponent - 1;
+ end loop;
+ else
+ Result := One;
+ while Exponent < 0 loop
+ Result := Result / Left;
+ Exponent := Exponent + 1;
+ end loop;
+ end if;
+ return Result;
+ end "**";
+
+
+
+ function "+"
+ (Right : in Imaginary)
+ return Imaginary is
+ begin
+ return Right;
+ end "+";
+
+ function "-"
+ (Right : in Imaginary)
+ return Imaginary is
+ begin
+ return Imaginary (-Real'Base (Right));
+ end "-";
+
+ function "abs"
+ (Right : in Imaginary)
+ return Real'Base is
+ begin
+ return abs (Real'Base (Right));
+ end "abs";
+
+
+
+ function "+"
+ (Left, Right : in Imaginary)
+ return Imaginary is
+ begin
+ return Imaginary (Real'Base (Left) + Real'Base (Right));
+ end "+";
+
+ function "-"
+ (Left, Right : in Imaginary)
+ return Imaginary is
+ begin
+ return Imaginary (Real'Base (Left) - Real'Base (Right));
+ end "-";
+
+ function "*"
+ (Left, Right : in Imaginary)
+ return Real'Base is
+ begin
+ return -(Real'Base (Left) * Real'Base (Right));
+ end "*";
+
+ function "/"
+ (Left, Right : in Imaginary)
+ return Real'Base is
+ begin
+ return Real'Base (Left) / Real'Base (Right);
+ end "/";
+
+
+
+ function "**"
+ (Left : in Imaginary;
+ Right : in Integer)
+ return Complex
+ is
+ Result : Complex;
+ Exponent : Integer := Right;
+ begin
+ if Exponent >= 1 then
+ Result := Cartesian (Left);
+ while Exponent > 1 loop
+ Result := Result * Left;
+ Exponent := Exponent - 1;
+ end loop;
+ else
+ Result := One;
+ while Exponent < 0 loop
+ Result := Result / Left;
+ Exponent := Exponent + 1;
+ end loop;
+ end if;
+ return Result;
+ end "**";
+
+
+
+ function "<"
+ (Left, Right : in Imaginary)
+ return Boolean is
+ begin
+ return Real'Base (Left) < Real'Base (Right);
+ end "<";
+
+ function "<="
+ (Left, Right : in Imaginary)
+ return Boolean is
+ begin
+ return Real'Base (Left) <= Real'Base (Right);
+ end "<=";
+
+ function ">"
+ (Left, Right : in Imaginary)
+ return Boolean is
+ begin
+ return Real'Base (Left) > Real'Base (Right);
+ end ">";
+
+ function ">="
+ (Left, Right : in Imaginary)
+ return Boolean is
+ begin
+ return Real'Base (Left) >= Real'Base (Right);
+ end ">=";
+
+
+
+ function "+"
+ (Left : in Complex;
+ Right : in Real'Base)
+ return Complex is
+ begin
+ return (Re => Left.Re + Right, Im => Left.Im);
+ end "+";
+
+ function "+"
+ (Left : in Real'Base;
+ Right : in Complex)
+ return Complex is
+ begin
+ return Right + Left;
+ end "+";
+
+ function "-"
+ (Left : in Complex;
+ Right : in Real'Base)
+ return Complex is
+ begin
+ return (Re => Left.Re - Right, Im => Left.Im);
+ end "-";
+
+ function "-"
+ (Left : in Real'Base;
+ Right : in Complex)
+ return Complex is
+ begin
+ return (Re => Left - Right.Re, Im => -Right.Im);
+ end "-";
+
+ function "*"
+ (Left : in Complex;
+ Right : in Real'Base)
+ return Complex is
+ begin
+ return (Re => Left.Re * Right, Im => Left.Im * Right);
+ end "*";
+
+ function "*"
+ (Left : in Real'Base;
+ Right : in Complex)
+ return Complex is
+ begin
+ return Right * Left;
+ end "*";
+
+ function "/"
+ (Left : in Complex;
+ Right : in Real'Base)
+ return Complex is
+ begin
+ return (Re => Left.Re / Right, Im => Left.Im / Right);
+ end "/";
+
+ function "/"
+ (Left : in Real'Base;
+ Right : in Complex)
+ return Complex is
+ begin
+ return Cartesian (Left) / Right;
+ end "/";
+
+
+
+ function "+"
+ (Left : in Complex;
+ Right : in Imaginary)
+ return Complex is
+ begin
+ return (Re => Left.Re, Im => Left.Im + Real'Base (Right));
+ end "+";
+
+ function "+"
+ (Left : in Imaginary;
+ Right : in Complex)
+ return Complex is
+ begin
+ return Right + Left;
+ end "+";
+
+ function "-"
+ (Left : in Complex;
+ Right : in Imaginary)
+ return Complex is
+ begin
+ return (Re => Left.Re, Im => Left.Im - Real'Base (Right));
+ end "-";
+
+ function "-"
+ (Left : in Imaginary;
+ Right : in Complex)
+ return Complex is
+ begin
+ return (Re => -Right.Re, Im => Real'Base (Left) - Right.Im);
+ end "-";
+
+ function "*"
+ (Left : in Complex;
+ Right : in Imaginary)
+ return Complex is
+ begin
+ return (Re => -Left.Im * Real'Base (Right),
+ Im => Left.Re * Real'Base (Right));
+ end "*";
+
+ function "*"
+ (Left : in Imaginary;
+ Right : in Complex)
+ return Complex is
+ begin
+ return Right * Left;
+ end "*";
+
+ function "/"
+ (Left : in Complex;
+ Right : in Imaginary)
+ return Complex is
+ begin
+ return (Re => Left.Im / Real'Base (Right),
+ Im => -Left.Re / Real'Base (Right));
+ end "/";
+
+ function "/"
+ (Left : in Imaginary;
+ Right : in Complex)
+ return Complex is
+ begin
+ return Cartesian (Left) / Right;
+ end "/";
+
+
+
+ function "+"
+ (Left : in Imaginary;
+ Right : in Real'Base)
+ return Complex is
+ begin
+ return (Re => Right, Im => Real'Base (Left));
+ end "+";
+
+ function "+"
+ (Left : in Real'Base;
+ Right : in Imaginary)
+ return Complex is
+ begin
+ return (Re => Left, Im => Real'Base (Right));
+ end "+";
+
+ function "-"
+ (Left : in Imaginary;
+ Right : in Real'Base)
+ return Complex is
+ begin
+ return (Re => -Right, Im => Real'Base (Left));
+ end "-";
+
+ function "-"
+ (Left : in Real'Base;
+ Right : in Imaginary)
+ return Complex is
+ begin
+ return (Re => Left, Im => -Real'Base (Right));
+ end "-";
+
+ function "*"
+ (Left : in Imaginary;
+ Right : in Real'Base)
+ return Imaginary is
+ begin
+ return Imaginary (Standard."*" (Right, Real'Base (Left)));
+ end "*";
+
+ function "*"
+ (Left : in Real'Base;
+ Right : in Imaginary)
+ return Imaginary is
+ begin
+ return Imaginary (Standard."*" (Left, Real'Base (Right)));
+ end "*";
+
+ function "/"
+ (Left : in Imaginary;
+ Right : in Real'Base)
+ return Imaginary is
+ begin
+ return Imaginary (Standard."/" (Real'Base (Left), Right));
+ end "/";
+
+ function "/"
+ (Left : in Real'Base;
+ Right : in Imaginary)
+ return Imaginary is
+ begin
+ return -Imaginary (Standard."/" (Left, Real'Base (Right)));
+ end "/";
+
+end Complex_Fixed_Points;
+
diff --git a/src/complex_fixed_points.ads b/src/complex_fixed_points.ads
new file mode 100644
index 0000000..97667b1
--- /dev/null
+++ b/src/complex_fixed_points.ads
@@ -0,0 +1,117 @@
+
+generic
+ type Real is delta <> digits <>;
+package Complex_Fixed_Points is
+
+ pragma Pure (Complex_Fixed_Points);
+
+ type Complex is private;
+ type Imaginary is private;
+ pragma Preelaborable_Initialization (Imaginary);
+
+ i : constant Imaginary;
+ j : constant Imaginary;
+
+ Zero : constant Complex;
+ One : constant Complex;
+
+ function Re (X : in Complex) return Real'Base;
+ function Im (X : in Complex) return Real'Base;
+ function Im (X : in Imaginary) return Real'Base;
+
+ procedure Set_Re
+ (X : in out Complex;
+ Re : in Real'Base);
+ procedure Set_Im
+ (X : in out Complex;
+ Im : in Real'Base);
+ procedure Set_Im
+ (X : in out Imaginary;
+ Im : in Real'Base);
+
+ function Cartesian (Re, Im : in Real'Base) return Complex;
+ function Cartesian (Re : in Real'Base) return Complex;
+ function Cartesian (Im : in Imaginary) return Complex;
+
+ function Modulus (X : in Complex) return Real'Base;
+ function "abs" (Right : in Complex) return Real'Base renames Modulus;
+
+ -- Like hell am I writing fixed point trig functions right now.
+
+ -- function Argument (X : in Complex) return Real'Base;
+ -- function Argument (X : in Complex; Cycle : in Real'Base) return Real'Base;
+
+ -- function Polar (Modulus, Argument : in Real'Base) return Complex;
+ -- function Polar (Modulus, Argument, Cycle : in Real'Base) return Complex;
+
+ function "+" (Right : in Complex) return Complex;
+ function "-" (Right : in Complex) return Complex;
+ function Conjugate (X : in Complex) return Complex;
+
+ function "+" (Left, Right : in Complex) return Complex;
+ function "-" (Left, Right : in Complex) return Complex;
+ function "*" (Left, Right : in Complex) return Complex;
+ function "/" (Left, Right : in Complex) return Complex;
+
+ function "**" (Left : in Complex; Right : in Integer) return Complex;
+
+ function "+" (Right : in Imaginary) return Imaginary;
+ function "-" (Right : in Imaginary) return Imaginary;
+ function Conjugate (X : in Imaginary) return Imaginary renames "-";
+ function "abs" (Right : in Imaginary) return Real'Base;
+
+ function "+" (Left, Right : in Imaginary) return Imaginary;
+ function "-" (Left, Right : in Imaginary) return Imaginary;
+ function "*" (Left, Right : in Imaginary) return Real'Base;
+ function "/" (Left, Right : in Imaginary) return Real'Base;
+
+ function "**" (Left : in Imaginary; Right : in Integer) return Complex;
+
+ function "<" (Left, Right : in Imaginary) return Boolean;
+ function "<=" (Left, Right : in Imaginary) return Boolean;
+ function ">" (Left, Right : in Imaginary) return Boolean;
+ function ">=" (Left, Right : in Imaginary) return Boolean;
+
+ function "+" (Left : in Complex; Right : in Real'Base) return Complex;
+ function "+" (Left : in Real'Base; Right : in Complex) return Complex;
+ function "-" (Left : in Complex; Right : in Real'Base) return Complex;
+ function "-" (Left : in Real'Base; Right : in Complex) return Complex;
+ function "*" (Left : in Complex; Right : in Real'Base) return Complex;
+ function "*" (Left : in Real'Base; Right : in Complex) return Complex;
+ function "/" (Left : in Complex; Right : in Real'Base) return Complex;
+ function "/" (Left : in Real'Base; Right : in Complex) return Complex;
+
+ function "+" (Left : in Complex; Right : in Imaginary) return Complex;
+ function "+" (Left : in Imaginary; Right : in Complex) return Complex;
+ function "-" (Left : in Complex; Right : in Imaginary) return Complex;
+ function "-" (Left : in Imaginary; Right : in Complex) return Complex;
+ function "*" (Left : in Complex; Right : in Imaginary) return Complex;
+ function "*" (Left : in Imaginary; Right : in Complex) return Complex;
+ function "/" (Left : in Complex; Right : in Imaginary) return Complex;
+ function "/" (Left : in Imaginary; Right : in Complex) return Complex;
+
+ function "+" (Left : in Imaginary; Right : in Real'Base) return Complex;
+ function "+" (Left : in Real'Base; Right : in Imaginary) return Complex;
+ function "-" (Left : in Imaginary; Right : in Real'Base) return Complex;
+ function "-" (Left : in Real'Base; Right : in Imaginary) return Complex;
+ function "*" (Left : in Imaginary; Right : in Real'Base) return Imaginary;
+ function "*" (Left : in Real'Base; Right : in Imaginary) return Imaginary;
+ function "/" (Left : in Imaginary; Right : in Real'Base) return Imaginary;
+ function "/" (Left : in Real'Base; Right : in Imaginary) return Imaginary;
+
+private
+
+ type Complex is record
+ Re, Im : Real'Base;
+ end record;
+
+ type Imaginary is new Real'Base;
+
+ i : constant Imaginary := 1.0;
+ j : constant Imaginary := 1.0;
+
+ Zero : constant Complex := (Re => 0.0, Im => 0.0);
+ One : constant Complex := (Re => 1.0, Im => 0.0);
+
+end Complex_Fixed_Points;
+
diff --git a/src/fluid_simulator.adb b/src/fluid_simulator.adb
new file mode 100644
index 0000000..01d3772
--- /dev/null
+++ b/src/fluid_simulator.adb
@@ -0,0 +1,197 @@
+
+with
+
+ Ada.Numerics.Generic_Complex_Types,
+ Ada.Containers.Vectors,
+ Ada.Characters.Latin_1,
+ Ada.Text_IO;
+
+procedure Fluid_Simulator is
+
+ package Latin renames Ada.Characters.Latin_1;
+ package IO renames Ada.Text_IO;
+
+ procedure Clear_Screen is
+ begin
+ -- ANSI control sequence Erase in Display
+ -- Variant to clear entire screen
+ IO.Put (Latin.ESC & "[2J");
+ end Clear_Screen;
+
+ procedure Reset_Cursor is
+ begin
+ -- ANSI control sequence Cursor Position
+ -- Parameters to move cursor to top left corner
+ IO.Put (Latin.ESC & "[1;1H");
+ end Reset_Cursor;
+
+
+
+ type Quantity is digits 18;
+
+ package Fixed is new Ada.Numerics.Generic_Complex_Types (Real => Quantity);
+ use type Fixed.Complex;
+
+ type Particle is record
+ Place : Fixed.Complex;
+ Solid : Boolean;
+ Density : Quantity := 0.0;
+ Acceleration : Fixed.Complex := Fixed.Compose_From_Cartesian (0.0, 0.0);
+ Velocity : Fixed.Complex := Fixed.Compose_From_Cartesian (0.0, 0.0);
+ end record;
+
+ function Create
+ (X, Y : Quantity;
+ Solid : in Boolean)
+ return Particle is
+ begin
+ return (Place => Fixed.Compose_From_Cartesian (X, Y), Solid => Solid, others => <>);
+ end Create;
+
+
+
+ -- Constant properties of particles
+ Particle_Radius : constant Quantity := 2.0;
+ Particle_Mass : constant Quantity := 1.0;
+
+ -- Constant used in calculating fluid interaction forces
+ P0 : constant Quantity := 1.5;
+
+ -- Other constant force factors
+ Gravity_Factor : constant Fixed.Complex := Fixed.Compose_From_Cartesian (0.0, 1.0);
+ Pressure_Factor : constant Quantity := 4.0;
+ Viscosity_Factor : constant Quantity := 8.0;
+
+
+
+ package Particle_Vectors is new Ada.Containers.Vectors
+ (Index_Type => Positive,
+ Element_Type => Particle);
+
+ Particles : Particle_Vectors.Vector := Particle_Vectors.Empty_Vector;
+
+
+
+ function Marching_Squares
+ (Input : in Particle_Vectors.Vector)
+ return String
+ is
+ Liquid_Chars : String (1 .. 16) := " .,_`/[/']\\-/\#";
+ Grid : array (0 .. 80, 0 .. 25) of Boolean := (others => (others => False));
+ Output : String (1 .. 2025);
+ X, Y : Integer;
+ begin
+ for P of Input loop
+ X := Integer (Fixed.Re (P.Place));
+ Y := Integer (Fixed.Im (P.Place));
+ if X >= 0 and X <= 80 and Y >= 0 and Y <= 25 then
+ Grid (X, Y) := True;
+ end if;
+ end loop;
+ for J in Integer range 1 .. 25 loop
+ for I in Integer range 1 .. 80 loop
+ Output ((J - 1) * 81 + I) := Liquid_Chars
+ (Boolean'Pos (Grid (I - 1, J)) + 1 +
+ Boolean'Pos (Grid (I, J)) * 2 +
+ Boolean'Pos (Grid (I, J - 1)) * 4 +
+ Boolean'Pos (Grid (I - 1, J - 1)) * 8);
+ end loop;
+ Output (J * 81) := Latin.LF;
+ end loop;
+ return Output;
+ end Marching_Squares;
+
+begin
+
+ declare
+ Input : Character;
+ X, Y : Quantity := 1.0;
+ begin
+ while not IO.End_Of_File (IO.Standard_Input) loop
+ IO.Get_Immediate (Input);
+ if Input = Latin.LF then
+ X := 1.0;
+ Y := Y + 1.0;
+ else
+ if Input > Latin.Space and Input < Latin.DEL then
+ Particles.Append (Create (X, Y, (Input = '#')));
+ end if;
+ X := X + 1.0;
+ end if;
+ end loop;
+ end;
+
+ loop
+ -- Calculate density
+ for P of Particles loop
+ if P.Solid then
+ P.Density := 9.0;
+ else
+ P.Density := 0.0;
+ end if;
+ for Q of Particles loop
+ declare
+ Rij : Quantity := Fixed.Modulus (P.Place - Q.Place);
+ W : Quantity := (Rij / Particle_Radius - 1.0) ** 2;
+ begin
+ if Rij < Particle_Radius then
+ P.Density := P.Density + Particle_Mass * W;
+ end if;
+ end;
+ end loop;
+ end loop;
+
+ -- Calculate interaction forces
+ for P of Particles loop
+ P.Acceleration := Gravity_Factor;
+ for Q of Particles loop
+ declare
+ Displacement : Fixed.Complex := P.Place - Q.Place;
+ Rij : Quantity := Fixed.Modulus (Displacement);
+ begin
+ if Rij < Particle_Radius then
+ declare
+ Pressure : Fixed.Complex :=
+ (P.Density + Q.Density - 2.0 * P0) * Pressure_Factor * Displacement;
+ Viscosity : Fixed.Complex :=
+ (P.Velocity - Q.Velocity) * Viscosity_Factor;
+ begin
+ P.Acceleration := P.Acceleration +
+ Fixed.Compose_From_Cartesian (1.0 - Rij / Particle_Radius) /
+ P.Density * (Pressure - Viscosity);
+ end;
+ end if;
+ end;
+ end loop;
+ end loop;
+
+ -- Render using marching squares
+ Clear_Screen;
+ Reset_Cursor;
+ IO.Put (Marching_Squares (Particles));
+
+ -- Update acceleration, velocity, and position
+ for P of Particles loop
+ if not P.Solid then
+ P.Velocity := P.Velocity + P.Acceleration / 10.0;
+ P.Place := P.Place + P.Velocity;
+ end if;
+ end loop;
+
+ -- Cull particles going out of bounds before any overflows happen
+ for C in reverse Particles.First_Index .. Particles.Last_Index loop
+ if Fixed.Re (Particles (C).Place) < -50.0 or
+ Fixed.Re (Particles (C).Place) > 130.0 or
+ Fixed.Im (Particles (C).Place) < -50.0 or
+ Fixed.Im (Particles (C).Place) > 75.0
+ then
+ Particles.Delete (C);
+ end if;
+ end loop;
+
+ -- Sleep for a bit
+ delay 0.012321;
+ end loop;
+
+end Fluid_Simulator;
+
generated by cgit (git 2.30.0) at 2024-12-22 08:34:49 +0000