How to properly apply NumericQ or s.th. similar here? - NDSolve inside NMaximize

Question

I wrote a function that numerically solves a given Schrodinger equation and returns replacement rules to obtain the time evolution operator. It looks and is used as follows:

solver[H_, time_] := 
    soln = Module[{d, init, eqs, vars, solargs, t, t0, tf},
    d = Dimensions[H][[1]];
    t0 = time[[2]];
    tf = time[[3]];
    t = time[[1]];
    u[t_] := Table[Subscript[u, i, j][t], {i, 1, d}, {j, 1, d}];
    init = Thread[Flatten /@ (u[t0] == IdentityMatrix[d])];
    eqs = Thread[Flatten /@ (I*u'[t] == H.u[t])];
    vars = Flatten[Table[Subscript[u, i, j], {i, 1, d}, {j, 1, d}]];
    solargs = Join[eqs, init];
    Return[NDSolve[solargs, vars, time, InterpolationOrder -> All, Method -> "ExplicitRungeKutta", AccuracyGoal -> 10]]
  ];

(* use like this, e.g. for 2D problem*)
ham[t_]:={{0,t},{t,0}};
solEvol=solver[ham[t],{t,0,20}];
evolOp[t_]=u[t]/.solEvol[[1]];

I do now need to work with the result obtained from my solver inside e.g. NMaximize to determine several parameters. To simplify things I reduced my original code/question to this snippet:

testFun[mat_] := Abs[mat];
uTemp2D[t_] := Table[Subscript[u, i, j], {i, 1, 2}, {j, 1, 2}];
hamParam[t_] = {{a, t}, {t, a}};
optimum = NMaximize[
    {testFun[uTemp2D[20] /. solver[hamParam[t], {t, 0, 20}]], 
    -10 < a < 10}, 
    a, 
    Method -> {"NelderMead"}
  ];

As far as I understand this now evaluates to

NDSolve::ndnum: Encountered non-numerical value for a derivative at t == 0.`.
ReplaceAll::reps: (* here come all entries of eqs in solver function*) is neither a list of replacement rules nor a valid dispatch table, and so cannot be used for replacing."

because NDSolve in solver does get non-numerical values (here a) as input. So the order of evaluation is incorrect. I thought of using something like _?NumericQ for matrices (https://mathematica.stackexchange.com/a/19600) but that is either working nor necessary accoring to a comment by MichaelE2.

Question

So in brief, what is the proper way to define my solver function or maybe the testFun to get results for my variable optimum?

---

Original code snippets (to which MichaelE2's answer refers)

Originally I was intending to execute the following (hParam[t] same as above):

partialTrace[states_, mat_] :=
  Module[{l = Length[states], trace},
   trace = Sum[{states[[i]]}\[Conjugate].mat.states[[i]], {i, 1, l}];
   Return[trace];
 ];

fidelPhase[evol_, Uideal_, states_] := 
  Module[{result}, 
   result = (1/(Length[states])^2)*(Abs[
   partialTrace[states, Uideal\[ConjugateTranspose].evol]])^2; 
   Return[result];];

solver[H_, time_, index_] := 
  soln = Module[{d, init, eqs, vars, solargs, t, t0, tf,meth = index}, 
  d = Dimensions[H][[1]];
  t0 = time[[2]];
  tf = time[[3]];
  t = time[[1]];
  u[t_] := Table[Subscript[u, i, j][t], {i, 1, d}, {j, 1, d}];
  init = Thread[Flatten /@ (u[t0] == IdentityMatrix[d])];
  eqs = Thread[Flatten /@ (I*u'[t] == H.u[t])];
  vars = Flatten[Table[Subscript[u, i, j], {i, 1, d}, {j, 1, d}]];
  solargs = Join[eqs, init];
  Return[
    Which[
      meth == 1, NDSolve[solargs, vars, time, InterpolationOrder -> All, AccuracyGoal -> 10],
      meth == 2, NDSolve[solargs, vars, time, InterpolationOrder -> All, Method -> {"FixedStep", Method -> {"ImplicitRungeKutta", "DifferenceOrder" -> 15}}, StartingStepSize -> 1/300, AccuracyGoal -> 10],
      meth == 3, NDSolve[solargs, vars, time, InterpolationOrder -> All, Method -> "ExplicitRungeKutta", AccuracyGoal -> 10]
     ]
   ];
 ];

optimizer[gateTime_, ham_, ideal_, vars_, range_, states_] :=
  Module[{params},
  params = Map[Flatten, Transpose[{vars, range}]];
  solsOpt = NMaximize[
    fidelPhase[uTemp3D[gateTime]/.solver[ham/.tg -> gateTime, {t,0,gateTime},3][[1]], ideal,states],
  params,
  Method -> {"NelderMead", "Tolerance" -> Sqrt[$MachineEpsilon]}
  ];
  Return[solsOpt];
 ];

uIdeal = {{0, 1, 0}, {1, 0, 0}, {0, 0, 1}};
test = optimizer[
    20, hParam[t], uIdeal, 
    {a, b}, {{-10, 10}, {-10, 10}}, {{1, 0, 0}, {0, 1, 0}}
  ];

which now works with MichaelE2's answer.

Answer 1

Primary fixes: Adding First to fidelPhase and adding t to Subscript[u, i, j][t]. Making an objective function obj that is not evaluated until a and b are numeric may or may not be important. I won't have time to check it out (I've lost track of the original optimizer). The following works and it's not crucial to be so precise in one's fixes.

Clear[optimizer];
optimizer[gateTime_, ham_, ideal_, vars_, range_, states_] := 
 Module[{params, obj},
  obj[v_ /; VectorQ[v, NumericQ]] := 
   First@fidelPhase[
     uTemp3D[gateTime] /. 
      solver[ham /. Thread[vars -> v] /. tg -> gateTime, {t, 0, 
         gateTime}, 3][[1]], ideal, states];
  params = Map[Flatten, Transpose[{vars, range}]];
  solsOpt = 
   NMaximize[obj[vars], params, 
    Method -> {"NelderMead", "Tolerance" -> Sqrt[$MachineEpsilon]}, 
    EvaluationMonitor :> Print];
  Return[solsOpt];]

uTemp3D[t_] = Table[Subscript[u, i, j][t], {i, 1, 3}, {j, 1, 3}];