Fix rotation bug (#93)

* implement correct distribution for isotropic 3D reorientations * remove `azimuthal` from api, rename `polar` to `angle`, default angles to 3D isotropic pdf on sphere * replace `polar` with `angle` * bump version - breaking change to api * no defaults for tumble * define `Isotropic` function to produce correct pdf for given dimensionality * update tests * update examples * update docs * update defaults of chemotaxis models * update api * remove old unused files * fix docs * fix motile states constructors * update benchmarks
mastrof · Oct 5, 2024 · c75eed9 · c75eed9
1 parent ded80e0
commit c75eed9
Show file tree

Hide file tree

Showing 34 changed files with 171 additions and 223 deletions.
diff --git a/Project.toml b/Project.toml
@@ -1,7 +1,7 @@
 name = "MicrobeAgents"
 uuid = "b17a4bac-8667-4a2d-84f7-1883ae0b8dbb"
 authors = ["Riccardo Foffi <[email protected]> and contributors"]
-version = "0.4.2"
+version = "0.5.0"
 
 [deps]
 Agents = "46ada45e-f475-11e8-01d0-f70cc89e6671"

diff --git a/benchmark/benchmarks.jl b/benchmark/benchmarks.jl
@@ -53,7 +53,7 @@ open(filename, "a") do io
                 motility = RunReverseFlick(0.0, [10.0], 0.0, [10.0])
             else
                 model = StandardABM(T{D,2}, space, dt; container=Vector)
-                motility = RunTumble(0.0, [10.0], 0)
+                motility = RunTumble(0.0, [10.0], Isotropic(D), 0)
             end
             add_agent!(model; motility)
             nsteps = 1000

diff --git a/benchmark/benchmarks_0.5.0.md b/benchmark/benchmarks_0.5.0.md
@@ -0,0 +1,20 @@
+# MicrobeAgents.jl benchmarks - 0.5.0 
+
+## Simple stepping
+Type             | Time (μs)        | Allocs           | Memory          
+ --- | --- | --- | --- 
+Microbe{1}       | 4.98e-01         | 1.10e+01         | 1.79e+02        
+Microbe{2}       | 6.31e-01         | 1.50e+01         | 3.07e+02        
+Microbe{3}       | 7.35e-01         | 1.50e+01         | 3.07e+02        
+BrownBerg{1}     | 6.01e-01         | 1.70e+01         | 7.55e+02        
+BrownBerg{2}     | 7.43e-01         | 2.10e+01         | 9.79e+02        
+BrownBerg{3}     | 8.57e-01         | 2.10e+01         | 1.08e+03        
+Brumley{1}       | 6.40e-01         | 1.65e+01         | 7.47e+02        
+Brumley{2}       | 7.58e-01         | 2.05e+01         | 9.71e+02        
+Brumley{3}       | 8.59e-01         | 2.05e+01         | 1.07e+03        
+Celani{1}        | 5.48e-01         | 1.30e+01         | 3.71e+02        
+Celani{2}        | 6.81e-01         | 1.70e+01         | 5.31e+02        
+Celani{3}        | 8.05e-01         | 1.70e+01         | 5.63e+02        
+SonMenolascina{1} | 6.08e-01         | 1.60e+01         | 7.39e+02        
+SonMenolascina{2} | 7.39e-01         | 2.00e+01         | 9.64e+02        
+SonMenolascina{3} | 8.51e-01         | 2.00e+01         | 1.06e+03        
diff --git a/docs/src/api.md b/docs/src/api.md
@@ -14,8 +14,7 @@ motilestate
 states
 transition_weights
 duration
-polar
-azimuthal
+angle
 distance
 distancevector
 random_speed

diff --git a/docs/src/examples/Chemotaxis/1_linear_ramp.md b/docs/src/examples/Chemotaxis/1_linear_ramp.md
@@ -88,7 +88,7 @@ for i in 1:n
     p0 = spacesize(model) ./ 2
     delta = rand(abmrng(model), SVector{2}) .* 10
     pos = p0 .+ delta # random scatter around center
-    motility = RunTumble(0.67, [30.0], 0.1)
+    motility = RunTumble(0.67, [30.0], Isotropic(2), 0.1)
     add_agent!(pos, model; motility)
 end
 model

diff --git a/docs/src/examples/Chemotaxis/2_celani_gauss2D.md b/docs/src/examples/Chemotaxis/2_celani_gauss2D.md
@@ -35,7 +35,7 @@ properties = Dict(
 )
 
 model = StandardABM(Celani{2,2}, space, timestep; properties)
-motility = RunTumble(0.67, [30.0], 0.1)
+motility = RunTumble(0.67, [30.0], Isotropic(2), 0.1)
 
 for _ in 1:300
     add_agent!(model; motility,

diff --git a/docs/src/examples/Chemotaxis/4_response_functions.md b/docs/src/examples/Chemotaxis/4_response_functions.md
@@ -61,11 +61,11 @@ properties = Dict(
 
 model = StandardABM(Union{BrownBerg{3},Celani{3}}, space, dt; properties)
 
-add_agent!(BrownBerg{3}, model; motility=RunTumble(1.0, [0]),
+add_agent!(BrownBerg{3}, model; motility=RunTumble(1.0, [0], Isotropic(3)),
     memory=1,
 )
-add_agent!(Celani{3}, model; motility=RunTumble(1.0, [0]), gain=4)
-add_agent!(Celani{3}, model; motility=RunTumble(1.0, [0]), gain=4,
+add_agent!(Celani{3}, model; motility=RunTumble(1.0, [0], Isotropic(3)), gain=4)
+add_agent!(Celani{3}, model; motility=RunTumble(1.0, [0], Isotropic(3)), gain=4,
     chemotactic_precision=50.0
 )
 

diff --git a/docs/src/examples/RandomWalks/1_randomwalk1D.md b/docs/src/examples/RandomWalks/1_randomwalk1D.md
@@ -37,18 +37,22 @@ from position `(0,)`.
 
 Moreover, it's always required to specify the motility of the microbes via
 the keyword argument `motility`.
-For example, `RunTumble` needs as inputs the average duration of runs
-and the distribution of run speeds. Passing a one-element vector `[U]`
+For example, `RunTumble` needs as inputs the average duration of runs,
+the distribution of run speeds and the distribution of reorientation angles.
+The latter is irrelevant in 1D, since the bacterium can only revert its
+direction along the line, but we still have to pass it.
+We can just use an empty vector.
+Passing a one-element vector `[U]` as the speed distribution
 means that all runs will have the same velocity `U`.
 We could have also specified the average duration of tumbles
-(equal to reversals in 1 dimension) via an optional 3rd argument (try it! e.g. 0.5);
+(equal to reversals in 1 dimension) via an optional 4th argument (try it! e.g. 0.5);
 when unspecified, the tumbles are taken to be instantaneous.
 
 ````@example 1_randomwalk1D
 n = 10 # number of microbes to add
 τ_run = 1.0 # average run duration in s
 U = 30.0 # swimming speed in μm/s
-motility = RunTumble(τ_run, [U], 0.5)
+motility = RunTumble(τ_run, [U], [], 0.5)
 foreach(_ -> add_agent!((0,), model; motility), 1:n)
 ````
 

diff --git a/docs/src/examples/RandomWalks/2_randomwalk2D_motilepatterns.md b/docs/src/examples/RandomWalks/2_randomwalk2D_motilepatterns.md
@@ -22,17 +22,23 @@ model = StandardABM(Microbe{2}, space, dt)
 We will now add microbes individually, choosing different properties for each.
 The motile pattern can be customized through the `motility` keyword.
 
-The `RunTumble` motility consists of straight runs interspersed with isotropic
+The `RunTumble` motility consists of straight runs interspersed with
 reorientations (tumbles).
 With the first argument, we set the average duration of tumbles, 1 s in this case.
 We can then define the speed to follow a a `Normal` distribution
 (from Distributions.jl) with mean 30 μm/s and standard deviation 6 μm/s.
 This means that, after every tumble, the microbe will change its speed following
 this distribution.
-Further, we can set tumbles to also have a finite duration, let's say 0.1 s.
+We should then specify the distribution of angular reorientations.
+For convenience, MicrobeAgents implements an `Isotropic` function
+which produces the appropriate distribution to have isotropic reorientations
+in a space with given dimensionality.
+In this case, using `Isotropic(2)` produces a `Uniform` distribution
+of angles between -π and +π.
+Finally, we can set tumbles to also have a finite duration, let's say 0.1 s.
 
 ````@example 2_randomwalk2D_motilepatterns
-add_agent!(model; motility=RunTumble(1.0, Normal(30,6), 0.1))
+add_agent!(model; motility=RunTumble(1.0, Normal(30,6), Isotropic(2), 0.1))
 ````
 
 The `RunReverse` motility consists of alternating straight runs and 180-degree reversals.

diff --git a/docs/src/examples/RandomWalks/3_randomwalk3D.md b/docs/src/examples/RandomWalks/3_randomwalk3D.md
@@ -17,7 +17,7 @@ dt = 0.1
 model = StandardABM(Microbe{3}, space, dt)
 
 add_agent!(model; motility=RunReverse(1.0, [55], 1.0, [55]), rotational_diffusivity=0.2)
-add_agent!(model; motility=RunTumble(2.0, Normal(30,6)))
+add_agent!(model; motility=RunTumble(2.0, Normal(30,6), Isotropic(3)))
 add_agent!(model; motility=RunReverseFlick(1.0, [30], 1.0, [6]), rotational_diffusivity=0.1)
 
 nsteps = 600

diff --git a/docs/src/introduction.md b/docs/src/introduction.md
@@ -69,16 +69,18 @@ model = StandardABM(Microbe{1}, space, dt)
 Now, with the `add_agent!` function we will populate the model with microbes of
 the specified type (`Microbe{1}`).
 The only argument we must *always* specify for `add_agent!` is the motility of the microbe, via the `motility` keyword. An overview of the motility interface is given later;
-For now we will just use a Run-Tumble motility with an average run duration of 1 second
-and a constant swimming speed of 20 micron / second.
+For now we will just use a Run-Tumble motility with an average run duration of 1 second and a constant swimming speed of 20 micron / second.
+The third required argument to the call is the distribution of reorientation
+angles, but it is irrelevant in 1-dimensional systems so we can just pass
+an empty array.
 If unspecified, position, direction and speed of the microbe will be assigned randomly;
 all the other fields will be assigned default values from the constructor (unless specified).
 To select a position, it can be passed as the first argument to the `add_agent!` call,
 and any other bacterial parameter can be defined via keyword arguments.
 All of the following are valid calls
 ```
-motility = RunTumble(1.0, [20.0])
-# a Microbe with large radius and low tumble rate
+motility = RunTumble(1.0, [20.0], [])
+# a Microbe with large radius
 add_agent!(model; motility, radius=10.0)
 # a Microbe with custom position and high coefficient of rotational diffusion
 add_agent!((53.2,), model; motility, rotational_diffusivity=0.5)

diff --git a/docs/src/randomwalks.md b/docs/src/randomwalks.md
diff --git a/docs/src/validation.md b/docs/src/validation.md
diff --git a/docs/src/validation/msd_runtumble.md b/docs/src/validation/msd_runtumble.md
@@ -15,7 +15,7 @@ where ``v`` is the microbe velocity, and ``\tilde{\tau} = \tau/(1-\cos\theta)``
 is the correlation-corrected effective run length.
 
 To test this, we will simulate run-tumble motility with different distributions
-of `polar` reorientation angles and evaluate their MSDs.
+of reorientation angles and evaluate their MSDs.
 For simplicity, since only the average angle matters rather than the entire
 distribution, we will only allow microbes to perform rotations of angle `θ`
 (and `-θ` for symmetry).
@@ -39,7 +39,7 @@ U = 30.0 # μm/s
 models = map(_ -> StandardABM(Microbe{3}, space, dt; container=Vector), θs)
 nmicrobes = 100
 for (i,θ) in enumerate(θs)
-    motility = RunTumble(τ, [U], 0.0, [θ,-θ])
+    motility = RunTumble(τ, [U], [θ,-θ])
     foreach(_ -> add_agent!(models[i]; motility), 1:nmicrobes)
 end
 

diff --git a/docs/src/validation/velocity_autocorrelations.md b/docs/src/validation/velocity_autocorrelations.md
@@ -35,7 +35,7 @@ model = StandardABM(Microbe{3}, space, Δt; container=Vector)
 n = 200
 for Mot in (RunTumble, RunReverse, RunReverseFlick), i in 1:n
     if Mot == RunTumble
-        motility = RunTumble(τ_run, [U], 0.0)
+        motility = RunTumble(τ_run, [U], Isotropic(3), 0.0)
     else
         motility = Mot(τ_run, [U], τ_run, [U])
     end

diff --git a/examples/Chemotaxis/1_linear_ramp.jl b/examples/Chemotaxis/1_linear_ramp.jl
@@ -84,7 +84,7 @@ for i in 1:n
     p0 = spacesize(model) ./ 2
     delta = rand(abmrng(model), SVector{2}) .* 10
     pos = p0 .+ delta # random scatter around center
-    motility = RunTumble(0.67, [30.0], 0.1)
+    motility = RunTumble(0.67, [30.0], Isotropic(2), 0.1)
     add_agent!(pos, model; motility)
 end
 model

diff --git a/examples/Chemotaxis/2_celani_gauss2D.jl b/examples/Chemotaxis/2_celani_gauss2D.jl
@@ -32,7 +32,7 @@ properties = Dict(
 )
 
 model = StandardABM(Celani{2,2}, space, timestep; properties)
-motility = RunTumble(0.67, [30.0], 0.1)
+motility = RunTumble(0.67, [30.0], Isotropic(2), 0.1)
 
 for _ in 1:300
     add_agent!(model; motility,

diff --git a/examples/Chemotaxis/4_response_functions.jl b/examples/Chemotaxis/4_response_functions.jl
@@ -57,11 +57,11 @@ properties = Dict(
 
 model = StandardABM(Union{BrownBerg{3},Celani{3}}, space, dt; properties)
 
-add_agent!(BrownBerg{3}, model; motility=RunTumble(1.0, [0]),
+add_agent!(BrownBerg{3}, model; motility=RunTumble(1.0, [0], Isotropic(3)),
     memory=1,
 )
-add_agent!(Celani{3}, model; motility=RunTumble(1.0, [0]), gain=4)
-add_agent!(Celani{3}, model; motility=RunTumble(1.0, [0]), gain=4,
+add_agent!(Celani{3}, model; motility=RunTumble(1.0, [0], Isotropic(3)), gain=4)
+add_agent!(Celani{3}, model; motility=RunTumble(1.0, [0], Isotropic(3)), gain=4,
     chemotactic_precision=50.0
 )
 

diff --git a/examples/RandomWalks/1_randomwalk1D.jl b/examples/RandomWalks/1_randomwalk1D.jl
@@ -34,18 +34,22 @@ from position `(0,)`.
 
 Moreover, it's always required to specify the motility of the microbes via
 the keyword argument `motility`.
-For example, `RunTumble` needs as inputs the average duration of runs
-and the distribution of run speeds. Passing a one-element vector `[U]`
+For example, `RunTumble` needs as inputs the average duration of runs,
+the distribution of run speeds and the distribution of reorientation angles.
+The latter is irrelevant in 1D, since the bacterium can only revert its
+direction along the line, but we still have to pass it.
+We can just use an empty vector.
+Passing a one-element vector `[U]` as the speed distribution
 means that all runs will have the same velocity `U`.
 We could have also specified the average duration of tumbles
-(equal to reversals in 1 dimension) via an optional 3rd argument (try it! e.g. 0.5);
+(equal to reversals in 1 dimension) via an optional 4th argument (try it! e.g. 0.5);
 when unspecified, the tumbles are taken to be instantaneous.
 =#
 
 n = 10 # number of microbes to add
 τ_run = 1.0 # average run duration in s
 U = 30.0 # swimming speed in μm/s
-motility = RunTumble(τ_run, [U], 0.5)
+motility = RunTumble(τ_run, [U], [], 0.5)
 foreach(_ -> add_agent!((0,), model; motility), 1:n)
 
 #=