5.5 How to plot a function on 2D using coordinates of grid locations?

Sometimes in finite difference we want to evaluate the force or the function in the RHS of the equation \(f(x,y)\) using the coordinates \(x,y\) on a 2D grid. The grid spacing can be \(h\) and the coordinates extend from some \(x\) value to some other \(x\) value and from some \(y\) value to some other \(y\) values. In Matlab, the coordinates are setup using meshgrid.

In Mathematica, there is really no need for meshgrid as the build-in command Map can be used to evaluate the function at each grid point. The coodinates physical values are generated using the Table command. Here is example evaluating and plotting \(f(x,y)= \sin (10 x y) e^{-x y}\) on a grid of spacing \(h=1/N\) where \(N\) is number of elements. The grid goes from \(x=-1\) to \(x=1\) by spacing of \(h\). Same for the \(y\) direction.

Mathematica

Clear[x, y, z]; 
nElem = 10.; 
h = 1/nElem 
N@Range[-1, 1, h] 
grid = N@Table[{i, j}, {i, -1, 1, h}, {j, -1, 1, h}]; 
f[x_, y_] := Sin[x*10 y] Exp[-x y]; 
force = Map[f[#[[1]], #[[2]]] &, grid, {2}]; 
ListPlot3D[force, AxesLabel -> {x, y, z}]
 

To use the actual physical coordinates on the axis above, then replace the plot above with this below. This also includes an implementation of meshgrid like function in Mathematica which makes it easier to work with for someone familiar with Matlab:

meshgrid[x_List,y_List]:={ConstantArray[x, Length[x]], 
  Transpose@ConstantArray[y, Length[y]]}; 
 
nElem = 10.; 
h = 1/nElem; 
{x, y} = meshgrid[Range[-1, 1, h], Range[-1, 1, h]]; 
f[x_, y_] := Sin[x*10 y] Exp[-x y]; 
force = f[x, y]; 
pts = Flatten[{x, y, force}, {2, 3}]; 
ListPlot3D[pts, PlotRange -> All, 
      AxesLabel -> Automatic, ImagePadding -> 20, 
      ColorFunction -> "Rainbow", Boxed -> False]
 

To see each point value, use  InterpolationOrder -> 0 in the plot command.

Matlab clear all; close all; nElem = 10; h = 1/nElem; [X,Y] = meshgrid(-1:h:1,-1:h:1); f = @(x,y) sin(x*10.*y) .* exp(-x.*y); force = f(X,Y); surf(X,Y,force)