Sigh. 2d Lagrange? You do realize that some time after developing Lagrange interpolation, Lagrange died? There must be some causal influence here. He is dead now, and alive before he developed it. Case closed. ;-) ;-)
Ok, that was said very much with tongue in cheek. But somewhat seriously, I want to impress upon you that Lagrange interpolation is actually not a good tool to use in the future for you, even though you have been tasked with this as a project. (I'll look at your code after this digression is done.) Why not?
We can see the problem even in 1-dimension. Consider a simple function of 1 variable, to be interpolated over a set of points using say a moving cubic Lagrange. (You REALLY don't want to consider the use of high order Lagrange, as high order polynomials oscillate terribly. Even cubic polynomials can be problematic. Worse, you cannot perform the computations in any reasonable amount of precision as the polynomial order grows at all large.)
x = linspace(0,pi,10);
y = sin(x);
plot(x,y,'o')
Now, imagine that between a pair of points, you will use the 2 points on either side of the point in question to interpolate. So as you follow along the curve between that pair of points, you get a nice, smoooth interpolant.
But once you cross a node of the function, suddnly, you are uing a DIFFERENT collection of 4 points from before. The result is the function you will see interpolated will be merely continuous, but NOT differentiable across that node of the interpolant.
The net result is a nice, smooth interpolant BETWEEN nodes, but first derivative singularities at every node. I don't want to get into too much depth here, beyond explaining why Lagrange should never be used in practice.
So you really don't have any reason to be using Lagrange interpolation beyond a first order interpolant, since splines exist, and can produce twice continuously differentiable results across the nodes. Get it? Lagrange is a bad idea, at least beyond the linear case. And high order Lagrange is REALLY a bad idea.
One caveat applies in this, in that FIRST order Lagrange interpolants are indeed a viable tool, since if you are willing to do connect the dots interpolation, then you already expect an interpolant that is piecewise linear, with derivative singularities at each break.
plot(x,y,'-o')
So when I did the connect-the-dots plot here, I effectively employed a first order moving Lagrange interpolant. As you can see, it has slope changes at each node, and this is expected.
Ok. Rant mode over. Yes, use Lagrange for your project. As long as you understand this is a tool to be used for your project, and then never again, then go for it.
So next, let me loo at what you wrote.... (give me a couple of minutes to read your code.)
Ugh, retching sound. You are building a Lagrange polynomial of order 201? And that in only ONE dimension? You need to understand that these computations cannot be performed in double precision? Do you accept just how nastily these polynomials will act, since you are effectively raising numbers to powers as large as 201?
A Lagrange interpolant makes sense only if it is low order. Thus over a span, you use nodes [1 2 3 4]. Then you swap to nodes [2 3 4 5],, then nodes [3 4 5 6], etc., moving along the curve. Seriously, if you really expect to perform tensor product high order Lagrange interpolation, where the order of the polynommial is 201 in each dimension, then your project is dead in the water immediately. You cannot do it. Tensor product LINEAR interpolation is fine. That is done all the time. Even a tensor product cubic is doable, though a bit complex. Sorry, but this is time to cut your losses. Talk to your instructor, and make sure you understand what is the goal of your task.
