Jak zamienić wyrażenie w gnieździe?

Question

Muszę obliczyć następującą funkcję z wymianą w cyklu

f[x0_, y0_]:= (s1 + s2)/t1 /. {s1 -> 
NDSolveValue[{x''[t] + x[t] == 0., y''[t] + x[t]^2 y[t] == 0., 
   x[0.] == x0, x'[0.] == 0., y[0.] == y0, y'[0.] == 0.}, 
  y, {t, 0, t1}][t1], s2 -> NDSolveValue[{x''[t] + x[t] == 0, y''[t] + x[t]^2 y[t] == 0.,
    x[0.] == x0, x'[0.] == 0., y[0.] == y0, y'[0.] == 1.}, 
  y, {t, 0, t1}][t1]} /. {t1 -> Take[Reap[
    NDSolve[{x''[t] + x[t] == 0., x[0.] == x0, x'[0.] == 0., 
      WhenEvent[x'[t] > 0., {Sow[t], "StopIntegration"}]}, 
     x, {t, 0., 100.}, 
     MaxStepSize -> 0.001]], {2, -1}][[1]][[1]][[1]]}

Ale napotykam błąd mówiący „NDSolveValue: Endpoint t1 in {t, 0., t1} nie jest liczbą rzeczywistą”.

Kluczowym wyrażeniem jest

NDSolveValue[{x''[t] + x[t] == 0, y''[t] + x[t]^2 y[t] == 0., 
 x[0.] == 1., x'[0.] == 0., y[0.] == 1., y'[0.] == 0.}, 
y, {t, 0, Re[t1]}][Re[t1]] /. {t1 -> Take[Reap[
    NDSolve[{x''[t] + x[t] == 0, x[0.] == 1., x'[0.] == 0., 
      WhenEvent[x'[t] > 0., {Sow[t], "StopIntegration"}]}, 
     x, {t, 0., 100.}, 
     MaxStepSize -> 0.001]], {2, -1}][[1]][[1]][[1]]}

który przychodzi błąd. Jak mogę rozwiązać ten problem?

Answer 1

Dokonaj t1zamiany w s1,s2podstawieniu, grupując w nawiasach:

ff[x0_, y0_] := s1 + s2 /. (
    {s1 :> 
       NDSolveValue[{x''[t] + x[t] == 0., y''[t] + x[t]^2 y[t] == 0., 
          x[0.] == x0, x'[0.] == 0., y[0.] == y0, y'[0.] == 0.}, 
         y, {t, 0, t1}][t1], 
      s2 :> NDSolveValue[{x''[t] + x[t] == 0, 
          y''[t] + x[t]^2 y[t] == 0., x[0.] == x0, x'[0.] == 0., 
          y[0.] == y0, y'[0.] == 1.}, y, {t, 0, t1}][t1]} /. {t1 ->
       NDSolveValue[{x''[t] + x[t] == 0., x[0.] == x0, x'[0.] == 0., 
         WhenEvent[x'[t] > 0., {"StopIntegration"}]}, 
        Indexed[x["Domain"], {1, -1}], {t, 0., 100.}, 
        MaxStepSize -> 0.001]}
    );

ff[-1, -3]

(*  -2.54137  *)

Alternatywny:

ff[x0_, y0_] := 
  Block[{t1 = 
     NDSolveValue[{x''[t] + x[t] == 0., x[0.] == x0, x'[0.] == 0., 
       WhenEvent[x'[t] > 0., {"StopIntegration"}]}, 
      Indexed[x["Domain"], {1, -1}], {t, 0., 100.}, 
      MaxStepSize -> 0.001]},
   s1 + s2 /.
    {s1 :> 
      NDSolveValue[{x''[t] + x[t] == 0., y''[t] + x[t]^2 y[t] == 0., 
         x[0.] == x0, x'[0.] == 0., y[0.] == y0, y'[0.] == 0.}, 
        y, {t, 0, t1}][t1], 
     s2 :> NDSolveValue[{x''[t] + x[t] == 0, 
         y''[t] + x[t]^2 y[t] == 0., x[0.] == x0, x'[0.] == 0., 
         y[0.] == y0, y'[0.] == 1.}, y, {t, 0, t1}][t1]}
   ];