Multivariable Gradients, Hessian Matrices & Newton-Raphson Optimization

In real-world computing, machine learning, physics simulation, and engineering optimization, systems rarely depend on a single variable. Neural networks optimize loss functions across billions of parameters, climate simulations track fluid equations across 3D grids, and autonomous vehicles plan trajectories across multidimensional state spaces. Multivariable calculus provides the geometric and analytical framework for navigating, differentiating, and optimizing multi-dimensional scalar and vector fields.

1. Differentiability in Multidimensional Euclidean Spaces#

Extending single-variable calculus to functions of several variables f: R^n → R requires moving beyond simple difference quotients. In single-variable functions, an input can approach a point from only two directions (left or right). In R^n, there exist uncountably many paths and curves approaching a target point.

Partial Derivatives: Axis-Aligned Rates of Change#

Let f: R^n → R be a scalar field. The partial derivative of f with respect to the coordinate variable x_i evaluates the instantaneous rate of change along the i-th standard basis vector while holding all other n - 1 variables fixed as constants:

∂f/∂x_i = lim_{h → 0} [ f(x1, ..., x_i + h, ..., xn) - f(x1, ..., x_i, ..., xn) ] / h.

The Differentiability Paradox: Why Partial Derivatives Are Not Enough#

In single-variable calculus, if f'(c) exists, f is guaranteed to be continuous at c. In multivariable calculus, this guarantee completely fails!

A function can have all partial derivatives exist at a point while failing to even be continuous at that point!

Canonical Counterexample:

Consider the function f: R^2 → R:

f(x, y) = (x · y) / (x^2 + y^2) for (x, y) ≠ (0, 0), and f(0, 0) = 0.

1. Along the coordinate axes:

∂f/∂x(0, 0) = lim_{h → 0} [ f(h, 0) - f(0, 0) ] / h = lim_{h → 0} [ 0 - 0 ] / h = 0.

∂f/∂y(0, 0) = lim_{h → 0} [ f(0, h) - f(0, 0) ] / h = lim_{h → 0} [ 0 - 0 ] / h = 0.

Both partial derivatives exist and equal zero at the origin.

1. Along the diagonal line y = m · x:

f(x, m x) = (x · (m x)) / (x^2 + (m x)^2) = (m · x^2) / (x^2 (1 + m^2)) = m / (1 + m^2).

As (x, y) → (0, 0) along y = x (m = 1), f approaches 1/2.

Along y = -x (m = -1), f approaches -1/2.

Because the limit depends on the direction of approach, the limit lim_{(x, y) → (0, 0)} f(x, y) does not exist! The function is discontinuous at the origin, despite having well-defined partial derivatives.

Total Differentiability (Fréchet Differentiability)#

To ensure true geometric smoothness, multivariable analysis requires total differentiability.

Definition: A scalar field f: R^n → R is differentiable at point x0 in R^n if there exists a linear transformation L: R^n → R (represented by a row vector) such that:

lim_{h → 0} [ |f(x0 + h) - f(x0) - L · h| ] / ||h|| = 0

where h = [ h1, h2, ..., hn ]^T and ||h|| is the Euclidean norm.

If such a linear operator L exists, its components are uniquely given by the partial derivatives:

L = [ ∂f/∂x1, ∂f/∂x2, ..., ∂f/∂xn ].

Theorem (Sufficient Condition for Differentiability):

If all first-order partial derivatives ∂f/∂x_i exist in an open neighborhood around x0 and are continuous at x0 (meaning f is of class C^1), then f is guaranteed to be totally differentiable at x0.

2. The Gradient Vector Field & Directional Derivatives#

When a function f is differentiable, its linear operator can be packaged as a vector in the domain space R^n.

The Gradient Operator#

The gradient vector of f at point x, denoted ∇f(x) (pronounced "del f" or "grad f"), gathers all first-order partial derivatives into a column vector:

∇f(x) = [ ∂f/∂x1, ∂f/∂x2, ..., ∂f/∂xn ]^T.

