Vector Calculus

Interactive version: https://gowda.ai/ml-refresh/vector-calculus.html

A scalar field gives a height at each point. Its gradient gives a slope in every input direction. The Hessian tells how those slopes change. A vector-valued function instead gives several outputs, and its Jacobian maps small input changes to output changes.

Function	Derivative	Shape
$f:\mathbb{R}^2\to\mathbb{R}$	$\nabla f$	2-vector
$F:\mathbb{R}^2\to\mathbb{R}^2$	$J_F$	$2\times2$ matrix
$\nabla f:\mathbb{R}^2\to\mathbb{R}^2$	$H_f=J_{\nabla f}$	$2\times2$ matrix

Python runs in your browser through Pyodide. The notebook is also executable in ordinary Jupyter. Numerical arrays use float32; derivatives here are analytic, not automatic differentiation.

import sys
if sys.platform == 'emscripten':
    import piplite
    await piplite.install(['plotly==6.3.1', 'nbformat==5.10.4'])

import numpy as np
import plotly.graph_objects as go
import plotly.io as pio
from plotly.subplots import make_subplots
from IPython.display import display, Math

pio.renderers.default = 'plotly_mimetype'
pio.templates.default = 'plotly_white'
plot_config = {'responsive': False, 'displaylogo': False, 'scrollZoom': True}
print(f'NumPy {np.__version__}; interactive Plotly ready')

1. A scalar function of a vector

The four examples distinguish different geometry:

bowl: $f=x^2+y^2$ . Constant positive curvature.
coupled bowl: $f=x^2+xy+y^2$ . Mixed partials tilt the level curves.
saddle: $f=x^2-y^2$ . Upward curvature in one direction and downward in another.
quartic: $f=x^4/4+y^2/2$ . Curvature changes with position.

The same polynomial coefficients define $f$ , $\nabla f$ , and $H_f$ . If you replace the function with another expression, update its derivatives too. The controls below are ordinary Python variables: change them and rerun the following cells.

example = 'bowl'
position = np.array([0.8, 0.5], dtype=np.float32)
direction_degrees = 30.0
step_size = 0.3

examples = {
    'bowl': (1.0, 0.0, 1.0, 0.0),
    'coupled bowl': (1.0, 1.0, 1.0, 0.0),
    'saddle': (1.0, 0.0, -1.0, 0.0),
    'quartic': (0.0, 0.0, 0.5, 0.25),
}
square_x, mixed, square_y, quartic = examples[example]

def field(point):
    horizontal, vertical = np.asarray(point, dtype=np.float32)
    return square_x * horizontal**2 + mixed * horizontal * vertical + square_y * vertical**2 + quartic * horizontal**4

def gradient(point):
    horizontal, vertical = np.asarray(point, dtype=np.float32)
    return np.array([2 * square_x * horizontal + mixed * vertical + 4 * quartic * horizontal**3,
                     mixed * horizontal + 2 * square_y * vertical], dtype=np.float32)

def hessian(point):
    horizontal, vertical = np.asarray(point, dtype=np.float32)
    return np.array([[2 * square_x + 12 * quartic * horizontal**2, mixed],
                     [mixed, 2 * square_y]], dtype=np.float32)

radians = np.deg2rad(np.float32(direction_degrees))
direction = np.array([np.cos(radians), np.sin(radians)], dtype=np.float32)
display(Math(r'f(x,y)=a x^2+bxy+c y^2+q x^4'))
display(Math(r'\nabla f=\begin{bmatrix}2ax+by+4qx^3\\ bx+2cy\end{bmatrix}'))
display(Math(r'H_f=\begin{bmatrix}2a+12qx^2&b\\b&2c\end{bmatrix}'))
print(f'{example}: a={square_x:g}, b={mixed:g}, c={square_y:g}, q={quartic:g}')
print('p =', position, '  f(p) =', field(position))
print('gradient =', gradient(position))
print('Hessian =\n', hessian(position))

2. Height and gradient

