49 views (last 30 days)

Show older comments

Hi,

I have 2D-surface functions z(x,y) that I want to sample efficiently. Does matlab include any such addaptive mesh-refinement functions and, if so, how to use it?

Backgrund: I want to find sample points on my functions, so that I later can intepolate (using e.g. linear 'scatteredInterpolant') between those sampled points and by this effectively recover my "full" surface function. I will probably need around 1e6 sample points on each surface -- in order to precisely enough capture the full shape of my 2D surfaces z(x,y) -- but still these points must be mainly concentrated in those regions where z have peaks or edges.

Let me give a simple example (my functions typically also have narrow edges, and higher peaks):

z = @(x,y) 3*(1-x).^2.*exp(-(x.^2) - (y+1).^2)- 10*(x/5 - x.^3 - y.^5).*exp(-x.^2-y.^2) - 1/3*exp(-(x+1).^2 - y.^2);

Now, I want to adaptaively sample x, y points such that it allows to capture the full behaviour of z.

The only function I could find that does something like this was 'fsurf':

h = fsurf(z, [-5 5 -5 5],'MeshDensity',10,'Visible','off');

I can then access fsurf's sampled points via:

x=h.XData; y=h.YData; z = h.ZData;

However, I want to control the adaptive refinments. fsurf gives no such options (excpet the initial MeshDensity). I accidentally noticed that I can also set a parameter called h.AdaptiveMeshDensity to larger values (e.g. to 9) after fsurf is executed (if I keep fplot's 'empty' window open!). This setting seems to allow to futher refine the sampling in the most relevant regions. Let me illustrate it by this code:

figure()

h.AdaptiveMeshDensity=9;

plot3(h.XData,h.YData,h.ZData,'.b')

hold on

h.AdaptiveMeshDensity=2;

plot3(h.XData,h.YData,h.ZData,'.r','MarkerSize',10)

So I do obtain sampling points adaptively in this way, but:

- I can not find any documentation about the setting`AdaptiveMeshDensity', nor is this atribute listed when I run get(h), but still this setting influences the adaptive sampling in fsurf?
- I don't want fsurf to open a figure, but I don't understand how to prevent it, if using fsurf?
- I imagine that there exists a more fundamental matlab function/operation (which fsurf might even be using) to adaptively find good mesh/sampling points. However, I could not find any such adaptive-mesh-refinement routine(s) within matlab?

Thanks for any useful help and inputs!

(p.s. I noticed there is a function generateMesh, but unclear if this can be of any use here)

KSSV
on 29 Mar 2021

z = @(x,y) 3*(1-x).^2.*exp(-(x.^2) - (y+1).^2)- 10*(x/5 - x.^3 - y.^5).*exp(-x.^2-y.^2) - 1/3*exp(-(x+1).^2 - y.^2);

m = 100 ; n = 100 ; % change these values

xi = linspace(-5,5,m) ;

yi = linspace(-5,5,n) ;

[x,y] = meshgrid(xi,yi) ;

z = z(x,y) ;

surf(x,y,z)

Bjorn Gustavsson
on 5 Apr 2021

Edited: Bjorn Gustavsson
on 5 Apr 2021

It is not entirely clear what your "problem situation" or "design desires" are. These might be very important for what type of solution to suggest. There is a sparse-grid-interpolation tool on the file exchange that is rather powerful and clever. However, you indicate that you have a general 2-D function you want to approximate efficiently on a completely unstructured grid - which sounds like quite a challenge. Once upon a time a coleague and I tried to make a "sparse image representer" by iteratively adding points to a Delaunay-triangulation where the absolute difference were largest:

I = imread('liftingbody.png');

I = double(I);

X = [1 1 512 512];

Y = [1 512 512 1];

F = scatteredInterpolant(X(:),Y(:),(diag(I(Y,X))));

Iapp = F(x,y);

[Imax,i1,i2] = max2D(abs(double(I)-Iapp)); % My own max2D-version

IMax = Imax;

subplot(1,3,1)

imagesc(abs(double(I))),colorbar

while Imax > 10

subplot(1,3,2)

imagesc(abs(double(Iapp))),colorbar

subplot(1,3,3)

imagesc(abs(double(I)-Iapp)),colorbar,hold on,plot(X,Y,'r.'),hold off,

X = [X,i2];

Y = [Y,i1];

F = scatteredInterpolant(X(:),Y(:),(diag(I(Y,X))),'natural');

% 'natural' in this case gives vibes of Turner, 'linear' reminds me of cubist-period Picasso

Iapp = F(x,y);

title(numel(X))

[Imax,i1,i2] = max2D(abs(double(I)-Iapp));

IMax(end+1) = Imax;

drawnow

end

This was mainly for laughs and giggles. One could also turn to a quad-tree-type parameterization (the attached function uses a general polynomial model-fiting-function orignating from John d'Ericco similar to: polyfitn).

Additional tools that might be of interest on the file exchange might be: griddatalsc, gaussian-interpolation-with-successive-corrections and barnes-interpolation-barnes-objective-analysis.

In the end it all comes down to how you plan to calculate the fit between your function and approximant. Is the function very expensive to calculate? Can you live with a regular grid or fixed number of test-points for evaluating the fit quality, are those points to be adapted? For more detailed advice more details about the problem is necessary.

HTH

Bjorn Gustavsson
on 7 Apr 2021

If you'd be happy to get the grid from the quadrature-functions then it sounds like you could just have a look at the source-code of quadgk and integral2 and see what you should extract from those functions. They are normal .m-files at least in 2020a.

HTH

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!