Directional Derivatives#

To evaluate how f changes along an arbitrary direction in space, let u in R^n be a unit vector (where ||u|| = 1).

The directional derivative of f at point x in direction u is defined as:

D_u f(x) = lim_{t → 0} [ f(x + t · u) - f(x) ] / t.

Theorem: If f is differentiable at x, the directional derivative equals the standard inner product of the gradient vector with the unit direction vector:

D_u f(x) = ∇f(x) · u = ||∇f(x)|| · ||u|| · cos(θ) = ||∇f(x)|| · cos(θ)

where θ is the angle between ∇f(x) and direction vector u.

+-----------------------------------------------------------------------------+
|                      GEOMETRY OF THE GRADIENT VECTOR                        |
+-----------------------------------------------------------------------------+
|                                                                             |
|            Level Curve S: f(x, y) = c                                       |
|                  .---------.                                                |
|                .'           `.   Tangent Vector T (nabla f · T == 0)        |
|               /       P       \  <----------------------------------        |
|              |        *================> Gradient Vector nabla f            |
|               \                / (Points along Steepest Ascent: theta = 0)  |
|                `.            .'                                             |
|                  `----------'                                               |
|                                                                             |
|  Key Invariants:                                                            |
|  1. theta = 0   ===> Direction of STEEPEST ASCENT (maximal increase)        |
|  2. theta = pi  ===> Direction of STEEPEST DESCENT (maximal decrease: -grad)|
|  3. theta = pi/2===> Tangent to Level Set (rate of change is zero: D_u f = 0)|
+-----------------------------------------------------------------------------+

The Three Fundamental Geometric Invariants of the Gradient#

1. Direction of Steepest Ascent:

Because cos(θ) achieves its absolute maximum of +1 when θ = 0, the directional derivative D_u f(x) is maximized when u is aligned with ∇f(x). The gradient points strictly in the direction of greatest instantaneous increase of f.

1. Magnitude as Maximal Rate of Change:

The maximum rate of increase equals the Euclidean magnitude ||∇f(x)||.

1. Orthogonality to Level Sets:

Let S = { x in R^n : f(x) = c } be a level hypersurface. For any smooth curve r(t) lying on S, f(r(t)) = c.

Differentiating both sides using the chain rule:

d/dt [ f(r(t)) ] = ∇f(r(t)) · r'(t) = 0.

Because r'(t) represents an arbitrary tangent vector to the level set, ∇f(x) is strictly perpendicular (orthogonal) to the level hypersurface!

The Tangent Hyperplane Equation#

For a surface defined implicitly by F(x, y, z) = c, the tangent plane at point (x0, y0, z0) has normal vector n = ∇F(x0, y0, z0). Its equation is:

∇F(x0, y0, z0) · [ x - x0, y - y0, z - z0 ]^T = 0

Fx(x0, y0, z0)(x - x0) + Fy(x0, y0, z0)(y - y0) + Fz(x0, y0, z0)(z - z0) = 0.

3. The Jacobian Matrix for Vector-Valued Transformations#

When functions map vector spaces to vector spaces, F: R^n → R^m, the derivative generalizes from a gradient vector into a matrix of partial derivatives.

The Jacobian Matrix#

Let F(x) = [ F1(x), F2(x), ..., Fm(x) ]^T be a differentiable mapping from R^n to R^m.

The Jacobian Matrix J_F(x) in R^(m×n) is defined as:

J_F(x) = [

[ ∂F1/∂x1, ∂F1/∂x2, ..., ∂F1/∂xn ],

[ ∂F2/∂x1, ∂F2/∂x2, ..., ∂F2/∂xn ],

[ ..., ..., ..., ... ],

[ ∂Fm/∂x1, ∂Fm/∂x2, ..., ∂Fm/∂xn ]

]

The i-th row of the Jacobian is the transposed gradient of the i-th component function: (∇F_i)^T.

The Jacobian represents the optimal local linear approximation of the mapping:

