WIP comparison tmp

Added
- apophis gravity and temperature examples
- read body data file
- simulation state file behavior and indices
- python plot routines

Started some work for comparison with other models, temporary prints and tests.
WIP on computation speed to be studied

.gitignore

@@ -9,6 +9,7 @@ target/
 
 # Python generics
 .env
+__pycache__/
 python/__pycache__/
 python/kalast/kalast.abi3.so
 python/kalast/__pycache__/__init__.cpython-311.pyc
@@ -17,6 +18,7 @@ python/kalast/__pycache__/__init__.cpython-311.pyc
 assets/perso/
 bin/releases/
 runs/
+out/
 kalast.log
 kalast.err
-kalast.out
+kalast.out

cfg/cfg.toml

@@ -1,23 +1,53 @@
 [window]
-ambient = [0.07, 0.07, 0.07]
-export_frames = true
-shadows = true
+colormap = { name = "inferno", vmax = 400.0, scalar = "temperature" }
+# export_frames = true
+# high_dpi = true
+# width = 1024
+# height = 768
 
-[scene]
-sun.position.cartesian = [100, 0, 0]
+[simulation]
+routines = "thermal"
+step = 198
+duration = 94780800
 
-[scene.camera]
-name = "test"
-position.cartesian = [10, 0, 0]
-projection.perspective = 0.5236
-# projection = "orthographic"
+[simulation.export]
+step = 198
+duration = 172800
+period = 172800
+cooldown_start = 94608000
 
-[[bodies]]
-name = "A"
-state.cartesian.position = [0, 0, 0]
-material.albedo = 0.1
+# 1728000 = 20 days
+# 1900800 = 22 days
+# 31536000 = 365 days
+# 31708800 = 367 days
+# 94608000 = 365*3 days
+# 94780800 = 365*3+2 days
+
+[scene]
+sun.position.cartesian = [149597870.7, 0, 0]
+camera.position.cartesian = [5, 0, 0]
 
 [[bodies]]
-name = "B"
-state.cartesian.position = [3, 0, 0]
-material.albedo = 0.1
+name = "Simple"
+color = "data"
+interior.grid1d.increasing = { skin = "one", m = 10, n = 5, b = 20 }
+record = { columns = [0], mesh = true, depth = true }
+spin = { period = 86400 }
+temperature.scalar = 1
+
+[bodies.mesh]
+shape.shape = "triangle"
+position = [0.0, -0.3333, 0.3333]
+
+# 1) Identity
+# 2 to 3) Rotation 90° around X and Y axis
+# orientation = [1, 0, 0, 0, 1, 0, 0, 0, 1]
+# orientation = [1, 0, 0, 0, 0, 1, 0, -1, 0]
+orientation = [0, 0, -1, 0, 1, 0, 1, 0, 0]
+
+[bodies.material]
+albedo = 0.1
+emissivity = 0.9
+thermal_inertia = 500.0
+density = 2000.0
+heat_capacity = 600.0

examples/apophis-gravity/cfg/cfg.toml

@@ -0,0 +1,27 @@
+[window]
+# colormap = { name = "cividis", vmax = 100, scalar = "file" }
+colormap = { name = "cividis", vmax = 50, scalar = "file" }
+export_frames = true
+high_dpi = true
+# width = 1024
+# height = 768
+
+[simulation]
+step = 60
+duration = 12000
+
+[simulation.export]
+step = 60
+duration = 12000
+period = 12000
+cooldown_start = 0
+
+[scene]
+camera.position.cartesian = [0, 1, 0]
+
+[[bodies]]
+name = "Apophis"
+mesh = { shape.path = "/Users/gregoireh/data/kalast/runs/Apophis flyby Hakan/input/apophis_new.obj", factor = [1e-3, 1e-3, 1e-3] }
+color = "data"
+# file_data = { path = "/Users/gregoireh/data/kalast/runs/Gravity Elisa/input/results_slope_change_percentage-2.csv" }
+file_data = { path = "/Users/gregoireh/data/kalast/runs/Gravity Elisa/input/results_slope_change-2.csv" }

