Skip to main content
added 2452 characters in body
Source Link
Alex Trounev
  • 48.8k
  • 3
  • 51
  • 115

Update 1. Second toy example can be solved with using "Do" loop as follows

Unprotect[HeavisideTheta]; 
HeavisideTheta[0.0] = 1; Protect[HeavisideTheta];

Legendre[n_, t_] = LegendreP[n - 1, (2 t - 1)]*Sqrt[2 (n - 1) + 1];
p[n_, t_] = Legendre[n, Exp[-t]]*HeavisideTheta[t];
\[Omega][t_] = Exp[-t];
Kk[n_, t_] = p[n, t]*\[Omega][t];
func[t_] = Sin[2 Pi*t];
d = 2;
T = 3;
f = Function[\[Tau], 
   Evaluate@
    Table[Integrate[
      func[t]*p[k, t - \[Tau]]*\[Omega][
        t - \[Tau]], {t, \[Tau], \[Infinity]}, 
      Assumptions -> \[Tau] >= 0], {k, d}]];


interactions = 
  Function[\[CapitalDelta], 
   Evaluate@
    Table[Integrate[
      Kk[i, t]*Kk[j, t + \[CapitalDelta]], {t, 0, \[Infinity]}, 
      Assumptions -> \[CapitalDelta] >= 0], {i, d}, {j, d}]];

