\documentclass[reqno]{amsart}
\usepackage{hyperref}
\usepackage{graphicx}

\AtBeginDocument{{\noindent\small
\emph{Electronic Journal of Differential Equations},
Vol. 2015 (2015), No. 206, pp. 1--17.\newline
ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu
\newline ftp ejde.math.txstate.edu}
\thanks{\copyright 2015 Texas State University - San Marcos.}
\vspace{9mm}}

\begin{document}
\title[\hfilneg EJDE-2015/206\hfil Piecewise uniform optimal design]
{Piecewise uniform optimal design of a bar with an attached mass}

\author[B. P. Belinskiy, J. W. Hiestand, J. V. Matthews \hfil EJDE-2015/206\hfilneg]
{Boris P. Belinskiy, James W. Hiestand, John V. Matthews}


\address{Boris P. Belinskiy \newline
 Department of Mathematics, University of Tennessee at Chattanooga,
615 Mccallie  Avenue, Chattanooga, TN 37403-2598, USA}
\email{Boris-Belinskiy@utc.edu}

\address{James W. Hiestand \newline
 College of Engineering,  University of Tennessee at Chattanooga,
 615 Mccallie  Avenue,  Chattanooga, TN 37403-2598, USA}
\email{James-Hiestand@utc.edu}

\address{John V. Matthews \newline
Department of Mathematics, University of Tennessee at Chattanooga,
615 Mccallie  Avenue, Chattanooga, TN 37403-2598, USA}
\email{Matt-Matthews@utc.edu}

\thanks{Submitted December 14, 2014. Published August 10, 2015.}
\subjclass[2010]{62K05, 80A20, 49R05, 35K05, 34B24, 65H10}
\keywords{Optimal design; heat transfer; heat equation;
least eigenvalue; \hfill\break\indent
 Sturm-Liouville problem; calculus of variations;
 transcendental equation; computer algebra}

\begin{abstract}
 We minimize, with respect to the cross sectional area, the mass of a
 bar given the rate of heat transfer. The bar enhances the heat transfer
 surface of a larger known mass to which the bar is attached.
 This article is an extension of a previous publication by two coauthors,
 where  heat transfer from the sides of the bar was neglected and only
 conduction through its length was considered. The rate of cooling is
 defined by the first eigenvalue of the corresponding Sturm-Liouville problem.
 We compare the mass of the computed variable cross-section bar with the mass
 of a bar with constant cross-sectional area and the same rate of heat transfer,
 and conclude that a fin design with constant, or near constant, cross-sectional
 area is best.
\end{abstract}

\maketitle
\numberwithin{equation}{section}
\allowdisplaybreaks

\section{Introduction} \label{intro_sect}

Removal of waste heat to another material or the environment by convection 
and radiation is important in everyday life and industrial applications.  
Heat must be removed from devices such as computers, refrigerators, 
and engines, large and small. Convective heat transfer utilizes the bulk 
motion of a fluid; heat transfer by radiation utilizes electromagnetic waves 
driven by the temperature difference between the source and the target.
	
The rate of heat transfer depends on the temperature difference, the area 
of the heat transfer surface, the heat transfer coefficient (for convection) 
and surface conditions and orientation (for radiation).  
Sometimes the temperature difference is fixed.  The convective heat 
transfer coefficient depends on the flow rate between the surface and the 
surrounding medium.  This may be enhanced by a fan, if available.  
Computers have fans, for example. We note that radiation is not considered in 
the analysis that follows. This effect is more important at high temperatures.

The condition most easily modified is the surface area between the source 
and the disposal medium for the heat transfer.  To increase the surface area
 without unduly increasing the weight, material used, and hence cost of the 
application, extended surfaces are often used. Such surface extensions for 
convective heat transfer frequently are called fins.  
Several examples of fins are 
(a) donuts placed at regular intervals around pipes, 
(b) protrusions like hairs from surfaces, and  
(c), layers of thin sheets of material, which have high surface-to-volume ratios.  
Human and animal bodies are also dependent on adequate heat transfer to prevent 
overheating.  Elephants, for example, increase their effective heat transfer 
surface areas through their large ears, which function as fins. 
On the other hand animals in cold climates need to retain body heat and 
hence tend to have larger bodies with lower surface-to-volume ratios.

\begin{figure}[ht]\begin{center}
\includegraphics[width=0.6\textwidth]{fig1} %diagram-cropped.pdf
\end{center}
\caption{Example of a fin attached to a larger mass, indicated by dotted boundary. 
Not to scale.}
\label{fig1}
\end{figure}


As shown in the tapered fin in Figure \ref{fig1}, in usual engineering practice 
the area $A(x)$ decreases in the direction away from the body to which the 
extended surface is attached. This practice is based on the steady-state 
analysis as justified in \cite[p. 79]{22} to minimize the area and hence the 
material used when a fixed amount of heat is to be continuously removed 
from the body to which the extended surface is attached. 
In contrast the present analysis is for transient cooling for a fixed amount 
of heat. This results in the opposite dependence of $A(x)$; it is slightly 
increasing in the direction away from the body. However, this result is 
consistent with Turner's and Taylor's results \cite{Turner,Taylor2} 
for a parallel problem in mechanical vibration 
(see also \cite{TaylorLiu,Taylor}).

This article is a continuation of a previous publication \cite{BeHMcC} 
where we discussed a particular case when convective heat transfer from 
the side of a fin is neglected. In the present work we do account for 
convective heat transfer but assume that the cross-section of the fin varies 
little along its length. To make the presentation here coherent, we sometimes 
repeat the results of \cite{BeHMcC} briefly.

Convective heat transfer is modeled by the equation
\begin{equation}
\dot Q=h A_s (T-T_\infty) \label{1.1}
\end{equation}
where $\dot Q$ is the heat transfer rate, $h$ is an empirical heat transfer 
coefficient, $A_s$ is the surface area, $T$ is the temperature of the surface 
and $T_\infty$ is the temperature of the surrounding medium. We want to maximize 
the surface-to-volume ratio since the heat transfer rate is proportional 
to the surface area. On the other hand, we want to minimize the volume of 
the heat transfer surface, in order to keep its weight and material cost as 
low as possible.  Hence fin design is an optimization problem.

Assume the extended surface is a surface of revolution with radius \(r(x)\) 
attached to a given base mass, $M_0$. The mass of the added extended surface
is small compared to the base mass $M_0$. Heat transfer along a bar is by 
conduction, with the rate given by
\begin{equation}
\dot Q=-\frac{kA \Delta T}{L}
\label{1.2}
\end{equation}
where $k$ is the thermal conductivity of the material, $A$ is the cross-sectional 
area, and $L$ is the length of the material. Moreover, $\Delta T$ is the 
temperature difference between the end points of the heat transfer. 
The equations above, along with the corresponding physical background, may be 
found in \cite{13}.

If an energy balance is performed for the region of the bar between $x$ and
 $x + \Delta x$, energy enters by conduction at $x$ and leaves by conduction
at $x + \Delta x$ and also from the side by convection (see Equation \eqref{1.1}). 
The difference is the rate of change of the energy content of that region of the bar. 
We find
\begin{quote}
 Rate of change of the region energy content equals
 conducted  energy in minus conducted energy out and convected energy out.