examples/apophis/cfg/cfg.toml

@@ -0,0 +1,47 @@
+[window]
+colormap = { name = "inferno", vmax = 400.0, scalar = "temperature" }
+export_frames = true
+high_dpi = true
+# width = 1024
+# height = 768
+
+[simulation]
+routines = "thermal"
+step = 600
+duration = 25920000
+file = { path = "/Users/gregoireh/Downloads/hakan/forkalast-2.csv", behavior = "loop" }
+
+[simulation.export]
+step = 60
+duration = 109800
+period = 1296000
+cooldown_start = 13399200
+
+# 1296000 = 360h
+# 12960000 = 3600h
+# 25920000 = 3600h * 2
+# 109800 = 30.5h
+# 439200 = 30.5h * 4
+# 549000 = 30.5h * 5
+# 1098000 = 30.5h * 10
+# 13399200 = 30.5h * 4 + 3600h * 1
+# 117079200 = 30.5h * 4 + 3600h * 9
+
+[scene]
+sun.position = "file"
+camera.position.cartesian = [-1, 0, 0]
+
+[[bodies]]
+name = "Apophis"
+color = "data"
+mesh = { shape.path = "/Users/gregoireh/Downloads/hakan/apophis_new.obj", factor = [1e-3, 1e-3, 1e-3] }
+interior.grid1d.linear = { size = 40, a = 2e-2 }
+state = "file"
+record = { rows = [0], mesh = true, depth = true }
+
+[bodies.material]
+albedo = 0.3
+emissivity = 0.9
+thermal_inertia = 500.0
+density = 2100.0
+heat_capacity = 600.0

examples/juventas/cfg/cfg.toml

@@ -1,10 +1,8 @@
 [window]
 ambient = [0.07, 0.07, 0.07]
 export_frames = true
-# width = 2048
-# height = 1536
-width = 1024
-height = 768
+width = 2048
+height = 1536
 shadows = true
 
 [spice]
@@ -13,12 +11,20 @@ frame = "JUVENTAS_NAVCAM"
 origin = "DIDYMOS"
 
 [simulation]
-start = "2027 MAR 23 16:02:00.0 TDB"
+# start = "2027 MAR 23 16:02:00.0 TDB"
+# start = "2027 MAR 26 16:16:10.0 TDB"
+# start = "2027 MAR 30 16:16:10.0 TDB"
+# start = "2027 APR 02 16:16:10.0"
+# start = "2027 APR 05 16:16:10.0"
+# start = "2027 MAY 03 16:16:10.0"
+start = "2027 MAY 06 06:08:47.0"
+
+
+duration = 2678400
 step = 60
-duration = 259080
-# 86400 * 60 = 5184000
-# 86400 * 3 = 259200
-# ^-120 = 259080
+export = { step = 60, duration = 2678400 }
+
+# 86400 * 31 = 2678400
 
 [scene]
 sun.position = "spice"
@@ -29,8 +35,7 @@ name = "JUVENTAS"
 position = "spice"
 up = [1, 0, 0]
 direction.vector = [0, 0, 1]
-# projection.perspective = 0.4974188  # 28.5°
-projection.perspective = 0.09
+projection.perspective = 0.4974188  # 28.5°
 
 [[bodies]]
 name = "DIDYMOS"

plot/apophis.py

