Got a function depending on single parameter "t" defined below which returns a 2 by 2 matrix

H[t_] := PauliMatrix[1] (3 - Cos[Pi t]) + PauliMatrix[3];
U[t_] := Total[Map[With[{x = Normalize[#]}, TensorProduct[x, x /. t -> 0]] &, Eigenvectors[H[t]]], 1];

When I call the function the result gets rounded to the identity

 {{1., 0.}, {0., 1.}}

But if I use a replacement I get the answer to a suitable precision

 U[t] /. t -> 0.2
 {{0.999843, -0.0177345}, {0.0177345, 0.999843}}

How do I change this behaviour?


This is not a precision problem. You can check that with

U[2/10] - U[0.2] // N

{{1.11022*10^-16, 5.55112*10^-17}, {5.55112*10^-17, 1.66533*10^-16}}

Eigenvectors simply does not (ortho-)normalize eigenvectors for symbolic input. Mostly because this would blow-up the complexity of the returned expressions even more.

  • $\begingroup$ Im sorry I don't see what your example illustrates. Why does U[2/10] give anything different to U[0.2]? I know Eigenvectors doesn't orthonormalise for symbolic input, which is why I'm forcing it to normalise the symbolic expression? If I evaluate the RHS side with "t' left symbolic, and then define the function from that expression it behaves as I want it to. How do you suggest I write the function? $\endgroup$
    – oweydd
    Nov 21 '18 at 13:10
  • $\begingroup$ Up to machine precision, U[2/10] does not give anything different from U[0.2]. $\endgroup$ Nov 21 '18 at 13:24
  • 2
    $\begingroup$ I don't see how you orthonormalize the eigenvectors. Did it occur to you that the replacement rule t -> 0 does nothing if you call U with anything else but t? Better use Orthogonalize. $\endgroup$ Nov 21 '18 at 13:24
  • $\begingroup$ The eigenvectors are automatically orthogonal, so I don't need to worry about that. The replacement rule does in fact work if I call U with any symbolic input. None of this answers my question. $\endgroup$
    – oweydd
    Nov 21 '18 at 13:37
  • $\begingroup$ Yes you are right: The replacement rule does something also for other symbols. But it does not what you expect for numerical input. Here is something to mull over: f[t_] := {1, 2, 3, t} /. t -> 0; f /@ Range[3] f /@ N@Range[3]. $\endgroup$ Nov 21 '18 at 13:48

The problem with the original code is an evaluation problem. @Henrik points this out in a comment.

The rule t -> 0 becomes 0.2 -> 0 in the call U[0.2]:

Hold[U[0.2]] /. DownValues@U
  Hold[Total[(With[{x$ = Normalize[#1]}, 
        x$ \[TensorProduct] (x$ /. 0.2 -> 0)] &) /@ Eigenvectors[H[0.2]], 1]]

However the rule is applied to the eigenvectors of H[0.2], which do not contain a 0.2 to replace:

(*  {{0.540734, -0.841193}, {-0.841193, -0.540734}}  *)

It's different in the call U[t], in which the rule is t -> 0 and the eigenvectors of H[t] contain a t to replace with zero.

One way around the difficulty is to use a different symbol for t in H[t]. Formal symbols cannot normally be assigned a value, so \[FormalT] is a good candidate:

UU[t_] := With[{v = Normalize /@ Eigenvectors[H[\[FormalT]]]},
     TensorProduct, {v /. \[FormalT] -> t, v /. \[FormalT] -> 0}], 1]

Since Eigenvectors[H[t]] is always the same, one might want to precompute the formulas for U[t] and use them to define the function U. It's a good idea to protect t with Block, in case it has been assigned a value. Also, it turns out formal symbols can be assigned a value by a user with Block[{\[FormalT] = 5}, <code>], so it would be even safer to use Module to create a temporary symbol for H[t]:

  UU[t_] = Module[{$t, v},
v = Normalize /@ Eigenvectors[H[$t]];
    Total[MapThread[TensorProduct, {v /. $t -> t, v /. $t -> 0}], 1]

Defining the function as

 U[t_] := Evaluate[Total[Map[With[{x = Normalize[#]}, TensorProduct[x,x /. t -> 0]] &, Eigenvectors[H[t]]], 1]];

gives the behaviour I want. But I still don't understand why the original definition did't work.

  • $\begingroup$ If this resolves your problem, then you should learn about the difference between Set (=) and SetDelayed (:=)... $\endgroup$ Nov 21 '18 at 13:45

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.