Re: Re: DiffMaps package

• To: mathgroup at smc.vnet.net
• Subject: [mg91420] Re: Re: DiffMaps package
• From: "Jean-Marc Gulliet" <jeanmarc.gulliet at gmail.com>
• Date: Thu, 21 Aug 2008 04:14:47 -0400 (EDT)
• References: <20080820143021.4EB9F70@mach.vub.ac.be>

```On Wed, Aug 20, 2008 at 4:30 PM, eabad wrote:

> a couple of months ago, you were kind enough to reply to my email on how to compute extrema of 2D functions (see below). Imagine now that I have a list of data for the values of a two variable function, from which I can get an Interpolating function f[x,y] in mathematica following the standard procedure. When trying to apply your code below to f[x,y], the Solve or Nsolve command fails to provide the solutions of the equation D[f[x,y],{{x,y}},1]==0 for the zeros of the gradient in a specified x-y range. I wonder if you know a getaround for this problem. I am also somewhat worried about the discontinuity of the derivatives of f[x,y], I wonder how to obtain a smoother interpolation than the one provided by the default procedure. I'd be very grateful if you could briefly answer to these questions.

[... Message crossposted to MathGroup ...]

(Do not forget to post your query to the list/newsgroup MathGroup
(comp.soft-sys.math.newsgroup) to get a broader on how to solve your
problem.)

With an interpolating function, you can use *FindRoot* to find zeros.
Pay particular attention to the various ways of entering the initial
value(s) for each variable for they influence the choice of a search
method and allow to restrict the search to a specific interval.

You can change the interpolating order with the option
*InterpolatingOrder->n* and also tell Mathematica whether the fuctiion
is periodic with the option *PeriodicInterpolation -> True*.

The following example illustrate some technics you could use and some
pitfall you may encounter.

f[x_, y_] := Sin[x ] Cos[y]

Plot3D[f[x, y], {x, 0, 2 Pi}, {y, 0, 2 Pi}]

Plot3D[f[x, y] + 1/10 RandomReal[], {x, 0, 2 Pi}, {y, 0, 2 Pi}]

data = Flatten[
Table[{{x, y}, f[x, y] + 1/10 RandomReal[]}, {x, 0, 2 Pi}, {y, 0,
2 Pi}], 1]

g = Interpolation[data, InterpolationOrder -> 4]

Plot3D[g[x, y], {x, 0, 2 Pi}, {y, 0, 2 Pi}]

ListPlot3D[Flatten /@ data]

grad = D[g[x, y], {{x, y}, 1}]

FindRoot[grad, {{x, 1, 0, 2 Pi}, {y, 1, 0, 2 Pi}}]

Table[FindRoot[grad, {{x, x0, 0, 2 Pi}, {y, y0, 0, 2 Pi}}], {x0, 0,
2 Pi}, {y0, 0, 2 Pi}]

Regards,
- Jean-Marc

>>Modeler wrote:
>>
>>> Hi everybody,
>>>
>>> I need to determine the critical points of a 3D surface parameterized by a function f[x,y]. I was told there is a package called diffMaps (apparently not installed by default) which allows one to easily compute local minima and  maxima as well as saddle points. Does anyone know how to get it or else can give any ideas about computation of critical points? Thanks a lot.
>>
>>
>>Say that the function we want to study is the following:
>>
>>     f[x_, y_] := 3*x*y - x^3 - y^3
>>
>>We compute the first derivatives of f w.r.t. to x and y:
>>
>>     D[f[x, y], {{x, y}, 1}]
>>
>>     {-3 x^2 + 3 y, 3 x - 3 y^2}
>>
>>We equate each derivative to zero and solve for x and y:
>>
>>     Solve[% == 0]
>>
>>     {{x -> 0, y -> 0}, {x -> 1, y -> 1}, {x -> -(-1)^(1/3),
>>       y -> (-1)^(2/3)}, {x -> (-1)^(2/3), y -> -(-1)^(1/3)}}
>>
>>We want only the real points, so we discard the complex solutions:
>>
>>     pts =
>>         Cases[%, {x -> vx_, y -> vy_} /; Im[vx] == 0 && Im[vy] == 0]
>>
>>     {{x -> 0, y -> 0}, {x -> 1, y -> 1}}
>>
>>So we have two critical points: (0, 0) and (1, 1).
>>
>>Now, we compute the Hessian matrix,
>>
>>     D[f[x, y], {{x, y}, 2}]
>>
>>     {{-6 x, 3}, {3, -6 y}}
>>
>>and take its determinant:
>>
>>     % // Det
>>
>>     -9 + 36 x y
>>
>>Then we compute the values of the Hessian determinant at the critical
>>points:
>>
>>     % /. pts
>>
>>     {-9, 27}
>>
>>Since the Hessian is negative at the point (0, 0), we conclude that this
>>
>>Since the Hessian is positive at the point (1, 1), we conclude that this
>>is a local maximum.
>>
>>We can "check" (or anticipate) these results by plotting the function:
>>
>>     Plot3D[f[x, y], {x, -2, 2}, {y, -2, 2}, PlotRange -> All]
>>
>>     ContourPlot[f[x, y], {x, -2, 2}, {y, -2, 2},
>>         ColorFunction -> "SunsetColors"]

```

• Prev by Date: Mathematica Special Interest Group (Washington DC Area) UPDATED
• Next by Date: Re: Re: NDSolve[] with nested If[] and Piecewise[] usage:
• Previous by thread: Mathematica Special Interest Group (Washington DC Area) UPDATED
• Next by thread: Hypergeometric1F1 polynomial