Phase Portraits of Linear Systems

The Linear System

A planar linear autonomous system has the form \(\dot{\mathbf{x}}=A\mathbf{x}\), where \(A\) is a constant \(2\times2\) matrix and \(\mathbf{x}=(x,y)^T\). The origin is always an equilibrium. The qualitative behaviour of all solutions — whether they spiral in, radiate outward, oscillate, or approach along straight lines — is determined entirely by the eigenvalues of \(A\).

Classification by Eigenvalues

Let \(\tau=\operatorname{tr}A=a+d\), \(\delta=\det A=ad-bc\), and \(\Delta=\tau^2-4\delta\).

Numerical Method

Trajectories are integrated forward and backward from the click point with a choice of one-step methods:

• RK4 — classical explicit Runge–Kutta, 4th order. Good general default for non-stiff systems.
• Tsit5 — Tsitouras 5(4) explicit Runge–Kutta, 6 stages, 5th-order accuracy. Smoother curves than RK4 at the same step size; same stability footprint (no benefit on stiff problems).
• DOP853 — Dormand–Prince 8(5,3) explicit Runge–Kutta, 12 stages, 8th order (Hairer, Nørsett & Wanner). The highest-accuracy explicit choice: at the same step it reaches near machine precision on this linear system, as the accuracy panel shows.
• SDIRK2 — singly-diagonally-implicit Runge–Kutta, 2 stages, 2nd order with \(\gamma = 1-\tfrac{1}{\sqrt{2}}\). The method is L-stable, so on stiff linear systems it stays bounded and usable at step sizes where explicit methods would blow up. It is only 2nd order, so its error is larger than the explicit methods on non-stiff problems, and a coarse step still leaves a visible error across the fast initial transient of a stiff mode — but the solution never diverges.

In Auto mode the page picks SDIRK2 whenever the linear system is "obviously stiff" — both eigenvalues real, the same sign, with a magnitude ratio \(|\lambda_\text{fast}/\lambda_\text{slow}|\ge 50\) — for stable and unstable nodes alike, and falls back to RK4 otherwise. Tracing a trajectory means integrating toward the fixed point (forward for a stable node, backward for an unstable one), and in that direction the fast mode is a stiff decay: a coarse explicit step oscillates or blows up on it — the backward branch of a stiff unstable node is where RK4/Tsit5/DOP853 visibly go haywire — while the L-stable SDIRK2 stays smooth. Centers are integrated separately with a fixed-step explicit scheme to keep closed orbits closed (an implicit method would introduce numerical damping).

Accuracy panel

Because the system is linear with constant coefficients, its solution is known exactly: \(\mathbf{x}(t)=e^{At}\mathbf{x}_0\), where the matrix exponential is evaluated in closed form from the eigenvalues (handling distinct-real, complex, and repeated cases). The Accuracy panel runs every method at a common, deliberately coarse fixed step and reports each one's global error against this exact reference in the normalised state 2-norm \(\max_{t_i\ge T/2}\lVert \mathbf{x}_\text{num}(t_i)-\mathbf{x}_\text{exact}(t_i)\rVert_2 \big/ \max_i\lVert \mathbf{x}_\text{exact}(t_i)\rVert_2\). The error is maximised over the second half of the interval, i.e. after the startup transient: on a stiff system the coarse step cannot resolve the fast mode (that is what makes it stiff), and an L-stable method needs a few steps to damp the resulting spurious component — so measuring from the start would report that unavoidable first-step overshoot rather than what a stiff solver is actually for, staying stable and tracking the slow manifold. On non-stiff systems the ordering of the errors reflects the methods' orders of accuracy — SDIRK2 (2) > RK4 (4) > Tsit5 (5) > DOP853 (8) — with DOP853 typically at machine precision. On a stiff node the picture inverts: the explicit methods cross their stability limit and diverge, while SDIRK2 stays bounded and, once the transient has settled, tracks the dynamics accurately — so the low-order implicit method becomes the only usable one. This holds for both stable and unstable stiff nodes, because the panel integrates toward the fixed point — backward in time for an unstable node, the contracting direction a phase portrait actually traces — so in both cases the fast mode it must resolve is a stiff decay, exactly the regime where the explicit methods oscillate and SDIRK2 does not. The initial state is deliberately chosen to excite both eigenmodes — otherwise a matrix whose slow eigenvector happens to align with \(\mathbf{x}_0\) would never activate its fast mode and would look non-stiff.

How to Use

• Enter the matrix entries \(a,b,c,d\) in the matrix cells. The direction field and classification update automatically after a short pause.
• Select an example from the dropdown to load a preset matrix with strategic trajectories illustrating the classification.
• Click anywhere on the phase plane to integrate a trajectory through that initial condition (forward and backward), or type an exact \((x_0, y_0)\) in the Initial condition box and press Add trajectory (keyboard-friendly). Trajectories are color-coded.
• Use Generate to redraw the field and clear trajectories; Clear trajectories removes curves while keeping the field.
• Set the viewing window with the \(x_\text{min}, x_\text{max}, y_\text{min}, y_\text{max}\) fields (numbers or expressions such as -pi or 2*pi) to zoom in or out on the phase plane.
• Adjust arrows per axis, arrow length, field opacity, and colors. Toggle grid and initial point markers using the checkboxes.
• Open the Time Series panel to plot \(x(t)\) and/or \(y(t)\) for all current trajectories.
• Download a PNG of the current portrait using Download PNG.