# 13.1. Tensor-product discretizations¶

As you learned when starting double integration in vector calculus, the simplest extension of an interval to two dimensions is a rectangle. We will use a particular notation for rectangles:

(13.1.1)$[a,b] \times [c,d] = \bigl\{ (x,y)\in\mathbb{R}^2 : a\le x \le b,\; c\le y \le d \bigr\}.$

The $$\times$$ in this notation is called a tensor product, and a rectangle is the fundamental example of a tensor-product domain. The implication of the tensor product is that each variable independently varies over a fixed set. The simplest three-dimensional tensor-product domain is the cuboid $$[a,b]\times[c,d]\times[e,f]$$. When the interval is the same in each dimension (that is, the region is a square or a cube), we may write $$[a,b]^2$$ or $$[a,b]^3$$. We will limit our discussion to two dimensions henceforth.

The discretization of a two-dimensional tensor-product domain is straightforward.

Definition 13.1.1 :  Tensor-product grid

Given discretizations of two intervals,

(13.1.2)$a= x_0< x_1 < \cdots < x_m = b, \qquad c = y_0 < y_1 < \cdots < y_n = d,$

then a tensor-product grid on $$[a,b]\times[c,d]$$ is the set

(13.1.3)$\bigl\{ (x_i,y_j): i=0,\ldots,m,\; j=0,\ldots,n \bigr\}.$

## Functions on grids¶

The double indexing of the grid set (13.1.3) implies an irresistible connection to matrices. Corresponding to any function $$f(x,y)$$ defined on the rectangle is an $$(m+1)\times(n+1)$$ matrix $$\mathbf{F}$$ defined by collecting the values of $$f$$ at the points in the grid. This transformation of a function to a matrix is so important that we give it a formal name:

(13.1.4)$\begin{split}\mathbf{F} = \mtx(f) = \Bigl[f(x_i,y_j)\Bigr]_{\substack{i=0,\ldots,m\\j=0,\ldots,n}}.\end{split}$

Caution

There is potential for confusion because the first dimension of a matrix varies in the vertical direction, while the first coordinate $$x$$ varies horizontally. In fact, the Julia plotting routines we use expect the transpose of this arrangement, so that $$x$$ varies along columns and $$y$$ along rows.

Example 13.1.2

Let the interval $$[0,2]$$ be divided into $$m=4$$ equally sized pieces, and let $$[1,3]$$ be discretized in $$n=2$$ equal pieces. Then the grid in the rectangle $$[0,2]\times[1,3]$$ is given by all points $$(i/2,1+j)$$ for all choices $$i=0,1,2,3,4$$ and $$j=0,1,2$$. If $$f(x,y)=\sin(\pi xy)$$, then

$\begin{split} \mtx(f) = \begin{bmatrix} \sin(\pi\cdot 0\cdot 1) & \sin(\pi\cdot0\cdot 2) & \sin(\pi\cdot0\cdot 3) \\[1mm] \sin\left(\pi\cdot\tfrac{1}{2} \cdot 1 \right) & \sin\left(\pi\cdot\tfrac{1}{2} \cdot 2 \right) & \sin\left(\pi\cdot\tfrac{1}{2} \cdot 3 \right) \\[1mm] \sin\left(\pi \cdot 1 \cdot 1 \right) & \sin\left(\pi \cdot 1 \cdot 2 \right) & \sin\left(\pi \cdot 1 \cdot 3 \right) \\[1mm] \sin\left(\pi\cdot \tfrac{3}{2} \cdot 1 \right) & \sin\left(\pi\cdot\tfrac{3}{2} \cdot 2 \right) & \sin\left(\pi\cdot\tfrac{3}{2} \cdot 3 \right) \\[1mm] \sin\left(\pi \cdot 2 \cdot 1 \right) & \sin\left(\pi \cdot 2 \cdot 2 \right) & \sin\left(\pi \cdot 2 \cdot 3 \right) \end{bmatrix} = \begin{bmatrix} 0 & 0 & 0 \\ 1 & 0 & -1 \\ 0 & 0 & 0 \\ -1 & 0 & 1 \\ 0 & 0 & 0 \end{bmatrix}. \end{split}$
Demo 13.1.3

Here is the grid from Example 13.1.2.

m = 4;   x = range(0,2,length=m+1);
n = 2;   y = range(1,3,length=n+1);


