theta[t_] := 0.5 (1 + Tanh[100 t])
f[t_] := theta[t]
solw = NDSolve[{w'''[t] + 6 w[t] w'[t] - w'[t] == f[t], w[0] == 0,
    w'[0] == 0, w''[0] == 0}, w, {t, 0, 20}];
wsol[t_] := Evaluate[w[t] /. solw]

solG = ParametricNDSolve[{G'''[t] + 6 G[t] G'[t] - G'[t] == 0,
    G[0] == 0, G'[0] == 0, G''[0] == s}, G, {t, 0, 20}, {s},
   MaxSteps -> Infinity, Method -> {"ParametricCaching" -> None}];
GGreen[t_, s_] := Evaluate[G[s][t] /. solG]

En = 50;
T[n_] := t/En n
GsolSimp[t_, s_, s0_] :=
 s0 t/(3 En) (GGreen[T[0], s] f[t - T[0]] + 2 Sum[GGreen[T[2  n], s] f[
        t - T[2  n]],{n, 1, En/2 - 1}] + 4 Sum[GGreen[T[2  n - 1], s] f[
        t - T[2  n - 1]],{n, 1, En/2}] + GGreen[T[En], s] f[t - T[En]])

NMinimize[{Max[
   Table[Abs[wsol[t] - GsolSimp[t, s, s0]], {t, 0, 5, 0.1}]],
  0 <= s <= 1, -10 <= s0 <= 10}, {s, s0},
 Method -> {"RandomSearch", "SearchPoints" -> 100}]