F(x + Δx) ≈ F(x) + J_F(x) · Δx.

The Multivariable Chain Rule in Matrix Form#

Let y = g(x) be a mapping from R^n to R^m, and let z = f(y) be a scalar function from R^m to R.

The composite function is h(x) = f(g(x)) from R^n to R.

By the chain rule:

∂h/∂x_j = sum(i=1 to m, (∂f/∂y_i) · (∂y_i/∂x_j)).

In compact matrix-vector notation:

∇_x h(x) = J_g(x)^T · ∇_y f(y).

This matrix-vector equation is the exact mathematical foundation of Reverse-Mode Automatic Differentiation (Backpropagation) in deep learning!

The Jacobian Determinant and Differential Volume Scaling#

When m = n, the mapping F: R^n → R^n is square. The absolute value of the determinant |det(J_F)| measures the local volumetric distortion factor under the transformation.

When changing coordinate systems in multiple integrals:

dx dy = |det(J)| dr dθ = r dr dθ (Polar Coordinates)

dx dy dz = |det(J)| dρ dφ dθ = ρ^2 sin(φ) dρ dφ dθ (Spherical Coordinates).

4. Second-Order Derivatives & The Hessian Curvature Matrix#

While the gradient indicates the slope of a multidimensional landscape, optimizing functions requires analyzing curvature to classify critical points.

The Hessian Matrix#

For a twice continuously differentiable scalar field f: R^n → R (class C^2), the Hessian Matrix H_f(x) in R^(n×n) is the matrix of second-order partial derivatives:

H_f(x) = [

[ ∂^2 f / ∂x1^2, ∂^2 f / (∂x1 ∂x2), ..., ∂^2 f / (∂x1 ∂xn) ],

[ ∂^2 f / (∂x2 ∂x1), ∂^2 f / ∂x2^2, ..., ∂^2 f / (∂x2 ∂xn) ],

[ ..., ..., ..., ... ],

[ ∂^2 f / (∂xn ∂x1), ∂^2 f / (∂xn ∂x2), ..., ∂^2 f / ∂xn^2 ]

]

Clairaut's / Schwarz's Theorem: Symmetry of Mixed Partials#

Theorem: If all second-order partial derivatives of f are continuous on an open domain, then mixed partial derivatives commute:

∂^2 f / (∂x_i ∂x_j) = ∂^2 f / (∂x_j ∂x_i).

Consequence: The Hessian matrix is guaranteed to be a real symmetric matrix (H_f = H_f^T)!

By the Spectral Theorem of linear algebra, all eigenvalues of H_f are purely real numbers, and its eigenvectors form an orthogonal coordinate system defining the principal axes of curvature.

Multivariable Second-Order Taylor Expansion#

Near a reference point x*, the second-order Taylor expansion of f is:

f(x + h) ≈ f(x) + ∇f(x)^T · h + (1/2) · h^T · H_f(x) · h.

At a critical point x where the gradient vanishes (∇f(x) = 0):

f(x + h) - f(x) ≈ (1/2) · h^T · H_f(x*) · h.

The local behavior is governed entirely by the quadratic form q(h) = h^T · H_f · h.

Critical Point Classification via Hessian Eigenvalues#

Let x be a critical point (∇f(x) = 0):

1. Local Minimum:

H_f(x) is Positive-Definite* (all eigenvalues λ_i > 0).

The quadratic term is strictly positive for all h ≠ 0. The surface curves upwards in all directions like a bowl.

1. Local Maximum:

H_f(x) is Negative-Definite* (all eigenvalues λ_i < 0).

The surface curves downwards in all directions.

1. Saddle Point:

H_f(x) is Indefinite* (Hessian has both positive and negative eigenvalues).

Traversing along eigenvectors corresponding to positive eigenvalues increases f, while traversing along eigenvectors with negative eigenvalues decreases f. (The landscape resembles a mountain pass or Pringles potato chip).

1. Inconclusive (Degenerate):

