One way is to create a function that uses Check
to see if the message was printed.
InDomainQ[ip_,x_]:=Quiet[Check[ip@@x;True,False,InterpolatingFunction::dmval]]
For the function in the question:
region = RegionPlot[InDomainQ[ip,{x,y}],
{x,ip[[1,1,1]],ip[[1,1,2]]},
{y,ip[[1,2,1]],ip[[1,2,2]]},
MaxRecursion->4,
PlotPoints->50,
Epilog->{Red,Point[{x,y}]}]

This seems to be the convex hull of the points:
<<ComputationalGeometry`
hullpts = ip["Grid"][[ConvexHull[ip["Grid"]]]]
hull = Graphics[{
{Red, Point[{x, y}]},
{Opacity[0.2], Green, Polygon[hullpts]}
}];
Show[region, hull]

But how often this is the case I don't know.
Edit: For functions $\mathbb{R}^3 \rightarrow \mathbb{R}$ it actually seems the idea of using the convex hull to determine the domain is more reliable than checking for the message
SeedRandom[1];
dat = RandomReal[1, {9, 4}];
paramdat = Map[{Most[#], Last[#]} &, dat];
ip = Interpolation[paramdat, InterpolationOrder -> 1];
ip[10000,0,0]
is clearly outside the domain but does not trigger a message, and checking that the argument is inside the rectangle given by ip["Domain"]
would only be partially helpful.
InDomainQ[ip_InterpolatingFunction, x_List] :=
Quiet[Check[ip @@ x; True, False, InterpolatingFunction::dmval],
InterpolatingFunction::dmval]
(* Plot InDomainQ region together with convex hull *)
d = 0.05;
pts = Select[
Flatten[Table[{x, y, z}, {x, 0, 1, d}, {y, 0, 1, d}, {z, 0, 1, d}],2],
(InDomainQ[ip, #] &)
];
(* ConvexHull3D is in MathWorldPackages.zip
available at http://library.wolfram.com/infocenter/MathSource/4775 *)
hull = ConvexHull3D[ip["Grid"]];
Graphics3D[{
{Green, Opacity[0.2], EdgeForm[Opacity[0.3]],
Cuboid[# - d/2, # + d/2] & /@ pts},
{Blue, Sphere[ip["Grid"], 0.05]},
{Opacity[0.3], hull}
}]

Green is where InDomainQ
is True
. The blue spheres is the domain part of the interpolated points, and the polygons connecting them form their convex hull.
Since ComplexHull3D
returns such nicely oriented Polygons
it's easy to construct a new InDomainQ
function simply by checking if the signed distance from a given point to each polygon has the same sign:
(* pt's side of the plane passing through the p1,p2,p3 *)
SideOfPlane[{p1_,p2_,p3_},pt_]:=Module[
{v1=p2-p1,v2=p3-p1,n},
n=Cross[v1,v2];
Sign[n.(pt-p1)]
]
GetDomainFunction[ip_]:=Module[
{ hull=ConvexHull3D[ip["Grid"]], planes},
planes=hull/.Polygon[a_]:>a;
(* Here is where the orientation from ConvexHull3D comes in
if the sign from SideOfPlane is same for all the polygons
it is in the convex hull *)
Function[{x,y,z},Evaluate[
Equal@@Function[plane,SideOfPlane[plane,{x,y,z}]]/@planes
]]
];
(* InDomainQ[x,y,z], only for this instance of ip *)
InDomainQ=GetDomainFunction[ip];
hull = ConvexHull3D[ip["Grid"]];
Show[{
Graphics3D[{Opacity[0.4],Green,hull}],
RegionPlot3D[InDomainQ[x,y,z],{x,0,1},{y,0,1},{z,0,1},
MaxRecursion->4,
PlotPoints->30,
PlotStyle->Directive[Opacity[0.4],Blue]]
}]

<< DifferentialEquations`InterpolatingFunctionAnatomy`
. It's documented in the advanced NDSolve docs, and if you look at the package source, you'll find all theseip["something"]
things. The package is simply a wrapper over these. $\endgroup$InterpolatingFunction
itself did not seem very well documented. That documentation wasn't linked either. So unless I were usingNDSolve
or looking here, that would have been hard to find. $\endgroup$