Checking Definiteness of a Matrix

Author

Silu Muduli

Published

August 18, 2025

Introduction

In linear algebra, the definiteness of a matrix is a property of symmetric matrices that determines the behavior of the quadratic form \(q(\mathbf{x}) = \mathbf{x}^T A \mathbf{x}\) , where \(( A )\) is a symmetric matrix and \(( \mathbf{x} )\) is a non-zero vector. A matrix can be:

  • Positive Definite: \(q(\mathbf{x}) > 0\) for all \(\mathbf{x} \neq 0\) .
  • Negative Definite: \(q(\mathbf{x}) < 0\) for all \(\mathbf{x} \neq 0\) .
  • Positive Semi-Definite: \(q(\mathbf{x}) \geq 0\) for all \(\mathbf{x}\) .
  • Negative Semi-Definite: \(q(\mathbf{x}) \leq 0\) for all \(\mathbf{x}\) .
  • Indefinite: \(q(\mathbf{x})\) takes both positive and negative values.

This document provides R code to check the definiteness of a matrix using eigenvalues and the quadratic form.

Methodology

To determine the definiteness of a matrix \(A\) :

  1. Ensure the matrix is symmetric: A matrix must be symmetric \(( A = A^T )\) for the quadratic form to be well-defined.
  2. Compute eigenvalues: The eigenvalues of ( A ) determine its definiteness:
    • All eigenvalues positive \(\implies\) Positive Definite.
    • All eigenvalues negative \(\implies\) Negative Definite.
    • All eigenvalues non-negative \(\implies\) Positive Semi-Definite.
    • All eigenvalues non-positive \(\implies\) Negative Semi-Definite.
    • Both positive and negative eigenvalues \(\implies\) Indefinite.
  3. Optional: Evaluate the quadratic form for sample vectors to confirm behavior.

R Code Implementation

We will use the matrixcalc package to check if a matrix is symmetric and compute eigenvalues using base R.

if (!require(matrixcalc)) install.packages("matrixcalc")
Loading required package: matrixcalc
library(matrixcalc)

Algorithm to Check Definiteness

# Function to check matrix definiteness
check_definiteness <- function(A) {
  # Check if matrix is symmetric
  if (!is.symmetric.matrix(A)) {
    return("Matrix is not symmetric")
  }
  
  # Compute eigenvalues
  eigenvalues <- eigen(A)$values
  
  # Remove small numerical errors (treat very small eigenvalues as zero)
  eigenvalues <- round(eigenvalues, 10)
  
  # Check definiteness
  if (all(eigenvalues > 0)) {
    return("Positive Definite")
  } else if (all(eigenvalues < 0)) {
    return("Negative Definite")
  } else if (all(eigenvalues >= 0)) {
    return("Positive Semi-Definite")
  } else if (all(eigenvalues <= 0)) {
    return("Negative Semi-Definite")
  } else {
    return("Indefinite")
  }
}

Examples

Positive Definite Matrix

A <- matrix(c(2, 1, 1, 2), nrow = 2)
cat("Matrix A:\n")
Matrix A:
print(A)
     [,1] [,2]
[1,]    2    1
[2,]    1    2
cat("Eigenvalues:", eigen(A)$values, "\n")
Eigenvalues: 3 1 
cat("Definiteness:", check_definiteness(A), "\n")
Definiteness: Positive Definite 

Matrix of \(n\) entries

n <- 100
A0 <- matrix(rnorm(n * n), nrow = n)

# Make it symmetric: A = (A + A^T)/2
A <- (A0 + t(A0))/2

cat("Eigenvalues:", eigen(A)$values, "\n")
Eigenvalues: 14.24295 13.4121 12.83995 12.07863 11.6172 11.10101 11.00185 10.90218 10.48453 10.02725 9.802502 9.555198 9.327673 9.020054 8.4252 8.237772 8.153419 7.951646 7.795358 7.529827 7.143569 6.955281 6.564925 6.477892 6.30818 5.858358 5.64292 5.52661 5.217162 4.750112 4.598369 4.480842 4.169449 3.93506 3.7726 3.551029 3.291429 3.089458 2.867726 2.735783 2.41935 1.916225 1.74884 1.5271 1.265522 0.8535213 0.561144 0.4338773 0.3803177 0.07963074 -0.0754434 -0.3033257 -0.4634509 -0.8478243 -1.048325 -1.259649 -1.471326 -1.570607 -1.724638 -2.001282 -2.299269 -2.595162 -2.849509 -3.009785 -3.385007 -3.50459 -3.704085 -3.910198 -4.113509 -4.507615 -4.670822 -4.771972 -5.272064 -5.430987 -5.622239 -5.859424 -6.0572 -6.238861 -6.410246 -6.688629 -7.187636 -7.612918 -8.067703 -8.1996 -8.333832 -8.615299 -9.037427 -9.177884 -9.666843 -9.725392 -9.923505 -10.20963 -10.79324 -10.93289 -11.53934 -11.5444 -12.24255 -12.59712 -12.79779 -14.21458 
cat("Definiteness:", check_definiteness(A), "\n")
Definiteness: Indefinite 