At least one eigenvalue is zero (det(H_f) = 0) and all other eigenvalues have matching signs. Higher-order terms must be examined.

+-----------------------------------------------------------------------------+
|                SECOND-ORDER CRITICAL POINT CLASSIFICATION                   |
+-----------------------------------------------------------------------------+
| Condition at grad(f)(x*) = 0 | Geometric Curvature | Classification         |
+------------------------------+---------------------+------------------------+
| All eigenvalues lambda_i > 0 | Bowl facing upwards | Strict Local Minimum   |
| All eigenvalues lambda_i < 0 | Dome facing down    | Strict Local Maximum   |
| Mixed signs (+ and -)        | Saddle surface      | Saddle Point           |
| At least one lambda_i == 0   | Ridge / Flat trough | Inconclusive Test      |
+------------------------------+---------------------+------------------------+

5. Multivariate Newton-Raphson Optimization#

In numerical optimization, finding the minimum of a cost function f(x) is paramount.

Gradient Descent (First-Order Optimization)#

First-order gradient descent updates parameters using negative gradient steps:

x_{k+1} = x_k - η · ∇f(x_k).

The Pathological Valley Problem:

When the Hessian has a large condition number κ(H) = λ_max / λ_min >> 1 (e.g., long, narrow curved canyons like the Rosenbrock function), gradient descent suffers severely. The gradient vector points almost perpendicular to the canyon floor, causing the optimizer to oscillate wildly back and forth between steep canyon walls while making microscopic progress along the floor.

Multivariate Newton-Raphson Method (Second-Order Optimization)#

Newton's method resolves this by incorporating local curvature. It approximates f around the current iterate x_k with its second-order Taylor polynomial:

m_k(p) = f(x_k) + ∇f(x_k)^T · p + (1/2) · p^T · H_f(x_k) · p.

To find the step p that minimizes this model, differentiate with respect to p and set to zero:

∇_p [ m_k(p) ] = ∇f(x_k) + H_f(x_k) · p = 0

H_f(x_k) · p = -∇f(x_k)

p = -[ H_f(x_k) ]^(-1) · ∇f(x_k).

The Multivariate Newton-Raphson Iteration:

x_{k+1} = x_k - [ H_f(x_k) ]^(-1) · ∇f(x_k).

Convergence Economics & Trade-Offs#

• Convergence Rate: Near a local minimum where H is positive-definite, Newton's method exhibits quadratic convergence:

||x_{k+1} - x|| <= C · ||x_k - x||^2.

The number of accurate decimal digits roughly doubles every iteration!

• Affine Invariance: Unlike gradient descent, which changes behavior under coordinate scaling or rotation, Newton's method is completely invariant to linear coordinate transformations.
• Computational Cost:

Computing the Hessian requires O(n^2) evaluations.

Solving the linear system H · p = -g requires O(n^3) floating-point operations.

For deep neural networks with n = 10^9 parameters, inverting H requires 10^(27) flops—impossible on any existing supercomputer.

Quasi-Newton Methods (BFGS and L-BFGS)#

To capture second-order curvature benefits without computing or inverting H, Quasi-Newton methods (such as BFGS: Broyden-Fletcher-Goldfarb-Shanno) iteratively construct a low-rank approximation B_k ≈ H^(-1) using observed changes in gradient positions:

s_k = x_{k+1} - x_k, y_k = ∇f(x_{k+1}) - ∇f(x_k).

Limited-memory BFGS (L-BFGS) stores only the last m gradient vectors (typically m = 10 to 20), achieving superlinear convergence with O(m · n) memory footprint.

6. Machine Learning, Computer Vision & Robotics Applications#

Multivariable calculus is the primary mathematical engine powering modern AI and robotics.

1. Neural Network Loss Landscapes & Saddle Point Proliferation#

In deep learning, loss functions L(W) exist in spaces with millions or billions of dimensions. A classical misconception was that optimization gets trapped in bad local minima.

