Merge main

Author: Sosnowsky

blobmodel/model.py

@@ -337,11 +337,12 @@ def _compute_start_stop(self, blob: Blob, speed_up: bool, error: float):
                     np.abs(blob.t_init * blob.v_x)
                     + blob.width_p * np.log(error * np.sqrt(np.pi))
                     + blob.pos_x
-                    - self._geometry.t[0]
+                    - self._geometry.t[0] * blob.v_x
                 )
                 / np.abs(self._geometry.dt * blob.v_x)
             ),
         )
+        start = np.minimum(self._geometry.t.size, start)
         # ignores t_drain when calculating stop time
         stop = np.minimum(
             self._geometry.t.size,

examples/compare_to_analytical_sol.py

@@ -12,38 +12,34 @@
 # with analytically derived results. See O. E. Garcia, et al.; Phys. Plasmas 1 May 2016; 23 (5): 052308. https://doi.org/10.1063/1.4951016
 
 # use DefaultBlobFactory to define distribution functions fo random variables
-bf = DefaultBlobFactory(A_dist=DistributionEnum.deg, vy_dist=DistributionEnum.zeros)
+bf = DefaultBlobFactory(
+    A_dist=DistributionEnum.deg, vy_dist=DistributionEnum.zeros, vy_parameter=10
+)
+t_drain = 1e10
 
 tmp = Model(
-    Nx=20,
-    Ny=1,
+    Nx=10,
+    Ny=10,
     Lx=10,
-    Ly=0,
-    dt=1,
-    T=1000,
-    blob_shape=BlobShapeImpl(BlobShapeEnum.exp, BlobShapeEnum.gaussian),
-    t_drain=2,
-    periodic_y=False,
+    Ly=10,
+    dt=0.1,
+    T=10000,
+    blob_shape=BlobShapeImpl(BlobShapeEnum.gaussian, BlobShapeEnum.gaussian),
+    t_drain=1e10,
+    periodic_y=True,
     num_blobs=10000,
     blob_factory=bf,
-    one_dimensional=True,
+    one_dimensional=False,
+    t_init=10,
 )
 
 ds = tmp.make_realization(speed_up=True, error=1e-10)
-x = np.linspace(0, 10, 20)
-t_p = 1
-t_w = 1 / 10
-amp = 1
-v_p = 1.0
-t_loss = 2.0
-t_d = t_loss * t_p / (t_loss + t_p)
-
-analytical_profile = t_d / t_w * amp * np.exp(-x / (v_p * t_loss))
-
-ds.n.isel(y=0).mean(dim=("t")).plot(label="blob_model")
-plt.yscale("log")
-plt.xlabel("x")
-plt.ylabel("Time average")
-plt.plot(x, analytical_profile, label="analytical solution")
-plt.legend()
+
+mean = ds.n.mean(dim="t").values
+
+fig, ax = plt.subplots()
+
+im = ax.imshow(mean, clim=(0, np.max(mean)))
+plt.colorbar(im, ax=ax)
+
 plt.show()

examples/point_process.py

@@ -13,7 +13,7 @@
 
 
 dt = 0.1
-T = int(1e5)
+T = int(1e4)
 num_blobs = int(1e5)
 start_index = int(100 / dt)