Another way

n <- 100
A0 <- matrix(rnorm(n * n), nrow = n)

# Make it symmetric: A =A* A^T
A <- A0%*%t(A0)
cat("Eigenvalues:", eigen(A)$values, "\n")
Eigenvalues: 373.2558 350.9676 323.0334 313.8147 308.0619 296.0021 288.5377 278.8562 270.2092 263.5634 254.53 241.6039 237.6003 229.3329 224.53 215.0039 209.5668 205.232 196.3533 190.8453 184.3748 182.3511 172.6199 168.3618 160.0649 158.5669 153.9689 149.1809 147.9439 138.922 135.5441 131.4529 127.5906 120.6281 115.2475 113.2519 108.2044 106.7045 102.215 101.3788 92.58315 89.0944 84.40527 81.52255 77.34372 76.12581 72.6896 70.64024 70.20734 68.34225 64.4203 61.41703 59.39161 57.17392 54.01327 50.63814 47.03186 45.57693 42.81771 41.2767 39.64845 38.72164 36.39983 34.56768 32.26043 30.45593 28.38024 25.39653 24.73548 22.69498 22.0677 20.17226 18.57292 17.40926 16.53421 14.10781 13.47348 13.01447 12.2304 10.35326 9.265709 8.200884 7.092448 6.160402 5.291743 5.039683 4.205902 3.7516 3.415017 2.960643 2.277607 1.86299 1.262528 0.9350024 0.5750869 0.3113302 0.2254543 0.1933613 0.117922 0.005682838 
cat("Definiteness:", check_definiteness(A), "\n")
Definiteness: Positive Definite 

Hessian Matrix

For a function \(f: \mathbb{R} ^2 \to \mathbb{R}\) defined by

\[ f(x,y)=x^3+yx^{100}+45+y^{23} \]

We can calculate the Hessian matrix at \((x,y)=(1,2)\)

if (!require(numDeriv)) install.packages("numDeriv")
Loading required package: numDeriv
library(numDeriv)

# Define the function
f <- function(x) x[1]^3 + x[1]*x[2]^(100)+45 + x[2]^(23)

# Compute the Hessian at point (1, 2)
point <- c(1, 2)
H <- numDeriv::hessian(f, point)
print(H)
             [,1]         [,2]
[1,] 1.067470e+17 5.522425e+31
[2,] 5.522425e+31 3.075567e+33
cat("Eigenvalues:", eigen(H)$values, "\n")
Eigenvalues: 3.076559e+33 -9.912758e+29 
cat("Definiteness:", check_definiteness(H), "\n")
Definiteness: Indefinite 

Optimal of a Function of Two Variables

For a function \(f: \mathbb{R} ^2 \to \mathbb{R}\) defined by

\[ f(x,y)=x^2+ y^2 -4x-2y+xy \] Its first derivatives are

\[ f_x=2x-4+y \ \ , \ \ f_y=2y-2+x \]

Partial derivatives attend 0 at \((x^*, y^*)=(2,0)\)

We can calculate the Hessian matrix as

if (!require(numDeriv)) install.packages("numDeriv")
if (!require(rootSolve)) install.packages("rootSolve")
Loading required package: rootSolve

Attaching package: 'rootSolve'
The following object is masked from 'package:numDeriv':

    hessian
library(numDeriv)
library(rootSolve)

# Define the function f(x, y) = x^2 + y^2 - 4x - 2y + xy
f <- function(x) x[1]^2 + x[2]^2 - 4*x[1] - 2*x[2] + x[1]*x[2]

# Define the gradient function
grad_f <- function(x) {
  c(2*x[1] - 4+x[2], 2*x[2] - 2+x[1])  # [df/dx, df/dy]
}

# Find critical points by solving gradient = 0
critical_points <- multiroot(f = grad_f, start = c(0, 0))

# Extract the critical point
cp <- critical_points$root
cat("Critical point: (", round(cp[1], 4), ", ", round(cp[2], 4), ")\n", sep = "")
Critical point: (2, 0)
# Compute the Hessian at the critical point
H <- numDeriv::hessian(f, cp)

print(H)
     [,1]     [,2]
[1,]    2 1.000000
[2,]    1 2.000003
cat("Eigenvalues of H:", eigen(H)$values, "\n")
Eigenvalues of H: 3.000001 1.000001 
cat("Definiteness of H:", check_definiteness(H), "\n")
Definiteness of H: Positive Definite 
# Evaluate function at critical point
f_value <- f(cp)
cat("Function value at critical point:", f_value, "\n")
Function value at critical point: -4