In 2014, Dauphin et al. proved mathematically using random matrix theory that in high-dimensional non-convex landscapes:

• Critical points with high loss are overwhelmingly saddle points, NOT local minima!
• For a critical point to be a local minimum, all n eigenvalues of the Hessian must be simultaneously positive. If eigenvalues have random signs with probability 0.5, the probability of encountering a true local minimum scales as (0.5)^n ≈ 10^(-300) for n = 1000!

Optimization algorithms like Adam (Adaptive Moment Estimation) and momentum SGD are specifically engineered to escape saddle points by accumulating velocity along negative curvature directions.

2. Robotics: Inverse Kinematics via the Manipulator Jacobian#

In robotic arms (such as industrial KUKA or ABB welding robots), the forward kinematics mapping relates joint motor angles θ in R^6 to the end-effector gripper 3D position and orientation x in R^6:

x = F(θ).

When the robot must move its gripper along a straight trajectory in physical 3D space, it must compute joint velocities dθ/dt. Differentiating with the chain rule:

dx/dt = J_F(θ) · (dθ/dt)

where J_F(θ) is the 6×6 robotic manipulator Jacobian.

Solving for required motor speeds:

dθ/dt = [ J_F(θ) ]^(-1) · (dx/dt).

At kinematic singular positions (where det(J_F) = 0, such as when the arm is fully extended), the robot loses a degree of freedom and [ J_F ]^(-1) blows up, requiring damped least-squares (Levenberg-Marquardt) regularization.

3. Computer Vision: Optical Flow & Spatial Intensity Gradients#

In autonomous driving and drone tracking, measuring object velocity across video frames uses the Lucas-Kanade Optical Flow algorithm.

Let I(x, y, t) be video pixel intensity. Assuming brightness constancy along motion trajectory (u = dx/dt, v = dy/dt):

dI/dt = (∂I/∂x) · u + (∂I/∂y) · v + (∂I/∂t) = 0.

This is the Optical Flow Constraint Equation:

I_x · u + I_y · v = -I_t.

The spatial image gradients ∇I = [ I_x, I_y ]^T act as the measurement matrix, allowing real-time motion tracking via least-squares over local image patches.

7. Common Architectural Pitfalls & Misconceptions#

• Pitfall 1: Assuming All Low-Gradient Regions Are Minima.

In multi-dimensional spaces, a vanishing gradient ||∇f(x)|| ≈ 0 frequently signifies a saddle point or flat plateau rather than an optimal solution. Always examine second-order curvature or apply exploratory stochastic noise.

• Pitfall 2: Deploying Raw Newton's Method on Non-Convex Surfaces.

If you apply x_{k+1} = x_k - H^(-1) ∇f in a region where the Hessian has negative eigenvalues, Newton's method will happily step uphill toward a local maximum or along a saddle ridge! Production optimization software (e.g., PyTorch, SciPy) always adds a diagonal damping term: (H + μ · I)^(-1) (Levenberg-Marquardt regularization) to guarantee positive-definiteness.

• Pitfall 3: Forgetting Unit Vector Normalization in Directional Derivatives.

The formula D_u f = ∇f · u strictly requires ||u|| = 1. If you compute the directional derivative along v = [ 3, 4 ]^T without dividing by ||v|| = 5, your computed rate of change will be 5 times larger than reality.

• Pitfall 4: Conflating the Direction of the Gradient with the Tangent Vector.

The gradient ∇f is perpendicular (normal) to the contour curve f(x, y) = c. It does NOT point along the contour line.

8. Summary & Checkpoint Exercises#

Solidify your mastery of multivariable gradients, Jacobians, Hessians, and Newton optimization with these four checkpoint problems. Full step-by-step mathematical solutions are provided.

Problem 1 (Gradient, Directional Derivative, and Tangent Hyperplane Construction):

Consider the scalar field f(x, y, z) = x^2 · y - y · z^3 + 2 · x · z.

