# Lagrange Interpolation over quadrilateral points

I've been working on a program which calculates, given a point and 4 surrounding points, the Lagrange polynomial, in order to interpolate a value. Consider that I'm not a mathematician and I better understand code than formulas.

I've came up with the following code, which works, but I really don't think is general (and it's pretty ugly to me).

using System;
using System.Collections.Generic;
using System.Linq;

namespace LagrangeInterpolation
{
public class Point : ICloneable
{
public double X { get; set; }

public double Y { get; set; }

public double Value { get; set; }

public Point Clone()
{
return (Point)this.MemberwiseClone();
}
}

public static class Lagrange
{
public static Point Interpolate(Point[] controlPoints, Point point)
{
var A = -controlPoints.X - controlPoints.X + controlPoints.X + controlPoints.X;
var B = -controlPoints.X + controlPoints.X + controlPoints.X - controlPoints.X;
var C = +controlPoints.X - controlPoints.X + controlPoints.X - controlPoints.X;
var X = 4 * point.X - controlPoints.Sum(x => x.X);

var D = -controlPoints.Y - controlPoints.Y + controlPoints.Y + controlPoints.Y;
var E = -controlPoints.Y + controlPoints.Y + controlPoints.Y - controlPoints.Y;
var F = +controlPoints.Y - controlPoints.Y + controlPoints.Y - controlPoints.Y;
var Y = 4 * point.Y - controlPoints.Sum(x => x.Y);

var r = (X / 4 - B * Y / 4 * E) / (1 - D / 4 * E);
var s = (Y - D * r) / E;

var prevR = 0d;
var prevS = 0d;
const double precision = 0.00000001;
while (!(prevR - r < precision && prevS - s < precision))
{
prevR = r;
prevS = s;
r = (X - B * s) / (A + C * s);
s = (Y - D * r) / (E + F * r);
}

// Interpolate value
var result = point.Clone();
result.Value = ((1 - r) * (1 - s) * controlPoints.Value + (1 - r) * (1 + s) * controlPoints.Value + (1 + r) * (1 + s) * controlPoints.Value + (1 + r) * (1 - s) * controlPoints.Value) / 4;

return result;
}
}
}


The input parameters are:

• controlPoints: the 4 points, each one with its Value.
• point: the point for which we want to calculate the interpolated value

The returned Point is a clone of the Input point with the Value property set. Every instance of Point have X and Y normalised within range -1..1 (I subtract the quadrilateral center from each point).

Example:

controlPoints = new [] {
new Point() { X = -0.033675000000000566, Y = -0.02564999999999884, Value = 1.2787 },
new Point() { X = -0.035524999999999807, Y = 0.024329999999999075, Value = 1.329 },
new Point() { X = 0.03370499999999943, Y = 0.02564999999999884, Value = 1.3376 },
new Point() { X = 0.035494999999999166, Y = -0.024329999999999075, Value = 1.302 }
}

point = new Point() { X = 0.018148174616284152, Y = -0.014201699949808244 }


Expected result is Point.Value = 1.3044829106888913

Can someone suggest a better way (formally and mathematically) to perform this calculation?

• A link or explanation of Lagrange interpolation would be nice. – paparazzo Jun 13 '18 at 16:08
• Done! I'm also looking for some plots but I can't find any... – Dan Jun 13 '18 at 16:41
• I have rolled back your last edit. Please don't change or add to the code in your question after you have received answers. See What should I do when someone answers my question? Thank you. – Phrancis Jun 13 '18 at 17:41
• What do you mean by "Point have x and y normalized within range -1..1"? Could you update the question with a data set and the desired result (It is allowed to update the question with additional information, but not to change or add code blocks to be reviewed)? – Henrik Hansen Jun 14 '18 at 8:59
• Done! Let me know if there is something else I could improve. – Dan Jun 18 '18 at 8:11

I would use a struct with copy constructor instead of a class with clone. Performance wise it would be better and would be more readable.

Extract function for your condition in the while.

• Even better would be the new readonly ref struct :-o or the usage of the new in parameter Reference semantics with value types – t3chb0t Jun 13 '18 at 17:35