The surface is the scalar value $z=f(x,y)$ . In the top-down view, color is height and each contour has constant height. The gradient is perpendicular to smooth level curves where it is nonzero, pointing toward steepest increase. Its length is scaled in the drawing; the printed vector gives its actual magnitude.

The dashed line passes through $p$ in the unit direction $d$ . The next section takes a vertical slice along that line.

axis = np.linspace(-2, 2, 81, dtype=np.float32)
grid_x, grid_y = np.meshgrid(axis, axis)
heights = field(np.array([grid_x, grid_y], dtype=np.float32))
surface = go.Figure(go.Surface(x=axis, y=axis, z=heights, colorscale='RdBu', colorbar={'title': 'f(x,y)', 'orientation': 'h', 'thickness': 12, 'len': 0.7, 'y': -0.12}, hovertemplate='x=%{x:.2f}<br>y=%{y:.2f}<br>f=%{z:.3f}<extra></extra>'))
surface.add_trace(go.Scatter3d(x=[position[0]], y=[position[1]], z=[field(position)], mode='markers', marker={'color': 'black', 'size': 5}, name='Point p'))
surface.update_layout(title={'text': 'Scalar height: z = f(x, y)', 'font': {'size': 16}}, height=400, margin={'l': 10, 'r': 10, 't': 55, 'b': 75}, scene={'xaxis_title': 'x', 'yaxis_title': 'y', 'zaxis_title': 'f(x, y)', 'aspectmode': 'cube', 'camera': {'eye': {'x': 1.7, 'y': 1.7, 'z': 1.7}}}, showlegend=False)
surface.show(config=plot_config)

contour = go.Figure(go.Contour(x=axis, y=axis, z=heights, colorscale='RdBu', colorbar={'title': 'f', 'thickness': 12}, contours={'showlabels': True}, hovertemplate='x=%{x:.2f}<br>y=%{y:.2f}<br>f=%{z:.3f}<extra></extra>'))
grad = gradient(position)
arrow = grad / max(1, np.linalg.norm(grad))
contour.add_annotation(x=float(position[0] + arrow[0]), y=float(position[1] + arrow[1]), ax=float(position[0]), ay=float(position[1]), axref='x', ayref='y', text='', showarrow=True, arrowhead=3, arrowwidth=3, arrowcolor='black')
cut = position[:, None] + direction[:, None] * np.array([-3, 3], dtype=np.float32)
contour.add_trace(go.Scatter(x=cut[0], y=cut[1], mode='lines', line={'color': '#a96c08', 'dash': 'dash'}, name='Slice direction d'))
contour.add_trace(go.Scatter(x=[position[0]], y=[position[1]], mode='markers', marker={'color': 'black', 'size': 9}, name='Point p'))
contour.update_layout(title={'text': 'Level curves and gradient', 'font': {'size': 16}}, height=460, margin={'l': 45, 'r': 20, 't': 65, 'b': 65}, legend={'orientation': 'h', 'y': -0.2}, xaxis={'title': 'x', 'range': [-2, 2], 'constrain': 'domain'}, yaxis={'title': 'y', 'range': [-2, 2], 'scaleanchor': 'x', 'scaleratio': 1}, dragmode='pan')
contour.show(config=plot_config)

3. Directional slope and curvature

Along a unit direction $d$ , set $g(t)=f(p+td)$ . The scalar $t$ is distance from $p$ along the dashed line.

$g'(0)=\nabla f(p)^T d$

$g''(0)=d^T H_f(p)d$

The tangent and quadratic approximations are

$L(t)=f(p)+t\nabla f(p)^Td$

$Q(t)=L(t)+\tfrac12 t^2d^TH_f(p)d.$

For a quadratic function, $Q$ and the exact slice coincide. For the quartic, they agree only locally. Rotating $d$ changes which slope and curvature you measure.