For a given $$f(x,y)$$ we can find $$\operatorname{mtx}(f)$$ by using a comprehension syntax.

f = (x,y) -> cos(π*x*y-y)
F = [ f(x,y) for x in x, y in y ]

5×3 Matrix{Float64}:
0.540302  -0.416147  -0.989992
0.841471   0.416147  -0.14112
-0.540302  -0.416147   0.989992
-0.841471   0.416147   0.14112
0.540302  -0.416147  -0.989992


The plots of this section look better using a different graphics engine on the back end:

plotlyjs();  # use better 3D renderer

ArgumentError: Package PlotlyJS not found in current path:
- Run import Pkg; Pkg.add("PlotlyJS") to install the PlotlyJS package.

Stacktrace:
[1] require(into::Module, mod::Symbol)
[2] top-level scope
@ ~/.julia/packages/Plots/8K4be/src/backends.jl:663
[3] eval
@ ./boot.jl:373 [inlined]
[4] _initialize_backend(pkg::Plots.PlotlyJSBackend)
@ Plots ~/.julia/packages/Plots/8K4be/src/backends.jl:662
[5] backend(pkg::Plots.PlotlyJSBackend)
@ Plots ~/.julia/packages/Plots/8K4be/src/backends.jl:176
[6] #plotlyjs#248
@ ~/.julia/packages/Plots/8K4be/src/backends.jl:31 [inlined]
[7] plotlyjs()
@ Plots ~/.julia/packages/Plots/8K4be/src/backends.jl:31
[8] top-level scope
@ In[4]:1
[9] eval
@ ./boot.jl:373 [inlined]
[10] include_string(mapexpr::typeof(REPL.softscope), mod::Module, code::String, filename::String)


We can make a nice plot of the function by first choosing a much finer grid. However, the contour and surface plotting functions expect the transpose of mtx($$f$$).

To emphasize departures from a zero level, use a colormap such as redsblues, and use clims to set balanced color differences.

m = 70;   x = range(0,2,length=m+1);
n = 50;   y = range(1,3,length=n+1);
F = [ f(x,y) for x in x, y in y ];

plot(x,y,F',levels=10,fill=true,aspect_ratio=1,
color=:redsblues,clims=(-1,1),
xlabel="x",ylabel="y")

surface(x,y,F',l=0,leg=:none,
color=:redsblues,clims=(-1,1),
xlabel="x",ylabel="y",zlabel="f(x,y)")


## Parameterized surfaces¶

We are not limited to rectangles by tensor products. Many regions and surfaces may be parameterized by means of $$x(u,v)$$, $$y(u,v)$$, and $$z(u,v)$$, where $$u$$ and $$v$$ lie in a rectangle. Such “logically rectangular” surfaces include the unit disk,

(13.1.5)\begin{split}\left\{ \begin{aligned} x &= u \cos v, \\ y &= u \sin v,\\ \end{aligned} \right. \qquad \qquad \left. \begin{aligned} 0 & \le u < 1, \\ 0 &\le v \le 2\pi, \end{aligned} \right.\end{split}

and the unit sphere,

(13.1.6)\begin{split}\left\{ \begin{aligned} x &= \cos u \sin v,\\ y &= \sin u \sin v,\\ z &= \cos v, \end{aligned} \right. \qquad \qquad \left. \begin{aligned} 0 & \le u < 2\pi, \\ 0 &\le v \le \pi. \end{aligned} \right.\end{split}
Demo 13.1.4

For a function given in polar form, such as $$f(r,\theta)=1-r^4$$, construction of a function over the unit disk is straightforward using a grid in $$(r,\theta)$$ space.

r = range(0,1,length=41)
θ = range(0,2π,length=81)
F = [ 1-r^4 for r in r, θ in θ ]

surface(r,θ,F',legend=:none,l=0,color=:viridis,
xlabel="r",ylabel="θ",title="A polar function")


Of course, we are used to seeing such plots over the $$(x,y)$$ plane, not the $$(r,\theta)$$ plane. For this we create matrices for the coordinate functions $$x$$ and $$y$$.

X = [ r*cos(θ) for r in r, θ in θ ]
Y = [ r*sin(θ) for r in r, θ in θ ]

surface(X',Y',F',legend=:none,l=0,color=:viridis,
xlabel="x",ylabel="y",title="Function on the unit disk")