@@ -0,0 +1,72 @@
+from pathlib import Path
+
+import numpy
+from natsort import natsorted
+from pyarrow import csv
+
+import plot
+
+d = dict()
+
+path_data = Path("/Users/gregoireh/data/kalast")
+path_runs = path_data / "runs"
+path_draft = path_data / "graphs/draft"
+
+path_run = path_runs / "Apophis flyby Hakan"
+
+path_cfg = path_run / "cfg/cfg.toml"
+path_simu = path_run / "simu"
+path_rec = path_simu / "rec"
+
+d["nf"] = 3996
+d["body"] = "Apophis"
+path_setup_body = path_simu / f"{d['body']}"
+path_progress = path_simu / "progress.csv"
+
+path_sph = path_setup_body / "mesh.csv"
+read_options = csv.ReadOptions(
+    column_names=[
+        "x",
+        "y",
+        "z",
+        "lon",
+        "lat",
+        "rad",
+    ],
+    skip_rows=1,
+)
+tab = csv.read_csv(path_sph, read_options)
+df = tab.to_pandas()
+d["centers"] = numpy.array(
+    [df["x"].to_numpy(), df["y"].to_numpy(), df["z"].to_numpy()]
+).T
+d["sph"] = numpy.array(
+    [df["lon"].to_numpy(), df["lat"].to_numpy(), df["rad"].to_numpy()]
+).T
+
+path_depth = path_setup_body / "depth.csv"
+read_options = csv.ReadOptions(
+    column_names=[
+        "depth",
+    ],
+    skip_rows=1,
+)
+tab = csv.read_csv(path_depth, read_options)
+df = tab.to_pandas()
+d["depth"] = df["depth"].to_numpy()
+
+list_path_date = [p for p in natsorted((path_rec).glob("*"), key=str) if p.is_dir()]
+list_elapsed = [int(p.name) for p in list_path_date]
+it_list_path_date = iter(list_path_date)
+
+path_date = next(it_list_path_date)
+
+p_csv = path_date / d["body"] / "temperatures/temperatures-all.csv"
+read_options = csv.ReadOptions()
+tab = csv.read_csv(p_csv, read_options)
+df = tab.to_pandas()
+d["tmp"] = df["tmp"].array.reshape((d["nf"], -1)).T
+print(d["tmp"].shape)
+print(d["tmp"][0].min(), d["tmp"][0].max(), d["tmp"][0].mean())
+
+plot.surface.plot(d, save=True)

plot/main.mplstyle

@@ -0,0 +1,26 @@
+figure.figsize: 4.0, 4.0
+figure.autolayout: True
+
+svg.fonttype: none
+font.family: Helvetica, Arial, sans-serif
+
+xtick.top: True
+xtick.bottom: True
+xtick.major.size: 4.0
+xtick.major.width: 1.0
+xtick.minor.size: 2.0
+xtick.minor.width: 1.0
+xtick.direction: in
+
+ytick.left: True
+ytick.right: True
+ytick.major.size: 4.0
+ytick.major.width: 1.0
+ytick.minor.size: 2.0
+ytick.minor.width: 1.0
+ytick.direction: in
+
+# grid.color: "#f8f8f8ff"
+# grid.linewidth: 0.25
+
+savefig.format: svg

plot/plot/__init__.py

@@ -0,0 +1 @@
+from plot import daily, depth, smap

plot/plot/daily.py

@@ -0,0 +1,48 @@
+#!/usr/bin/env python
+
+from pathlib import Path
+
+import matplotlib
+import numpy
+import scipy
+import util
+from matplotlib import pyplot, ticker
+from mpl_toolkits.basemap import Basemap
+
+Data = dict[str, numpy.array]
+Data_m = dict[str, str]
+
+
+def plot(d: dict[str, Data], save: bool = False):
+    pyplot.style.use("plot/main.mplstyle")
+    fig, ax = pyplot.subplots(figsize=(15, 7.3))
+
+    ax.set_xlabel("Time elapsed (h)")
+    ax.set_ylabel("Temperature (K)")
+
+    dt = 200
+    xmax = d["tmp-cols"].shape[0] * dt
+    x = numpy.arange(0, xmax, dt) / 3600
+    y = d["tmp-cols"][:, 0]
+
+    ax.plot(x, y)
+
+    ax.grid(True)
+
+    ax.set_xlim(0, x.max())
+    ax.set_ylim(200, 400)
+
+    # loc = ticker.MultipleLocator(base=30.0)
+    # ax.xaxis.set_major_locator(loc)
+
+    loc = ticker.MultipleLocator(base=50)
+    ax.yaxis.set_major_locator(loc)
+
+    if save:
+        out = Path("out/surface")
+        if not out.exists():
+            out.mkdir(parents=True)
+
+        fig.savefig(out / "daily.png", bbox_inches="tight", dpi=300)
+
+    pyplot.show()