offsets = np.linspace(-1, 1, 201, dtype=np.float32)
value = field(position)
slope = gradient(position) @ direction
curvature = direction @ hessian(position) @ direction
exact = field(position[:, None] + direction[:, None] * offsets)
linear = value + offsets * slope
quadratic = linear + np.float32(0.5) * offsets**2 * curvature
slice_plot = go.Figure()
for values, name, color, dash in [(quadratic, 'Q(t): quadratic', '#ae426a', 'dot'), (exact, 'g(t): exact', '#087f72', 'solid'), (linear, 'L(t): tangent', '#a96c08', 'dash')]:
    slice_plot.add_trace(go.Scatter(x=offsets, y=values, mode='lines', name=name, line={'color': color, 'dash': dash, 'width': 3}))
slice_plot.add_trace(go.Scatter(x=[0], y=[value], mode='markers', marker={'color': 'black', 'size': 9}, name='Point p'))
slice_plot.update_layout(title={'text': f'Slope = {slope:.3f}<br>Curvature = {curvature:.3f}', 'font': {'size': 16}}, height=460, margin={'l': 55, 'r': 15, 't': 80, 'b': 110}, xaxis_title='t: displacement along d', yaxis_title='Height', hovermode='x unified', legend={'orientation': 'h', 'y': -0.3})
slice_plot.show(config=plot_config)

4. What the Hessian tells us

Writing $f_{xy}$ for a mixed second partial derivative,

$H_f=\begin{bmatrix}f_{xx}&f_{xy}\\f_{yx}&f_{yy}\end{bmatrix}.$

The first-order change in the gradient is

$\begin{aligned}\nabla f(p+\Delta p)-\nabla f(p)\\\approx H_f(p)\Delta p.\end{aligned}$

The Hessian is the Jacobian of the gradient. For continuous second partial derivatives it is symmetric. Off-diagonal entries mean that moving one coordinate changes the slope in another.

At a stationary point ( $\nabla f=0$ ): positive eigenvalues imply a strict local minimum; negative eigenvalues imply a strict local maximum; mixed signs imply a saddle. A zero eigenvalue makes this test inconclusive. Away from a stationary point, the signs describe curvature, not an optimum.

At the saddle's origin the gradient is zero, but curvature along x is +2 and along y is -2. For the quartic at the origin, the Hessian has eigenvalues 0 and 1, yet the higher-order term gives a strict minimum.

print('Hessian at p:\n', hessian(position))
print('Eigenvalues:', np.linalg.eigvalsh(hessian(position)))
print('Gradient at p:', gradient(position))
for name, unit_direction in [('x', [1, 0]), ('y', [0, 1]), ('diagonal', [1, 1])]:
    unit_direction = np.array(unit_direction, dtype=np.float32)
    unit_direction /= np.linalg.norm(unit_direction)
    print(name, 'directional curvature:', unit_direction @ hessian(position) @ unit_direction)

5. A vector-valued function and its Jacobian

This is a different function, with two outputs rather than one height:

$F(x,y)=\begin{bmatrix}x^2-y\\xy\end{bmatrix}$

$J_F(x,y)=\begin{bmatrix}2x&-1\\y&x\end{bmatrix}.$

Each row is the gradient of one output; each column gives the output response to one input coordinate. The Jacobian need not be symmetric.

Actual output change:

$\Delta F=F(p+\Delta p)-F(p).$

Linear prediction:

$\Delta F\approx J_F(p)\Delta p.$

A circle of input displacements becomes an ellipse under a linear map, possibly collapsing to a line or point. The nonlinear image deviates from that ellipse. Both panels below use the same coordinate scale; arrows track the selected displacement $\Delta p=rd$ . Decrease $r$ to see the relative error shrink.

For a composition $G(F(p))$ , the chain rule is $J_{G\circ F}=J_G(F(p))J_F(p)$ . For a scalar loss, backpropagation applies $J_F^T$ to the output gradient.

def vector_map(point):
    horizontal, vertical = np.asarray(point, dtype=np.float32)
    return np.array([horizontal**2 - vertical, horizontal * vertical], dtype=np.float32)

def jacobian(point):
    horizontal, vertical = np.asarray(point, dtype=np.float32)
    return np.array([[2 * horizontal, -1], [vertical, horizontal]], dtype=np.float32)

