[Date Index] [Thread Index] [Author Index]

Re: ContourFindRoot

• To: mathgroup at smc.vnet.net
• Subject: [mg71327] Re: ContourFindRoot
• From: "dimitris" <dimmechan at yahoo.com>
• Date: Wed, 15 Nov 2006 06:43:23 -0500 (EST)
• References: <ejc4rs$6nr$1@smc.vnet.net>

In the previous post I forrgot to mention any kind of references.

So here they are

http://groups.google.com/group/comp.soft-sys.math.mathematica/browse_thread/thread/f042e78929ddb078/0a9c893ced104f6f?lnk=st&q=FindRoot&rnum=3#0a9c893ced104f6f

http://groups.google.com/group/comp.soft-sys.math.mathematica/browse_thread/thread/2587d839a40223bb/1a2c1c290d3c4b26?lnk=gst&q=FindRoot+David+Park&rnum=1#1a2c1c290d3c4b26

http://groups.google.com/group/comp.soft-sys.math.mathematica/browse_thread/thread/9836d5644ef5b0a1/49569b37ba055a99?lnk=gst&q=FindRoot&rnum=7#49569b37ba055a99

http://groups.google.com/group/comp.soft-sys.math.mathematica/browse_thread/thread/7612c1825cd938da/f5c53fdd2c7b737e?lnk=gst&q=Find+crossings&rnum=3#f5c53fdd2c7b737e

http://www.amazon.com/Mathematica-Action-Stan-Wagon/dp/0387986847/sr=1-2/qid=1162830266/ref=sr_1_2/103-5414680-6091022?ie=UTF8&s=books

http://www.amazon.com/Complex-Analysis-MATHEMATICA-William-Shaw/dp/0521836263/sr=8-1/qid=1163459323/ref=pd_bbs_sr_1/103-5414680-6091022?ie=UTF8&s=books

Dimitris





