"Scattered Interpolant" produces invalid values

Hi there, i am trying to interpolate the attached file and in some cases it the "scatteredinterpolant" produces invalid values.

Here is my code:

Tb1   = importdata('EngineModelMap.csv',';');
Tb    = Tb1.data;
Tb    = sortrows(Tb,[1 2]);
TET   = Tb(:,3); 
TH    = Tb(:,4); 
SFC   = Tb(:,5); 
ALT   = Tb(:,1); 
MACH  = Tb(:,2);
f_TET = scatteredInterpolant(ALT,MACH,TH,TET,'natural','none');
f_sfc = scatteredInterpolant(ALT,MACH,TET,SFC,'natural','none');
f_Thrust = scatteredInterpolant(ALT,MACH,TET,TH,'natural','none');
Alt  = 7000; 
M    = 0.6; 
Tr   = 13.21;
tet  = f_TET(Alt,M,Tr);
sfc  = f_sfc(Alt,M,tet);
ff   = Tr*sfc;
% SEP,fs (max Thrust)
tetM = 1500; %max(TET);
sfcM = f_sfc(Alt,M,tetM);
Tmax = f_Thrust(Alt,M,tetM);

The issue i encounter is that some combinations that should be outside the engine map, are treated as they are inside it. For example if you set as initial values Alt=7000,M=0.6,Tr=13.21 and take a look at the imported file you will notice that Tr is out of range and thus when the tet=f_TET(Alt,M,T) is calculated should return tet=NaN. In this case though, the interpolation returns tet=1500, which is the maximum value for tet. A way to validate this miscalculation is that later on when Tmax is calculated and tet is set to 1500 (maximum tet for max thrust) the returning Tmax is ~10.

It seems like it is extrapolating, although i have disabled extrapolation.

How can i fix this issue?

[darova]

Please show some images

plot3(x,y,z,'.r')       % original data
surface(X,Y,Z)          % interpolated data

[Laz-os]

I uploaded the file

[darova]

What is this?

[Laz-os]

Just a typo...also corrected!

[darova]

Is it possible to see the problem visually? Something like slice or isosurface?

[Laz-os]