1. Compute the gradient vector field ∇f(x, y, z).
2. Evaluate ∇f at point P = (1, 2, -1).
3. Compute the directional derivative of f at P in the direction of vector v = [ 2, -1, 2 ]^T.
4. Derive the equation of the tangent plane to the level surface f(x, y, z) = c passing through point P.

Solution:

Part 1: Compute Partial Derivatives.

• ∂f/∂x = 2 · x · y + 2 · z
• ∂f/∂y = x^2 - z^3
• ∂f/∂z = -3 · y · z^2 + 2 · x

The gradient vector field is:

∇f(x, y, z) = [ 2 x y + 2 z, x^2 - z^3, -3 y z^2 + 2 x ]^T.

Part 2: Evaluate at P = (1, 2, -1).

• fx(1, 2, -1) = 2(1)(2) + 2(-1) = 4 - 2 = 2
• fy(1, 2, -1) = (1)^2 - (-1)^3 = 1 - (-1) = 2
• fz(1, 2, -1) = -3(2)(-1)^2 + 2(1) = -6(1) + 2 = -4

Gradient at P:

∇f(1, 2, -1) = [ 2, 2, -4 ]^T.

Part 3: Compute Directional Derivative.

Vector v = [ 2, -1, 2 ]^T.

Compute its Euclidean norm:

||v|| = sqrt(2^2 + (-1)^2 + 2^2) = sqrt(4 + 1 + 4) = sqrt(9) = 3.

Normalize v to a unit vector u:

u = v / ||v|| = [ 2/3, -1/3, 2/3 ]^T.

Compute D_u f(P) = ∇f(P) · u:

D_u f(P) = [ 2, 2, -4 ] · [ 2/3, -1/3, 2/3 ]^T

= (2)(2/3) + (2)(-1/3) + (-4)(2/3)

= 4/3 - 2/3 - 8/3

= -6/3 = -2.

The function decreases at an instantaneous rate of 2 units per distance along direction v.

Part 4: Tangent Plane Equation at P.

The normal vector to the level surface is n = ∇f(1, 2, -1) = [ 2, 2, -4 ]^T.

Equation: n · [ x - 1, y - 2, z - (-1) ]^T = 0

2(x - 1) + 2(y - 2) - 4(z + 1) = 0

2x - 2 + 2y - 4 - 4z - 4 = 0

2x + 2y - 4z - 10 = 0

Divide through by 2:

x + y - 2z = 5.

Problem 2 (Jacobian Matrix and Multivariable Chain Rule Application):

Let coordinate transformation g: R^2 → R^2 be defined by:

y1 = g1(x1, x2) = x1^2 + x2

y2 = g2(x1, x2) = 2 · x1 · x2

Let loss function f: R^2 → R be defined on the output space by:

f(y1, y2) = 3 · y1^2 + 4 · y2.

1. Compute the Jacobian matrix J_g(x).
2. Compute the gradient ∇_y f(y).
3. Using the multivariable chain rule ∇_x h = J_g(x)^T · ∇_y f, evaluate ∇_x h at point x = [ 1, 2 ]^T.

Solution:

Part 1: Compute Jacobian Matrix J_g(x).

g1(x1, x2) = x1^2 + x2

• ∂g1/∂x1 = 2 x1, ∂g1/∂x2 = 1

g2(x1, x2) = 2 x1 x2

• ∂g2/∂x1 = 2 x2, ∂g2/∂x2 = 2 x1

Jacobian matrix J_g(x):

J_g(x) = [

[ 2 x1, 1 ],

[ 2 x2, 2 x1 ]

]

At x = [ 1, 2 ]^T:

J_g(1, 2) = [

[ 2(1), 1 ],

[ 2(2), 2(1) ]

] = [

[ 2, 1 ],

[ 4, 2 ]

]

Part 2: Compute ∇_y f(y).

f(y1, y2) = 3 y1^2 + 4 y2.

• ∂f/∂y1 = 6 y1
• ∂f/∂y2 = 4