dimitris wrote:
> Below we construct a routine, which will take two equations {f[x,y]=0,
> g[x,y]=0} and return all the simultaneous solutions (only the
> crossings; not tangencies) in a rectangle.
>
> The basic steps of the routine are the following:
>
> First, we generate a contour plot of f[x,y]=0, go into the graphics
> objects that results, and gather up the data the data making up the
> curves. Then, on each curve, we evaluate g. The places where g changes
> sign correspond to solutions; perhaps not very accurate. We next send
> them to FindRoot in order to get exact solution.
>
> ----------------------------------------------------------------------------------------------------------------------
>
> Clear["Global`*"]
>
> Off[General::spell1]
>
> CurveFindRoot::usage =
> "CurveFindRoot[{f,g},{x,xmin,xmax},{y,ymin,ymax},curve\
> ,points,opts]\nfinds all nontangential solutions to {f=0, g=0}, on a
> contour \
> line.";
>
> ContourFindRoot::usage =
> "ContourFindRoot[{f,g},{x,xmin,xmax},{y,ymin,ymax},p\
> oints,opts]\nfinds all nontangential solutions to {f=0, g=0}. in a
> rectangle \
> in the 2D domain";
>
> ContourPlotFunctionRoots::usage =
> "ContourPlotFunctionRoots[{f,g},{x,xmin,xma\
> x},{y,ymin,ymax},points,roots,optsplot]\ngives the graph of curves and
> \
> zeros.";
>
> CurveFindRoot[{f_, g_}, {x_, x0_, x1_}, {y_, y0_, y1_}, curve_,
> points_,
>    opts___] := Block[{contourdata, pos, seeds, gg},
>    contourdata = Cases[Graphics[ContourPlot[f, {x, x0, x1}, {y, y0,
> y1},
>         Contours -> {0}, ContourShading -> False, PlotPoints -> points,
>
>         DisplayFunction -> Identity]], Line[a_] :> a, Infinity];
>     gg = Function[{x, y}, Evaluate[g]];
>     pos = Position[Partition[Apply[gg, contourdata[[curve]], {1}], 2,
> 1],
>       _?(Times @@ #1 < 0 & ), {1}]; seeds = contourdata[[curve,
>       Flatten[pos]]]; (FindRoot[{f == 0, g == 0}, {x, #1[[1]]}, {y,
> #1[[2]]},
>        opts] & ) /@ seeds]
>
> ContourFindRoot[{f_, g_}, {x_, x0_, x1_}, {y_, y0_, y1_}, points_,
>    opts___] := Block[{lengthcontourdata, roots},
>    lengthcontourdata = Length[Cases[Graphics[ContourPlot[f, {x, x0,
> x1},
>          {y, y0, y1}, Contours -> {0}, ContourShading -> False,
>          PlotPoints -> points, DisplayFunction -> Identity]], Line[a_]
> :> a,
>        Infinity]]; roots = (CurveFindRoot[{f, g}, {x, x0, x1}, {y, y0,
> y1},
>         #1, points, opts] & ) /@ Range[lengthcontourdata];
>     {x, y} /. Flatten[DeleteCases[roots, {}], 1]]
>
> ContourPlotFunctionRoots[{f_, g_}, {x_, x0_, x1_}, {y_, y0_, y1_},
> points_,
>    roots_, optsplot___] :=
>   Show[MapThread[ContourPlot[#1, {x, x0, x1}, {y, y0, y1}, optsplot,
>       PlotPoints -> points, Contours -> {0}, ContourShading -> False,
>       DisplayFunction -> Identity, ContourStyle -> {Hue[#2/3]}] & ,
>     {{f, g}, {1, 2}}], Graphics[
>     ({Red, AbsolutePointSize[6], Point[#1]} & ) /@ roots],
>    DisplayFunction -> $DisplayFunction]
>
> ---------------------------------------------------------------------------------------------------------------------
>
> Here are some applications on finding the roots of a function in a
> rectangle
> in the complex domain.
>
> Clear[f]
>
> f[z_] := Sin[z^3] + Cos[z^2] - z
>
> ref[x_, y_] = Simplify[Re[ComplexExpand[f[z] /. z -> x + I*y]],
>    Element[{x,y},Reals]]
> imf[x_, y_] = Simplify[Im[ComplexExpand[f[z] /. z -> x + I*y]],
>    Element[{x,y},Reals]]
> -x + Cos[x^2 - y^2]*Cosh[2*x*y] + Cosh[3*x^2*y - y^3]*Sin[x^3 -
> 3*x*y^2]
> -y - Sin[x^2 - y^2]*Sinh[2*x*y] + Cos[x^3 - 3*x*y^2]*Sinh[3*x^2*y -
> y^3]
>
> sols = ContourFindRoot[{ref[x, y], imf[x, y]}, {x, -2, 2}, {y, -2, 2},
> 100]
>
> ContourPlotFunctionRoots[{ref[x, y], imf[x, y]}, {x, -2, 2}, {y, -2,
> 2}, 100, sols];
>
> Here is another example adopted from a private communication with David
> Park
>
> Clear[f]
>
> f[z_] := Sin[Cos[z]] - z
>
> ref[x_, y_] = Simplify[Re[ComplexExpand[f[z] /. z -> x + I*y]],
>   Element[{x,y},Reals]]
> imf[x_, y_] = Simplify[Im[ComplexExpand[f[z] /. z -> x + I*y]],
>    Element[{x,y},Reals]]
> -x + Cosh[Sin[x]*Sinh[y]]*Sin[Cos[x]*Cosh[y]]
> -y - Cos[Cos[x]*Cosh[y]]*Sinh[Sin[x]*Sinh[y]]
>
> Timing[sols = ContourFindRoot[{ref[x, y], imf[x, y]}, {x, -Pi, Pi},
>     {y, -Pi, Pi}, 200]]
>
> ContourPlotFunctionRoots[{ref[x, y], imf[x, y]}, {x, -Pi, Pi}, {y, -Pi,
> Pi},
>    200, sols];
>
> Timing[sols = ContourFindRoot[{ref[x, y], imf[x, y]}, {x, -Pi, Pi},
>     {y, Pi, 5*(Pi/3)}, 200]]
>
> ContourPlotFunctionRoots[{ref[x, y], imf[x, y]}, {
>     x, -Pi, Pi}, {y, Pi, 5Pi/3}, 200, sols, ImageSize -> {600, 400}];
>
> Here is a system of two equations f[x,y]==0&&g[x,y]==0
>
> f[x_, y_] := Cosh[x + y] + x^2*y^2 - 7*y*x^3
> g[x_, y_] := Sinh[x + y] + x*y - 2*x
>
> Timing[sols = ContourFindRoot[{f[x, y], g[x, y]}, {x, -10, 10}, {y,
> -10, 10}, 200]]
>
> ContourPlotFunctionRoots[{f[x, y], g[x, y]}, {x, -10, 10}, {y, -10,
> 10}, 200, sols]
>
> I would appreciate a lot any kinds of comments for my code, possible
> modifications,
> suggestions for extension to treat more cases etc.
>
> I would further appreciate if you experiment with it with functions of
> your own
> and write down the results.
> 
> Best Regards
> Dimitris