10 Introduction
A classical problem in minimization is the brachistrochrone, i.e., find the curve (meaning determining the fastest trajectory) to move from a point A to a (lower) point B considering that the only action is the one of the gravity.
Code
# Parameters for the cycloid (optimal curve)
r <- 1
theta <- seq(0, pi, length.out = 500)
# Parametric coordinates of the cycloid
x_cycloid <- r * (theta - sin(theta))
y_cycloid <- -r * (1 - cos(theta))
# Coordinates of the cycloid's final point
x_end <- max(x_cycloid)
y_end <- min(y_cycloid)
# Straight line between (0,0) and (x_end, y_end)
x_line <- seq(0, x_end, length.out = 500)
y_line <- seq(0, y_end, length.out = 500)
# Parabola passing through (0,0) and (x_end, y_end)
a_parabola <- y_end / (x_end^2) # Coefficient for the parabola
x_parabola <- seq(0, x_end, length.out = 500)
y_parabola <- a_parabola * x_parabola^2 # Correct parabola equation
# Quadratic Bézier curve passing through (0,0) and (x_end, y_end)
# We choose a control point at 1.5 times the height of the endpoint
control_x <- x_end / 2
control_y <- 1.5 * y_end
# Parameters for the Bézier curve
t <- seq(0, 1, length.out = 500)
x_bezier <- (1 - t)^2 * 0 + 2 * (1 - t) * t * control_x + t^2 * x_end
y_bezier <- (1 - t)^2 * 0 + 2 * (1 - t) * t * control_y + t^2 * y_end
# Create data frames for ggplot
data_cycloid <- data.frame(x = x_cycloid, y = y_cycloid, curve = "Cycloid (optimal)")
data_line <- data.frame(x = x_line, y = y_line, curve = "Straight line")
data_parabola <- data.frame(x = x_parabola, y = y_parabola, curve = "Parabola")
data_bezier <- data.frame(x = x_bezier, y = y_bezier, curve = "Bézier spline")
# Combine all data into a single data frame
data_all <- rbind(data_cycloid, data_line, data_parabola, data_bezier)
# Create the plot
ggplot(data_all, aes(x = x, y = y, color = curve)) +
geom_path(linewidth = 1.2) +
labs(title = "Brachistochrone Problem: Curve Comparison",
subtitle = "Cycloid (optimal curve) vs Straight Line vs Parabola vs Bézier Spline",
x = "x-axis",
y = "y-axis") +
theme_minimal() +
scale_color_manual(values = c("Cycloid (optimal)" = "blue",
"Straight line" = "red",
"Parabola" = "green",
"Bézier spline" = "purple")) +
coord_fixed() + # Maintains correct aspect ratio
# Add points using annotate
annotate("point", x = 0, y = 0, color = "black", size = 4) + # Start point
annotate("point", x = x_end, y = y_end, color = "black", size = 4) # End pointMathematically the problem is:
Definition 10.1 (Brachistrochrone problem) Determine the curve
\mathcal C: [a,b] \to\Bbb{R}
such that
\mathcal C(a) = y_a \quad\textrm{and}\quad \mathcal C(b) = y_b
that minimize the function T(\mathcal{C}) which represent the time to travel from A=(a,y_a) to B=(b,y_b). Thus, given T(\mathcal{C}) the problem become:
\textrm{minimize } T(\mathcal{C}) \qquad \textrm{for all possible curves } \mathcal{C}
A function takes a number as input and returns a number, like f(x) = x \mathrm{e}^x. A functional is instead something that takes as input a function and returns a number and an example is
\mathcal{F}(x) = \int_a^b x(t)\,\mathrm{d}t
where x is a function. Considering for example x(t) = t^2 the previous functional becomes
\mathcal{F}(x) = \int_a^b t^2\,\mathrm{d}t = \frac{t^3}{3}\Big|_a^b= \frac{b^3-a^3}{3}
Another example of functional can be
\mathcal{G}(z) = z'(0) \int_0^1 z^2(t)\,\mathrm{d}t
and this expression can be evaluated for every generic function z(t).
The question connected to Definition 10.1 is how to define a minimum for a functional; in order to do so we have to firstly understand how to minimize a simple function and in particular given a function
f: A\subseteq\Bbb{R}\mapsto\Bbb{R}
has a minimum in the point x^\star if
f(x^\star) \leq f(x) \qquad \forall x \in A
Similarly given a function \mathcal{F}(x) a point x^\star (that in reality is a function) is a minimum if
\mathcal{F}(x^\star) \leq \mathcal{F}(x) \qquad \forall x \in ?
In this case we have to specify the class (functional space) of functions we are considering, as example x can be a function is continuous in the domain [a,b] and so we can define
x\in C([a,b]) = \big\{ g:[a,b]\mapsto R \ | \ g \textrm{ is continuous} \big\}
In general changing the domain of the functional \mathcal{F} may change the problem.
Example 10.1 (domain change) Considering the function f(x) = x^2+2, it’s roots can be computed if we allow complex solutions z\in\Bbb{C} (and in fact z = \pm \mathrm{i}\sqrt 2), while if the domain of the solution is the real set z\in\Bbb{R} no solutions exists.
Considering now the functional \mathcal{F} defined as
\mathcal{F}(x) = \int_{-1}^1 \big(x(t) - |t|\big)^2\, \mathrm{d}t
The |t| introduce a cuspid in t=0 that cannot be derived. If we minimize the functional in the continuous interval x\in C([-1,1]), the solution is the function x^\star(t) = |t|, in fact
\mathcal{F}(x^\star) = \int_{-1}^1 \big(|t|-|t|\big)^2 \, \mathrm{d}t = \int_{-1}^1 0\, \mathrm{d}t = 0
Choosing any other continuous function will result in a functional with a positive value.
Considering now to minimize the function in the domain of functions with continuous derivatives, and so x\in C^1([-1,1]). In this case |t|\notin C^1([-1,1]) (due to the cuspid). An approximation of the function |t| that is continuous with also continuous first derivative is the function
x_\varepsilon(t) = \begin{cases} t \qquad & t \geq \varepsilon \\ \frac 1 2 \left(\dfrac{t^2}\varepsilon + \varepsilon\right) \qquad & -\varepsilon < t< \varepsilon \\ -t & t \leq- \varepsilon \end{cases}
By pushing the limit \varepsilon\to 0 we can have the function that minimize the functional \mathcal{F} in the domain C^1([-1,1]). In fact we can see that the functional of x_\varepsilon becomes
\mathcal{F}(x_\varepsilon) = \int_{-\varepsilon}^\varepsilon \left( \frac{t^2}{2\varepsilon} + \frac \varepsilon 2 - |t| \right)^2\, \mathrm{d}t = \frac{\varepsilon^3}{10}
Some examples of space of functions
- Functions f:[a,b]\mapsto\Bbb{R} continuous
- Functions f:[a,b]\mapsto\Bbb{R} which are C^k (continuous with derivative up to order k).
- Functions f:[a,b]\mapsto\Bbb{R} piecewise continuos
Less trivial space of functions examples
Space L_1([a,b]) of functions f:[a,b]\mapsto\Bbb{R} which are module integrable
\int_a^b | f(x) |\,\mathrm{d}t < \infty
[the integral must exists and is finite]
Space L_2([a,b]) of functions f:[a,b]\mapsto\Bbb{R} which are square integrable
\int_a^b | f(x) |^2\,\mathrm{d}t < \infty
[the integral must exists and is finite]
Example 10.2 Let f:[0,1]\mapsto\Bbb{R}, f(x)=1/x is not L_1([0,1]) in fact
\int_0^1 \dfrac{1}{x}\mathrm{d}x \to \infty
The function g:[0,1]\mapsto\Bbb{R}, g(x)=1/\sqrt{x}, is in L_1([0,1]) in fact
\int_0^1 \dfrac{1}{\sqrt{x}}\mathrm{d}x = \left[2\sqrt{|x|}\right]_{x=0}^{x=1} =2<\infty
however this function g(x) is not in L_2([0,1]), in fact
\int_0^1 g(x)^2\mathrm{d}x \int_0^1 \dfrac{1}{(\sqrt{x})^2}\mathrm{d}x= \int_0^1 \dfrac{1}{|x|}\mathrm{d}x \to \infty
the function h(x)=x^{-1/4} is in L_2([0,1]), in fact
\begin{aligned} \int_0^1 h(x)^2\mathrm{d}x &= \int_0^1 (x^{-1/4})2\mathrm{d}x \\ &= \int_0^1 \dfrac{1}{\sqrt{x}}\mathrm{d}x\\ &= \left[2\sqrt{|x|}\right]_{x=0}^{x=1} \\ &=2<\infty \end{aligned}
10.1 Analogies with linear algebra
To solve the problem of minimizing a functional we can see some relations with the linear algebra. The domain of the functional can be in fact see as a functional space \Bbb{V} analogous to the vectorial one which the vectors are represented by the functions and the scalars are represented by real values. With this definition example of functional spaces might be
\begin{aligned} \Bbb{V}_1 & = \big\{ f : [0,1] \mapsto\Bbb{R}\textrm{ continuous}\big\} \\ \Bbb{V}_2 & = \big\{ f: [a,b]\mapsto\Bbb{R}\textrm{ such that } f\in C^k([a,b]) \big\} \end{aligned}
with a,b \in \Bbb{R} and k\in \Bbb{N}.
Given in fact two function f,g \in \Bbb{V} that are member of the same functional space, we can see that each linear combination of the functions determine a function h that’s still in the vectorial space:
h(x) = \alpha f(x) + \beta g(x) \in \Bbb{V} \qquad \forall \alpha,\beta \in \Bbb{R}, \ f,g\in\Bbb{V}
As example let’s consider the two continuous function
f(x) = x^2 \quad\textrm{and}\quad g(x) = \sin x
then we can clearly see that the function
h(x) = 2f(x) + \frac 13 g(x) = 2x^2 +\frac 13\sin x
is still continuous. This in general means that the functional space is closed respect to the operations of function summation and multiplication by a scalar.
Scalar product
In linear algebra given two vector \bm{v}_1, \bm{v}_2 it exists the bilinear operator scalar product \langle \bm{v}_1,\bm{v}_2 \rangle = \bm{v}_1 \cdot \bm{v}_2 that satisfy the following rules:
\langle \alpha \bm{v} + \beta \bm{w}, \bm{z}\rangle = \langle \bm{z} , \alpha \bm{v} + \beta \bm{w}\rangle = \alpha \langle \bm{v},\bm{z}\rangle + \beta \langle \bm{w},\bm{z}\rangle
for all \alpha,\beta \in \Bbb{R}, \bm{v},\bm{w},\bm{z} \in \Bbb{V} and
\langle \bm{v},\bm{v}\rangle \geq 0 \quad \textrm{ and } \quad \langle \bm{v},\bm{v} \rangle = 0 \quad \Leftrightarrow \quad \bm{v} = \bm{0}
Also in the functional space \Bbb{V} can exists definitions of product scalar such the one here presented:
\langle f,g\rangle = \int_0^1 f(x)g(x)\, \mathrm{d}x \tag{10.1}
Note
This is not the lonely function that can serve as product scalar of function and in this case the relation holds for the functional space
\Bbb{V} = \{ f:[0,1] \mapsto \Bbb{R} \textrm{ integrable} \}
In this case we can prove that this definitions meets the requirements stated for the scalar product considering the linear properties of the integrals as follows:
\begin{aligned} \langle \alpha f(\cdot) + \beta g(\cdot),h \rangle & = \int_0^1 \big(\alpha f(x) + \beta g(x)\big)h(x) \, \mathrm{d}x \\ & = \alpha \int_0^1 f(x)h(g) \, \mathrm{d}x + \beta \int_0^2 g(x)h(x) \, \mathrm{d}x \\ & = \alpha \langle f,h \rangle + \beta \langle g,h \rangle \end{aligned}
and also
\langle f,f \rangle = \int_0^1 f(x) f(x) \, \mathrm{d}x = \int_0^1 f^2(x) \geq 0
and in particular we can see that the scalar product \langle f,f\rangle will give as result 0 if and only if the function f is identically null, so such that f(x)=0 for all x in it’s domain (in this case [0,1]).
Norm
In linear algebra it’s also defined the norm of a vector as
\|\bm{v}\| := \sqrt{\bm{v} \cdot \bm{v}}
and this expression is used to compute, from a vector, a single positive value (and in particular \|\bm{v}\| = 0 if and only if \bm{v} = \bm{0}). This definition also holds for the functional space and considering the scalar product defined (as example) in equation Equation 10.1 we can see that one definition of norm for function can be the one
\|f(x)\| = \sqrt{\langle f(x),f(x)\rangle} = \sqrt{\int_0^1f^2(x) \, \mathrm{d}x}
In this case we can clearly see that the norm \|f\| \geq 0 is always positive defined (and is zero only in the case on which the function f is identically null).
Example 10.3 (Spaces of Functions) A space of functions can be thought of as a set of functions satisfying certain properties.
For example, one such space is the set of continuous functions f:[a,b]\to \mathbb{R}, \quad \text{with } f \in C([a,b]).
More generally, we may consider the space of k-times continuously differentiable functions f \in C^k([a,b]), \quad k \in \mathbb{N}.
Similarly, we can work with spaces of piecewise continuous functions.
A less trivial example is the space of integrable functions, namely f:[a,b]\to \mathbb{R}, \quad \text{such that } \int_a^b |f(x)|\, dx < \infty. This is precisely the space of functions that are absolutely integrable on [a,b].
In general we define as L^p([a,b]) the space of p-integrable functions, so such that
f\in L^p([a,b]) \qquad \iff \qquad \int_a^b |f(x)|^p\, \mathrm{d}x < \infty
10.2 First Variation
Let \mathcal{J}:\mathbb{V}\to\mathbb{R} be a functional defined on a function space \mathbb{V}.
To study conditions for \mathcal{J} to attain a minimum, we introduce the notion of the derivative of a functional, called the first variation (or directional derivative).
Given a function x \in \mathbb{V} and a perturbation \eta \in \mathbb{V}, we consider the variation
x_\alpha(\cdot) = x(\cdot) + \alpha \eta(\cdot), \qquad \alpha \in \mathbb{R}.
Defining the auxiliary function
g(\alpha) = \mathcal{J}\!\big(x(\cdot) + \alpha \eta(\cdot)\big).
then the first variation of \mathcal{J} at x in the direction \eta is
\delta \mathcal{J}(x;\eta) = g'(0). \tag{10.2}
Equivalent Forms
From the definition of the derivative, we can write \delta \mathcal{J}(x;\eta) = \lim_{\alpha\to 0} \frac{\mathcal{J}(x + \alpha \eta) - \mathcal{J}(x)}{\alpha}.
Equivalently, using the small-o notation, \mathcal{J}(x + \alpha \eta) = \mathcal{J}(x) + \alpha\,\delta\mathcal{J}(x;\eta) + o(\alpha).
Example 10.4 (directional derivative of a functional — first variation) Given the functional
\mathcal{F}(x) = \int_0^1 \big(x^2(t) + 1\big) \, \mathrm{d}t \ + x(1)
it’s first variation respect to a generic direction d(t) can be calculated by firstly determining the associated function g:\Bbb{R}\mapsto \Bbb{R} defined as
\begin{aligned} g(\alpha) & = \mathcal{F}\big(x(\cdot) + \alpha\, d(\cdot)\big) \\ & = \int_0^1 \big(x(t) + \alpha\, d(t)\big)^2 + 1\, \mathrm{d}t + x(1) + \alpha\, d(1) \end{aligned}
As reported in equation Equation 10.2, the first variation of the functional \mathcal{F}(x) can be regarded as the derivative of g(\alpha) respect to \alpha evaluated for \alpha = 0; the first step is so determining
\begin{aligned} g'(\alpha ) &= \frac{\mathrm{d}}{\mathrm{d}\alpha}g(\alpha) \\ &= \frac{\mathrm{d}}{\mathrm{d}\alpha } \int_0^1 \big(x(t) + \alpha\, d(t)\big)^2 + 1\, \mathrm{d}t \\ &+ \frac{\mathrm{d}}{\mathrm{d}\alpha} \big(x(1) + \alpha\, d(1)\big) \end{aligned}
Assuming that \frac{\mathrm{d}}{\mathrm{d}\alpha}\int = \int \frac{\mathrm{d}}{\mathrm{d}\alpha} (operation that cannot always be performed) then the derivative of g(\alpha) can be regarded as
\begin{aligned} g'(\alpha) &= \int_0^1 \frac{\mathrm{d}}{\mathrm{d}\alpha} \Big(x(t) + \alpha\, d(t)\Big)^2\, \mathrm{d}t + \frac{\mathrm{d}}{\mathrm{d}\alpha} \Big(x(1) + \alpha \, d(1)\Big) \\ & = \int_0^1 2\big(x(t) + \alpha\, d(t)\big) d(t)\, \mathrm{d}t + d(1) \end{aligned}
Evaluating this expression for \alpha = 0 then the first variation of \mathcal{F}(x) with direction d(t) becomes
\delta \mathcal{F}(x;d) = g'(0) = \int_0^1 2xd\, \mathrm{d}t + d(1)
10.3 Integrals and Limits
In the computation of variations it is often assumed that \dfrac{\mathrm{d}}{\mathrm{d}\alpha}\circ\int_a^b \quad = \quad \int_a^b \circ \dfrac{\mathrm{d}}{\mathrm{d}\alpha}.
However, this property does not hold in general: some regularity conditions on the functions involved are required.
To better understand the issue, let us first examine a simpler problem.
Commuting Limit and Integral
Since a derivative is defined as a limit of a quotient, the question reduces to whether the operations of limit and integration can always be interchanged: \lim\circ\int \quad \stackrel{?}{=}\quad \int\circ\lim.
Concretely, let \{f_n\}_{n=1}^\infty be a sequence of functions such that f(x) = \lim_{n\to\infty} f_n(x), \qquad \forall x \in \mathbb{R}.
We ask whether \lim_{n\to\infty}\int_{-\infty}^{\infty} f_n(t)\,dt \stackrel{?}{=} \int_{-\infty}^{\infty}\lim_{n\to\infty} f_n(t)\,dt = \int_{-\infty}^{\infty} f(t)\,dt.
A Counterexample
Consider the sequence f_n(x) = \begin{cases} 0, & x \leq 0, \\[6pt] n x, & 0 < x \leq \tfrac{1}{n}, \\[6pt] 2 - n x, & \tfrac{1}{n} < x \leq \tfrac{2}{n}, \\[6pt] 0, & x > \tfrac{2}{n}. \end{cases}
Code
# f_n(x) a tratti
f_n <- function(x, n) {
ifelse(x <= 0, 0,
ifelse(x <= 1/n, n * x,
ifelse(x <= 2/n, 2 - n * x, 0)))
}
# valori di n e griglia x
n_seq <- c(1, 5, 10, 20, 25)
x_vals <- seq(-0.25, 2, length.out = 1000)
# dati per il plot: etichette legenda come "n==1" ecc. (parse-abili)
plot_data <- do.call(rbind, lapply(n_seq, function(n) {
data.frame(
x = x_vals,
y = sapply(x_vals, f_n, n = n),
func = factor(paste0("n==", n), levels = paste0("n==", n_seq))
)
}))
# Plot
ggplot(plot_data, aes(x = x, y = y, color = func)) +
geom_line(linewidth = 1) +
labs(
title = expression("Plots of functions " * f[n](x) * " for various " * n),
x = expression(x),
y = expression(f[n](x))
) +
theme_minimal() +
theme(legend.position = "right") +
scale_color_manual(
values = scales::hue_pal()(length(n_seq)), # palette senza library(scales)
labels = scales::label_parse(), # converte "n==1" -> espressione n==1
name = expression(n) # titolo della legenda in math
)
Each f_n is supported in the interval (0,2/n) and takes the shape of a triangular “bump” of base 2/n and height 1.
For every fixed x \in \mathbb{R} we have \lim_{n\to\infty} f_n(x) = 0, so the pointwise limit is the zero function: f(x) \equiv 0.
However, \int_{-\infty}^{\infty} f_n(t)\,dt = \text{area of the triangle } = \frac{1}{2}\cdot\frac{2}{n}\cdot 1 = \frac{1}{n}.
or in more detail \begin{aligned} \int_{-\infty}^\infty f_n(x)\, \mathrm{d}x &= \int_{0}^{1/n} nx \,\mathrm{d}x+ \int_{1/n}^{2/n} 2-nx \,\mathrm{d}x \\ & = \left[n\frac{x^2}{2}\right]_{x=0}^{x=1/n}+ \left[2x-n\frac{x^2}{2}\right]_{x=1/n}^{x=2/n} \\ & = \left[ n\frac{(1/n)^2}{2}-n\frac{(0)^2}{2}\right]+ \left[ 2\frac{2}{n}-n\frac{(2/n)^2}{2} -2\frac{1}{n}+n\frac{(1/n)^2}{2} \right] \\ & = \frac{1}{2n}+ \frac{4}{n}-\frac{2}{n} -\frac{2}{n}+\frac{1}{2n}=\frac{1}{n} \end{aligned}
Thus
\lim_{n\to\infty}\int_{-\infty}^{\infty} f_n(t)\,dt = 0, \qquad \int_{-\infty}^{\infty} f(t)\,dt = 0,
so in this case the interchange works.
Consider now the functions
g_n(x) = n f_n(t) = \begin{cases} 0 & x \leq 0 \\ n^2x & 0 < x \leq \frac{1}{n} \\ 2n-n^2x & \frac{1}{n} < x \leq \frac{2}{n} \\ 0 & x > \frac{2}{n} \end{cases}
Code
library(ggplot2)
# Define function g_n(x) for a given n
g_n <- function(x, n) {
ifelse(x <= 0, 0,
ifelse(x <= 1/n, n^2 * x,
ifelse(x <= 2/n, 2 * n - n^2 * x, 0)))
}
# Values of n and x
n_seq <- c(1, 5, 10, 20, 25)
x_vals <- seq(-0.25, 2, length.out = 1000)
# Create data frame for plotting
plot_data <- do.call(rbind, lapply(n_seq, function(n) {
data.frame(
x = x_vals,
y = sapply(x_vals, g_n, n = n),
func = factor(paste0("n==", n), levels = paste0("n==", n_seq))
)
}))
# Plot g_n(x) for different n
ggplot(plot_data, aes(x = x, y = y, color = func)) +
geom_line(linewidth = 1) +
labs(
title = expression("Plots of function " * g[n](x) * " for various " * n),
x = expression(x),
y = expression(g[n](x))
) +
theme_minimal() +
theme(legend.position = "right") +
scale_color_manual(
values = scales::hue_pal()(length(n_seq)), # CORRETTO: hue_pal da scales
labels = scales::label_parse(), # convert "n==1" -> n==1 in legend
name = expression(n) # legend title in math
)
it is clear that
\lim_{n\to\infty} g_n(x) = 0,\qquad \forall x\in\Bbb{R}
and
\begin{aligned} \int_{-\infty}^\infty g_n(x)\,\mathrm{d}x &= \int_{0}^{1/n} n^2x\,\mathrm{d}x+ \int_{1/n}^{2/n} 2n-n^2x\,\mathrm{d}x \\ & = \frac{1}{2}+\frac{1}{2} = 1 \end{aligned}
and thus
\begin{aligned} \lim_{n\to\infty}\int_{-\infty}^\infty g_n(x) \mathrm{d}x & = \lim_{n\to\infty}1 = 1 \\ \int_{-\infty}^\infty \lim_{n\to\infty}g_n(x) \mathrm{d}x & = \int_{-\infty}^\infty 0 \mathrm{d}x = 0 \end{aligned}
and in this case limit and integral cannot be exchanged.
Thus in general, when limits, derivative and integral are exhanged we shuold check if the opearation is possibile.