assert np.isfinite(step_size) and step_size > 0, 'step_size must be positive and finite'
angles = np.linspace(0, 2 * np.pi, 129, dtype=np.float32)
circle = np.float32(step_size) * np.array([np.cos(angles), np.sin(angles)], dtype=np.float32)
mapped = vector_map(position[:, None] + circle) - vector_map(position)[:, None]
predicted = jacobian(position) @ circle
displacement = np.float32(step_size) * direction
actual_change = vector_map(position + displacement) - vector_map(position)
predicted_change = jacobian(position) @ displacement
bound = float(1.15 * max(np.abs(circle).max(), np.abs(mapped).max(), np.abs(predicted).max()))
mapping = make_subplots(rows=2, cols=1, subplot_titles=['Input displacement', 'Output displacement'], vertical_spacing=0.18)
for points, name, color, dash, row in [(circle, 'Input circle', '#087f72', 'solid', 1), (mapped, 'Exact change', '#087f72', 'solid', 2), (predicted, 'J times displacement', '#ae426a', 'dash', 2)]:
    mapping.add_trace(go.Scatter(x=points[0], y=points[1], mode='lines', name=name, line={'color': color, 'dash': dash}), row=row, col=1)
for change, color, row in [(displacement, '#087f72', 1), (actual_change, '#087f72', 2), (predicted_change, '#ae426a', 2)]:
    mapping.add_annotation(x=float(change[0]), y=float(change[1]), ax=0, ay=0, axref='x' if row == 1 else 'x2', ayref='y' if row == 1 else 'y2', text='', showarrow=True, arrowhead=3, arrowwidth=2, arrowcolor=color, row=row, col=1)
for row, labels in [(1, ('dx', 'dy')), (2, ('dF1', 'dF2'))]:
    mapping.update_xaxes(title_text=labels[0], range=[-bound, bound], constrain='domain', matches='x' if row == 2 else None, row=row, col=1)
    mapping.update_yaxes(title_text=labels[1], range=[-bound, bound], scaleanchor='x' if row == 1 else 'x2', scaleratio=1, matches='y' if row == 2 else None, row=row, col=1)
mapping.update_layout(height=720, margin={'l': 50, 'r': 15, 't': 55, 'b': 100}, legend={'orientation': 'h', 'y': -0.15}, dragmode='pan')
mapping.show(config=plot_config)
print('J(p) =\n', jacobian(position))
print('Actual change:', actual_change)
print('Linear prediction:', predicted_change)
print('Error norm:', np.linalg.norm(actual_change - predicted_change))

6. Experiments and derivative checks

Set the saddle's position to the origin. Compare directions 0, 45, and 90 degrees. Why can directional curvature be zero while the Hessian is nonzero?
Choose the coupled bowl. Which directions have the largest and smallest curvature? Compare with its Hessian eigenvalues, 1 and 3.
Choose the quartic at the origin. Why does the second-order approximation miss the growth along x?
Halve step_size in the vector example. The leading error is quadratic in the step: expect roughly one quarter of the previous absolute error.

Central differences check the selected analytic derivatives below. These are approximate numerical checks, not a symbolic proof.

epsilon = np.float32(0.002)
for point in [position, np.zeros(2, dtype=np.float32)]:
    for coordinate in range(2):
        delta = np.eye(2, dtype=np.float32)[coordinate] * epsilon
        np.testing.assert_allclose(gradient(point)[coordinate], (field(point + delta) - field(point - delta)) / (2 * epsilon), atol=2e-4, rtol=5e-4)
        np.testing.assert_allclose(hessian(point)[:, coordinate], (gradient(point + delta) - gradient(point - delta)) / (2 * epsilon), atol=2e-4, rtol=5e-4)
        np.testing.assert_allclose(jacobian(point)[:, coordinate], (vector_map(point + delta) - vector_map(point - delta)) / (2 * epsilon), atol=2e-4, rtol=5e-4)
print('Gradient, Hessian, and Jacobian checks passed.')