The mathematical underpinnings of vector analysis are reviewed as they are applied in the development of the ensemble of numeric statements for subsequent matrix solution. With the continued advances in computational power, there is increased interest in the field of atmospheric modelling to decrease the computational scale to a micro‐scale. This interest is partially motivated by the ability to solve large scale matrix systems in the number of occasions to enable a small‐scale time advancement to be approximated in a finite‐difference scheme. Solving entire large scale matrix systems several times a modelling second is now computationally feasible. Hence the motivation to increase computational detail by reducing modelling scale.
- vector fields
- Divergence Theorem
Computational Engineering Mathematics as applied to topics in Computational Geosciences and Computational Fluid Dynamics, among other themes and topics, is the foundation of computational modelling processes involved with Climate and Atmospheric processes. This chapter reviews the key mathematical underpinnings of computational mathematics used to represent climate, atmospheric and hydrologic related processes. The commonly used numerical and computational approaches of finite difference, finite element, and finite volume have common roots stemming from the mathematical analysis of partial differential equations, ordinary differential equations, vector differential and integral equations, among other topics including solution of large scale matrix systems involving tens of thousands of unknown variables, time‐stepped numerous times a modelling second, resulting in incredible numbers of computations. This chapter focuses on important aspects of matrix solutions and on the vector analysis towards use in development of finite volume and nodal domain numeric solutions of the governing Navier‐Stokes equations. Due to the continuing increasing computational capability forecast according to Moore’s Law, the use of Navier‐Stokes solvers in routine problems is becoming more prevalent, and will possibly be common practice in a few short years.
Methods for developing large scale matrix systems and their solution is a fundamental area of knowledge needed for the successful modelling of atmospheric processes and climate. This particular topic can be investigated in other references focused towards solving large scale matrix systems. Methods for efficiently constructing the matrix system for subsequent solution is also a subject for other references. Rather, in the current chapter, the focus is on the topics of matrix system solution and the existence and uniqueness of solutions. The reduction of the governing partial differential equations into numeric form or into finite element or finite volume form, and then further developed into numeric statements, is another subject dealt with in detail in other references.
In the current chapter, the mathematical underpinnings of vector analysis are reviewed as they are applied in the development of the ensemble of numeric statements for subsequent matrix solution. With the continued advances in computational power, there is increased interest in the field of atmospheric modelling to decrease the computational scale to a micro‐scale. This interest is partially motivated by the ability to solve large scale matrix systems in the number of occasions to enable a small‐scale time advancement to be approximated in a finite‐difference scheme. Solving entire large scale matrix systems several times a modelling second is now computationally feasible. Hence the motivation to increase computational detail by reducing modelling scale.
2. The Navier‐Stokes equations in computational engineering mathematics
The mathematical description of fluid processes are embodied in the well‐known Navier‐Stokes equations. These equations relate Newton’s second law of motion to a control volume of fluid. Their origins can be found in the work dated between 1827 and 1845, including the advances made by Saint‐Venant. Although all fluids possess viscous effects, in many problem situations the effects of viscosity can be neglected, reducing the governing equations considerably. The assumption of the fluid being incompressible in the range of temperatures and conditions under study further reduces the governing equations to being applicable to a fluid that is incompressible and having zero viscosity. That is, an ideal fluid. The analysis of ideal fluid flow is considerably simpler to undertake than the full Navier‐Stokes equations, and provides generally a good framework and basis for the analysis of the complete fluid flow mathematical regime. Ideal fluid flow ignores shear stresses in fluids and normal stresses (that is, thermodynamic pressure forces) are the internal force components carried forward into the analysis.
The resulting Euler equations are well known in fluid mechanics and are particularly amenable to solution when applied to flow on a streamline. For steady flow, the Euler equations indicate that flow velocity and pressure sum to as constant value along a horizontal streamline (or neglecting weight of the fluid control volume). Integration of the Euler equation along a streamline for steady flow conditions results in the Bernoulli equation, which is a commonly used relationship in the analysis of fluid flow. The Bernoulli equation is often misapplied by failure to adhere to the assumptions applied in the derivation of that equation. Namely, steady flow conditions, frictionless or shear‐free flow, incompressible fluid flow, and flow on a streamline. A common misapplication is to apply Bernoulli’s equation at two points within a fluid domain that do not lie on the same streamline. The sum of the Bernoulli equation components is called the Bernoulli constant, and has different values depending on the streamline under study. If the Bernoulli constant is constant throughout the entire fluid domain, regardless of streamline, then the flow is said to be irrotational. Irrotational flow can be assessed within a vector field by use of the vector operator Curl, defined in the development below. In particular, the streamline can be mathematically described in terms of vector notation with respect to usual measure of distance along the streamline, .
When dealing with fluid flow, particularly atmospheric fluid flow, thermodynamic pressure is a key component. Thermodynamic pressure is also called static pressure. The stagnation pressure is measured under an ideal frictionless process where the fluid is decelerated to zero speed (note that velocity is a vector whereas speed is a scalar term.) The dynamic pressure quantifies the magnitude of the flow velocity. The flow speed is determined given the stagnation pressure and the static pressure at a point. These relationships become impacted by the speed, particularly with higher Mach numbers greater than about 0.3.
Application of the Bernoulli equation finds good use in many problems of high interest involving fluid flow regimes involving contraction or expansion of flows. In all of these applications, however, care is required to preserve the fundamental assumptions employed in the derivation of the Bernoulli equation.
The above discussion focuses on the mass continuity and flow momentum equations embodied in the Navier‐Stokes equations. A similar resulting formulation arises using the energy equation. In this second situation, a control volume is determined that conforms to the stream lines of the flow regime, and steady incompressible flow is considered that is inviscid. The Bernoulli equation then arrives even though different concepts and boundary conditions are employed. Therefore, the Laws of Thermodynamics for steady, incompressible flow along a streamline simplifies to the well‐known Bernoulli equation which is also employed for analysis of mechanical energy fluid flow problems.
In civil engineering hydraulics problems, the energy equation formulation of the Bernoulli equation is typically shown graphically, as the Energy Grade Line or EGL displaying the sum of the three energy head components of pressure head plus velocity head plus gravity head, all in units of length. The sum of these three components is the Total Energy Head, which represents the height the fluid would rise up to in a tube open for the fluid to access. The Hydraulic Grade Line or HGL is plotted the velocity head below the EGL.
3. Irrotational flow and the vector operator curl
Another important vector calculus operator is divergence. This operator provides several measures of fluid flow and fluid properties including compressibility, and also a measure of the strength of source or sink conditions. A divergence value that is positive indicates outwards flux from the target point, whereas a negative value indicates flow trends indicative of a sink.
The vector operators of curl and divergence are key vector calculus tools used in the analysis of fluid flows of both compressible and incompressible flow regimes, and form the basis of continuity equations.
Modeling topics currently tend to divide along lines of numerical methods, particularly in leveraging parallelism versus increased processor performance. The partial differential equations describing the atmospheric processes involved still are heavily influenced by processes such as rainfall, convective processes, and others. Many models use grid‐based finite‐difference analogs, whereas other approaches are based upon finite element and finite volume schemes. Because many of these schemes involves grid density problems at the poles, further attention has been paid towards use of discretization's of the globe into different geometric finite volumes and tessellations. Finite volume methods, like finite element models, are attractive in conserving energy, mass, and momentum. Furthermore, the Divergence Theorem has a more direct application with use of finite volume or finite element methods.
In this chapter, the fundamental vector calculus principles are reviewed that are most relevant to fluid flow modeling such as involved with compressible and incompressible flow problems in water and the atmosphere. Details regarding these fields far surpass the scope of this chapter, but the fundamental elements are generally based in the computational mathematics involved in formulation of and solution of large scale matrix equations.
4.1. Definition and representation
A vector is used to describe a quantity such as displacement, velocity, or force that has both magnitude and direction. In two dimensions or three dimensions, we can represent the magnitude and direction of a vector with an arrow (Figure 1) .
By convention, a vector is usually written in bold, or with an arrow, such as or . We may also represent a vector as a displacement vector. For example, the displacement from point A to point B can be written as the vector (Figure 2).
It is important to note that a vector is not anchored to any particular place in the coordinate system. Rather, it represents only a
Most often when dealing with vectors of more than two components, it is easier to represent them algebraically. We write a vector as a list of its components, which are equivalent to the change in each of the coordinate directions (Figure 4). This is fairly easy to represent in two or even three dimensions, but also extends to dimensions.
A vector that denotes a specific point in the coordinate system is called a position vector (meaning a vector that, when its tail is placed at the origin, points to a specific point) Thus, a position vector for a point woud be . (Figure 5)
4.2. Vector addition and scalar multiplication
Two basic operations that can be performed on vectors are addition and scalar multiplication. In other words, we can add and subtract vectors, as well as multiply them by a constant which will lengthen or shorten the vector.
Algebraically, vector addition results in adding the vectors component‐wise:
In addition to representing a vector by its components, we can represent a vector by a linear combination of the basis vectors, a basis vector being a vector with a magnitude of one in the direction of only one of the coordinate axes.
Vector addition results in factoring each basis vector:
Algebraically, scalar multiplication results in multiplying each component of the vector by the scalar constant:
Vector subtraction results from applying both vector addition and scalar multiplication by Furthermore, from this we can derive from geometric or algebraic arguments the following properties of vectors:
Given , , are vectors and are scalars: 
The length of a vector can be found by utilizing the Pythagorean Theorem. If a vector has the components , then the length of
4.3. The dot product
Multiplication between vectors is not the same as between scalars. With vectors we have two “products”, the dot product and the cross product. The dot product between vectors will play a vital part in vector calculus. We can think of the dot product as being a metric of how much one vector is influencing another. An application for the dot product is work, where work is defined as . If we think of the force and distance as vectors, then the amount of work to move a sled between two points as one pulls on a rope is the dot product between the force and direction or distance vectors (Figure 7).
The dot product of two vectors is therefore a scalar (number), by definition, and not a vector. To find the dot product we multiply the corresponding components and add.
4.4. The cross product
The cross product is the other “product” between two vectors. The caveat is that the cross product is only defined for vectors in three dimensions. An application for the cross product is torque. This is similar to work which we saw with the dot product, but the cross product results in a vector orthogonal to the two vectors that created it.
The cross product of two vectors is therefore a vector (unlike the dot product). The cross product of two vectors and is a non‐zero vector that is
That is, the magnitude of the resultant vector from the cross‐product depends on the angle between them. The direction (while remaining orthogonal) depends on the order of the cross product, or
5. Parametric equations
An important concept within vector calculus is parameterization. In this we will create a vector function that can build a line, a curve (called a space curve) or a surface in space. The idea is that each component of the vector is a function of a single independent variable called the parameter. The parameter dictates the value of the component for each value in the domain of the parameter. The terminal point of this vector will “draw out” the line, curve or surface as the parameter varies through its domain.
5.1. Equation of a line
To parameterize a line we will create a vector equation in which is a position vector for all points on the line. The addition of two vectors creates a vector that will “draw” the line with its terminal point based upon values of a parameter. The equation is:
The parameter is , is a vector parallel to the line, and is a position vector to any single point on the line. A sketch of the mechanics of this equation is shown in Figure 9.
As varies, the terminal point of
If we wanted to parameterize a straight line between the points (0,1) to (2,1) then we could use as the parameter, and since the value is constant across the entire line, our vector equation becomes .
5.2. Equation of a plane
The vector equation for a plane is: . In this case
6. The gradient
An important calculation in vector calculus is the gradient. We will see that the gradient of a function will return a vector that points in the direction of greatest increase from any point in the function. First we will inspect where the gradient comes from.
6.1. Partial derivatives
As we move into three dimensions we see our functions represented as . Now our functions are of two variables. Remembering that a derivative is the rate of change of a function at a specific point, the partial derivative with respect to x or y becomes the rate of change at a point in either the x or y direction. There are several notations for partial derivatives .
The partial derivatives are found, largely, the same way as single derivatives. The only new aspect is the other variable, which we treat as a constant when taking the derivative.
Find , where
6.2. Directional derivative
Now that we know how to find the rate of change in either the x or y direction, what about the other infinite directions we may be interested in? We will utilize a directional derivative to find the rate of change from a point in a particular direction.
The directional derivative appears to be a dot product between two vectors. In fact if we assume that it is and write the function as a dot product between two vectors the resulting equation reveals:
We immediately see our vector
6.3. The gradient
When talking about the gradient we will introduce an operator called the dell operator, which is similar to an upside down triangle. The del operator is a vector.
The gradient is the product of the del operator and the function
Thus another way to annotate a directional derivative would be:
It is important to note that the gradient of a function will point in the direction of greatest increase of a function at a specific point. The magnitude of the gradient at that point will be the greatest rate of change. When viewing a function as its level curves. The gradient will be orthogonal to the curve at any point.
Find the maximum rate of change of at (2,1)
6.4. Tangent planes
In a similar way to finding a line tangent to a curve at a point, we can find a plane tangent to a surface at a point. The idea of a tangent plane will be applied later in the chapter.
Suppose that we have a surface and we want to find the plane tangent to the surface at the point . The equation is:
If we think of the differences between the values of , , as the change in the , , and coordinates, then as that difference approaches zero, we can use the differentials for x,y,z.
7. Vector fields
We saw that a vector has both a magnitude and a direction. When we assign a vector to every point in space we build a vector field. A practical example in the real world would be wind speeds. It is easy to see that at every point the wind has both a velocity and a direction.
Perhaps, the most important attribute of a vector field is the curl. The curl will determine if there is some sort of twist or spin within the vector field. This will become very important when we start to discuss work through a vector field. The curl is determined by taking the cross product between the del operator and the vector field.
To find the actual curl vector we will find the determinate of the matrix :
If the curl results in a zero vector, then the vector field is determined to have zero curl and the vector field is conservative.
When we compare the plots of the two similar vector fields in Figure 10, we notice a subtle difference between them. In the vector field on the left there are vectors pointing towards the center as well as vectors pointing away from the center. This indicates curl. In the vector field on the right, all points are pointing away from the center which indicates that there is no curl in the vector field.
The divergence of a vector field can be thought of as a quantification of the rate that the vector field is expanding or contracting per unit volume at any point in the field . If the divergence is equal to zero the vector field is said to be incompressible. Divergence is calculated by taking the dot product between the del operator and the vector field .
If the divergence at a point is positive the point is considered a source. If the divergence is negative the point is a sink. If the divergence is zero, then the point is neither a sink nor a source.
8. Line integrals
If we remember back to the dot product, we defined work as . We now think of a particle moving through a vector field along a curve between two points described by the vector The work done by the vector field on the particle as it is moving between two points can be thought of as the dot product between the vector field and the direction vector of the particle. If we inspect a small segment of the curve, the work between two points can be approximated by taking the dot product of
If is the length of the curve between two points, then the approximation of work over the entire curve would be the sum of work in all subsections
As goes to zero we have an infinite sum, or integral .
8.1. The fundamental theorem for line integrals
We talked earlier about conservative vector fields, where the curl was zero. There is an important implication to the curl being zero; there exists a function such that the gradient of that function is equal to the vector field. This function is called the potential function.
Clearly, if the path is closed, where the start point is the same as the end point, the work will be zero. Another implication to the conservative vector field is path independence. This means that the work done by a conservative vector field on an object moving between two points is the same regardless of the path.
9. Green's Theorem
In essence, Green's Theorem allows us to calculate the work done around a closed path through a non‐conservative vector field by looking at the area enclosed by the path rather than the path itself.
We can recognize in this integral that we are translating the line integral into a double integral of the curl in the x‐y plane over the area enclosed by the curve. The positive orientation indicates that we are moving in a counter clockwise direction around the curve. Green's Theorem can also be applied to the plane and the plane if the curve or path is in either of those two planes. Green's Theorem is especially useful if the path encloses an area with a hole in it. In this case we are still integrating the curl over the area.
Evaluate the given line integral around the path defined by the circles:
We are going to utilize the area between the two circles to solve this line integral.
By Green's Theorem:
Converting to polar coordinates:
10. Surface area
Previously we saw how to parameterize a line in which the
Finding the surface area of a flat or geometric surface is relatively easy as there are known formulas for such cases. However, it is more difficult to find the surface area of a more complexly defined surface. In this case we will see that the relationship between surface area and surface integrals is similar to the relationship between arc length and line integrals.
If we have a space curve defined by the vector function where , the length of the of the arc between and is found by :
If the curve is in the x‐y plane where then
In a similar sense to the arc length, where we utilized a secant line to approximate the length of a curve, we can utilize a tangent plane to a surface to approximate the surface area.
The primary idea is that to get as close to the true surface area as possible we must use an ever increasing number of planes approximating our surface. Then as the number of planes goes towards infinity the sum of their areas approaches the area of the surface.
It is important to remember how to define a plane. Two ways that will be important to surface area are finding two vectors within the plane and by finding a point with a normal vector to the plane. In order to find two vectors within the plane we can utilize the resulting vectors from the directional derivatives of our surface at our point within the surface. With the parameterization of the surface the vectors would become and . To find the normal vector to the plane containing the two directional derivative we will take their cross product.
In order to find the surface area of the plane we can scale each vector by the change of the parameter of the subinterval. We then take the cross product of the scaled vectors and by the rules for the cross product we can factor out the change in variables such that:
If we then sum the area of each plane, the approximation to the surface becomes:
Then the surface area for our surface S is: 
If the surface is defined in terms of
then the formula becomes: 
11. Surface integrals
In this section we will see that the surface integral is very similar to what a line integral is to curves in space.
It is clear that each coordinate is a function of our two parameters. The parameterized surface, , is considered smooth if both of the partial derivatives are continuous and the cross product, , is never zero on the interior of the parameter domain. This cross product is called the standard normal to . .
The first step is to assume the parameter domain,
Our surface is similarly divided into corresponding subsections, each defined by . If we form the Riemann sum of the product of a function evaluated at some point in each subregion of our surface with the area of each subregion, . The resulting summation is: 
As the area goes to zero and
Let's take a second to remember what the
If we remember back to how surface area was calculated utilizing an infinite number of tangent planes, then the change in the surface area clearly becomes
Remember that is the length of a vector normal to the surface. This allows us to compute a surface integral by converting it into a double integral over the parameter domain
Then the resulting surface integral becomes: 
Similar formulas can be easily derived should the surface be defined in a way other than by .
12. Oriented surfaces
It is important to talk about the “orientation” of a surface. It is easy to see that there are always two vectors that are normal to any surface as they point in exactly opposite directions. We can think about a “positive” orientation being one defined by a normal vector pointing away from the origin and a “negative” orientation being defined by a normal vector pointing towards the origin. In vector calculus, the only surfaces that are used are those that are orientable. A type of surface that is not orientable, such as the Mobius Strip, is not applicable within vector calculus.
13. Surface integrals and flux
Suppose that we have an oriented surface,
The rate of flow (mass per unit time) per unit area is . If we divide our surface into sub regions defined by where each subregion is small and nearly rectangular we can approximate the mass or amount of fluid crossing in the direction of
If we add this mass for each sub section as the area of each subsection goes to 0 it results in an infinite sum or integral which calculates the rate of flow through
This type of surface integral occurs frequently in physics. If the density function is constant then we are only left with the velocity vector field.
More simply, if
We are then able to translate this integral into our previously shown double integral over a region such that: 
With regards to the actual mechanics of how this integral is calculated it is important to have an understanding of the relationship between the parts. It is also important to pay attention to how the surface is defined. Let's assume is our surface. We must express the surface as a function . Then, . The resulting normal vector of the surface becomes:
The differential for surface area is:
The resulting integral then becomes:
The denominator within the dot product and the radical within the surface differential will cancel which leaves us with the final integral: 
14. Stokes Theorem
Stokes’ Theorem is very similar to Green's Theorem. In fact it can be thought of that Green's Theorem is a special case of Stokes’ Theorem, or that Stokes’ Theorem is a higher‐dimensional version of Green's Theorem. What Stokes is going to do is relate a line integral with a surface integral. Stokes’ Theorem will find the work done around a closed path by computing a surface integral of the curl of a vector field over a “surface cap ” on top of the path.
In other words the circulation of a vector field around the boundary of an oriented surface in space in the direction counter clockwise with respect to the surface's unit normal vector equals the integral of the curl for vector field over the normal component of the surface. The surface must also be piecewise smooth.
We can think of this integral as relating a line integral around a closed path to a surface integral over a “capping surface” of the path. The interesting idea is that it doesn't matter what the surface is over the path, once the surface integral is calculated it reaches the same value!
While this result is difficult to prove in the general case, an easier intuition is gained through inspection of Green's Theorem demonstrated earlier. This can then be extended to Stokes, which is the same concept in three dimensions.
Evaluate where is the curve of intersection of the paraboloid and the plane , oriented up, through .
First check the curl to see if the vector field is conservative.
, not conservative
We are going to utilize Stokes’ Theorem to solve this problem. Remembering that we can define our surface over the path as anything, let's use the plane .
If we project the surface into the x‐y plane it will form a circle. By making the substitution 5 for then , and we see our area is a circle of radius 2. We will convert into polar coordinates.
Further substituting 5 for becomes
15. Divergence Theorem
In this final section we move away from line integrals and again visit flux. In this case we will be looking at how a vector field moves through a solid shape in space.
The divergence theorem states: Let be a simple solid region and let
In other words the Divergence Theorem states that the flux of
The Divergence Theorem is to surface integrals what Green's Theorem is to line integrals in that it allows us to convert a surface integral over a closed surface into a triple integral over a closed region.
While the equation looks to be a direct equivalence between the double and triple integrals, it is not quite that simple. If we were to try and utilize the surface integral definition for flux on a closed surface, we would have to compute the flux through each sub surface individually then add the result. In contrast, the divergence theorem allows us to find the flux through the closed surface by utilizing the divergence through the volume.
Find the outward flux of
If we project our region into the x‐y plane, we see that it is a quarter circle with a radius of 2. In this problem we will utilize cylindrical coordinates to solve the integral.