At point x = [ 1, 2 ]^T:

y1 = 1^2 + 2 = 3.

y2 = 2(1)(2) = 4.

Evaluating gradient at y = [ 3, 4 ]^T:

∇_y f = [ 6(3), 4 ]^T = [ 18, 4 ]^T.

Part 3: Apply Multivariable Chain Rule.

∇_x h = J_g(1, 2)^T · ∇_y f.

Transpose the Jacobian:

J_g(1, 2)^T = [

[ 2, 4 ],

[ 1, 2 ]

]

Multiply:

∇_x h = [

[ 2, 4 ],

[ 1, 2 ]

] · [

18,

4

]

• Component 1: 2(18) + 4(4) = 36 + 16 = 52.
• Component 2: 1(18) + 2(4) = 18 + 8 = 26.

Result:

∇_x h(1, 2) = [ 52, 26 ]^T.

(Verification: h(x1, x2) = 3(x1^2 + x2)^2 + 4(2 x1 x2) = 3(x1^4 + 2 x1^2 x2 + x2^2) + 8 x1 x2.

∂h/∂x1 = 12 x1^3 + 12 x1 x2 + 8 x2 = 12(1) + 12(2) + 8(2) = 12 + 24 + 16 = 52.

∂h/∂x2 = 6 x1^2 + 6 x2 + 8 x1 = 6(1) + 6(2) + 8(1) = 6 + 12 + 8 = 26. Exact match!).

Problem 3 (Constructing the Hessian Matrix and Classifying Critical Points):

Given the multivariable function:

f(x, y) = 2 x^3 + y^3 - 3 x^2 - 12 x - 3 y.

1. Find all critical points of f(x, y).
2. Construct the general Hessian matrix H_f(x, y).
3. Evaluate the Hessian at each critical point, compute its eigenvalues, and classify the critical point as a local minimum, local maximum, or saddle point.

Solution:

Part 1: Find Critical Points.

Set partial derivatives to zero:

∂f/∂x = 6 x^2 - 6 x - 12 = 0

Divide by 6: x^2 - x - 2 = 0 => (x - 2)(x + 1) = 0 => x = 2, x = -1.

∂f/∂y = 3 y^2 - 3 = 0

Divide by 3: y^2 - 1 = 0 => y = 1, y = -1.

Combining the independent roots yields 4 distinct critical points:

P1 = (2, 1), P2 = (2, -1), P3 = (-1, 1), P4 = (-1, -1).

Part 2: Construct General Hessian Matrix.

Compute second partial derivatives:

• ∂^2 f / ∂x^2 = d/dx [ 6 x^2 - 6 x - 12 ] = 12 x - 6.
• ∂^2 f / (∂x ∂y) = 0.
• ∂^2 f / ∂y^2 = d/dy [ 3 y^2 - 3 ] = 6 y.

H_f(x, y) = [

[ 12 x - 6, 0 ],

[ 0, 6 y ]

]

Because the off-diagonal entries are 0, the eigenvalues of H_f are simply the diagonal entries:

λ1 = 12 x - 6, λ2 = 6 y.

Part 3: Evaluate and Classify Each Critical Point.

1. Point P1 = (2, 1):

H_f(2, 1) = [ [ 12(2) - 6, 0 ], [ 0, 6(1) ] ] = [ [ 18, 0 ], [ 0, 6 ] ].

Eigenvalues: λ1 = 18 > 0, λ2 = 6 > 0.

Hessian is Positive-Definite.

Classification: Strict Local Minimum.

1. Point P2 = (2, -1):

H_f(2, -1) = [ [ 12(2) - 6, 0 ], [ 0, 6(-1) ] ] = [ [ 18, 0 ], [ 0, -6 ] ].

Eigenvalues: λ1 = 18 > 0, λ2 = -6 < 0.

Hessian is Indefinite (mixed signs).

Classification: Saddle Point.

1. Point P3 = (-1, 1):

