Ah, if I understand your question well, then you’ve run into one of the more elegant but challenging bits of math in finite element analysis: flux extraction.
Essentially, one should not evaluate the fluxes (i.e., traction in the case of elasticity/linear momentum equations) as simply the stress operator acting on the discrete solution. While this operator is “permitted” in volumes, regularity issues do not “permit” this operation on surfaces. While it may seem like the absolutely logical thing to do, in fact, it is not the flux that is balanced in the finite element scheme.
Which fluxes do balance in FE is not straightforward to see. It relates to the question in fluid mechanics “is FEM conservative?”, which has been debated heavily (and misconceptions still persist). The answer is yes; but only w.r.t. the appropriate definition of the flux (see Redirecting)
What you should not do
Evaluate the traction as \boldsymbol{\sigma}(\boldsymbol{u}^h)\cdot \boldsymbol{n}
What you should do
Post process the solution \boldsymbol{u}^h to extract the corresponding flux:
Define a new variational form, find t^h \in V^h s.t. \forall \, v \in V^h :
Where V^h is your discrete space, without enforcement of the Dirichlet conditions. Once solved for t (which nodal valued should only be nonzero for nodes on the boundary), you can use this field in your subsequent computation in place of your t_N
I quite like this paper on this concept: https://doi.org/10.1115/1.4005187
p.s.to lighten the above computation, you can also constrain all nodes that do not live on the boundary with a dirichletbc