buildV = 
  Function[spkList, 
   Function[\[Tau], 
    Evaluate[
     f[\[Tau]] - 
      Table[Total[
        HeavisideTheta[\[Tau] - #[[2]]]*#[[3]]*
           interactions[\[Tau] - #[[2]]][[m, #[[1]]]] & /@ 
         spkList], {m, d}]]]];


tStart = 0;
spikes = Table[{i, tStart - eps, f[tStart - eps][[i]]}, {i, d}]; V = 
 buildV[spikes];


tS[0] = tStart; spike[0] = spikes; 
spk[i_] := {{{1, tS[i], 1}}, {{1, tS[i], -1}}, {{2, tS[i], 1}}, {{2, 
    tS[i], -1}}}; Do[{V1[i], V2[i]} = 
  NDSolveValue[{{v1'[t], v2'[t]} == V'[t], 
    WhenEvent[v1[t] == 1/4, {tS[i] = t, j = 1}; "StopIntegration"], 
    WhenEvent[v1[t] == -1/4, {tS[i] = t, j = 2}; "StopIntegration"], 
    WhenEvent[v2[t] == 1/4, {tS[i] = t, j = 3}; "StopIntegration"], 
    WhenEvent[v2[t] == -1/4, {tS[i] = t, j = 4}; "StopIntegration"],
    WhenEvent[Abs[v2[t]] == 1/4, tS[i] = t; 
     "StopIntegration"], {v1[tS[i - 1]], v2[tS[i - 1]]} == 
     V[tS[i - 1]]}, {v1, v2}, {t, tS[i - 1], tS[i - 1] + 1}]; 
 spike[i] = Join[spike[i - 1], spk[i][[j]]]; 
 V = Evaluate[buildV[spike[i]]];, {i, 10}];

Visualization

VV1 = Piecewise[
  Table[{V1[i], tS[i - 1] <= t < tS[i]}, {i, 10}]]; VV2 = 
 Piecewise[Table[{V2[i], tS[i - 1] <= t < tS[i]}, {i, 10}]];

Plot[{VV1[t], VV2[t]}, {t, tS[0], tS[10]}, Frame -> True, 
 PlotLegends -> {"v1", "v2"}]

Fugure 2

Update 1. Second toy example can be solved with using "Do" loop as follows

Unprotect[HeavisideTheta]; 
HeavisideTheta[0.0] = 1; Protect[HeavisideTheta];

Legendre[n_, t_] = LegendreP[n - 1, (2 t - 1)]*Sqrt[2 (n - 1) + 1];
p[n_, t_] = Legendre[n, Exp[-t]]*HeavisideTheta[t];
\[Omega][t_] = Exp[-t];
Kk[n_, t_] = p[n, t]*\[Omega][t];
func[t_] = Sin[2 Pi*t];
d = 2;
T = 3;
f = Function[\[Tau], 
   Evaluate@
    Table[Integrate[
      func[t]*p[k, t - \[Tau]]*\[Omega][
        t - \[Tau]], {t, \[Tau], \[Infinity]}, 
      Assumptions -> \[Tau] >= 0], {k, d}]];


interactions = 
  Function[\[CapitalDelta], 
   Evaluate@
    Table[Integrate[
      Kk[i, t]*Kk[j, t + \[CapitalDelta]], {t, 0, \[Infinity]}, 
      Assumptions -> \[CapitalDelta] >= 0], {i, d}, {j, d}]];

buildV = 
  Function[spkList, 
   Function[\[Tau], 
    Evaluate[
     f[\[Tau]] - 
      Table[Total[
        HeavisideTheta[\[Tau] - #[[2]]]*#[[3]]*
           interactions[\[Tau] - #[[2]]][[m, #[[1]]]] & /@ 
         spkList], {m, d}]]]];


tStart = 0;
spikes = Table[{i, tStart - eps, f[tStart - eps][[i]]}, {i, d}]; V = 
 buildV[spikes];


tS[0] = tStart; spike[0] = spikes; 
spk[i_] := {{{1, tS[i], 1}}, {{1, tS[i], -1}}, {{2, tS[i], 1}}, {{2, 
    tS[i], -1}}}; Do[{V1[i], V2[i]} = 
  NDSolveValue[{{v1'[t], v2'[t]} == V'[t], 
    WhenEvent[v1[t] == 1/4, {tS[i] = t, j = 1}; "StopIntegration"], 
    WhenEvent[v1[t] == -1/4, {tS[i] = t, j = 2}; "StopIntegration"], 
    WhenEvent[v2[t] == 1/4, {tS[i] = t, j = 3}; "StopIntegration"], 
    WhenEvent[v2[t] == -1/4, {tS[i] = t, j = 4}; "StopIntegration"],
    WhenEvent[Abs[v2[t]] == 1/4, tS[i] = t; 
     "StopIntegration"], {v1[tS[i - 1]], v2[tS[i - 1]]} == 
     V[tS[i - 1]]}, {v1, v2}, {t, tS[i - 1], tS[i - 1] + 1}]; 
 spike[i] = Join[spike[i - 1], spk[i][[j]]]; 
 V = Evaluate[buildV[spike[i]]];, {i, 10}];

Visualization

VV1 = Piecewise[
  Table[{V1[i], tS[i - 1] <= t < tS[i]}, {i, 10}]]; VV2 = 
 Piecewise[Table[{V2[i], tS[i - 1] <= t < tS[i]}, {i, 10}]];

Plot[{VV1[t], VV2[t]}, {t, tS[0], tS[10]}, Frame -> True, 
 PlotLegends -> {"v1", "v2"}]

Fugure 2

added 181 characters in body
Source Link
Alex Trounev
  • 48.8k
  • 3
  • 51
  • 115

Do you mean something like that?We can use discrete variable to control w as follows

buildV = Function[w, Function[t, Sin[w*t]]];
V = buildV[5*Random[]];
sol = NDSolve[{{v'[t] == V'[t]buildV[w[t]]'[t], v[0] == 0, 
w[0] == 5 RandomReal[]},  


WhenEvent[v[t] == 1/2, V{w[t] =-> buildV[5*Random[]];5 RandomReal[],
     v[t] -> v[t] + 0.05]05}]}, v[0] == V[0]}{v, v[t]w}, {t, 0, 10}][[1]], 
  DiscreteVariables -> w]

Visualization

Plot[v[t]Plot[Evaluate[{v[t], w[t]} /. solsol[[1]]], {t, 0, 10}, 
 PlotPoints -> 200]200, PlotRange -> All, PlotLegends -> {"v", "w"}, 
 AxesLabel -> Automatic] 

Figure 1Figure 1

Do you mean something like that?

buildV = Function[w, Function[t, Sin[w*t]]];
V = buildV[5*Random[]];
sol = NDSolve[{v'[t] == V'[t], 
    WhenEvent[v[t] == 1/2, V = buildV[5*Random[]];
     v[t] -> v[t] + 0.05], v[0] == V[0]}, v[t], {t, 0, 10}][[1]]

Visualization

Plot[v[t] /. sol, {t, 0, 10}, PlotPoints -> 200]

Figure 1

We can use discrete variable to control w as follows

buildV = Function[w, Function[t, Sin[w*t]]];
sol = NDSolve[{{v'[t] == buildV[w[t]]'[t], v[0] == 0, 
w[0] == 5 RandomReal[]},  


WhenEvent[v[t] == 1/2, {w[t] -> 5 RandomReal[],
     v[t] -> v[t] + 0.05}]}, {v, w}, {t, 0, 10}, 
  DiscreteVariables -> w]

Visualization

Plot[Evaluate[{v[t], w[t]} /. sol[[1]]], {t, 0, 10}, 
 PlotPoints -> 200, PlotRange -> All, PlotLegends -> {"v", "w"}, 
 AxesLabel -> Automatic] 

Figure 1

Source Link
Alex Trounev
  • 48.8k
  • 3
  • 51
  • 115

Do you mean something like that?

buildV = Function[w, Function[t, Sin[w*t]]];
V = buildV[5*Random[]];
sol = NDSolve[{v'[t] == V'[t], 
    WhenEvent[v[t] == 1/2, V = buildV[5*Random[]];
     v[t] -> v[t] + 0.05], v[0] == V[0]}, v[t], {t, 0, 10}][[1]]

Visualization

Plot[v[t] /. sol, {t, 0, 10}, PlotPoints -> 200]

Figure 1