%{
This scripts runs through the steps to create the de-nuovo spiral used as
IC for all simulations in the publication.

Parameter references
[1] Sgro et al, From intracellular signaling to population oscillations: bridging size- and time-scales in collective behavior
[2] Hoefer et al, Dictyostelium discoideum: Cellular Self-Organization in an Excitable Biological Medium
[3] Vidal-Henriquez et at, Spontaneous center formation in Dictyostelium discoideum
[HF] From experiments by Hugh Ford

Written 12/2021 by Angelika Manhart

%}

%% Set base parameters for spiral creation

p.Lx=4; % size of domain in x in mm
p.Ly=4; % size of domain in y in mm
p.maxTime=100; % in min
p.nX=round(p.Lx*50); % nr of grid points in x
p.nY=round(p.Ly*50); % nr of grid points in y

% cell movement
cellDens=8000; % per mm^2: 4000-5000 [HF]
p.N=round(p.Lx*p.Ly*cellDens); % nr of cells
p.movNoise=0;
p.repRad=5e-3; % cell repulsion range in mm(to do with cell size) in mm [HF]
p.repStrength=0; % cell repulsion strength

% chemotaxis parameters
p.chemDesens=2.67*1e-5; % [HF] & numerics
p.chemTimeScale=6.7*1e-3; % [HF] & numerics
p.turnRate=3*1e-6; % [HF] & numerics
p.speedModelP1=0; % [HF] & numerics, if set to zero no cell movement
p.speedModelP2=0.6; % [HF] & numerics

% internal cell dynamics
p.epsi=0.1;   %[1] ratio between activator and inhibitor timescales
p.preA=0.058;  % [1] magnitude of response of pre-processing (Pp) unit
p.preKd=1e-5;  % [1] threshold for response of pre-processing (Pp) unit
p.inhDecay=0.5; % [1] inhibitor decay
p.inhBasal=1.2; % [1] inhibitor steady state value w/o signal (not really true)
p.intTimeScale=1/6; % (estimated from data) 1/n means roughly n spike every 30mins in a constant ext cAMP scenario

% signal dynamics
p.sigDiff=2.4e-2; % signal diffusion constant in mm^2/min: 1.8e-3 [2], 2.4e-2 [3]
p.sigDecay=5; % signal decay rate in /min: 1-10 [3]
p.sigSecret=1e3; % % signal secretion rate per active cell in molecules/min
p.sigSpreadRad=0; % the radius over which the secreted signal is initially spread in mm;

p.dt=0.05; % time step in minutes

% simulation setting
p.timeDepParas=0;
p.optBC=0; % 0 .. PBC, 1.. no flux
p.optSigDisc=1; % 0 .. explicit, 1.. implcit
p.optICsig=0;
p.optICcells=0;
p.saveEvery=0;
p.plotEvery=5;
p.isVideo=0;

%% Creat linear wave and let it equilibrate

addpath('Program_Files')
p.maxTime=20;
results=main_Dicty(p);
ICs=results(2);

%% Delete half of external cAMP and de-activate cells in that half. Spiral forms.
%{
In main_Dicty.m activate these lines (in "Prelude"):
sigConc(1:p.nX/2,:)=0;
inhConc(yPos<0)=2.5;
%}

p.optBC=1; % set BCs to no-flux
p.maxTime=10;
results=main_Dicty(p,ICs);
ICs=results(2);

% save('ICs\ICs_naturalSpiral.mat', 'ICs')