H_f(-1, 1) = [ [ 12(-1) - 6, 0 ], [ 0, 6(1) ] ] = [ [ -18, 0 ], [ 0, 6 ] ].

Eigenvalues: λ1 = -18 < 0, λ2 = 6 > 0.

Hessian is Indefinite.

Classification: Saddle Point.

1. Point P4 = (-1, -1):

H_f(-1, -1) = [ [ 12(-1) - 6, 0 ], [ 0, 6(-1) ] ] = [ [ -18, 0 ], [ 0, -6 ] ].

Eigenvalues: λ1 = -18 < 0, λ2 = -6 < 0.

Hessian is Negative-Definite.

Classification: Strict Local Maximum.

Problem 4 (Executing Multivariate Newton-Raphson Optimization by Hand):

Consider the quadratic bowl function:

f(x, y) = 2 x^2 + 2 x y + 3 y^2 - 6 x - 8 y.

Starting from initial guess (x0, y0) = (0, 0):

1. Compute the gradient vector ∇f(0, 0).
2. Compute the Hessian matrix H_f(0, 0) and verify that it is positive-definite.
3. Compute the inverse Hessian [ H_f ]^(-1).
4. Perform one multivariate Newton-Raphson step (x1, y1) = (x0, y0) - [ H_f ]^(-1) · ∇f(0, 0).
5. Prove that this single step lands exactly at the global minimum of the function.

Solution:

Part 1: Compute Gradient.

∂f/∂x = 4 x + 2 y - 6

∂f/∂y = 2 x + 6 y - 8

At (x0, y0) = (0, 0):

∇f(0, 0) = [ -6, -8 ]^T.

Part 2: Compute Hessian and Check Definiteness.

∂^2 f / ∂x^2 = 4

∂^2 f / (∂x ∂y) = 2

∂^2 f / ∂y^2 = 6

H_f = [

[ 4, 2 ],

[ 2, 6 ]

]

Note: Because f is a quadratic polynomial, H_f is constant everywhere!

Check definiteness:

• det(H_f) = (4)(6) - (2)(2) = 24 - 4 = 20 > 0.
• Leading principal minor: 4 > 0.
• Trace: 4 + 6 = 10 > 0.

Both eigenvalues are strictly positive. H_f is strictly positive-definite.

Part 3: Compute [ H_f ]^(-1).

The inverse of [ [ a, b ], [ c, d ] ] is (1 / det) · [ [ d, -b ], [ -c, a ] ]:

[ H_f ]^(-1) = (1 / 20) · [

[ 6, -2 ],

[ -2, 4 ]

] = [

[ 0.3, -0.1 ],

[ -0.1, 0.2 ]

]

Part 4: Execute Newton-Raphson Step.

(x1, y1) = (x0, y0) - [ H_f ]^(-1) · ∇f(x0, y0)

[ x1, y1 ]^T = [ 0, 0 ]^T - (1 / 20) · [

[ 6, -2 ],

[ -2, 4 ]

] · [

-6,

-8

]

Multiply matrix by gradient:

• Row 1: 6(-6) + (-2)(-8) = -36 + 16 = -20
• Row 2: -2(-6) + 4(-8) = 12 - 32 = -20

Substitute:

[ x1, y1 ]^T = [ 0, 0 ]^T - (1 / 20) · [ -20, -20 ]^T

= [ 0, 0 ]^T - [ -1, -1 ]^T

= [ 1, 1 ]^T.

Part 5: Verification of Global Minimum.

Evaluate gradient at (x1, y1) = (1, 1):

∂f/∂x(1, 1) = 4(1) + 2(1) - 6 = 4 + 2 - 6 = 0.

∂f/∂y(1, 1) = 2(1) + 6(1) - 8 = 2 + 6 - 8 = 0.

∇f(1, 1) = [ 0, 0 ]^T.

Because H_f is positive-definite everywhere, f is strictly convex.

The point (1, 1) is the unique global minimum of the function.

Newton's method reached the exact global minimum in a single iteration!