Sorry but i havent tried any visual representation! The best derscription is the one i wrote in the main post :-(. What should i do to help you? The case i give as an example is just a sample that i observed the anomaly...

Answer 1

[John D'Errico]

I can shed light on your problem. What you think of as extrapolation is different from the way scatteredInterpolant sees the concept. The problem is, your data lives in a domain that is not a convex set.

And, while you THINK that the point you show is not inside the data, it IS inside the convex hull of the data. (Not difficult to prove. I've done so.)

So what does extrpolation mean to a tool like scatteredInterpolant? That means predicting a value that lies OUTSIDE of the convex hull of the data. These tools work via triangulations of the domain - Delaunay triangulations, which result in convex things. The outer boundary surface of a Delaunay triangulation is in fact the convex hull of the data. So even though your data happens to look non-convex, scatteredInterpolant does not care in the least. It creates a triangulation that contains everything that lies inside the convex hull of your data. Anything inside that hull is interior to the data, as far as it cares.

Extrapolation is only needed when a point falls outside of that convex hull. That is where your assumptions fall apart.

[Laz-os]

Though the outcome of this colculation is not physicaly plausible....it is unnatural to have greater thrust than the one that that you get with maximun engine temperature! I have also tried to normalize my data before entering them in scatteredInterpolant, but i get the same results!

I get that it not actually extrapolating, i just described it that way because that is how i perceive it, but through the validation i mentioned above https://www.mathworks.com/matlabcentral/answers/536515-scattered-interpolant-produces-invalid-values#comment_871619 the value f_TET(7000,0.6,13.21)=1500 is invalid.

Hope i helped you understand my point...

[John D'Errico]

I understand your point. However, your point is not relevant. Software is written as it is written, to do what it is designed to do. That you want it to perform in a different way is not going to get you far. If you use the wrong software, or if you use that software in the wrong way, what do you expect?

scatteredInterpolant is designed to work WITHIN the convex hull of the data. It does not know the physics of your problem, nor could it possibly care. It is just software. You gave it some numbers. Effectively the code follows a simple flow chart:

Is the point within the convex hull?
  yes?: Do the interpolation as designed.
  No?: Extrapolation! BAD! DANGER!

And you need to understand that what scatteredInterpolant defines as an extrapolation is not what you think extrapolation means. They are different things, but when you use scatteredInterpolant, you play by its rules.

Now to your comment. "it is unnatural to have greater thrust than the one that that you get with maximun engine temperature".

You have used software to interpolate data that is not convex. It produces an approximation based on interpolation as it is designed to do, where the point happens to lie in a region that is outside of the valid parameter space. But since the point is inside the convex hull of the data, garbage arises because the interpolant was still happy. And that should not be a surprise.

Remember that interpolation is just an approximation! When used to predict a value at some location that is between points in your data, there is no presumption it will produce the theoretically correct value, because the interpolation scheme does not understand the physics of your problem.

You are thinking physics with all of its glorious and sutble nonlinearities, but forgetting that interpolation is not physics.

[Laz-os]

I didnt mean to change the software's principles, nor did i expected to make it work in a way that its not designed to. I preseneted my case to find out if there is any workaround. I get that this is they way scatteredInterpollant works but is there any other way i can make my estimations? Maybe by reforming my data or by using any other workflow or function ?!

Answer 2

[darova]

Here are some of your data

scatter3(ALT,MACH,TH,20,TET,'fill')

I used slice to see how the interpoalted data looks

[x,y,z] = meshgrid(0:1e3:15e3,0:0.1:1,0:1:20);
v = f_TET(x,y,z);
hold on
for iso = linspace(min(z(:)),max(z(:)),5)
%     isosurface(x,y,z,v,iso)
    slice(x,y,z,v,[],[],iso)
end
hold off

And it looks pretty logical. Isn't it?

I tried isosurface

[x,y,z] = meshgrid(0:1e3:15e3,0:0.1:1,0:1:20);
v = f_TET(x,y,z);
hold on
for iso = linspace(min(v(:)),max(v(:)),5)
    isosurface(x,y,z,v,iso)
%     slice(x,y,z,v,[],[],iso)
end
hold off

And it looks ok for me as well

Can explain more about the problem?

[darova]

Almost identical

[x,y,z] = meshgrid(0:1e3:15e3,0:0.1:1,0:1:20);
v = f_TET(x,y,z);
hold on
for iso = linspace(min(v(:))*1.02,max(v(:))*0.98,4)
    isosurface(x,y,z,v,iso)
%     slice(x,y,z,v,[],[],iso)
end
hold off

[Laz-os]

Let me explain...These values are data from a jet engine and they represent Altitude, Mach number, Turbine Engine Temperature, Thrust, Specific Fuel Consumption. f_TET aims to find the TET for which the combination of Altitude, Mach number, Thrust would be possible, but for the example we are talking about the combination should not be feasible because at that specific Alt and Mach the maximum thrust would be produced for max TET(=1500). So if we check f_Thrust(7000,0.6,1500), which is the combination for maximum Thrust output gives Tmax=~10. This contradicts the calculation of our example : f_TET(7000,0.6,13.21)=1500. For the engine to achive Thrust=13.21 the TET should be bigger which is not possible for our engine.

Dont know if helped you ?!

[darova]

I'm afraid not. There is too much explanations and theory. Is it possible to make a sketch or something?

[Laz-os]

Check this out...

f_TET(7000,0.6,13.21)=1500

while

f_Thrust(7000,0.6,1500)=~10. These two relationships are in contradicton. Either the first should result to tet>1500 or the second Tr=13.21. You see that for the same engine flight parameters (Alt,M) they have different output.

[darova]

I understand. Maybe it's because of interpolation method. Did you try something different? griddata? It has 'cubic' method, maybe it will be more precise

[darova]

'natural' method is the worst

x = 1:10;
y = sin(x);
x1 = 1:0.1:10;
y1 = interp1(x,y,x1,'natural');
y2 = interp1(x,y,x1,'cubic');

plot(x,y,'o-b')
hold on
plot(x1,y1,'.-r')
plot(x1,y2,'.-g')
hold off

legend('original data','natural interpolation','cubic interpolation')

[Laz-os]

I tried griddata and it produces the same output :-(.

[darova]

Hm. I don't have any other ideas. Maybe John D'Errico can help you, see answer below

[Laz-os]

Anyway thanks for your effort! :-D