\[Theta] = \[Pi]/4; \[Epsilon]r = FullSimplify@RotationMatrix[\[Theta], {1, 0, 0}].DiagonalMatrix[{\[Epsilon]1, \[Epsilon]1, \[Epsilon]2}].RotationMatrix[-\[Theta], {1, 0, 0}]; eq = Table[    (D[{Ex[x], Ey[x], Ez[x]}, {x, 2}])[[j]] == -(\[Epsilon]r.{Ex[x], Ey[x], Ez[x]})[[j]], {j, 1, 3, 1}]; sol = FullSimplify[DSolveValue[Flatten[{eq, Ex[0] == 0, Ex'[0] == 0, Ey[0] == 0, Ey'[0] == 0, Ez[0] == 1, Ez'[0] == 0}], {Ex[x], Ey[x], Ez[x]}, x] ]; p1 = Table[ Show[     Graphics3D[{Lighting -> "Neutral", Gray, Opacity[0.3], Cuboid[{0, -1, -1}, {62.8, 1, 1}], Opacity[1], Thick, Black, Arrowheads[0.002], Arrow[{{0, -1, -1}, {0, 1, 1}}]}]     ,     ParametricPlot3D[Piecewise[{{{x, 0, Cos[x - \[Tau]]} /. {\[Epsilon]1 -> 1, \[Epsilon]2 -> 1.1}, x < 0}, {{x, sol2[[2]], sol2[[3]]} /. {\[Epsilon]1 -> 1, \[Epsilon]2 -> 1.1}, 0 < x < 62.8}, {{x, Cos[x - \[Tau]], 0} /. {\[Epsilon]1 -> 1, \[Epsilon]2 -> 1.1}, x > 62.8}}], {x, -20, 80}, BoxRatios -> {3, 1, 1}, PlotStyle -> {Orange}, PlotPoints -> 50]     , PlotRange -> {{-20, 80}, {-1, 1}, {-1, 1}}, BoxRatios -> {3, 1, 1}, Boxed -> False     ]    , {\[Tau], 0, 2 \[Pi], (2 \[Pi])/30}]; ListAnimate[p1]