Browse Source

[nucleus] recompute coordinates once in a while

Gaspard Jankowiak 11 months ago
parent
commit
0cb10429f2
2 changed files with 36 additions and 6 deletions
  1. 5
    1
      src/CellSim.jl
  2. 31
    5
      src/Nucleus.jl

+ 5
- 1
src/CellSim.jl View File

@@ -195,6 +195,10 @@ function read_config(config_filename::String)
195 195
         config["output_prefix"] = "runs/"
196 196
     end
197 197
 
198
+    if !haskey(config, "recompute_nucleus_coords")
199
+        config["recompute_nucleus_coords"] = 5
200
+    end
201
+
198 202
     if (config["metrics"]["start_iteration"] <= 0) && !F.nucleus
199 203
         println("[WARNING] metrics.start_iteration is non-positive and nucleus is deactivated, metrics will never start!")
200 204
     end
@@ -457,7 +461,7 @@ function launch(P::CSC.Params, F::CSC.Flags, config)
457 461
             if F.centrosome
458 462
                 Nucleus.compute_centronuclear_force(potentials, coords, nucleus_coords, P, F)
459 463
             end
460
-            Nucleus.update_coords(old_nucleus_coords, nucleus_coords, potentials, P, F, temparrays)
464
+            Nucleus.update_coords(old_nucleus_coords, nucleus_coords, potentials, P, F, temparrays, k)
461 465
 
462 466
         end
463 467
 

+ 31
- 5
src/Nucleus.jl View File

@@ -36,7 +36,7 @@ mutable struct NucleusCoords
36 36
     # normal vector
37 37
     n::Matrix{Float64}
38 38
 
39
-    # local length element
39
+    # local length element, i.e. log(r)
40 40
     η::Vector{Float64}
41 41
 
42 42
     # nuclear membrane length
@@ -70,6 +70,28 @@ function g_p(x::Vector{Float64}, α::Float64)
70 70
     return -(2α*min.(α*x.-1, 0.0).*log.(α*x) .+ min.(α*x.-1, 0.0).^2 ./x)
71 71
 end
72 72
 
73
+function recompute_nucleus_coords(c::NucleusCoords)
74
+    delta_Y = c.Y - c.Y[c.circ_idx.m1,:]
75
+    c.r[:] = sqrt.(sum(abs2, delta_Y; dims=2))
76
+    c.q[:] = 0.5*(c.r + c.r[c.circ_idx.p1])
77
+
78
+    c.θ[:] = angle.(delta_Y[:,1] + delta_Y[:,2]*im)
79
+    delta_θ = rem.(c.θ - c.θ[c.circ_idx.m1], 2π, RoundNearest)
80
+    c.θ[:] = cumsum([c.θ[1]; delta_θ[1:end-1]])
81
+
82
+    # http://cral-labo.univ-lyon1.fr/labo/fc/Ateliers_archives/ateliers_2005-06/cercle_3pts.pdf
83
+    xc = (0.5*(c.Y[c.circ_idx.p1,1].^2 - c.Y[:,1].^2 + c.Y[c.circ_idx.p1,2].^2 - c.Y[:,2].^2) ./ (c.Y[c.circ_idx.p1,2]-c.Y[:,2])
84
+          - 0.5*(c.Y[:,1].^2 - c.Y[c.circ_idx.m1,1].^2 + c.Y[:,2].^2 - c.Y[c.circ_idx.m1,2].^2) ./ (c.Y[:,2]-c.Y[c.circ_idx.m1,2])) ./ ((c.Y[:,1] - c.Y[c.circ_idx.m1,1])./(c.Y[:,2] - c.Y[c.circ_idx.m1,2]) - (c.Y[c.circ_idx.p1,1] - c.Y[:,1])./(c.Y[c.circ_idx.p1,2] - c.Y[:,2]))
85
+
86
+    yc = -(c.Y[:,1] - c.Y[c.circ_idx.m1,1])./(c.Y[:,2] - c.Y[c.circ_idx.m1,2]) .* xc + 0.5*(c.Y[:,1].^2 - c.Y[c.circ_idx.m1,1].^2 + c.Y[:,2].^2 - c.Y[c.circ_idx.m1,2].^2) ./ (c.Y[:,2]-c.Y[c.circ_idx.m1,2])
87
+    radius = sqrt.((xc - c.Y[:,1]).^2 + (yc - c.Y[:,2]).^2)
88
+    c.k[:] = 1 ./radius
89
+    c.k[:] = 0.5*(c.k + c.k[c.circ_idx.m1])
90
+
91
+    c.η[:] = log.(c.r)
92
+    c.L = sum(c.r)
93
+end
94
+
73 95
 function compute_contact_force(pots::CSC.InteractionPotentials,
74 96
                                cor_coords::PointCoords, c::NucleusCoords,
75 97
                                P::Params, F::Flags)
@@ -352,10 +374,10 @@ function initialize_coords(P::Params, F::Flags, cortex_c::PointCoords; fill_wall
352 374
     # The nucleus is initialized as a circle centered on the
353 375
     # barycenter of the cell, with radius P.N_r_init
354 376
 
355
-    t = collect(range(0, stop=2pi, length=P.Nnuc+1))[1:P.Nnuc]
377
+    t = collect(range(-0.5pi, stop=1.5pi, length=P.Nnuc+1))[1:P.Nnuc]
356 378
 
357
-    θ0 = π/P.Nnuc
358
-    θ = collect(π/2 + θ0 .+ (-1:(P.Nnuc-2))*2*θ0)
379
+    θ0 = 2π/P.Nnuc
380
+    θ = collect(-0.5*θ0 .+ (0:(P.Nnuc-1))*θ0)
359 381
 
360 382
     if fill_wall
361 383
         bc = [0.0 0.5π/P.f_ω0]
@@ -384,7 +406,7 @@ end
384 406
 
385 407
 function update_coords(c::NucleusCoords, new_c::NucleusCoords,
386 408
                        potentials::CSC.InteractionPotentials,
387
-                       P::Params, F::Flags, temparrays::CSC.TempArrays6)
409
+                       P::Params, F::Flags, temparrays::CSC.TempArrays6, iter::Int64)
388 410
 
389 411
     N_W = potentials.N_W
390 412
     N_∇W = potentials.N_∇W
@@ -402,6 +424,10 @@ function update_coords(c::NucleusCoords, new_c::NucleusCoords,
402 424
     update_θ(c, new_c, N_W, N_∇W, P, F, temparrays)
403 425
     update_Y(c, new_c, N_W, N_∇W, P, F, temparrays)
404 426
 
427
+    if iter%5 == 0
428
+        recompute_nucleus_coords(new_c)
429
+    end
430
+
405 431
     if F.DEBUG
406 432
         println("POST alpha, beta, r, q, K, θ")
407 433
         # display([new_c.α new_c.β new_c.r new_c.q new_c.k new_c.θ])

Loading…
Cancel
Save