Re: Fitting complicate functions

• To: mathgroup at smc.vnet.net
• Subject: [mg111070] Re: Fitting complicate functions
• From: janos <janostothmeister at gmail.com>
• Date: Mon, 19 Jul 2010 02:09:48 -0400 (EDT)
• References: <i1u21r\$7sd\$1@smc.vnet.net>

```On j=FAl. 18, 07:04, Marco Masi <marco.m... at ymail.com> wrote:
> I have a very complicate function which has to be inserted into a loop that calculates it zillions of times. This slows the simulation speed and I'm wondering if it is possible to fit the function with a simpler one by adjusting its fitting parameters. I used "Fit" of Mathematica and was partially succesful: however, I tried several functions (polynomial, logarithmic, etc.) and parameters but the difference between the original and the fitting function always tends to diverge in some regions. I'm wondering if this can be done better? The program is appended, let me know if you have a better idea and/or parameterizing fitting function or method.
>
> Regards, Mark.
>
> ----------------------------------------------------------------------
>
> (*Fundamental parameter setting*)
> pc = 3600*24*365*299792458*3.2638; M0 = 1.989 10^30; Myr =
>  3600*24*365*10^6;
> \[Alpha] = (M0*Myr^2)/pc^3; G = \[Alpha]*6.672 10^-11; Rmax = 10000=
0;
>
> Subscript[M, BH] = 10^6;
> Subscript[M, BG] = 1.6*10^10; rc = 420;
> \[CapitalSigma] = 492; rd = 3500;(*Central density and scale radius \
> of the Freeman disk*)\[Rho] = 0.05; rh = 4600;(*Central density and \
> scale radius of the Psiso halo*)\[Rho]mpsiso = 0.014; rmpsiso = \
> 13000;(*Central density and scale radius of MPSISO*)\[Rho]nfw = \
> 0.012; rnfw = 15000;(*Central density and scale radius of NFW*)(*The \
> Miyamoto-Nagai potential parameters.*)M1 = 2.4 10^9; M2 =
>  18.7 10^9; M3 = 40.2 10^9; M4 = 3.4 10^9;
> a1 = 80; a2 = 500; a3 = 4600; a4 = 9000; b1 = 95; b2 = 200; b=
3 =
> = 400; \
> b4 = 200;
>
> Subscript[\[Phi], BH][X_,
>   Y_] := -((G Subscript[M, BH])/
>     Sqrt[X^2 + Y^2]);(*Potential of the Black Hole*)
> Subscript[\[Phi], BG][X_,
>   Y_] := -((G Subscript[M, BG])/Sqrt[X^2 + Y^2 + rc^2]);
> Subscript[\[Phi], Free][X_,
>    Y_] := -\[Pi] G \[CapitalSigma] Sqrt[
>     X^2 + Y^2] (BesselI[0, Sqrt[X^2 + Y^2]/(2 rd)] BesselK[1,
>        Sqrt[X^2 + Y^2]/(2 rd)] -
>      BesselI[1, Sqrt[X^2 + Y^2]/(2 rd)] BesselK[0,
>        Sqrt[X^2 + Y^2]/(2 rd)]);
> Subscript[\[Phi], PSISOMOD][X_,
>    Y_] := -4 \[Pi] G \[Rho]mpsiso rmpsiso^2 (rmpsiso/
>        Sqrt[X^2 + Y^2] ArcSinh[Sqrt[X^2 + Y^2]/rmpsiso] -
>      1/Sqrt[1 + Rmax^2/rmpsiso^2]);
> (*Subscript[\[Phi],tot][X_,Y_]=Subscript[\[Phi],BH][X,Y]+Subscript[\
> \[Phi],BG][X,Y]+Subscript[\[Phi],Free][X,Y]+Subscript[\[Phi],PSISOMOD]\
> [X,Y];*)
> Subscript[\[Phi], tot][X_, Y_] = Subscript[\[Phi], Free][X, Y];
>
> Rrange = 15000; PotRange = -50000; Rstep = 301;
>
> Plot3D[Subscript[\[Phi], tot][X, Y], {X, -Rrange,
>   Rrange}, {Y, -Rrange, Rrange}, PlotRange -> {PotRange, 0},
>  PlotStyle -> RGBColor[1, 0, 0]]
> Print["The above figure: plot of the Galactic potential. Wait..."];
>
> GalPotPoints =
>   Flatten[Table[{X, Y, Subscript[\[Phi], tot][X, Y]}, {X, -Rrange,
>      Rrange, Rstep}, {Y, -Rrange, Rrange, Rstep}], 1];
> MinPot = Min[
>    Table[GalPotPoints[[i, 3]], {i, 1, Length[GalPotPoints]}]];
> ListPointPlot3D[GalPotPoints, PlotRange -> {MinPot, 0},
>  PlotStyle -> RGBColor[0, 0, 1]]
> Print["The above figure: plot of the Galactic potential points mesh. \
> Wait..."];
>
> (* Here "quad" is the fitting function *)
> quad = Fit[GalPotPoints, {1, Log[Sqrt[X^2 + Y^2]]}, {X, Y}];
>
> Show[Plot3D[quad, {X, -Rrange, Rrange}, {Y, -Rrange, Rrange},
>   PlotStyle -> Directive[PointSize[Medium], Blue],
>   PlotRange -> {MinPot, 0}],
>  Graphics3D[{Red, PointSize[0.005], Map[Point, GalPotPoints]}]]
> Print["The above figure: plot of the Galactic potential fit function",
>    quad, ".  Wait..."];
>
> DiffFitGalPotPoints =
>   Table[Subscript[\[Phi], tot][GalPotPoints[[nr, 1]],
>      GalPotPoints[[nr, 2]]] - (quad /. X -> GalPotPoints[[nr, 1]] /=
=
> .
>       Y -> GalPotPoints[[nr, 2]]), {nr, 1, Length[GalPotPoints]}]=
=
> ;
> MaxDiffPot =
>   Max[Table[DiffFitGalPotPoints[[nr]], {nr, 1, Length[GalPotPoints]}]];
> MinDiffPot =
>   Min[Table[DiffFitGalPotPoints[[nr]], {nr, 1, Length[GalPotPoints]}]];
> DiffFitGalPotPointsPlot =
>   Table[{GalPotPoints[[nr, 1]], GalPotPoints[[nr, 2]],
>     DiffFitGalPotPoints[[nr]]}, {nr, 1, Length[GalPotPoints]}];
> ListPointPlot3D[DiffFitGalPotPointsPlot,
>  PlotRange -> {MinDiffPot, MaxDiffPot}, PlotStyle -> RGBColor[1, 0, 1]]
> Print["The above figure: plot of the difference between the Galactic \
> potential and the Galactic potential fit function."];
> Print["Minimum Galactic potential of ", MinPot,
>   "    Maximun positive difference= ", MaxDiffPot,
>   "     Minimum negative difference= ", MinDiffPot];
> PotErrRange = 100*(MaxDiffPot + Abs[MinDiffPot])/Abs[MinPot];
> Print["Maximum error of fit over Galactic potential= ", PotErrRange,
>   " %"];

Something is wrong here:
DiffFitGalPotPoints =
Table[Subscript[\[Phi], tot][GalPotPoints[[nr, 1]],
GalPotPoints[[nr, 2]]] - (quad /. X -> GalPotPoints[[nr, 1]] /=

```

• Prev by Date: Re: Show left hand side
• Next by Date: Re: Transform differential equation by tranformation rule
• Previous by thread: Fitting complicate functions
• Next by thread: Re: Fitting complicate functions