Engee documentation
Notebook

Finding the tangent plane

We will show how to approximately calculate the gradient of a function using finite differences and construct a tangent plane to the surface at a given point.

In [ ]: 
Pkg.add( "Zygote" )

Let's set the initial function

Let's create a surface function for further analysis:

In [ ]: 
f( xy ) = xy[1].^2 .+ xy[2].^2;

In this example, it is important for us that all parameters are passed to the function using a vector

To get the parameters of the tangent plane at the point of interest, we need to take the derivative of this function.

In Julia, there are several ways to get the derivative of functions, for example, using auto-differentiation and the library Zygote.

In [ ]: 
using Zygote
fx, fy = f'( [-1, 3] )

This library can differentiate a very wide range of functions in the Julia language: for example, their differentiable code can contain complex numbers, loops, recursion, conditional operations, and command calls from standard libraries (for example, FFTW).

Equation of the tangent plane at a point you can set it as follows:

Coefficients and we got by using f'( 1,2 ). Setting the function ∇f(x,y), which will allow us to find the values of the tangent plane at any point:

In [ ]: 
x₀,y₀ = -1,3
fx,fy = f'( [x₀,y₀] )
∇f(x,y) = f( [x₀,y₀] ) .+ fx.*(x.-x₀) .+  fy.*(y.-y₀)
Out[0]:
∇f (generic function with 1 method)

Output the original function f(x,y) and the tangent plane at the point P = [1,2]:

In [ ]: 
gr()

x_range = -5:0.25:5
y_range = -5:0.25:5
xx,yy = repeat( x_range, inner = length(y_range)), repeat( y_range, outer = length(x_range));
Nx = Ny = round( Int32, sqrt(length(xx)) );
xm = reshape( xx, Nx, Nx );
ym = reshape( yy, Ny, Ny );

wireframe( x_range, y_range, ∇f(xm,ym), cbar=false )
surface!( x_range, y_range, f([xm, ym]), c=:viridis, cbar=false )
scatter!( [x₀], [y₀], [f([x₀, y₀])], c=:white, leg=false )
Out[0]:

Let's look at this graph from a different perspective.

In [ ]: 
plot!( camera=(-86, 9) )
Out[0]:

Numerical solution (diff function)

Now let's assume that we only have data, not a function, and we need partial derivatives in both coordinates. We obtain an approximate value of the derivatives at each point of a certain coordinate grid:

In [ ]: 
dfx = diff( f( [xm, ym] ), dims=2 );
dfy = diff( f( [xm, ym] ), dims=1 );

plot( surface(dfx), surface(dfy), cbar=false, size=(500,250) )
Out[0]:

To build a tangent at a point (x₀, y₀) – one of the points on the coordinate grid – find the index of this point:

In [ ]: 
t = (xm .== x₀) .& (ym .== y₀)
indt = findfirst( t )
Out[0]:
CartesianIndex(33, 17)
In [ ]: 
dT = yy[2]-yy[1] # Расстояние между всеми точками: 0.25
fx0 = dfx[ indt ] ./ dT;
fy0 = dfy[ indt ] ./ dT;

We get the same tangent graph again:

In [ ]: 
∇f2( x,y ) = f([x₀,y₀]) .+ fx0.*(x.-x₀) .+  fy0.*(y.-y₀)

wireframe( x_range, y_range, ∇f2(xm,ym) )
surface!( x_range, y_range, f([xm, ym]), c=:viridis, cbar=false )
plot!( camera=(-88, 8.5) )
Out[0]:

Conclusion

We have constructed the tangent plane to the function in two different ways:

  • by differentiating the function code in Julia language (at any point),
  • using numerical differentiation (on a certain coordinate grid).