Take the 2-minute tour ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.

In this paper about Interpolation Methods, I am trying to learn Minimum curvature method. I have not done partial differential equations before; hence I am finding it tough to penetrate through this article.

Here is an extract:

"This method and namely its computer implementation was developed by W.H.F. Smith and P.Wessel in1990. The interpolated surface by the Minimum Curvature method is analogous to a thin, linearly elastic plate passing through each of the data values with a minimum amount of bending. The algorithm of the Minimum curvature method is based on the numerical solution of the modified biharmonic differential equation $$(1-T)\nabla^{4}f(x,y)-(T)\nabla^{2}f(x,y)=0$$ with three boundary conditions: $$ (1-T)\frac{\partial^2f}{\partial n^2}+(T)\frac{\partial f}{\partial n}=0$$ $$\frac{\partial \nabla^2 f}{\partial n}=0 \text {on the edges}$$ $$\frac{\partial^2 f}{\partial x \partial y}=0 \text{at the corners}$$ where T $\in \langle 0,1 \rangle$ is a tensioning parameter, $\nabla^2$ is the Laplacian operator: $\nabla^2 f$=$\frac{\partial^2 f}{\partial x^2} + \frac{\partial^2 f}{\partial y^2},$ $\nabla^4=(\nabla^2)^2$ is the biharmonic operator, and $\nabla^4 f=\frac{\partial^4 f}{\partial x^4}+\frac{\partial^4 f}{\partial y^4}+2\frac{\partial^4 f}{\partial x^2 \partial y^2}$ and n is the boundary normal. If $T=0$, the biharmonic differential equation is solved; if $T=1$, the Laplace differential equation is solved - in this case the resulting surface may have local extremes only at points $XYZ.$"

Can anyone explain to me in simpler terms about this method. I would be highly grateful to you, as always. Thanks.

share|improve this question

1 Answer 1

up vote 4 down vote accepted

To start, let's say that T=1, so we are solving the Laplace equation. In my work, I had to solve this equation over a "hole" in a surface, which surface was specified by height values over a grid of pixels. In this case, the boundary conditions are the heights of the surface at the pixels bordering the hole. (Without explanation, I'll simply say that we can ignore the other boundary conditions given above regarding curvature at the edges, because the Laplace equation doesn't care about this.) This means that I was trying to find a surface, defined over the hole (including the border pixels of data), where the Laplacian is 0 everywhere (i.e. at each pixel). If you don't know what the Laplacian operator is and what it does, read the short, excellent book by H.M. Schey titled "Div, Grad, Curl and All That." Or, suffice it to say for now that the Laplacian of a surface at a certain point is the curvature as you travel in the x direction plus the curvature as you travel in the y-direction. Saying that the Laplacian is zero means that, at every point, the x-curvature should always be negative the y-curvature. That is, at every point, the surface locally looks like a "saddle" or a plane, and never like a peak or a pit.

This can be solved multiple ways, but it turns out that there is a simple way given "Numerical Recipes" called the method of "relaxation", which can be improved to a similar method called "successive over-relaxation". "Relaxation" is a commonly used term among people who solve diff eqs numerically, and so is "SOR", so you can google these terms. The method makes use of the property that every point on a Laplacian surface is the average of its neighbors. I'm speaking roughly now, because this property has to be defined a bit more formal for continuous surfaces, but it is in fact perfectly adequate for pixellated surfaces. "Relaxation" in this case means you iterate over all the points, setting each equal to the average of its neighbors, again and again, and the surface converges to the Laplacian surface. SOR means you overshoot a bit, which can speed things up considerably, but in order to calculate how much you should overshoot for maximum speed (for an arbitrary given hole shape), you need to solve a giant sparse matrix, and I haven't gotten around to figuring out how to do that quickly. So I just use a hard-coded "relaxation" value, and it's pretty fast.

(Believe it or not, there's a patch for the above method if the hole touches the edge of your window and so you don't have boundary data everywhere -- you have to use "von Neumann stencils", and I think I read about them in a textbook on math methods with Matlab.)

Incidentally, a Laplacian surface for a given boundary is the same as (or, very very close to?) the shape a film of soap takes if you dip the boundary in soapy water. I'm going to wave my hands here and only half-lie to you by saying that nature as usual is finding its "minimum energy" state as the soap "relaxes" to the shape it finaly takes, and that that is why you may have heard Laplacian surfaces referred to as "minimum energy" surfaces. Think about the fact that you have to add energy (say, by pushing with your finger) in order to make the soap film have a peak or a pit anywhere -- you have to stretch the soap, and it naturally doesn't want to be stretched.

Anyway, now let's say that T=0, so that we are solving purely the biharmonic equation. Just like the Laplacian operator talks about how much a surface is like a "peak" or a "pit" (since, when the Laplacian is 0, we are like a "saddle" instead -- there is no surplus stretching-energy stored anywhere) -- I am going to go out on a limb and hazard to say that biharmonic operator behaves the same way, but regarding "bending" energy instead of "stretching" energy.

Imagine if you had a surface shaped like a mountain, and you chopped off the top. If you filled the hole with a Laplacian surface, it would look flattish -- it wouldn't continue the upward slope and reconstruct a mountain, because soap only cares about being stretched, not bent. If you wanted to fill with a biharmonic surface, however, it would be like nailing a thin plastic plate down all around the edge -- with two rows of nails, so that the plate is flush with the upward slope of the mountain all around. The plate doesn't want to bend, but you are forcing it to have a certain position and slope all around the edge. So in the middle where it isn't nailed down, it is going to assume the shape that requires it to bend the least overall, which means not having a giant crease all around the edge of the mountain. The result is that you're going to reconstruct basically a mountain shape.

If T is between 0 and 1, then you have a surface that doesn't want to stretch OR bend.

I won't go much further, because I'm not extremely strong in this domain, but I hope this is a helpful start.

share|improve this answer
    
appreciate your answer. I had completely forgot this question as no one replied. But you did. Thank you. –  math Oct 25 '13 at 23:16

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.