Derivation of Divergence of a Vector Field Formula

Question

Assume a vector field $\boldsymbol{V}=M(x,y)\boldsymbol{i}+N(x,y)\boldsymbol{j}$ The divergence is then: $$\frac{\partial M}{\partial x}+\frac{\partial N}{\partial y}$$

My book(Thomas's calculus) explains this using a square in the x-y plane whose area goes to zero which is similar to this answer How would one arrive at the formulas for divergence and curl? Following the same argument, why is the change in the x For N (i.e $N(x+\Delta x,y))$) as we move along the x-axis not included, similarly for the change in y for M.

i.e why are these not included in the derivation of the divergence formula?

$$(-N(x,y)+N(x+\Delta x,y))\Delta x$$ and $$(-M(x,y)+M(x,y+\Delta y))\Delta y$$

Answer 1

So, as you go around the little rectangle, the flux is (approximately) the sum of the $\mathbf V\cdot\mathbf n \Delta s$ contributions. Starting with the bottom edge and proceeding counterclockwise, we have \begin{align*} \text{FLUX} &\approx\mathbf V\cdot (-\mathbf j) \Delta x + \mathbf V\cdot\mathbf i \Delta y + \mathbf V\cdot\mathbf j \Delta x + \mathbf V\cdot(-\mathbf i)\Delta y\\ &= -N(x,y)\Delta x + M(x+\Delta x,y)\Delta y + N(x,y+\Delta y)\Delta x + M(x,y)\Delta y \\&= \left(M(x+\Delta x,y)-M(x,y)\right)\Delta y + \left(N(x,y+\Delta y)-N(x,y)\right)\Delta x \\ &\approx \left(\frac{\partial M}{\partial x} + \frac{\partial N}{\partial y}\right)\Delta x\Delta y. \end{align*}

The terms to which you refer don't show up because they would be higher order in $\Delta x$. As we move along the bottom edge, of course $x$ is varying, but this variation will disappear in the limit. By the intermediate value theorem for integrals, the bottom edge will actually be $-N(\xi,y)\Delta x$ for some $\xi$ satisfying $x\le\xi\le x+\Delta x$. We then approximate this as $$-N(\xi,y)\Delta x \approx -\big(N(x,y)+\left(\frac{\partial N}{\partial x}\right)(\xi-x)\big)\Delta x.$$ Since $\xi-x \le \Delta x$, this term will contribute something bounded by $(\Delta x)^2$. When we combine this with the corresponding term from the upper edge, we'll have something that looks like $$\frac{\partial}{\partial y}\left(\frac{\partial N}{\partial x}\right) (\Delta x)^2\Delta y,$$ and such terms will disappear in the limit when we convert to a double integral.

Answer 2

I'd like to offer an alternative explanation as to why the divergence is equal to the sum of the partial derivatives of $M(x,y)$ and $N(x,y)$ in the $x$ and $y$ directions respectively. This will be done alongside with a pictorial representation of the same arbitrarily "small" rectangle depicted below.

Intuitively (and also formally provable), the net flux can be split into four parts: $$\text{net flux} = \text{right side flux}+\text{left side flux}+\text{top side flux}+\text{bottom side flux}$$

This is because the boundary of the rectangle $\text{Rect}$ (denoted $\partial \text{Rect}$) consists of all the points on the "outline" of this rectangle.

More formally,

$$\text{net flux}=\oint_{\partial \text{Rect}} \vec{V}\cdot \hat{n}\;\mathbb{d}s=\int_{R} \vec{V}\cdot \hat{n}\;\mathbb{d}s+\int_{L} \vec{V}\cdot \hat{n}\;\mathbb{d}s+\int_{T} \vec{V}\cdot \hat{n}\;\mathbb{d}s+\int_{B} \vec{V}\cdot \hat{n}\;\mathbb{d}s$$

Note that this "decomposition" of $\partial \text{Rect}$ makes easy work of using $\vec n$ on each subset $R,L,T$ and $B$ (assuming $\vec n$ points outwards from $\text{Rect}$).

Following the coordinate system provided, every line integral is expressible as a definite integral over the intervals representing each side.

$$F_R=\int_{R} \vec{V}\cdot \hat{n}\;\mathbb{d}s=\int_{y}^{y+\Delta y} \left(M(x+\Delta x,Y)\hat i+N(x+\Delta x,Y)\hat j\right)\cdot \hat i\;\mathbb{d}Y=\int_{y}^{y+\Delta y}M(x+\Delta x,Y)\;\mathbb{d}Y$$ Similarly for the left side L (but now noting that $\hat n=-\hat i$) $$\int_{L} \vec{V}\cdot \vec{n}\;\mathbb{d}s=-\int_{y}^{y+\Delta y}M(x,Y)\;\mathbb{d}Y$$

Likewise for the top and bottom sides $T$ and $B$, $$\int_{T} \vec{V}\cdot \hat{n}\;\mathbb{d}s=\int_{x}^{x+\Delta x} \left(M(X,y+\Delta y)\hat i+N(X,y+\Delta y)\hat j\right)\cdot \hat j\;\mathbb{d}X=\int_{x}^{x+\Delta x}N(X,y+\Delta y)\;\mathbb{d}X$$ Similarly for the bottom side B (but now noting that $\hat n=-\hat j$) $$\int_{B} \vec{V}\cdot \hat{n}\;\mathbb{d}s=-\int_{x}^{x+\Delta x}N(X,y)\;\mathbb{d}X$$

Therefore, we have that the net flux on $\partial R$ is $$\text{net flux}=\int_{y}^{y+\Delta y}M(x+\Delta x,Y)\;\mathbb{d}Y-\int_{y}^{y+\Delta y}M(x,Y)\;\mathbb{d}Y+\\\int_{x}^{x+\Delta x}N(X,y+\Delta y)\;\mathbb{d}X-\int_{x}^{x+\Delta x}N(X,y)\;\mathbb{d}X\tag{1}$$

Relationship to the book's approximation

Note that at this point it can be seen that the (approximation of the) flux on the right side R used in the book is the first term of the Taylor expansion of the exact flux, $F_R(\Delta y)$ (can be calculated using the Leibniz integral rule). Written using Big-O notation, the flux on the right side R is

$$\int_{y}^{y+\Delta y}M(x+\Delta x,Y)\;\mathbb{d}Y=M(x+\Delta x,y)\Delta y+\mathcal{O}(\Delta y^2)\tag{2} \quad \quad \text{as } \Delta y \to 0$$

At this point we could say that they are approximately equal to one another. However, there is no fine line for a small enough $\Delta y$ besides one's judgement or application: $$\int_{y}^{y+\Delta y}M(x+\Delta x,Y)\;\mathbb{d}Y\approx M(x+\Delta x,y)\Delta y$$

Below is the corresponding picture:

Note that they would be equal if $M(x+\Delta x,y)$ was constant for all $y\in [y,y+\Delta y]$, but that is not necessarily the case.

Relationship between this rectangle construction and the divergence of V at a point (x,y)

Combining integrals with "like" bounds in (1) you obtain: $$\oint_{\partial \text{Rect}} \vec{V}\cdot \hat{n}\;\mathbb{d}s=\int_{y}^{y+\Delta y}M(x+\Delta x,Y)-M(x,Y)\;\mathbb{d}Y+\int_{x}^{x+\Delta x}N(X,y+\Delta y)-N(X,y)\;\mathbb{d}X$$

Applying the definition of the divergence for this rectangle construction we have:

$$\text{div}\;\vec V\;\Biggr|_{(x,y)}=\lim_{(\Delta x,\Delta y)\to \vec 0}\frac{1}{\Delta x\Delta y}\oint_{\partial \text{Rect}} \vec{V}\cdot \hat{n}\;\mathbb{d}s$$

Substituting in the expression for the net flux for the rectangle:

$$\small{\begin{align}\text{div}\;\vec V\;\Biggr|_{(x,y)}&=\lim_{(\Delta x,\Delta y)\to \vec 0}\frac{1}{\Delta x\Delta y}\left(\int_{y}^{y+\Delta y}M(x+\Delta x,Y)-M(x,Y)\;\mathbb{d}Y+\int_{x}^{x+\Delta x}N(X,y+\Delta y)-N(X,y)\;\mathbb{d}X\right)\\&=\underbrace{\lim_{(\Delta x,\Delta y)\to \vec 0}\frac{1}{\Delta x\Delta y}\left(\int_{y}^{y+\Delta y}M(x+\Delta x,Y)-M(x,Y)\;\mathbb{d}Y\right)}_{A}\\&+\underbrace{\lim_{(\Delta x,\Delta y)\to \vec 0}\frac{1}{\Delta x\Delta y}\left(\int_{x}^{x+\Delta x}N(X,y+\Delta y)-N(X,y)\;\mathbb{d}X\right)}_B\end{align}}$$

Using iterated limits, term $A$ can be evaluated: $$\begin{align}A&=\lim_{(\Delta x,\Delta y)\to \vec 0}\frac{1}{\Delta x\Delta y}\left(\int_{y}^{y+\Delta y}M(x+\Delta x,Y)-M(x,Y)\;\mathbb{d}Y\right)\\&=\lim_{\Delta x\to 0}\lim_{\Delta y\to 0}\frac{1}{\Delta x\Delta y}\left(\int_{y}^{y+\Delta y}M(x+\Delta x,Y)-M(x,Y)\;\mathbb{d}Y\right)\\&=\lim_{\Delta x\to 0}\frac{1}{\Delta x}\lim_{\Delta y\to 0}\frac{1}{\Delta y}\left(\int_{y}^{y+\Delta y}M(x+\Delta x,Y)-M(x,Y)\;\mathbb{d}Y\right)\\&=\lim_{\Delta x\to 0}\frac{1}{\Delta x}\left(M(x+\Delta x,y)-M(x,y)\right)\end{align}$$

Going from line 3 to 4 can be achieved in at least two ways I'm aware of. Either by using the Taylor expansion as discussed in (2) and taking the limit (all higher order terms disappear), or by realising that:

$$\lim_{h\to 0}\frac{1}{h}\int_{a}^{a+h}f(x)\;\mathbb{d}x=\lim_{h\to 0}\left(\text{average value of } f \text{ on } [a,a+h]\right)=f(a)$$ when $f(x)$ is continuous. In our case the inside of the integral is assumed continuous.

Notice that this is precisely the definition of the partial derivative of $M$ at $(x,y)$, $$\therefore A=\frac{\partial M}{\partial x}$$

Applying the same process for term $B$, but in the opposite order for iterated limits, you obtain that

$$B=\frac{\partial N}{\partial y}$$

Therefore, the result is obtained: $$\text{div} \vec V=\frac{\partial M}{\partial x}+\frac{\partial N}{\partial y}\quad \blacksquare$$

This is a result that holds for all $(x,y)$ in the domain of $\vec V(x,y)$.

Remark on transitioning from integrals over boundary (line integral) to definite integrals with bounds:

This is a standard result that takes into consideration the paramaterisation of the curve. Assuming the curve R (or L) is paramaterised by arclength, then the transition from a line integral to definite integral (with bounds) is straightforward (simply add the bounds and change the differential to $\mathbb{d}Y$ as done here).