\end{quote}
Since conduction in the positive direction requires a negative temperature gradient, 
we find
\begin{equation}
\rho c A \Delta x \frac{\Delta T}{\Delta t}
=-kA\frac{\Delta T}{\Delta x}\Big |_x+kA\frac{\Delta T}{\Delta x}
\Big |_{x+\Delta x}-hA_s(T-T_{\infty}).
\label{1.3}
\end{equation}	
Here $T(x,t)$ is the temperature distribution. The surface area of a differential 
slice through the fin is $A_s(x) = P(x)\Delta x$, where
\[
P(x) = 2\pi r(x) \sqrt{ 1 + \left(r'(x)\right)^2}.
\]
The ratio $\frac{\Delta T}{\Delta t}$ is the rate of change of temperature 
with time and $\frac{\Delta T}{\Delta x}$ is the local temperature gradient. 
The bar material parameters are the material density, $\rho$, the specific
heat capacity, $c$, and the thermal conductivity, $k$. The convective heat
transfer coefficient is $h$. It is assumed that $\rho,c,k,\,{\rm and}\,h$ are positive constants.

Dividing by $\Delta x$, and taking the limit as $\Delta x$ and $\Delta t \to 0$, 
yields the partial differential equation
\begin{equation} \label{1.4}
A{\frac{\partial T}{\partial t}}
=\frac{k}{\rho c}\frac{\partial}{\partial x}
\big(A\frac{\partial T}{\partial x}\big)
-\frac{hP}{\rho c}(T-T_{\infty}),\quad (x,t)\in (0,l)\times (0,\infty).
\end{equation}

In a previous paper \cite{BeHMcC} we discussed a particular case of this general 
equation when convective heat transfer from the side of the bar is neglected, 
i.e., the limiting case $h \to 0$ was considered. 
Now we shall remove this restriction and allow $h$ to have a constant and positive 
value along the bar.  The presence of positive $h$ is very realistic, 
for the purpose of the fin is to increase the surface area involved in convection.  
Our purpose is to find the optimal distribution of the cross-sectional area, 
$A$, of a surface of revolution of a given heat transfer rate for a given $h$. 
This will produce minimum mass and may be considered the optimum solution.

When studying the limiting case $h\to 0$ in \cite{BeHMcC}, i.e. convective 
heat transfer from the side of the bar is neglected, it was possible to use 
the methods of the Calculus of Variations to find explicitly the optimal 
form of the cross-sectional area, $A(x)$, that maximizes the cooling rate 
subject to the given mass of the bar. The same form minimizes the mass given 
the cooling rate.

The methods of the Calculus of Variations are usually applied to optimal design 
problems of this sort. Many structural problems have been considered using these 
methods, including the maximization of a column's buckling 
load \cite{Taylor,TaylorLiu}, the minimization of the mass of an oscillating 
bar \cite{Taylor2,Turner}, the maximization of a column's 
height \cite{KellerNiordson} and the minimization of the moment of inertia of 
an oscillating turbine \cite{BeMcCW}. See also \cite{Cox,CoxMcCarthy,CoxOverton}. 
Only a few of the design problems mentioned, for example \cite{Turner} 
and \cite{BeHMcC}, have boundary conditions that contain a spectral parameter, 
as we will have here. A more detailed review of results in the area can be 
found in \cite{Dym} and the references therein.

In general, convective heat transfer may not be neglected. 
It is the purpose of this paper to find the optimal distribution of the 
cross-sectional area $A(x)$ of a surface of revolution with the given heat 
transfer rate such that the mass of the bar is minimal. It appears, though, 
that the application of the Calculus of Variation methods does not lead in 
this case to an explicitly solvable equation as in the case of $h=0$.

The optimal shape determined below was obtained by dividing the fin into a 
series of constant cross area segments.  This was mathematically necessary 
to solve the formulated equations. In practice these segments could be fit 
by a smooth curve. This is facilitated by the small variation in area among 
adjacent segments found below.

Hence, we find a piecewise uniform design for the bar under consideration. 
The corresponding mechanical problem, i.e. piecewise uniform optimum design 
of a bar with the minimal mass and of a given first eigenfrequency, was 
considered in \cite{Turner,Sheu,Sippel}. The exact solution for two uniform 
regions was found. The next step was made in \cite{Cardou}, where a closed 
form of the exact solution was found for an arbitrary number of uniform regions. 
The author then used it to show, at the numerical level, the convergence of 
discrete optimization to the continuous one, which is known explicitly 
in this case, see \cite{Turner,Taylor}.

Computationally we use two different approaches on the two different problems. 
Specifically, in the problem without convection for a bar of \(n\) equal-length 
pieces each with its own constant cross-sectional area, we may write the 
optimization in terms of quantities known from the problem with \((n-1)\) pieces, 
eliminating all but a single unknown. As a result, we may obtain exact solutions 
for an arbitrary number of pieces, and demonstrate this for several small values 
of \(n\).

For the problem with convection and a bar of \(n\) equal-length pieces each with 
its own constant cross-sectional area, the resulting equations permit no such 
exact solution. Rather, we formulate the optimization problem with the help 
of a Lagrange multiplier. With the help of a computer algebra system we obtain 
an exact reformulation of the optimization problem as a large system of nonlinear 
equations and then use numerical software to obtain an approximation of the exact 
solution. We verify that as the convective term approaches zero, our numerical 
solutions converge to the value of an exact solution of the system without 
convection.

The paper is organized as follows. In Section \ref{sect2} we give the mathematical 
description of the model, apply separation of variables, and describe some 
spectral properties of the corresponding Sturm--Liouville problem. 
We also give the solution for an elementary case of the problem when the 
cross-sectional area is constant, which we need later for comparison 
with the solution in case of variable cross-sectional area. 
In Section \ref{sect3} we formulate the design problem and we derive the 
necessary conditions of optimality in the form of a nonlinear differential 
equation for the first eigenfunction of the Sturm--Liouville problem. 
We briefly summarize the analysis from \cite{BeHMcC} in Section \ref{sect5} 
for a bar with no convective heat transfer and find an optimal form of the bar. 
In Section \ref{sect6} we develop another approach based on the discretization 
of the bar, i.e. we assume the cross-section $A(x)$ to be piece--wise constant.
Our algorithm allows us to find an optimal form of the bar (within the class 
of the piece-wise constant cross-sections). In this Section we proceed for 
the case without convective heat transfer. In Section \ref{sect7} we develop 
the similar approach for a bar with convective heat transfer. 
In Section \ref{sect8} we give the numerical comparison of the bar of optimal 
shape and the bar having the same cooling properties but with the constant 
cross-sectional area, and discuss the results of the study from an engineering 
point of view. We conclude that, from this perspective, our results suggest a 
fin design with constant, or near constant, cross-sectional area.

\section{Heat transfer of a bar of a variable cross-sectional area: 
separation of variables and the Sturm-Liouville problem}
\label{sect2}

We consider the heat transfer in a cylindrical bar $\{0<x<l\}$ with a base mass 
$M_0$ attached at the end point $x=0$. The temperature distribution 
$T:[0,l]\times [0,\infty)\to\mathbb{R}$ satisfies the transient one-dimensional 
conduction equation
\begin{equation} \label{2.1}
A{\frac{\partial T}{\partial t}}
=\frac{k}{\rho c} \frac{\partial}{\partial x}
\Big(A\frac{\partial T}{\partial x}\Big)-\frac{h P}{\rho c}(T-T_{\infty}),
\quad (x,t)\in (0,l)\times (0,\infty).
\end{equation}
Here $A, P$ are continuous differentiable positive functions from \([0,l]\) 
to \(\mathbb{R}_+\).
 As was mentioned in the Introduction, parameters $k, \rho, c$ are positive 
constants. The end point $x=l$ is kept at the (constant) temperature of the 
surrounding medium,
\begin{equation} \label{2.2}
T(l,t)=T_{\infty},\quad t\in [0, \infty).
\end{equation}
The rate of change of the energy content of the base is given by the difference 
between the energy flow into and out of it, as in the derivation of  \eqref{1.3}.
 For energy flow only by conduction outward at $x = 0$ this becomes
\begin{equation} \label{2.3-prel}
c M_0 \frac{\Delta T}{\Delta t}=k A \frac{\Delta T}{\Delta x}\big|_{x=0}.
\end{equation}
In the limit as $\Delta t$ and $\Delta x \to 0$ this becomes
\begin{equation} \label{2.3}
cM_0\frac{\partial T}{\partial t}\big(0,t\big)
=kA(0)\frac{\partial T}{\partial x}\big(0,t\big)\quad t\in [0, \infty).
\end{equation}
The initial distribution $T_0:[0,l]\to\mathbb{R}$ of the temperature is given,
\begin{equation} \label{2.4}
T(x,0)=T_0(x).
\end{equation}
It is well known that the initial boundary value problem \eqref{2.1}--\ref{2.4} 
has a unique solution \cite{Powers,Lad}. It is convenient to extract the 
term $T_{\infty}$ from the solution,
\begin{equation}\tau(x,t):= T(x,t)-T_{\infty}. \label{2.5}
\end{equation}
The new unknown function $\tau:[0,l]\times [0,\infty)\to\mathbb{R}$
 is the unique solution of the initial boundary value problem
\begin{gather} \label{2.6}
A\frac{\partial\tau}{\partial t}=\frac{k}{\rho c}
\frac{\partial}{\partial x}\Big(A\frac{\partial\tau}{\partial x}\Big)
-\frac{h P}{\rho c}\tau,\quad (x,t)\in (0,l)\times (0,\infty),
\\
\tau(l,t)=0,\quad t\in [0, \infty), \label{2.7}
\\
cM_0\frac{\partial\tau}{\partial t}\big(0,t\big)=k A(0)\frac{\partial\tau}{\partial x}\big(0,t\big),
\quad t\in [0, \infty),\label{2.8}
\\
\tau(x,0)=\tau_0(x),\quad x\in [0, l];\text{ here }\tau_0(x):= T_0(x)-T_{\infty}.
\label{2.9}
\end{gather}
If we use the standard procedure of the separation of variables
\begin{equation}
\tau(x,t):= e^{-\sigma t}\,u(x)
\label{2.10}
\end{equation}
and introduce the notation
\begin{equation} \label{2.11}
\frac{\rho c}{k}\sigma:= \lambda\quad\text{so that }
\frac{c \sigma}{k}=\frac{\lambda}{\rho},
\end{equation}
then the function $u:[0,l]\to\mathbb{R}$ satisfies the  
Sturm-Liouville problem
\begin{gather} \label{2.12}
(Au')'+\lambda Au-\frac{h P}{k}u=0,\quad x\in (0,l); \\
\label{2.12-1}
u(l)=0,\quad A(0)u'(0)+\frac{M_0}{\rho}\lambda u(0)=0.
\end{gather}
Although the spectral parameter $\lambda$ appears in the second boundary 
condition, a general theory for Sturm-Liouville problems of this type developed 
in \cite{Walter,Fulton,Hinton,Bel-Dau} may be used.
 It can be verified that the conditions of the corresponding theorems are satisfied. 
In particular, the eigenparameter dependent Sturm-Liouville problem 
\eqref{2.12}--\eqref{2.12-1} has a pure discrete positive real spectrum with 
the only point of accumulation at $+\infty$. The set of eigenfunctions 
satisfies two orthogonality relations for \(\lambda_n \neq \lambda_j\):
\begin{equation}
\begin{gathered}
\int_0^l Au_n u_j dx+\frac{M_0}{\rho} u_n(0)u_j(0)=0,\\
\int_0^l \Big(Au_n' u_j'+\frac{h P}{k}u_n u_j\Big)\,dx=0.
\end{gathered} \label{2.13}
\end{equation}
The Rayleigh quotient
\begin{equation} \label{2.14}
\lambda_n=\frac{\int_0^l Au_n'^2 dx+\frac{h}{k}\int_0^l P u^2 dx}
{\int_0^l Au_n^2 dx+{M_0 u_n^2(0)}/{\rho}}
\end{equation}
immediately follows.

The unique solution of the initial boundary value problem \eqref{2.6}-\eqref{2.9} 
may be represented by
\begin{equation}
\tau(x,t):= \sum_{n\ge 1} b_n e^{-\sigma_n\,t} u_n(x) \label{3.1}
\end{equation}
where
$$
\sigma_n:= \frac{k}{\rho c}\lambda_n
$$
and
\begin{equation} \label{3.2}
b_n=\frac{\rho \int_0^l A\tau_0 u_n dx+M_0 \tau_0(0)u_n(0)}
{\rho \int_0^l Au_n^2 dx+M_0 u_n^2(0)}.
\end{equation}

We note that if the cross-sectional area is constant, then $A(x)=A$ 
and \(P(x) = \alpha \sqrt{A}\) where the parameter 
$\alpha:=2\sqrt{\pi}$ connects the perimeter and cross-sectional area of 
a cylindrical bar.
The mass of the bar is $M=\rho Al$, and
the exact solution of the problem \eqref{2.6}-\eqref{2.9} is given by
$$
\tau(x,t):= \sum_{n\ge 1} b_n e^{-\sigma_nt}\sin
\Big(\Big(\lambda_n - \frac{\alpha h}{k \sqrt{A}} \Big)^{1/2}(x-l)\Big).
$$
Here $\lambda_n$ are the positive solutions of the transcendental equation
\begin{equation}
\frac{\lambda_n l^2 \tan\Big( \sqrt{ \lambda_n - \frac{\alpha h}{k \sqrt{A}}}\,
 l \Big)}
{ \sqrt{ \lambda_n - \frac{\alpha h}{k \sqrt{A}} } \, l }
= \frac{M}{M_0}, \quad  n = 1, 2, \ldots
\label{4.7}
\end{equation}

\section{Design problem. 
Necessary condition of optimality (smooth bar surface)}
\label{sect3}
The representation \eqref{3.1}-\eqref{3.2}, and \eqref{2.5} for the solution 
shows that the temperature $T(x,t)$ approaches the level $T_{\infty}$ 
exponentially fast, and the rate of approach is determined by the first 
eigenvalue $\lambda_1$. We now formulate the problem of optimal design.

\subsection*{Design Problem}
 Given the first eigenvalue $\lambda_1$ of the bar find 
the shape $A(x)>0$ of the bar such that its mass
\begin{equation} \label{M}
M=\rho\int_0^l A(x)\,dx
\end{equation}
is a minimum.

Our admissible class of designs is given by
\begin{equation}
ad = \Big\{ r : 0<r(x)<\infty,\; x\in [0,l];\; 
\frac{\int_0^l Au_1'^2 dx+\frac{h}{k}\int_0^l\;P\,u_1^2 dx}
{\int_0^l Au_1^2 dx+M_0 u_1^2(0)/\rho}=\lambda_1\Big\}
\label{5.1}
\end{equation}
where $u_1(x)$ is the first eigenfunction.

We seek the necessary conditions of optimality. It appears that the 
corresponding differential equation is too complex to be solved explicitly. 
For the purposes of illustration, we first consider the specific case \(h = 0,\) 
following the derivation originally given in \cite{BeHMcC}.

Since we have already known that the spectrum is discrete, we may use standard 
Calculus of Variations techniques to derive an optimality condition in 
the form of a differential equation. Consider the functional
\begin{equation} \label{6.0}
\begin{aligned}
F(A)&:=\int_0^l \rho\,A\,dx
+\int_0^l{\Lambda_1\Big((Au')'+\lambda_1 Au\Big) }dx \\ 
&\quad+ \Lambda_2\Big(A(0)u'(0)+\frac{M_0}{\rho}\lambda_1 u(0)\Big)
\end{aligned}
\end{equation}
where $\Lambda_j$, $j=1,2$ are Lagrange multipliers, and equate the first 
variation to zero. Note that while \(\Lambda_1\) is a function of \(x\), 
the multiplier \(\Lambda_2\) is merely a constant. The next few equations 
are quite cumbersome. To simplify them, we omit the limits of integration 
since they are the same in all integrals. We find
\begin{align*}
\delta F &=  \int \rho \delta A\,dx \\
&\quad +  \int \underline{\Lambda_1} \Big(\underline{(A\delta u')'
+(u'\delta A)'}+\lambda_1(A\delta u+u\delta A) \Big) \\
 & \quad +  \Lambda_2\Big(u'(0)\delta A(0)+A(0)\delta u'(0)
+\frac{M_0}{\rho}\lambda_1 \delta u(0)\Big)=0.
 \end{align*}

The underlined terms are integrated by parts and produce the following terms
\begin{equation} 
\Lambda_1 (A\delta u'+u'\delta A)\big|_{x=0}^{x=l}
-\Lambda_1'A\delta u\big|_{x=0}^{x=l}
+\int\big((A\Lambda_1')'\delta u-\Lambda_1'u'\delta A\big)\,dx
\label{6.4}
\end{equation}
We further consider the equality $\delta F=0$ and equate to zero the 
coefficients in front of all variations, $\delta A(x)$, $\delta u(x)$, etc. 
thus deriving the system of necessary conditions for optimality
\begin{gather} \label{var-A}
\delta A(x): \rho-\Lambda_1'u'+\lambda_1 \Lambda_1 u=0; \\
\label{var-u}
\delta u(x): (A\Lambda_1')'+\lambda_1 A\Lambda_1=0; \\
 \label{var-A(0)}
\delta A(0): \Lambda_2u'(0)-\Lambda_1(0)u'(0)=0; \\
 \label{var-u'(0)}
\delta u'(0): \Lambda_2A(0)-\Lambda_1(0)A(0)=0; \\
 \label{var-u(0)}
\delta u(0): \Lambda_2\frac{M_0}{\rho}\lambda_1 +A(0)\Lambda_1'(0)=0; \\
 \label{var-A(l)}
\delta A(l): \Lambda_1(l)u'(l)=0; \\
 \label{var-u(l)}
\delta u(l): \Lambda_1'(l)A(l)=0; \\
 \label{var-u'(l)}
\delta u'(l): \Lambda_1(l)A(l)=0.
\end{gather}
Recall that we chose $\lambda=\lambda_1$, which is positive 
(see the Rayleigh quotient \eqref{5.1}). Hence, the boundary condition 
at the end of the bar $x=0$ (see \eqref{2.12-1}) implies that $u'(0)\ne 0$. 
The equality \eqref{var-A(0)} then implies
\begin{equation} \label{Lambda=Lambda}
\Lambda_2=\Lambda_1(0).
\end{equation}
The very same equality appears in \eqref{var-u'(0)} since the admissible 
class of designs \eqref{5.1} requires $A(0)>0$. With \eqref{Lambda=Lambda} 
in mind, we conclude that the equality \eqref{var-u(0)} is similar to the 
boundary condition at the end $x=0$ but for $\Lambda(x)$,
\begin{equation} \label{b-cond-0-Lambda}
A(0)\Lambda_1'(0)+\frac{M_0}{\rho}\lambda_1 \Lambda_1(0)=0.
\end{equation}
Similarly, since the admissible class of designs \eqref{5.1} requires $A(l)>0$, 
the equalities \eqref{var-u(l)}--\eqref{var-u'(l)} imply the boundary condition 
at the end $x=l$ of the bar,
\begin{equation} \label{b-cond-l-Lambda}
\Lambda_1(l)=0\quad\text{or}\quad \Lambda_1'(l)=0.
\end{equation}
Following \cite{Turner,Taylor2}, we fix $u(l)$ to be any non-zero number, 
since the admissible class remains the same for any factor included in $u(x)$. 
After that, the term \(\delta u(l) \Lambda_1'(l)A(l)\) vanishes and we obtain 
only $\Lambda_1(l)=0$.

Equations \eqref{var-u} and \eqref{2.12} now imply that the functions 
$\Lambda_1(x)$ and $u(x)$ are proportional. After that, the equation \eqref{var-A} 
leads to the following equation for the optimal eigenfunction $u(x)$:
\[
u'^2-\lambda_1 u^2-\rho=0.
\]

For the case where \(h\neq0\) and an arbitrary function $A(x)$, we probably 
have no hope to solve the corresponding differential equation analytically. 
For example, for \(|r'| \ll 1\), the equation has the form
\begin{equation} \label{difficult-eqn}
u'^2+\Big(-\lambda_1 +\frac{\alpha h}{2k\sqrt{A(x)}}\Big)u^2-\rho=0.
\end{equation}

\section{A bar with no convective heat transfer ($h=0$) - 
analytic approach}
\label{sect5}

As shown above, we begin with the simplified version of the optimization 
problem, $h=0$. In this case, the equation has the form
\begin{equation} \label{nec-cond-2}
u'^2-\lambda_1 u^2-\rho=0
\end{equation}
and, along with the boundary condition at the end $x=l$ of the bar 
(see \eqref{2.12-1}), may be solved explicitly,
\begin{equation} \label{opt-u-1}
u(x)=\sqrt{R}\sinh\sqrt{\lambda_1}(x-l),
\end{equation}
for an arbitrary constant \(\sqrt{R}\).
The original Sturm--Liouville problem \eqref{2.12} may now be considered as 
a differential equation for the cross-sectional area $A(x)$,
\begin{equation} \label{opt-A-1}
(A\sqrt{\lambda_1}\cosh\sqrt{\lambda_1}(x-l))'+\lambda_1 A\sinh 
\sqrt{\lambda_1} (x-l)=0.
\end{equation}
After integration of this equation we arrive at
\begin{equation} \label{opt-A-2}
A(x)=\frac{C}{\cosh^2 \sqrt{\lambda_1} (x-l)}.
\end{equation}
The arbitrary constant $C$, as well as the connection between the parameter 
$\lambda_1$ and parameters of the model, may be found if we use the boundary 
condition at $x=0$ from \eqref{2.12-1} and evaluate the total mass of the bar. 
We omit the details since they were discussed in \cite{BeHMcC}, and only 
give the results. The optimal form of the bar is
\begin{equation} \label{opt-A-3}
A(x)=\frac{M_0\,\sqrt{\lambda_1} \sinh{\sqrt{\lambda_1}l} 
\cosh{\sqrt{\lambda_1}l}}{\rho\,\cosh^2{\sqrt{\lambda_1}(x-l)}}
\end{equation}
and the connection holds that
\begin{equation} \label{6.16}
M=M_0\int_0^l\frac{\sqrt{\lambda_1}\sinh{\sqrt{\lambda_1}l} 
\cosh{\sqrt{\lambda_1}l}} {\cosh^2{\sqrt{\lambda_1}(x-l)}}\,dx
=M_0\sinh^2{\sqrt{\lambda_1}l}.
\end{equation}
We note that the optimal form of the bar is the same as that found in  
\cite{Turner,Taylor2}, see also \cite{BeMcCW,BeHMcC}.
Equation \eqref{6.16} may be solved for $\lambda_1$ to produce the explicit 
expression for the optimal rate of cooling for the bar with the given mass
\begin{equation} \label{6.17}
\lambda_1= \Big(\frac{1}{l}\ln\Big(\sqrt{\frac{M}{M_0}}
 +\sqrt{\frac{M}{M_0}+1}\,\Big) \Big)^2.
\end{equation}
We introduce the dimensionless cooling rate, as
\begin{equation} \label{6.18}
z^{\rm opt}:= \sqrt{\lambda_1}\, l
\end{equation}
and find, from the representation \eqref{6.17}
\begin{equation} \label{6.18-1}
z^{\rm opt}:= \sqrt{\lambda_1} l
=\ln\Big(\sqrt{\frac{M}{M_0}}+\sqrt{\frac{M}{M_0}+1}\,\Big).
\end{equation}
In particular, we may compare the dimensionless cooling rate with the 
similar parameter $z$ for the constant cross-section. The last parameter 
satisfies 
\begin{equation} \label{6.19}
z \tan z= \frac{M}{M_0}
\end{equation}
which is just \eqref{4.7} for \(h=0\).
In \cite{BeHMcC}, we show that the inequality
\begin{equation} \label{6.20}
z^{\rm opt}\,>z
\end{equation}
holds for any (positive) ratio $M/M_0$ and confirm it by the numerical results. 
Hence, indeed, the optimal cooling rate for the bar with the optimal 
cross-section $A(x)$ given by \eqref{opt-A-3} is higher than for the 
bar with a constant cross-section.

We now describe the strategy for the remaining part of the paper. 
Since we see no possibility to proceed with the equation \eqref{difficult-eqn} 
for the optimal eigenfunction analytically, we use a numerical algorithm. 
The natural idea is to use piecewise uniform design, as it was suggested for 
axial vibration problems (see \cite{Cardou}, \cite{Sheu}, \cite{Sippel} and 
the references therein). To accumulate numerical experience, we first study 
the optimization problem from the current section, i.e. for a bar without 
convective heat transfer.

\section{A bar without convective heat transfer - discretized approach}
\label{sect6}

We consider the differential equation \eqref{2.12} with $h=0$ subject to 
the boundary conditions \eqref{2.12-1} and discretize this problem. 
Let the bar be split into $n$ equal pieces, 
$x_j-x_{j-1}=\frac{l}{n}:=\Delta$, $j=1,2,\dots,n$, so that $x_1=0$, $x_{n+1}=l$,
 and the cross-section $A(x)$ be piecewise constant,
 $A(x)=A_j$  for $x\in (x_j,x_{j+1})$. It is more convenient to introduce 
local coordinates, so that $0<x<\Delta$ on each piece and $u(x):=u_j(x)$ on 
$x\in (x_j,x_{j+1})$. Along with the boundary conditions, we need to introduce 
the continuity conditions, i.e.
\begin{equation}\label{T-cont}
\begin{gathered}
 u_j(\Delta-0)=u_{j+1}(0+), \quad j=1,2,\dots,n-1\\
 \frac{du_j}{dx}(\Delta-0)=\frac{du_{j+1}}{dx}(+0),
\quad j=1,2,\dots,n-1 \end{gathered} 
\end{equation}
if there are at least two pieces ($n>1$). The boundary conditions \eqref{2.12-1}
 result in
\begin{equation} \label{T-bc}
A_1 u_1'(0)+\frac{M_0}{\rho}\lambda u_1(0)=0,\quad u_n(\Delta)=0.
\end{equation}
Since the cross-section $A(x)$ is piecewise constant, the differential 
equation \eqref{2.12} may be solved explicitly. We find
\begin{equation} \label{discr-1}
u_j(x)=C_j \sin (\sqrt{\lambda}\,x+\varphi_j),\quad j=0,1,\dots,n.
\end{equation}
Substituting the form $u_j(x)$ of the solution in the boundary conditions 
\eqref{T-cont} and \eqref{T-bc} leads to a system of transcendental equations. 
To formulate it in a relatively simple form, we introduce the parameters
\begin{equation} \label{mu}
\mu:=\frac{M_0\sqrt{\lambda}}{\rho},\quad 
\nu:=\cot \sqrt{\lambda} \Delta.
\end{equation}
The system of equations has the form
\begin{gather}
\label{trans-syst-1a} A_1\cot \varphi_1+\mu =  0; \\
\label{trans-syst-1b} C_j\sin(\sqrt{\lambda} \Delta+\varphi_j) =  C_{j+1}\sin \varphi_{j+1},\\
\label{trans-syst-1c} A_j C_j\sqrt{\lambda}\cos(\sqrt{\lambda} \Delta+\varphi_j) 
 =  A_{j+1} C_{j+1}\sqrt{\lambda}\cos \varphi_{j+1},\\
\label{trans-syst-1d} \sqrt{\lambda}\Delta+\varphi_n  = \pi m
\end{gather}
for \( j=1,\ldots,n-1\) with an arbitrary integer $m$, so that
\begin{equation} \label{trans-syst-n}
\cot \varphi_n=-\nu.
\end{equation}
For each $j=1,\ldots,n-1$, dividing \eqref{trans-syst-1c} by \eqref{trans-syst-1b} 
yields
\begin{equation} \label{trans-syst-j}
A_j \cot(\sqrt{\lambda} \Delta+\varphi_j)=A_{j+1} \cot \varphi_{j+1}.
\end{equation}
Now with the help of an angle sum identity for cotangent,
\begin{equation*}
\cot(\theta_1+\theta_2) = \frac{\cot(\theta_1)\cot(\theta_2)-1}{\cot(\theta_1) 
+ \cot(\theta_2)},
\end{equation*}
we may rewrite the system of equations \eqref{trans-syst-1a}--\eqref{trans-syst-j} 
in the form
\begin{gather} \label{trans-syst-2-a}
A_1\cot \varphi_1+\mu=0; \\
 \label{trans-syst-2-b}
A_j \frac{\nu\cot \varphi_j-1}{\nu+\cot \varphi_j}=A_{j+1}\cot \varphi_{j+1},
\quad j=1,\dots,n-1;\\
\label{trans-syst-2-c}
\cot \varphi_n=-\nu.
\end{gather}
We can express $\cot \varphi_1$ from the first of these equations, substitute 
into the second equation, solve it for $\cot \varphi_2$, etc. We see that 
the transcendental functions may be excluded from the system 
\eqref{trans-syst-2-a}-\eqref{trans-syst-2-c} finally yielding the algebraic 
connection between the $A_j$'s.

We then evaluate the volume of the bar (with a piece-wise constant cross section) 
as follows
\begin{equation} \label{V-1}
V(A_1,\dots,A_n)=\Delta \sum_{j=1}^n\,A_j.
\end{equation}
Minimization of the volume may be finally performed by the standard calculus method; 
the algebraic connection between the cross-sections $A_j$ is a constraint for 
this minimization problem.

We now make the observation that the system of equations given by 
\eqref{trans-syst-2-a}-\eqref{trans-syst-2-c} for a bar subdivided into \(n\) 
pieces of equal length is a subset of the corresponding system of equations 
given by \eqref{trans-syst-2-a}-\eqref{trans-syst-2-c} for a bar subdivided 
into \((n+1)\) pieces. Given the pairwise dependence of the variables, 
we can use the solution of the optimization problem \eqref{V-1} with \(n\) 
pieces to reduce the optimization problem \eqref{V-1} with \((n+1)\) pieces 
to a problem in a single unknown.

By way of explanation, when \(n=2\) the optimization problem \eqref{V-1}
 yields optimal \(A_1, A_2\) with the ratio
\begin{equation} \label{a1a2ratio}
\frac{A_1}{A_2} = \frac{\nu^2}{\nu^2+2}.
\end{equation}
To see how one arrives at this value, the three relations  
\eqref{trans-syst-2-a}-\eqref{trans-syst-2-c} can be reduced to the 
single relation
\[
\frac{A_1}{A_2} = \frac{\nu (\nu A_1 - \mu)}{\nu\mu + A_1}
\]
and from that point one may write \(A_2\) in terms of \(A_1\) as
\[
A_2 = \frac{A_1(\nu\mu + A_1)}{\nu (\nu A_1 - \mu)}.
\]
To minimize \eqref{V-1}, we may instead simply minimize
\[
A_1 + A_2 = A_1 + \frac{A_1(\nu\mu + A_1)}{\nu (\nu A_1 - \mu)}.
\]
This explicit expression for \(A_1\) has its minimum at
$A_1 = 2\mu/\nu$ 
and then the corresponding value of \(A_2\) is
\begin{equation} \label{optimalA2}
A_2 = \frac{A_1 (2+\nu^2)}{\nu^2}.
\end{equation}
Immediately one arrives at the ratio \eqref{a1a2ratio} indicated above.


Our observation above then leads to the conclusion that when \(n=3\)
 we will find an optimal solution to \eqref{V-1} for which
\[
\frac{A_2}{A_3} = \frac{\nu^2}{\nu^2+2}.
\]
With this additional relationship, we may find an expression for
 \eqref{V-1} in which we have replaced \(A_2\) and \(A_3\) with functions 
of \(A_1,\) and we may then optimize in the single variable \(A_1.\)

In a similar way, once that optimal solution is obtained in this way for
 \(n=3\), the ratios amongst the three pieces \(A_1, A_2, A_3\) may then be
 used to reduce the optimization problem for \(n=4\) to a single variable problem.

We give only a few examples and compare the results with the absolute minimum 
of the volume given by \eqref{6.16}.
\begin{itemize}
\item \(n=1\) We find
\begin{equation} \label{V-n=1}
A_1=\frac{\mu}{\cot \sqrt{\lambda} l};\quad
\frac{V_{1,\rm min}}{\mu l}=\frac{1}{\cot \sqrt{\lambda} \Delta}
\end{equation}
where we assume $\sqrt{\lambda} \Delta<\pi/2$.

\item \(n=2\) See the process outlined above, starting from 
\eqref{a1a2ratio} and arriving at \eqref{optimalA2}, which leads to
\begin{equation} \label{V-n=2}
{\frac{V_{2,\rm min}}{\mu l}=\frac{2\big(1+\cot^2\big(\sqrt{\lambda} l/2\big)\big)}
{\cot^3 (\sqrt{\lambda} l/2)}}
\end{equation}
provided $\sqrt{\lambda} l<\pi$.

\item \(n=3\) Following an optimization process analogous to the above, we find
\begin{equation}
\frac{V_{3,\rm min}}{\mu l} = \frac{1}{3} 
\frac{\big(3 \cot^2\big({\sqrt{\lambda}l}/{3}\big) + 4\big)^2}
{\cot^5\big({\sqrt{\lambda}l}/{3}\big)}.
\end{equation}

\item \(n=4\) Again, using the same optimization process, we find
\begin{equation}
\frac{V_{4,\min}}{\mu l} = 4 
\frac{\big(\cot^2\big({\sqrt{\lambda}l}/{4}\big) + 2\big)^2 
\big(\cot^2\big({\sqrt{\lambda}l}/{4}\big) + 1\big)}
{\cot^7\big({\sqrt{\lambda}l}/{4}\big)}.
\end{equation}

\item \(n=\infty\) It easily follows from the representation \eqref{6.16} 
for the optimal mass (with no convective heat transfer and for the smooth surface 
bar $A(x)$) that
\begin{equation} \label{V-n=inf}
{\frac{V_{\infty,\min}}{\mu l}
=\frac{\sinh^2\big( \sqrt{\lambda} l \big)}{\sqrt{\lambda} l}.}
\end{equation}
\end{itemize}
	Normalization of the volume in the formulas \eqref{V-n=1} through 
\eqref{V-n=inf}, and below is based on the fact that $\mu l$ has the dimension 
of volume and $\sqrt{\lambda} \, l$ is dimensionless. It is easy to check that
\begin{equation} \label{ineq-1}
V_{1}>V_{2}>V_{3}>V_{4}>V_{\infty}.
\end{equation}

\section{A bar with convective heat transfer - discretized numerical approach}
\label{sect7}

The algorithm in Section 5, allows us to minimize the mass of the bar with 
the given cooling rate and without convective heat transfer. 
For a general discretization, the solution may be found exactly. 
In this section, we generalize that same discretized approach to the
 optimization problem for a bar with convective heat transfer. 
Specifically, we consider the differential equation \eqref{2.12} 
with $h\ne 0$ subject to the boundary conditions \eqref{2.12-1} and 
then discretize this problem as we did in Section 5. We represent the general 
solution to the equation \eqref{2.12} in the local coordinates as follows
\begin{equation} \label{h-ne-0-1}
\begin{gathered}
{u_j(x)=C_j\sin(\gamma_j x+\varphi_j),\quad j=1,\dots,n}\,,\\
\text{where } {\gamma_j:=\sqrt{\lambda-\frac{\alpha h}{k\sqrt{A_j}}}}
\end{gathered}
\end{equation}
The boundary conditions \eqref{T-cont} and \eqref{T-bc} result in the 
following system of equations (we skip the derivation that is similar 
to one in Section 5)
\begin{equation} \label{trans-syst-h}
\begin{gathered}
A_1\gamma_1\cot\varphi_1+\mu\sqrt{\lambda}=0;\\
A_j\gamma_j\cot(\gamma_j \Delta+\varphi_j)
=A_{j+1}\gamma_{j+1}\cot \varphi_{j+1},\quad j=1,\dots,n-1;\\
\cot \varphi_n=-\cot \gamma_n \Delta.
\end{gathered}
\end{equation}
Since $\gamma_j$ is dependent on $A_j$ (see \eqref{h-ne-0-1}), we may 
not expect an explicit solution, as in the case of a bar without convective 
heat transfer ($h=0$). Yet, having a system of equations of the same 
structure as in \eqref{trans-syst-2-a}-\eqref{trans-syst-2-c} we may use 
a numerical approach.

The formulation of the problem which is then solved numerically proceeds as follows:
\begin{itemize}
\item Starting with the first equation in \eqref{trans-syst-h}, 
isolate \(\cot\varphi_1\).

\item In the next equation, which includes only \(\cot\varphi_1\) 
and \(\cot\varphi_2\), isolate \(\cot\varphi_2\).

\item Continuing in the natural way through each subsequent equation 
 in \eqref{trans-syst-h}, eliminate \(\cot\varphi_j\) and leave an equation 
 with only \(\cot\varphi_{j+1}\).

\item The final equation in \eqref{trans-syst-h} allows the elimination 
 of \(\cot\varphi_n\), and what remains is a nonlinear equation only in 
 terms of \(A_1, \ldots , A_n\). We write this in the form
\[ 
f( A_1, \ldots, A_n\ ) = 0. 
\]
\end{itemize}
Finally, we utilize a Lagrange multiplier, \(\beta\), and minimize an expression 
of the form
\begin{equation}
\label{lagrange-mult-system}
F( A_1, \ldots , A_n , \beta ) = \Big(\sum_{j=1}^n A_j \Big) 
+ \beta f( A_1, \ldots , A_n )
\end{equation}
by finding roots of \(\partial F/\partial A_j = 0 \) and 
\(\partial F/\partial \beta = 0\). Compare \eqref{lagrange-mult-system} 
with  \eqref{V-1}.

All of the algebraic manipulation and computation of partial derivatives 
is handled algorithmically through the computer algebra system. 
The computer algebra system also automatically generates expressions for 
the partial derivatives, which then make up the function whose roots are 
found with a nonlinear solver based on the Newton-Raphson method.

The solutions found are indeed minimizers of the \(F\) above and the 
corresponding values for \(A_1, \ldots , A_n\) are roots of \(f\).

\section{Discussion and comparison of cooling properties for $h=0$ and $h>0$}
\label{sect8}

For our numerical experiments we chose the following physical and material 
parameters. We consider a steel fin of length \(l = 0.02\,m\). 
The fin is assumed to have \(\rho = 7800\,kg/m^3\) and thermal conductivity 
of \(k = 40.0\, W/(m \cdot K)\). For the convective heat transfer coefficient 
\(h\) we considered a wide range of values, from \(h=0\) through \(h=24\). 
Typical values of $h$ for natural convection (density driven rather than 
forced by a fan) are \(6 - 30\,W/(m^2 \cdot K)\) (see \cite[p. 17]{KrBo}).

We first summarize the numerical results for \(h=0\) within the context 
of the prior results in \cite{BeHMcC} where it is assumed that \(h=0\) 
and \(\lambda\) is maximized subject to the given mass \(M\). 
The dimensionless product $z=\sqrt{\lambda} l$ is a function of the ratio 
$M/M_0$ in both the constant area case \eqref{6.19} and the optimal 
case \eqref{6.17}.
In \cite{BeHMcC} the equation \eqref{6.19} was solved numerically.

Numerical results in \cite{BeHMcC} and in this work show that the advantage 
of the optimum cross-section over the constant cross-section is small 
and becomes less so as the base mass, $M_0$, increases. This is 
physically reasonable. Indeed, recall that convective heat transfer 
from the side of the area has been neglected. Hence addition of the extended 
surface does nothing but move the boundary condition at $x=l$ that distance 
from $M_0$. Furthermore, as $M$ becomes small compared to $M_0$, its very 
presence becomes negligible and hence its shape does not matter.

Though in agreement with the general results of Section \ref{sect5}, 
\(z^{opt}\ge z\), the effect is not large because of the physical reasons 
explained above. Moreover, for $M/M_0\,\to 0$ the optimum and constant area 
results merge.

We now discuss the numerical results for \(h\geq0\) and a given value 
of \(\lambda\) chosen according to \eqref{6.17} with \(M=0.2 M_0\). 
We use the numerical approach of Section \ref{sect7} to minimize the 
mass of the fin. A range of values for \(n\) is considered here, up to \(n=5\).
 For a larger value, say \(n = 10\), the nonlinear system of 11 equations 
in 11 unknowns can be generated procedurally by a computer algebra system 
and imported into a numerical suite such as Matlab. However, just the 
representations of these equations (in the millions of characters for 
\(n=10\)) taxes the system to the point of impracticality. 
An alternative would be to code the program in a compiled language (like C) 
and then apply a similar Newton-type solver, but given the guidance 
in \cite{Cardou}, where an engineering problem with a similar discrete 
structure was considered, we anticipate that the relative gains made would 
be quite small, as the solution for \(n\) pieces is nearly optimal even 
for \(n\) small.

For the quantities \(A_j\) found exactly for the special case \(h = 0\) 
in Section \ref{sect6} and the corresponding quantities found numerically 
for \(h>0\) in Section \ref{sect7}, we have computed the shape of the 
corresponding bar. Specifically, for \(n\) pieces, with \(n\) from 2 to 5, 
and for \(h=0, 0.01, 10\), and \(24\), the results are displayed in 
Figure \ref{fig2}. The \(A_j\) visualized in this figure arise from 
nondimensionalized versions of the corresponding equations which result 
in \(0 \leq A_j \leq 1\). Consequently, we consider the shape of the 
resulting bar to be the important feature, not the specific values of 
the cross-sectional areas.

\begin{figure}[ht]
\begin{center}
\includegraphics[width=0.8\textwidth]{fig2} %bars-n-2-to-5-h-0-to-24-cropped.pdf}
\end{center}
\caption{Shapes of bars showing influence of the number of segments \(n\) 
and the parameter \(h\).}
\label{fig2}
\end{figure}

As we observe, the cross-sectional area $A(x)$ is increasing slightly as we 
move toward the tip of the extended surface. This is physically reasonable 
for the transient problem. Indeed, as we proceed in this direction 
(away from the heated mass) with $h>0$, heat is convected away from the 
extended surface. This lowers the temperature within the material as 
\(x\) increases. Hence the temperature difference between the material and 
the surrounding ambient medium decreases, and convective heat flow to the 
surrounding medium would decrease for a constant surface area 
(see  \eqref{1.1}). However, if the surface area, $A_s$, increases the 
reduction in the temperature difference may be offset by increasing $A_s$. 
A uniform heat flow along the external surface is desirable; our results 
show this.

We also observe that the results presented herein were for a range of values 
of $h$. Over this range the areas computed changed little. Hence the 
trends observed are essentially independent of this parameter. However, 
for the refined case \(n=5\), we can compare the masses of the fins for 
the values \(h=0.0, 0.01, 10,\) and \(24\) \(W/(m^2\cdot K)\) to the 
fin with homogeneous cross-section. Table 1 shows the computed masses
of these fins and the relative decrease in mass for a given value of \(h\).

\begin{table}[ht]
\caption{Relative improvements in mass for a range of values of \(h\)}
\label{table:data}
\begin{center}
\begin{tabular}{cccc}
\hline
\(h\) & Mass (kg)  & Mass (kg) & Relative \\
\(W/(m^2\cdot K)\)       & homogeneous & piecewise \(n=5\) & decrease \\
\hline
  0.00 & \(3.12 \times 10^{-3}\) & \(3.1157 \times 10^{-3}\) & \(0.14\%\)  \\
  0.01 & \(3.12 \times 10^{-3}\) & \(3.1156 \times 10^{-3}\) & \(0.14\%\)  \\
  10.0 & \(3.12 \times 10^{-3}\) & \(3.0732 \times 10^{-3}\) & \(1.50\%\)  \\
  24.0 & \(3.12 \times 10^{-3}\) & \(3.0107 \times 10^{-3}\) & \(3.50\%\)  \\
\hline
\end{tabular}
\end{center}
\end{table}

However in usual engineering practice, see \cite{22}, extended surfaces 
are not built this way for the following reasons.  It would be structurally 
undesirable to have the mass increase away from the body. This would result 
in the thinnest part of the surface being adjacent to the larger body.  
The resulting rapid transition would be an area most likely to break 
if subjected to an external load. Furthermore the flow area between 
adjacent extended surfaces would decrease in the direction of the tip. 
This would restrict the flow of the external medium towards the tip and 
thus reduce the ability of the medium to carry away heat. 
As further noted by \cite{22}, (see p. 79), since the heat transfer 
rate $\dot Q$ conducted through the extended surface decreases in the 
direction of flow because of heat loss from the surface, reducing the 
cross-sectional area in this direction tends to equalize the heat flux 
(heat/cross-sectional area) in the direction of flow, a desirable result.

\section{Conclusion} \label{sect11}

We have found the optimal distribution of the cross-sectional area of a bar 
in the form of a surface of revolution for a given heat transfer rate, 
such that the mass of the bar is a minimum. 
We use a variational principle in conjunction with numerical 
analysis, which is based on a piece--wise approximation of the surface. 
It may be expected by the Duality Principle that the very same form of 
the bar produces the maximum cooling rate if the total mass is given. 
For \(h=0\), the optimal distribution coincides with one found by 
 Taylor \cite{Taylor} and  Turner \cite{Turner} for the design of a bar 
having a maximum lowest eigenfrequency with the given mass. 
The influence of convection on the optimal form of the bar is studied. 
We emphasize that no analytic approach seems to be available in this case. 
In particular, the comparison of the optimal design with the constant 
cross-section bar is studied.

Numerical experiments demonstrate that for a range of \(h\) values the areas 
computed change very little. Moreover, the problem is generally not sensitive 
to the number of segments in the bar. If one assumes a bar with constant 
cross-section (i.e. \(n=1\) and the area \(A\) from \eqref{4.7}) and 
a range of \(h\) values, the new designs we find in Section \ref{sect7}, 
say for \(n=5\), differ in mass by less than 3.5\%.

In practical engineering the design of the extended surfaces is based on a 
steady-state analysis. However, the current work has opened our way forward 
to the corresponding analysis for the steady-state optimization problem. 
Our approach to that problem will be the focus of a future work.

\subsection*{Acknowledgments}
We wish to thank the anonymous referees who worked  diligently to improve 
and elucidate our work, offering many helpful and insightful suggestions.

\begin{thebibliography}{99}

\bibitem{Bel-Dau}
    \newblock B.~P. Belinskiy, J.~P. Dauer;
    \newblock \emph{Eigenoscillations of mechanical systems with boundary 
conditions containing the frequency},
    \newblock Quarterly of Applied Mathematics, \textbf{56} (1998), pp.~521--541.

\bibitem{BeHMcC}
    \newblock B.~P. Belinskiy, J.~W. Hiestand,  M.~L. McCarthy;
    \newblock \emph{Optimal design of a bar with an attached mass for maximizing 
the heat transfer},
    \newblock Electron. J. Diff. Equations, \textbf{2012} (2012), pp.~1--13.

\bibitem{BeMcCW}
    \newblock B.~P. Belinskiy, C.~M. McCarthy,  T.~J. Walters;
    \newblock \emph{The optimal design of a turbine},
    \newblock European Series in Applied and Industrial Mathematics: 
Control, Optimization and Calculus of Variations (ESAIM: COCV), 
\textbf{9} (2003), pp.~217--230.

\bibitem{Cardou}
    \newblock A. Cardou;
    \newblock \emph{Piecewise uniform optimum design for axial vibration requirement},
    \newblock AIAA J., \textbf{11} (1973), pp.~1760--1761.

\bibitem{Cox}
   \newblock S.~J. Cox;
   \newblock \emph{The two phase drum with the deepest base note},
   \newblock Japan J. of Industrial and Applied Mathematics, \textbf{8} (1991), pp.~345--355.

\bibitem{CoxMcCarthy}
   \newblock {S.~J. Cox, C.~M. McCarthy};
    \newblock \emph{The shape of the tallest column},
    \newblock SIAM J. on Mathematical Analysis, \textbf{29} (1998), pp.~1--8.

\bibitem{CoxOverton}
    \newblock {S.~Cox, M.~Overton};
   \newblock \emph{On the optimal design of columns against buckling},
   \newblock  SIAM J. on Math. Anal., \textbf{23} (1992), pp.~287-325.

\bibitem{Dym}
   \newblock C.~L. Dym;
   \newblock \emph{On some recent approaches to structural optimization},
   \newblock J. Sounds \& Vibration, \textbf{32} (1974), pp.~49--70.

\bibitem{Fulton}
    \newblock C.~T. Fulton;
    \newblock \emph{Two-point boundary value problems with eigenvalue parameter 
contained in the boundary conditions},
    \newblock Proc. Royal Soc. Edin., \textbf{77A} (1977), pp.~293--308.

\bibitem{Hinton}
    \newblock D. Hinton;
    \newblock \emph{Eigenfunction expansions for a singular eigenvalue problem 
with eigenparameter in the boundary condition},
    \newblock SIAM J. Math. Anal., \textbf{12} (1981), pp.~572--584.

\bibitem{13}
    \newblock F.~P. Incropera, D.~P. Dewitt;
    \newblock Fundamentals of Heat and Mass Transfer
     \newblock John Wiley $\&$ Sons, New York, 4th edition, 1996.

\bibitem{KellerNiordson}
   \newblock J.~B. Keller, F.~I. Niordson;
   \newblock \emph{The tallest column},
    \newblock J. Math. Mech., \textbf{16} (1966), pp.~433--446.

\bibitem{KrBo}
    \newblock F.  Kreith, M.~S. Bohn;
    \newblock Principles of Heat Transfer
    \newblock Harper $\&$ Row, 4th edition, 1986.

\bibitem{Lad}
    \newblock O.~A. Ladyzhenskaia;
    \newblock The Boundary Value Problems of Mathematical Physics
    \newblock Springer-Verlag, New York, 1985.

\bibitem{Powers}
    \newblock {D.~L. Powers};
    \newblock {Boundary Value Problems},
    \newblock Academic Press, New York, 2nd edition, 1979.

\bibitem{22}
    \newblock P.~J. Schneider;
    \newblock Conduction Heat Transfer
    \newblock Addison-Wesley, Reading, Mass., 1955.

\bibitem{Sheu}
    \newblock C.~Y. Sheu;
    \newblock \emph{Minimum-weight design for specified fundamental frequency},
    \newblock International J. of Solids and Structures, \textbf{4} (1968), 
pp.~953--958.

\bibitem{Sippel}
    \newblock D.~L. Sippel, W.~H. Warner;
    \newblock \emph{Minimum-weight design of sandwich structures under a 
frequency constraint},
    \newblock AIAA J., \textbf{11} (1973), pp.~483--489.

\bibitem{Taylor}
   \newblock J.~E. Taylor;
   \newblock \emph{The strongest column: the energy approach},
   \newblock J. of Appl. Mech., \textbf{34} (1967), pp.~486--487.

\bibitem{Taylor2}
   \newblock J.~E. Taylor;
   \newblock \emph{Minimum mass bar for axial vibrations at specified natural 
frequency},
   \newblock AIAA J., \textbf{5} (1967), pp.~1911--1913.

\bibitem{TaylorLiu}
   \newblock J.~E. Taylor, C.~Y. Liu;
   \newblock \emph{On the optimal design of columns},
   \newblock{AIAA J.}, \textbf{6} (1968), pp.~1497--1502.

\bibitem{Turner}
   \newblock M.~J. Turner;
   \newblock  \emph{Design of minimum mass structures with specified natural frequencies},
   \newblock AIAA J., \textbf{5} (1963), pp.~406--412.

\bibitem{Walter}
    \newblock J. Walter;
    \newblock \emph{Regular eigenvalue problems with eigenvalue parameter in the boundary condition},
    \newblock Math. Z., \textbf{133} (1973), pp.~301--312.

\end{thebibliography}

\end{document}