Plotting the energy function found by summation of the series

Mathematica program

(* Plotting the energy function found by summation of the series *)

sect = "function";
entity = "figure";

(*Options*)

mindlist = 10; (*number of approximants*)
ncrop = 20; (*disregard coefficients with N > ncrop*)
pltpnts = 999;
If[iftest, pltpnts = 111];

f0 = func[[1]];
(*
  e0true = 0/10 f0;
  tr[e_, z_] := z e + (1 - z)e0true;
  *)
tr[e_, z_] := e;
cf = If[nm > ncrop, Take[func, ncrop], func];
nmax = Length[cf];
(*Create table of indexes of approximants*)
indlist = {};
Do[Do[ind = {n1, n2, n3};
      nuse = Plus @@ ind + 2;
      prior = 1.1^nuse;
      If[n1 < 0 || n2 < 0 || n3 < 1 || nuse > nmax, Continue[]];
      If[Abs[n1 - n2] > 1 || Abs[n2 - n3] > 1 || Abs[n3 - n1] > 1, 
        Continue[]];
      If[! (n1 <= n2 <= n3 <= n1 + 1), Continue[]];
      indlist = Append[indlist, {ind, prior, nuse}],
      {n2, n1, n1 + 1}, {n3, n1, n1 + 1}], {n1, 0, nm}];
indlist = Sort[indlist, #1[[2]] > #2[[2]] &];
mi = Length[indlist];
If[mindlist < mi, indlist = Take[indlist, mindlist]; mind = mindlist, 
    mind = mi];
indlist = Sort[indlist, #1[[3]] < #2[[3]] &];
tind = Transpose[indlist];
pind = {ToString /@ tind[[1]], tind[[3]]};
Print["Indexes: ", pind // TableForm];

(* Choosing range of z and E *)
If[NumberQ[zmin] && NumberQ[zmax], ifzbounds = True; {zmin0, zmax0} = {zmin, zmax}, ifzbounds = False];
If[NumberQ[fmin] && NumberQ[fmax], iffbounds = True; {fmin0, fmax0} = {fmin, fmax}, iffbounds = False];
(* Approx. radius of convergence *)
zr = 1/asave[name];
(* Approx. behavior of the function *)
nf = 20;
{zmin, zmax} = 0.85{-zr, -zr/3};
ztab = Table[zmin + (n - 1)/(nf - 1)(zmax - zmin), {n, nf}];
ftab = sumser[#, nm] & /@ ztab;
zftab = Transpose[{ztab, ftab}];
ffitm = Fit[zftab, {1, z}, z];
{zmin, zmax} = 0.85{zr/3, zr};
ztab = Table[zmin + (n - 1)/(nf - 1)(zmax - zmin), {n, nf}];
ftab = sumser[#, nm] & /@ ztab;
zftab = Transpose[{ztab, ftab}];
ffitp = Fit[zftab, {1, z}, z];
(* z - range *)
slistr = Select[slist0, Head[#] === Real &];
slistrm = Select[slistr, # < 0 &];
slistrp = Select[slistr, # > 0 &];
slistc = Select[slist0, Head[#] === Complex &];
If[ifzbounds, {zmin, zmax} = {zmin0, zmax0}, 
zmin = Min[2.4 slistrm, 1.8 Re[slistc]];
If[Abs[zmin]<3 zr, zmin = 1.4 zmin];
If[zmin > 0, zmin = -5 zr];
zmax = Max[1.5 slistrp, 1.5 Re[slistc]];
If[Abs[zmax]<2 zr, zmax = 1.4 zmax];
If[zmax < 0, zmax = 5 zr];
];
(* f - range *)
If[iffbounds, {fmin, fmax} = {fmin0, fmax0},
{zmin0, zmax0} = {zmin, zmax}; 
{zmin, zmax} = {zmin0, 0}; 
ztab = Table[zmin + (n - 1)/(nf - 1)(zmax - zmin), {n, nf}];
ftab = (ffitm /. z -> #) & /@ ztab;
fmax = Max[ftab];
{zmin, zmax} = {0, zmax0}; 
ztab = Table[zmin + (n - 1)/(nf - 1)(zmax - zmin), {n, nf}];
ftab = (ffitp /. z -> #) & /@ ztab;
fmin = Min[ftab];
{zmin, zmax} = {zmin0, zmax0}; 
(* Streching f - range *)
{fmin, fmax} = f0 + {1.4, 1.0} ({fmin, fmax} - f0);
(* Including en + engap into f - range *)
enex2 = enex + engap;
fd = 0.07(fmax - fmin);
If[NumberQ[engap] && enex2 + fd > fmax, fmax = enex2 + fd];
];
If[iprint, Print[{{"Min","Max"},{zmin, zmax}, {fmin, fmax}} // TableForm]];

If[newfigure == True,
scale = 1;
dat = If[NumberQ[engap], {{{1, enex}}, {{1, enex2}}}, {{{1, enex}}}];
pltex = MultipleListPlot[Sequence @@ dat,
      SymbolShape -> {shape[0], shape[1]},
      SymbolStyle -> {{Thickness[0.004 scale], RGBColor[1, 0, 0]}},
      DisplayFunction -> Identity];
grat = 1/GoldenRatio // N;
   pltseq = {};
   Clear[pol, c, funs];
   (*Cycle over indexes*)
Do[
    {{n1, n2, n3}, prior, nc} = indlist[[nindlist]];
    plts = Plot[tr[sumser[z, nc], z], {z, zmin, zmax},
        PlotRange -> {{zmin, zmax}, {fmin, fmax}},
        PlotStyle -> {Thickness[0.005], RGBColor[1, 0, 0]},
        AxesOrigin -> {0, f0},
        PlotPoints -> pltpnts,
        DisplayFunction -> Identity];
    padefrac = padesum[z, nc];
    padefun[t_] := (
        t1 = SetPrecision[t, 64];
        ft = padefrac /. z -> t1;
        If[ft > fmin && ft < fmax, ft, "near Pole"]);
    Off[Plot::"plnr"];
    pltp = Plot[tr[padefun[z], z], {z, zmin, zmax},
        PlotRange -> {{zmin, zmax}, {fmin, fmax}},
        PlotStyle -> {Thickness[0.005], RGBColor[0, 0, 1]},
        AxesOrigin -> {0, f0},
        PlotPoints -> pltpnts,
        Compiled -> False,
        DisplayFunction -> Identity];
    c[0] = 1;
    pol[0] = Sum[c[n]z^n, {n, 0, n1}];
    pol[1] = Sum[c[n + n1 + 1]z^n, {n, 0, n2}];
    pol[2] = Sum[c[n + n1 + n2 + 2]z^n, {n, 0, n3}];
    funs[z_] = Sum[cf[[n]]z^(n - 1), {n, nc}];
    (*Quadratic approximant*)
    pattern = pol[0] + pol[1]f[z] + pol[2]f[z]^2;
    eqs = Take[CoefficientList[pattern /. f -> funs, z], nc];
    s = Solve[eqs == Table[0, {nc}], Table[c[n], {n, nc}]];
    ls = Length[s];
    If[ls == 0, Print["No solutions of linear equations found!"]];
    If[ls > 1, Print["More than one solution of linear equations!"]];
    s1 = s[[1]];
    polfz = pattern /. s1 /. f[z] -> f;
    pola = pol[2] /. s1;
    dpola = D[pola, z];
    Clear[root];
    dim = 10.^-32;
    s = Solve[polfz == 0, f];
    dselect = 0.1(zmax - zmin)/pltpnts;
    root[t_, n_] := (t1 = SetPrecision[t, 64];
        fnt = (f /. s[[n]] /. z -> t1);
        {refnt, imfnt} = {Re[fnt], Im[fnt]};
        If[Abs[imfnt/fnt] < dim, 
          If[Abs[(pola/dpola) /. z -> t1] > dselect, refnt, 
            "almost Infinity"], "complex-valued"]);
    reroot[t_, n_] := (t1 = SetPrecision[t, 64];
        fnt = (f /. s[[n]] /. z -> t1);
        {refnt, imfnt} = {Re[fnt], Im[fnt]};
        If[Abs[imfnt/fnt] > dim, 
          If[Abs[(pola/dpola) /. z -> t1] > dselect, refnt, 
            "almost Infinity"], "real-valued"]);
    pltq = Plot[{tr[root[z1, 1], z1], tr[root[z1, 2], z1]}, {z1, zmin, zmax},
        PlotRange -> {{zmin, zmax}, {fmin, fmax}},
        PlotPoints -> pltpnts,
        PlotStyle -> {{Thickness[0.008], RGBColor[0, 1, 0]}},
        AxesOrigin -> {0, f0},
        DisplayFunction -> Identity
        ];
    pltqre = 
      Plot[{tr[reroot[z1, 1], z1], tr[reroot[z1, 2], z1]}, {z1, zmin, zmax},
        PlotRange -> {{zmin, zmax}, {fmin, fmax}},
        PlotPoints -> pltpnts,
        PlotStyle -> {{Thickness[0.008], RGBColor[0, 1, 0], 
              Dashing[{0.005, 0.015}]}},
        AxesOrigin -> {0, f0},
        DisplayFunction -> Identity
        ];
    On[Plot::"plnr"];
    pltg = Show[pltqre, pltq, pltp, plts, pltex,
      PlotRange -> {{zmin, zmax}, {fmin, fmax}},
      AxesOrigin -> {0, f0},
      TextStyle -> {FontFamily -> "Times New Roman", FontSize -> 12, 
                FontSubstitutions -> {"Math1" -> "Symbol"}},
      AspectRatio -> grat,
      Frame -> True,
FrameLabel -> {"\!\(\* StyleBox[\"z\",\nFontFamily->\"Times\",\nFontSize->16,\nFontSlant->\"Italic\"]\)",
    "\!\(\* StyleBox[\"E\",\nFontFamily->\"Times\",\nFontSize->16,\nFontSlant->\"Italic\"]\)"},
      ImageSize -> 400,
      DisplayFunction -> Identity];

        (*Show together*)
wr1 = 0.04;
wr2 = 0.06;
wr3 = 0.087;
wd = 0.006;
{ybot, ytop} = {0.122, 0.987} + 0.0{1, -1};
pltr = Rectangle[{0, 0} // Scaled, {1 - wr1 - wr2 - wr3 - 3 wd, 1} // Scaled, 
      pltg];
t1 = Text["Sum", {1/2, 1/2},
      TextStyle -> {FontFamily -> "Times New Roman",
          FontSize -> 12, FontWeight -> "Bold",
          FontColor -> RGBColor[1, 0, 0]}];
t2 = Text["Pade", {1/2, 1/2},
      TextStyle -> {FontFamily -> "Times New Roman",
          FontSize -> 12, FontWeight -> "Bold",
          FontColor -> RGBColor[0, 0, 1]}];
t3 = Text["Quadr.", {1/2, 1/2},
      TextStyle -> {FontFamily -> "Times New Roman",
          FontSize -> 12, FontWeight -> "Bold",
          FontColor -> RGBColor[0, 1, 0]}];
g = Graphics[#, Background -> RGBColor[0.9, 0.9, 1]] & /@ {t1, t2, t3};
mindlist1 = mindlist + 1;
y1 = ybot + (ytop - ybot)(mindlist1 - 1)/mindlist1;
y2 = ybot + (ytop - ybot)(mindlist1)/mindlist1;
x2 = 1;
x1 = x2 - wr3;
r3 = Rectangle[{x1, y1} // Scaled, {x2, y2} // Scaled, g[[3]]];
x2 = x1 - wd;
x1 = x2 - wr2;
r2 = Rectangle[{x1, y1} // Scaled, {x2, y2} // Scaled, g[[2]]];
x2 = x1 - wd;
x1 = x2 - wr1;
r1 = Rectangle[{x1, y1} // Scaled, {x2, y2} // Scaled, g[[1]]];
pltindh = {r1, r2, r3};
pltind = Table[{{n1, n2, n3}, prior, nc} = indlist[[n]];
      np1 = Floor[(nc - 1)/2];
      np2 = nc - 1 - np1;
      txt = 
        If[n <= mind, {ToString[nc], 
            "[" <> ToString[np1] <> "/" <> ToString[np2] <> "]", 
            "[" <> ToString[n1] <> ", " <> ToString[n2] <> ", " <> 
              ToString[n3] <> "]"}, {"", "", ""}];
      t = 
        Text[#, {1/2, 1/2}, 
              TextStyle -> {FontFamily -> "Times New Roman", FontSize -> 14, 
                  FontWeight -> If[n == nindlist, "Bold", "Plain"], 
                  FontColor -> RGBColor[0, 0, 0.5]}] & /@ txt;
      col = 
        If[n > mind, RGBColor[1, 1, 1], 
          If[n == nindlist, RGBColor[1, 1, 0], RGBColor[3/4, 3/4, 3/4]]];
      g = Graphics[#, Background -> col] & /@ t;
      y1 = ybot + (ytop - ybot)(mindlist1 - (n + 1))/mindlist1;
      y2 = ybot + (ytop - ybot)(mindlist1 - (n + 1) + 1)/mindlist1;
      x2 = 1;
      x1 = x2 - wr3;
      r3 = Rectangle[{x1, y1} // Scaled, {x2, y2} // Scaled, g[[3]]];
      x2 = x1 - wd;
      x1 = x2 - wr2;
      r2 = Rectangle[{x1, y1} // Scaled, {x2, y2} // Scaled, g[[2]]];
      x2 = x1 - wd;
      x1 = x2 - wr1;
      r1 = Rectangle[{x1, y1} // Scaled, {x2, y2} // Scaled, g[[1]]];
      {r1, r2, r3}, {n, mindlist}];
plt = Show[{pltr, pltindh, pltind} // Graphics,
      AspectRatio -> (1 - wr1 - wr2 - wr3 - 3 wd)grat,
      ImageSize -> 660, 
      DisplayFunction -> If[ifexport, Identity, $DisplayFunction]];
If[ifexport,
      outfile = ToFileName[dir, sect <>ToString[nindlist] <> ".gif"];
      Export[outfile, plt, "GIF", ImageSize -> 660, ImageResolution -> 72, 
    ConversionOptions -> {"AnimationDisplayTime" -> 1, "Loop" -> True, 
        "ColorReductionPalette" -> 32, "ColorReductionDither" -> False}];
        pltseq = Append[pltseq, plt] ];
    , {nindlist, mind}];
If[ifexport,
      outfile = ToFileName[dir, sect <> ".gif"];
      Export[outfile, pltseq, "GIF", ImageSize -> 660, ImageResolution -> 72, 
    ConversionOptions -> {"AnimationDisplayTime" -> 1, "Loop" -> True, 
        "ColorReductionPalette" -> 32, "ColorReductionDither" -> False}];
];
      ];

    "<A NAME=\"" <> sect <> 
        "\"></A><TABLE BORDER=1 CELLPADDING=8 CELLSPACING=0
 BGCOLOR=\"#FFFFFF\"><TR><TD>" // p;
    "<TABLE BORDER=0 CELLPADDING=8 BGCOLOR=\"#FFFFFF\">" // p;
    "<TR ALIGN=\"CENTER\" VALIGN=\"MIDDLE\" BGCOLOR=\"#FEF8CB\"><TH>
The function <I>E</I>(<I>z</I>) found by summation of its power series.
<BR>Dashed line indicates that the approximant is complex valued.
<BR>Red dot marks exact physical energy at <I>z</I>&nbsp;=&nbsp;1." // p;
If[NumberQ[engap],
"<BR>Red circle marks the lowest excited energy level at <I>z</I>&nbsp;=&nbsp;1." // p ];
"<BR>To view results of summation of a specific number of terms of the series, click on the right bar.
</TH></TR>" // p;
    "<TR ALIGN=\"CENTER\" VALIGN=\"MIDDLE\"><TD>" <> "<IMG SRC=\"" <> sect <> 
        ".gif?" <> rndm <> 
        "\" WIDTH=660 HEIGHT=324 BORDER=0 ALT=\"Partial sums, " <>
        "Pade and quadratic approximants\" USEMAP=\"#indbar1\"></TD></TR>" <> "</TABLE></TD></TR></TABLE>
<MAP name=\"indbar1\">" // p;
x0 = 527;
x1 = 658;
ytop = 28;
ybot = 283;
Do[y1 = ytop + n/mindlist (ybot-ytop);
   y0 = ytop + (n-1)/mindlist (ybot-ytop);
  {x0p, x1p, y0p, y1p} = ToString/@Round/@{x0, x1, y0, y1};
  ref = sect <> ToString[n] <> ".gif?" <> rndm;
  {{n1, n2, n3}, prior, nc} = indlist[[n]];
  app = "[" <> ToString[n1] <> ", " <> ToString[n2] <> ", " <> ToString[n3] <> "]";
"<AREA shape=\"rect\" coords=\"" <> x0p <> ","<> y0p <> ","<> x1p <> ","<> y1p <>
 "\" href=\"" <> ref <> "\" onclick=\"vappr=window.open('" <> ref <>
 "','Summation" <> (*ToString[n]*)"" <> "','toolbar=no,status=no,scrollbars=no,
location=no,menubar=no,directories=no,width=670,height=334'); vappr.focus(); return false\">" // p;
   , {n, mind}];
"<AREA shape=\"default\" nohref></MAP>" // p;
printprogr[sect]; 
printnav;

Molecule - icon for Allen-dataBlankExamples of MP seriesBlankSource code of Mathematica programBlankMathematica programsBlankWork in UMass DartmouthWork in UMassDBlankWaste iconUnpublished reports

Designed by A. Sergeev.