// Stephen C. Ferguson
// 11-9-00
//
// Compiled in Delphi 5

unit Formula;

interface

uses Windows, Messages, SysUtils,
     Classes, Controls,
     Dialogs, Menus, ExtCtrls, ExtDlgs,
     complex, Math, MainUnit;

  procedure DrawFractal(pData: Pointer);
  procedure render_method;
  procedure average;
  procedure standard_deviation;
  procedure min_max;
  procedure fractal_dimension;
  procedure gaussian;
  procedure submenu1;
  procedure submenu2;
  procedure submenu3;
  procedure submenu4;
  procedure transform1;
  procedure gaussian2;
  procedure gaussian3;
  procedure dimension4;
  procedure dimension5;
  procedure dimension6;
  procedure dimension8;
  procedure dimension9;
  procedure dimension10;
  procedure dimension11;

implementation

const
  constPI: Double = 3.14159265359;
  const2PI: Double = 6.28318530718;

  constZeroTol: Double = 1e-24;
  constDMax: Double = 1e24;
  constOverFlow = 1e2;

  constMANDEL: Integer = 0;
  constNEWTON: Integer = 1;
  constWARPED_MANDEL: Integer = 2;
  constWARPED_NEWTON: Integer = 3;

var

  //v: pDataPointer;

  xrange, yrange: double;
  xdata, ydata: double;
  xmean, ymean: double;
  x_std, y_std: double;
  dFDx, dFDy: double;
  numerator, numerator_x, numerator_y, denominator: double;
  i: Integer;
  xmax, xmin, ymax, ymin: double;
  temp, xtmp, ytmp: double;
  cx_tmp, cy_tmp: double;
  r, rmin, rmax: double;
  xrad, yrad: double;
  sum, sumsum: double;
  xsum, xsumsum: double;
  ysum, ysumsum: double;
  length_z, length_z1, length_zx, length_zy: double;
  mx, my, zz: double;
  P, K, LSum, LVal, rs: double;
  LSum_x, LSum_y: double;
  XVal, YVal: double;
  length_c: double;
  LTmp, LTmp_x, LTmp_y: double;

  min, max: double;
  ctemp: double;
  ztmp: tcomplex;
  X, Y, T: double;

  zri, zr, zi, cr, ci: double;
  zr1, zi1: double;
  zr2, zi2: double;
  zr3, zi3: double;
  zr4, zi4: double;
  tzr, tzi: double;
  dzr, dzi: double;
  zdr, zdi: double;

  nr, ni: double;  // numerator for complex divide
  dr, di: double;  // denominator for complex divide

  // dmj triangle variables
  sum2, ac, il, lp, az2, lowbound, f, ipower: Double;
  first: Bool;
  xdif, ydif: Double;
  t1, t2, t3, a, b, aa3: tcomplex;

procedure DrawFractal(pData: Pointer);
begin
  //MainForm.v := pData;

  MainForm.v.pXData := @xArray;
  MainForm.v.pYData := @yArray;

  //MainForm.v.StrFormula := 'fractal formula';
  //MainForm.v.StrRender := 'render formula';

  min := constOverFlow;
  max := 0;

  length_z := 0;
  LSum := 1;
  LSum_x := 1;
  LSum_y := 1;
  LVal := constOverFlow;
  XVal := constOverFlow;
  YVal := constOverFlow;
  LTmp := 0;
  LTmp_x := 0;
  LTmp_y := 0;
  MainForm.v.Iteration := -1;
  MainForm.v.xtot := 0;
  MainForm.v.ytot := 0;
  z.x := 0;
  z.y := 0;
  c.x := 0;
  c.y := 0;
  xrange := 0;
  yrange := 0;

  length_z := 0;
  length_c := (constZeroTol+sqrt(MainForm.v.cx*MainForm.v.cx))+(MainForm.v.cy*MainForm.v.cy);
  LSum := 1;
  LVal := 0;

  /////////////////////

  i := 0;
  while i < MainForm.v.maxit do
  begin
    xArray[i] := 0;
    yArray[i] := 0;
    inc(i);
  end;

  cx_tmp := MainForm.v.cx;
  cy_tmp := MainForm.v.cy;

  if MainForm.v.nRotation = 1 then  // Rotate 90
  begin
    xtmp := MainForm.v.cx;
    ytmp := MainForm.v.cy;
    MainForm.v.cx := ytmp;
    MainForm.v.cy := xtmp;

    xtmp := MainForm.v.zx;
    ytmp := MainForm.v.zy;
    MainForm.v.zx := ytmp;
    MainForm.v.zy := xtmp;
  end
  else
  if MainForm.v.nRotation = 2 then // Mirror and rotate 90
  begin
    xtmp := MainForm.v.cx*(-1.0);
    ytmp := MainForm.v.cy*(-1.0);
    MainForm.v.cx := ytmp;
    MainForm.v.cy := xtmp;

    xtmp := MainForm.v.zx*(-1.0);
    ytmp := MainForm.v.zy*(-1.0);
    MainForm.v.zx := ytmp;
    MainForm.v.zy := xtmp;
  end
  else
  if MainForm.v.nRotation = 3 then  // Mirror
  begin
    MainForm.v.cx := MainForm.v.cx*(-1.0);
    MainForm.v.cy := MainForm.v.cy*(-1.0);
    MainForm.v.zx := MainForm.v.zx*(-1.0);
    MainForm.v.zy := MainForm.v.zy*(-1.0);
  end;

  length_c := (constZeroTol+(MainForm.v.cx*MainForm.v.cx))+(MainForm.v.cy*MainForm.v.cy);

  MainForm.draw_fractal;

  MainForm.v.cx := cx_tmp;
  MainForm.v.cy := cy_tmp;
  render_method;

  //Application.ProcessMessages;
end;

//////////////////////////////////////////////////////////

procedure render_method;
begin
  // initialize
  MainForm.v.xtot := 0;
  MainForm.v.ytot := 0;

  xmean := 0;
  ymean := 0;

  xdata := 0;
  ydata := 0;

  x_std := 0;
  y_std := 0;

  dFDx  := 0;
  dFDy  := 0;

  xmax := -constDMAX;
  ymax := -constDMAX;

  xmin := constDMAX;
  ymin := constDMAX;

  rmin := constDMAX;
  rmax := 0;

  xrad := 0;
  yrad := 0;

  sum := 0;
  sumsum := 0;

  xsum := 0;
  xsumsum := 0;

  ysum := 0;
  ysumsum := 0;

  denominator := 0;
  numerator := 0;

  mx := 0;
  my := 0;

  XVal := 0;
  YVal := 0;

  // transform
  case MainForm.v.nTransform of
    1: begin
        i := 0;
        while i < MainForm.v.Iteration do
        begin
          //xArray[i] := 1e-15+(sin((xArray[i])-((1e-15+sin(yArray[i])))));
          //yArray[i] := 1e-15+(cos((yArray[i])-((1e-15+sin(xArray[i])))));

          // Swirly

          xtmp := xarray[i]*MainForm.v.dFactor3;
          ytmp := yarray[i]*MainForm.v.dFactor3;

          if (abs(xtmp) < overflow) and
             (abs(ytmp) < overflow) then
            r := (1e-12+xtmp*xtmp+ytmp*ytmp)
          else
            r := constDMax;

          zr1 := sin(r);
          zr2 := sin(r);

          xdata := zr1 * xtmp - zr2 * ytmp;
          ydata := zr2 * x + zr1 * ytmp;

          xarray[i] := xtmp+xdata*MainForm.v.dFactor5;
          yarray[i] := ytmp+ydata*MainForm.v.dFactor5;

          inc(i);
         end;
       end;
    2: begin
         i := 0;
         while i < MainForm.v.Iteration do
         begin
           //xArray[i] := 1e-15+(sin((xArray[i])+((1e-15+sin(yArray[i])))));
           //yArray[i] := 1e-15+(cos((yArray[i])+((1e-15+sin(xArray[i])))));

           (*
           xnew := sin(x);
           ynew := sin(y);
           fx := fx + MainForm.v.dFactor2 * xnew;
           fy := fy + MainForm.v.dFactor2 * ynew;
           *)

           xtmp := xarray[i]*MainForm.v.dFactor4;
           ytmp := yarray[i]*MainForm.v.dFactor4;

           xarray[i] := xtmp + sin(xtmp)*MainForm.v.dFactor5;
           yarray[i] := ytmp + sin(ytmp)*MainForm.v.dFactor5;

           inc(i);
         end;
       end;

    3: begin
         i := 0;
         while i < MainForm.v.Iteration do
         begin
           //xArray[i] := 1e-15+(sin((xArray[i])+((1e-15+sin(yArray[i])))));
           //yArray[i] := 1e-15+(cos((yArray[i])+((1e-15+sin(xArray[i])))));

           (*
           xnew := sin(x);
           ynew := sin(y);
           fx := fx + MainForm.v.dFactor2 * xnew;
           fy := fy + MainForm.v.dFactor2 * ynew;
           *)

           //xtmp := xarray[i]*MainForm.v.dFactor4;
           //ytmp := yarray[i]*MainForm.v.dFactor4;

           z.x := xarray[i]*MainForm.v.dFactor4;
           z.y := yarray[i]*MainForm.v.dFactor4;

           // complex sine
           z := csin(z);
           xtmp := z.x;
           ytmp := z.y;

           xarray[i] := xtmp*MainForm.v.dFactor5;
           yarray[i] := ytmp*MainForm.v.dFactor5;

           inc(i);
         end;
       end;

    4: begin
        i := 0;
        while i < MainForm.v.Iteration do
        begin
          //xArray[i] := 1e-15+(sin((xArray[i])*((1e-15+sin(yArray[i])))));
          //yArray[i] := 1e-15+(cos((yArray[i])*((1e-15+sin(xArray[i])))));

          (* ------ Inversion
	        xnew := x/(x*x+y*y);
          ynew := y/(x*x+y*y);
          *)

          xtmp := abs(xarray[i]);
          ytmp := abs(yarray[i]);

          if (xtmp = 0) or (ytmp = 0) then
            break;

          xdata := xtmp/(1e-15+xtmp*xtmp+ytmp*ytmp);
          ydata := ytmp/(1e-15+xtmp*xtmp+ytmp*ytmp);

          //xdata := zr1 * xtmp - zr2 * ytmp;
          //ydata := zr2 * x + zr1 * ytmp;

          xarray[i] := xtmp+xdata*MainForm.v.dFactor5;
          yarray[i] := ytmp+ydata*MainForm.v.dFactor5;

          inc(i);
         end;
       end;

    end;

case MainForm.v.render of
1:  begin  // iteraton escape count
      if (MainForm.v.Iteration > 20) then
      begin
        i := MainForm.v.Iteration-15;

        xtmp := abs(xArray[Round(i)]);
        ytmp := abs(yArray[Round(i)]);

        //xtmp := xtmp*MainForm.v.dFactor4;
        //ytmp := ytmp*MainForm.v.dFactor4;

        MainForm.v.xtot := ((i))/MainForm.v.Maxit;
        MainForm.v.ytot := ((i))/MainForm.v.Maxit;

        //MainForm.v.xtot := ln(i) + (ArcTan(xtmp/(constZeroTol+ytmp)))/const2PI;
        //MainForm.v.ytot := ln(i) + (ArcTan(ytmp/(constZeroTol+xtmp)))/const2PI;

        //MainForm.v.xtot := ln(i) - ArcTan(xtmp/(constZeroTol+ytmp));
        //MainForm.v.ytot := ln(i) - ArcTan(ytmp/(constZeroTol+xtmp));

        //MainForm.v.xtot := ln(i/MainForm.v.Maxit);
        //MainForm.v.ytot := ln(i/MainForm.v.Maxit);

        //MainForm.v.xtot := MainForm.v.xtot * (ArcTan(xtmp/(constZeroTol+ytmp)))/const2PI;
        //MainForm.v.ytot := MainForm.v.ytot * (ArcTan(ytmp/(constZeroTol+xtmp)))/const2PI;

        //MainForm.v.xtot := sqrt(1+xtmp);
        //MainForm.v.ytot := sqrt(1+ytmp);

        //MainForm.v.xtot := sqrt(1+i);
        //MainForm.v.ytot := sqrt(1+i);

        //MainForm.v.xtot := Sin(i);
        //MainForm.v.ytot := Sin(i);


      end
      else
      begin
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;
      end;
    end;

2:  begin
      //if MainForm.v.Iteration > 0 then
      begin
        i := 0;
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;

        while i < MainForm.v.Iteration do
        begin
          xtmp := abs(xArray[i]);
          ytmp := abs(yArray[i]);

          xtmp := xtmp*MainForm.v.dFactor4;
          ytmp := ytmp*MainForm.v.dFactor4;

          if (xtmp < constZeroTol) then
            xtmp := constZeroTol;

          if (ytmp < constZeroTol) then
            ytmp := constZeroTol;

          if (xtmp < MainForm.v.dStrands) then
            MainForm.v.xtot := MainForm.v.xtot+(1-xtmp/MainForm.v.dStrands);

          if (ytmp < MainForm.v.dStrands) then
            MainForm.v.ytot := MainForm.v.ytot+(1-ytmp/MainForm.v.dStrands);

          inc(i);
        end;

        if MainForm.v.xtot > constZeroTol then
          MainForm.v.xtot := MainForm.v.xtot/MainForm.v.Maxit;

        if MainForm.v.ytot > constZeroTol then
          MainForm.v.ytot := MainForm.v.ytot/MainForm.v.Maxit;

      end;
    end;

3:  begin
      average;
      standard_deviation;
      min_max;
      fractal_dimension;

      MainForm.v.xtot := abs(dFDx);
      MainForm.v.ytot := abs(dFDy);

    end;

4:  begin // variation of triangle inequality averaging by Damian Jones
(*
dmj-Triangle {
;
; This coloring method implements a variation of Kerry
; Mitchell's Triangle Inequality Average coloring
; method.  Because the smoothing used in this formula
; is based on the dmj-Smooth formula, which only works
; for z^n+c and derivates, the smoothing here will only
; work for those types as well.
;
init:
  float sum = 0.0
  float sum2 = 0.0
  float ac = cabs(#pixel)
  float il = 1/log(@power)
  float lp = log(log(@bailout)/2.0)
  float az2 = 0.0
  float lowbound = 0.0
  float f = 0.0
  BOOL first = true
  float ipower = 1/@apower

loop:
  sum2 = sum
  IF (!first)
    az2 = cabs(#z - #pixel)
    lowbound = abs(az2 - ac)
    IF (@aflavor == 0)
      sum = sum + ((cabs(#z) - lowbound) / (az2+ac - lowbound))^@apower
    ELSEIF (@aflavor == 1)
      sum = sum + 1-(1-(cabs(#z) - lowbound) / (az2+ac - lowbound))^ipower
    ENDIF
  ELSE
    first = false
  ENDIF

final:
  sum = sum / (#numiter)
  sum2 = sum2 / (#numiter-1)
  f = il*lp - il*log(log(cabs(#z)))
  #index = sum2 + (sum-sum2) * (f+1)

default:
  title = "Triangle Inequality Average"

  param apower
    caption = "Average Exponent"
    default = 1.0
    hint = "This skews the values averaged by raising them to \
            this power. Use 1.0 for the classic coloring."
  endparam
  param aflavor
    caption = "Average Flavor"
    default = 0
    enum = "normal" "reversed"
    hint = "Controls whether values are reversed before being \
            raised to a power. Has no effect if Average Exponent \
	    is 1.0."
  endparam
  param power
    caption = "Exponent"
    default = 2.0
    hint = "This should be set to match the exponent of the \
            formula you are using.  For Mandelbrot, this is 2."
  endparam
  param bailout
    caption = "Bailout"
    default = 1e20
    min = 1
    hint = "This should be set to match the bailout value in \
            the Formula tab.  Use a very high bailout!"
  endparam
*)
      // DMJ triangle variables
      sum := 0.0;
      sum2 := 0.0;
      il := 1/ln(1e-15+MainForm.v.ti_power);            // il = 1/log(@power)
      lp := ln(ln(constDMax)/MainForm.v.ti_power);   // lp = log(log(@bailout)/2.0)
      az2 := 0.0;
      lowbound := 0.0;
      f := 0.0;
      first := true;
      ipower := 1/MainForm.v.ti_power;          // ipower = 1/@apower

      cr := MainForm.v.cx;
      ci := MainForm.v.cy;

      ///////////////////////////

      if MainForm.v.Iteration > 2 then
      begin
        //i := 2;
        i := 0;
        while i < MainForm.v.Iteration do
        begin
          // ac is the pixel magnitude
          // az2 is difference between the pixel magnitude and the z magnitude
          // lowbound is the difference between ac and az2
          // sum is the current running average
          // sum2 is the previous running average
          //sum = sum + ((cabs(#z) - lowbound) / (az2+ac - lowbound))^@apower

          sum2 := sum;

          xtmp := xArray[i];
          ytmp := yArray[i];

          xtmp := xtmp*MainForm.v.dFactor4;
          ytmp := ytmp*MainForm.v.dFactor4;

          length_z := sqrt(xtmp*xtmp+ytmp*ytmp);  // z magnitude

          ac := sqrt(cr*cr+ci*ci);  // pixel magnitude
          az2 := sqrt((cr-xtmp)*(cr-xtmp)+(ci-ytmp)*(ci-ytmp));

          if (abs(az2) < constZeroTol) then
            az2 := 1e-24;

          if (abs(az2) > constDMax) then
            az2 := constDMax;

          lowbound := abs(az2-ac);

          numerator := (length_z - lowbound);
          denominator := (az2 + ac - lowbound);

          if (abs(numerator) < constZeroTol) then
            numerator := constZeroTol;

          if (abs(denominator) < constZeroTol) then
            denominator := constZeroTol;

          if (abs(numerator) > constDMax) then
            numerator := constDMax;

          if (abs(denominator) > constDMax) then
            denominator := constDMax;

          LTmp := numerator / denominator;
          sum := sum + LTmp*LTmp;

          inc(i);

        end;

        // final
        if (MainForm.v.Iteration > 2) then
        begin
          xtmp := (xArray[MainForm.v.iteration-1]);
          ytmp := (yArray[MainForm.v.iteration-1]);
          length_z := (constZeroTol+xtmp*xtmp+ytmp*ytmp);

          sum := sum / (MainForm.v.Iteration);
          sum2 := sum2 / (MainForm.v.Iteration-1);

          if (abs(1-length_z) < constZeroTol) then
             length_z := 1.0+constZeroTol;

          numerator := il*ln(abs(ln(length_z)));
          f := il*lp - numerator;
          MainForm.v.xtot := sum2 + (sum-sum2) * (f+1);    // (f+1) is default

          //if abs(MainForm.v.xtot) < constOverFlow then
            //MainForm.v.xtot := abs(MainForm.v.xtot*MainForm.v.xtot);

          //if abs(MainForm.v.ytot) < constOverFlow then

          MainForm.v.ytot := MainForm.v.xtot;
        end
        else
        begin
          MainForm.v.xtot := 0;
          MainForm.v.ytot := MainForm.v.xtot;
        end;

        MainForm.v.xtot := abs(MainForm.v.xtot);
        MainForm.v.ytot := abs(MainForm.v.ytot);

      end;
    end;

5:  begin
      //if MainForm.v.Iteration > 0 then
      begin
        i := 0;
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;

        while i < MainForm.v.Iteration do
        begin
          xtmp := abs(xArray[i]);
          ytmp := abs(yArray[i]);

          xtmp := (xtmp*xtmp+ytmp*ytmp)*MainForm.v.dFactor4;
          ytmp := xtmp;

          if (xtmp < constZeroTol) then
            xtmp := constZeroTol;

          if (ytmp < constZeroTol) then
            ytmp := constZeroTol;

          if (xtmp < MainForm.v.dStrands) then
            MainForm.v.xtot := MainForm.v.xtot+(1-xtmp/MainForm.v.dStrands);

          if (ytmp < MainForm.v.dStrands) then
            MainForm.v.ytot := MainForm.v.ytot+(1-ytmp/MainForm.v.dStrands);

          inc(i);
        end;

        if MainForm.v.xtot > constZeroTol then
          MainForm.v.xtot := MainForm.v.xtot/(MainForm.v.Maxit*MainForm.v.dStrands);

        if MainForm.v.ytot > constZeroTol then
          MainForm.v.ytot := MainForm.v.ytot/(MainForm.v.Maxit*MainForm.v.dStrands);

      end;
    end;

6:  begin
      //if MainForm.v.Iteration > 0 then
      begin
        i := 0;
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;

        while i < MainForm.v.Iteration do
        begin
          xtmp := abs(xArray[i]);
          ytmp := abs(yArray[i]);

          xtmp := (xtmp*xtmp+ytmp*ytmp)*MainForm.v.dFactor4;
          ytmp := xtmp;

          if (xtmp < constZeroTol) then
            xtmp := constZeroTol;

          if (ytmp < constZeroTol) then
            ytmp := constZeroTol;

          if (xtmp < MainForm.v.dStrands) then
            MainForm.v.xtot := MainForm.v.xtot/Mainform.v.Iteration+(1-xtmp/MainForm.v.dStrands);

          if (ytmp < MainForm.v.dStrands) then
            MainForm.v.ytot := MainForm.v.ytot/Mainform.v.Iteration+(1-ytmp/MainForm.v.dStrands);

          inc(i);
        end;

        if MainForm.v.xtot > constZeroTol then
          MainForm.v.xtot := MainForm.v.xtot/(MainForm.v.Maxit*MainForm.v.dStrands);

        if MainForm.v.ytot > constZeroTol then
          MainForm.v.ytot := MainForm.v.ytot/(MainForm.v.Maxit*MainForm.v.dStrands);

      end;
    end;

7:  begin
      //if MainForm.v.Iteration > 0 then
      begin
        i := 0;
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;

        while i < MainForm.v.Iteration do
        begin
          xtmp := abs(xArray[i]);
          ytmp := abs(yArray[i]);

          xtmp := xtmp*MainForm.v.dFactor4;
          ytmp := ytmp*MainForm.v.dFactor4;

          if (xtmp < constZeroTol) then
            xtmp := constZeroTol;

          if (ytmp < constZeroTol) then
            ytmp := constZeroTol;

          (*
          z.x := xtmp;
          z.y := ytmp;
          z := csin(z);
          xdata := z.x*z.x;
          ydata := z.y*z.y;
          *)

  		    //xdata := xtmp-(sin(ytmp-sin(xtmp)));
		      //ydata := ytmp-(sin(xtmp-sin(ytmp)));

  		    //xdata := abs(xtmp*xtmp-abs(ytmp-sin(xtmp)));
		      //ydata := abs(ytmp*ytmp-abs(xtmp-sin(ytmp)));

  		    //xdata := xtmp+(sin(const2PI*ytmp+cos(const2PI*xtmp)));
		      //ydata := ytmp+(sin(const2PI*xtmp+cos(const2PI*ytmp)));

		      xdata := abs(xtmp-ytmp-sin(ytmp));
		      ydata := abs(ytmp-xtmp-sin(xtmp));

          xtmp := abs(xdata);
          ytmp := abs(ydata);

          if (xtmp < MainForm.v.dStrands) then
            MainForm.v.xtot := MainForm.v.xtot+(1-xtmp/MainForm.v.dStrands);

          if (ytmp < MainForm.v.dStrands) then
            MainForm.v.ytot := MainForm.v.ytot+(1-ytmp/MainForm.v.dStrands);

          inc(i);
        end;

        if MainForm.v.xtot > constZeroTol then
          MainForm.v.xtot := MainForm.v.xtot/MainForm.v.Maxit;

        if MainForm.v.ytot > constZeroTol then
          MainForm.v.ytot := MainForm.v.ytot/MainForm.v.Maxit;

      end;
    end;

/////////////////////////////////////////

8:  begin
      //if MainForm.v.Iteration > 1 then
      begin
        i := 1;
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;

        while i < MainForm.v.Iteration do
        begin
          xtmp := abs(xArray[i]);
          ytmp := abs(yArray[i]);

          xtmp := xtmp*MainForm.v.dFactor4;
          ytmp := ytmp*MainForm.v.dFactor4;

          if (xtmp < constZeroTol) then
            xtmp := constZeroTol;

          if (ytmp < constZeroTol) then
            ytmp := constZeroTol;

  		    //xdata := abs(xtmp*xtmp-abs(ytmp-sin(xtmp)));
		      //ydata := abs(ytmp*ytmp-abs(xtmp-sin(ytmp)));

  		    //xdata := xtmp*xtmp-(sin(ytmp)-sin(xtmp));
		      //ydata := ytmp*ytmp-(sin(xtmp)-sin(ytmp));

  		    xdata := xtmp+(sin(const2PI*ytmp+cos(const2PI*xtmp)));
		      ydata := ytmp+(sin(const2PI*xtmp+cos(const2PI*ytmp)));

          xtmp := (xdata);
          ytmp := (ydata);

          if (xtmp < MainForm.v.dStrands) then
            MainForm.v.xtot := MainForm.v.xtot+(1-xtmp/MainForm.v.dStrands);

          if (ytmp < MainForm.v.dStrands) then
            MainForm.v.ytot := MainForm.v.ytot+(1-ytmp/MainForm.v.dStrands);

          inc(i);
        end;

        //MainForm.v.xtot := MainForm.v.xtot));
        //MainForm.v.ytot := (abs(MainForm.v.ytot));

        MainForm.v.xtot := 3/(1+MainForm.v.xtot/MainForm.v.Maxit*MainForm.v.dStrands);
        MainForm.v.ytot := 3/(1+MainForm.v.ytot/MainForm.v.Maxit*MainForm.v.dStrands);

      end;
    end;

9:  begin
      //if MainForm.v.Iteration > 0 then
      begin
        i := 0;
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;

        while i < MainForm.v.Iteration do
        begin
          xtmp := abs(xArray[i]);
          ytmp := abs(yArray[i]);

          xtmp := xtmp*MainForm.v.dFactor4;
          ytmp := ytmp*MainForm.v.dFactor4;

          if (xtmp < constZeroTol) then
            xtmp := constZeroTol;

          if (ytmp < constZeroTol) then
            ytmp := constZeroTol;

          xtmp := xtmp * (xtmp * xtmp + ytmp * ytmp);
          ytmp := ytmp * (ytmp * ytmp + xtmp * xtmp);

          //v.xtot := ArcTan(xtmp/ytmp);
          //v.ytot := ArcTan(ytmp/xtmp);

          // emberJ 1.)
		      //xdata := abs(abs(2*xtmp)-abs(ytmp)-sin(ytmp*ytmp));
		      //ydata := abs(abs(2*ytmp)-abs(xtmp)-sin(xtmp*xtmp));

          //xtmp := xdata;
          //ytmp := ydata;

          if (xtmp < MainForm.v.dStrands) then
            MainForm.v.xtot := MainForm.v.xtot+(1-xtmp/MainForm.v.dStrands);

          if (ytmp < MainForm.v.dStrands) then
            MainForm.v.ytot := MainForm.v.ytot+(1-ytmp/MainForm.v.dStrands);

          inc(i);
        end;

        MainForm.v.xtot := abs(MainForm.v.xtot);
        MainForm.v.ytot := abs(MainForm.v.ytot);

        MainForm.v.xtot := 3/(1+MainForm.v.xtot/MainForm.v.Maxit*MainForm.v.dStrands);
        MainForm.v.ytot := 3/(1+MainForm.v.ytot/MainForm.v.Maxit*MainForm.v.dStrands);

      end;
    end;

10:  begin
      //if MainForm.v.Iteration > 0 then
      begin
        i := 0;
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;

        while i < MainForm.v.Iteration do
        begin
          xtmp := abs(xArray[i]);
          ytmp := abs(yArray[i]);

          xtmp := xtmp*MainForm.v.dFactor4;
          ytmp := ytmp*MainForm.v.dFactor4;

          //xtmp := xtmp * xtmp;
          //ytmp := ytmp * ytmp;

          if (xtmp < constZeroTol) then
            xtmp := constZeroTol;

          if (ytmp < constZeroTol) then
            ytmp := constZeroTol;

          // emberJ 1.)
		      xdata := abs(xtmp-ytmp-sin(ytmp));
		      ydata := abs(ytmp-xtmp-sin(xtmp));

          //xdata := xtmp * ((sin(ytmp) - sin(xtmp)));
          //ydata := ytmp * ((sin(xtmp) - sin(ytmp)));

          xtmp := abs(xdata);
          ytmp := abs(ydata);

          if (xtmp < MainForm.v.dStrands) then
            MainForm.v.xtot := 0.2*MainForm.v.xtot+(1-xtmp/MainForm.v.dStrands);
            //MainForm.v.xtot := 0.2*MainForm.v.xtot+sin(1-xtmp/MainForm.v.dStrands);

          if (ytmp < MainForm.v.dStrands) then
            MainForm.v.ytot := 0.2*MainForm.v.ytot+(1-ytmp/MainForm.v.dStrands);
            //MainForm.v.ytot := 0.2*MainForm.v.ytot+sin(1-ytmp/MainForm.v.dStrands);

          inc(i);
        end;

        MainForm.v.xtot := 0.2*sqrt(abs(MainForm.v.xtot));
        MainForm.v.ytot := 0.2*sqrt(abs(MainForm.v.ytot));

        //MainForm.v.xtot := 10/(1+MainForm.v.xtot/MainForm.v.Maxit*MainForm.v.dStrands);
        //MainForm.v.ytot := 10/(1+MainForm.v.ytot/MainForm.v.Maxit*MainForm.v.dStrands);
      end;
    end;

11:  begin
      if MainForm.v.Iteration > 0 then
      begin
        i := 0;
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;

        while i < MainForm.v.Iteration do
        begin
          xtmp := abs(xArray[i]);
          ytmp := abs(yArray[i]);

          xtmp := xtmp*MainForm.v.dFactor4;
          ytmp := ytmp*MainForm.v.dFactor4;

          if (xtmp = 0) or (ytmp = 0) then
            break;

          xdata := xtmp*xtmp/ytmp*ytmp;
          ydata := ytmp*ytmp/xtmp*xtmp;

          xtmp := xdata;
          ytmp := ydata;

  		    xdata := abs(xtmp-sin(const2PI*(ytmp)));
	  	    ydata := abs(ytmp-sin(const2PI*(xtmp)));

          xtmp := abs(xdata);
          ytmp := abs(ydata);

          if (xtmp < MainForm.v.dStrands) then
            MainForm.v.xtot := 0.2*MainForm.v.xtot+(1-xtmp/MainForm.v.dStrands);

          if (ytmp < MainForm.v.dStrands) then
            MainForm.v.ytot := 0.2*MainForm.v.ytot+(1-ytmp/MainForm.v.dStrands);

          inc(i);
        end;

        MainForm.v.xtot := 0.2*sqrt(abs(MainForm.v.xtot));
        MainForm.v.ytot := 0.2*sqrt(abs(MainForm.v.ytot));

      end;
     end;

12:  begin
      if MainForm.v.Iteration > 0 then
      begin
        i := 0;
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;

        z1.x := 0;
        z1.y := 0;
        z2.x := 0;
        z2.y := 0;

        c.x := MainForm.v.cx;
        c.y := MainForm.v.cy;

        while i < MainForm.v.Iteration do
        begin
          z.x := xArray[i];
          z.y := yArray[i];

          //xtmp := abs(xArray[i]);
          //ytmp := abs(yArray[i]);

          //z := cmulti(z,z);
          z := csin(z);
          //z := ccos(cmulti(z,c));

          xtmp := (z.x);
          ytmp := (z.y);

          xtmp := xtmp*MainForm.v.dFactor4;
          ytmp := ytmp*MainForm.v.dFactor4;

          //xdata := xtmp;
          //ydata := ytmp;

          //xdata := (constZeroTol+xtmp*xtmp/(constZeroTol+ytmp*ytmp));
          //ydata := (constZeroTol+ytmp*ytmp/(constZeroTol+xtmp*xtmp));

          //xdata := (constZeroTol+xtmp+ytmp);
          //ydata := (constZeroTol+ytmp+xtmp);

          //xtmp := xdata;
          //ytmp := ydata;

          //if (xtmp < constZeroTol) then
            //xtmp := constZeroTol;

          //if (ytmp < constZeroTol) then
            //ytmp := constZeroTol;

          // emberJ 1.)
		      //xdata := abs(abs(2*xtmp)-abs(ytmp)-sin(ytmp*ytmp));
		      //ydata := abs(abs(2*ytmp)-abs(xtmp)-sin(xtmp*xtmp));

          xdata := ytmp-sin(ytmp) - cos(xtmp);
          ydata := xtmp-sin(xtmp) - cos(ytmp);

          xtmp := abs(xdata);
          ytmp := abs(ydata);

          //xtmp := abs(xtmp);
          //ytmp := abs(ytmp);

          if ((xtmp) < MainForm.v.dStrands) then
            MainForm.v.xtot := (MainForm.v.xtot*xtmp+(1-xtmp/MainForm.v.dStrands));

          if ((ytmp) < MainForm.v.dStrands) then
            MainForm.v.ytot := (MainForm.v.ytot*ytmp+(1-ytmp/MainForm.v.dStrands));

          (*
          if ((xtmp) < MainForm.v.dStrands) then
            MainForm.v.xtot := (MainForm.v.xtot+sin(1-xtmp/MainForm.v.dStrands));

          if ((ytmp) < MainForm.v.dStrands) then
            MainForm.v.ytot := (MainForm.v.ytot+sin(1-ytmp/MainForm.v.dStrands));
          *)

          (*
          if ((xtmp) < MainForm.v.dStrands) then
            z1 := cadd_CR(z1,(1-xtmp/MainForm.v.dStrands));

          if ((ytmp) < MainForm.v.dStrands) then
            z2 := cadd_CR(z2,(1-ytmp/MainForm.v.dStrands));

          MainForm.v.xtot := cmagnitude(z1);
          MainForm.v.ytot := cmagnitude(z2);
          *)

          inc(i);
        end;

        MainForm.v.xtot := 0.25*sqrt(constZeroTol+abs((MainForm.v.xtot)));
        MainForm.v.ytot := 0.25*sqrt(constZeroTol+abs((MainForm.v.ytot)));

        //MainForm.v.xtot := 3/(1+MainForm.v.xtot/MainForm.v.Maxit*MainForm.v.dStrands);
        //MainForm.v.ytot := 3/(1+MainForm.v.ytot/MainForm.v.Maxit*MainForm.v.dStrands);

      end;
    end;

13:  begin
      if MainForm.v.Iteration > 0 then
      begin
        i := 0;
        MainForm.v.xtot := 0;
        MainForm.v.ytot := 0;

        z1.x := 0;
        z1.y := 0;
        z2.x := 0;
        z2.y := 0;

        c.x := MainForm.v.cx;
        c.y := MainForm.v.cy;

        while i < MainForm.v.Iteration do
        begin
          z.x := xArray[i];
          z.y := yArray[i];

          //xtmp := abs(xArray[i]);
          //ytmp := abs(yArray[i]);

          //z := cmulti(z,z);
          //z := csin(z);
          z := ccos(cmulti(z,c));

          xtmp := (z.x);
          ytmp := (z.y);

          xtmp := xtmp*MainForm.v.dFactor4;
          ytmp := ytmp*MainForm.v.dFactor4;

          //xdata := xtmp;
          //ydata := ytmp;

          //xdata := (constZeroTol+xtmp*xtmp/(constZeroTol+ytmp*ytmp));
          //ydata := (constZeroTol+ytmp*ytmp/(constZeroTol+xtmp*xtmp));

          //xdata := (constZeroTol+xtmp+ytmp);
          //ydata := (constZeroTol+ytmp+xtmp);

          //xtmp := xdata;
          //ytmp := ydata;

          //if (xtmp < constZeroTol) then
            //xtmp := constZeroTol;

          //if (ytmp < constZeroTol) then
            //ytmp := constZeroTol;

          // emberJ 1.)
		      //xdata := abs(abs(2*xtmp)-abs(ytmp)-sin(ytmp*ytmp));
		      //ydata := abs(abs(2*ytmp)-abs(xtmp)-sin(xtmp*xtmp));

          xdata := ytmp-sin(ytmp) - cos(xtmp);
          ydata := xtmp-sin(xtmp) - cos(ytmp);

          xtmp := abs(xdata);
          ytmp := abs(ydata);

          //xtmp := abs(xtmp);
          //ytmp := abs(ytmp);

          //if ((xtmp) < MainForm.v.dStrands) then
            //MainForm.v.xtot := (MainForm.v.xtot*xtmp+(1-xtmp/MainForm.v.dStrands));

          //if ((ytmp) < MainForm.v.dStrands) then
            //MainForm.v.ytot := (MainForm.v.ytot*ytmp+(1-ytmp/MainForm.v.dStrands));

          if ((xtmp) < MainForm.v.dStrands) then
            MainForm.v.xtot := 0.25*(MainForm.v.xtot+(1-xtmp/MainForm.v.dStrands));

          if ((ytmp) < MainForm.v.dStrands) then
            MainForm.v.ytot := (0.25*MainForm.v.ytot+(1-ytmp/MainForm.v.dStrands));

          (*
          if ((xtmp) < MainForm.v.dStrands) then
            z1 := cadd_CR(z1,(1-xtmp/MainForm.v.dStrands));

          if ((ytmp) < MainForm.v.dStrands) then
            z2 := cadd_CR(z2,(1-ytmp/MainForm.v.dStrands));

          MainForm.v.xtot := cmagnitude(z1);
          MainForm.v.ytot := cmagnitude(z2);
          *)

          inc(i);
        end;

        MainForm.v.xtot := 0.25*sqrt(constZeroTol+abs((MainForm.v.xtot)));
        MainForm.v.ytot := 0.25*sqrt(constZeroTol+abs((MainForm.v.ytot)));

        //MainForm.v.xtot := 3/(1+MainForm.v.xtot/MainForm.v.Maxit*MainForm.v.dStrands);
        //MainForm.v.ytot := 3/(1+MainForm.v.ytot/MainForm.v.Maxit*MainForm.v.dStrands);

      end;
    end;

//////////////////////////////////

14:  begin
      gaussian;
      MainForm.v.xtot := xmean;
      MainForm.v.ytot := ymean;

      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);
    end;

15:  begin
      submenu1;
      MainForm.v.xtot := xmean;
      MainForm.v.ytot := ymean;

    end;

16:  begin
      submenu2;
      MainForm.v.xtot := xmean;
      MainForm.v.ytot := ymean;
      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);
    end;

17:  begin
      submenu3;
      MainForm.v.xtot := xmean;
      MainForm.v.ytot := ymean;
      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);
    end;

18:  begin
      i := 1;
      LSum_x := 0;
      LSum_y := 0;

      while i < MainForm.v.Iteration do
      begin
        xtmp := Abs(xArray[i]);
        ytmp := Abs(yArray[i]);

        xtmp := xtmp*MainForm.v.dFactor4;
        ytmp := ytmp*MainForm.v.dFactor4;

        length_z := xtmp*ytmp;

        tmp := xtmp;
        xtmp := ((tmp*ytmp)-tmp+sin(xtmp)-sin(xtmp*xtmp));
        ytmp := ((ytmp*tmp)-tmp+sin(ytmp)-sin(ytmp*ytmp));

        if (xtmp) < MainForm.v.dStrands then
          LSum_x := sqrt(LSum_x + (1 - (xtmp)/MainForm.v.dStrands));

        if (ytmp) < MainForm.v.dStrands then
          LSum_y := sqrt(LSum_y + (1 - (ytmp)/MainForm.v.dStrands));

        inc(i);
      end;

      MainForm.v.xtot := LSum_x/MainForm.v.Maxit;
      MainForm.v.ytot := LSum_y/MainForm.v.Maxit;

    end;

19:  begin
      i := 1;
      LSum_x := 0;
      LSum_y := 0;

      while i < MainForm.v.Iteration do
      begin
        xtmp := Abs(xArray[i]);
        ytmp := Abs(yArray[i]);

        xtmp := xtmp*MainForm.v.dFactor4;
        ytmp := ytmp*MainForm.v.dFactor4;

        length_z := xtmp*ytmp;

        xtmp := 2-1/(0.5+ytmp*xtmp-ytmp-sin(ytmp*ytmp));
        ytmp := 2-1/(0.5+ytmp*xtmp-xtmp-sin(xtmp*xtmp));

        if (xtmp) < MainForm.v.dStrands then
          LSum_x := sqrt(LSum_x + (1 - (xtmp)/MainForm.v.dStrands));

        if (ytmp) < MainForm.v.dStrands then
          LSum_y := sqrt(LSum_y + (1 - (ytmp)/MainForm.v.dStrands));

        inc(i);
      end;

      MainForm.v.xtot := LSum_x/MainForm.v.Maxit;
      MainForm.v.ytot := LSum_y/MainForm.v.Maxit;

    end;

20:  begin
      i := 1;
      LSum_x := 0;
      LSum_y := 0;

      XVal := 0;
      YVal := 0;

      while i < MainForm.v.Iteration do
      begin
        X := Abs(xArray[i]);
        Y := Abs(yArray[i]);

        X := X*MainForm.v.dFactor4;
        Y := Y*MainForm.v.dFactor4;

        length_z := X*X+Y*Y;

        tmp := X;
        xrad := (constZeroTol+sin(X*X));
        yrad := (constZeroTol+sin(Y*Y));

        T := X;
        X := 2-1/(0.5+Y*X-Y-sin(Y*Y));
        Y := 2-1/(0.5+Y*T-T-sin(T*T));

        if (X) < xrad then
          XVal := (1 - (X)/xrad);

        if (Y) < yrad then
          YVal := (1 - (Y)/yrad);

        LSum_x := LSum_x+(1/(1+abs(XVal)));
        LSum_y := LSum_y+(1/(1+abs(YVal)));

        inc(i);
      end;

      MainForm.v.xtot := LSum_x/MainForm.v.Maxit;
      MainForm.v.ytot := LSum_y/MainForm.v.Maxit;

    end;

21:  begin
      i := 0;
      LSum_x := 0;
      LSum_y := 0;

      XVal := 0;
      YVal := 0;

      while i < MainForm.v.Iteration do
      begin
        X := Abs(xArray[i]);
        Y := Abs(yArray[i]);

        if (X > constdMAX) or (Y > constdMAX) then
        begin
          //ShowMessage('break');
          break;
        end;

        X := X*MainForm.v.dFactor4;
        Y := Y*MainForm.v.dFactor4;

        length_z := X*X+Y*Y;

        tmp := X;
        xrad := (constZeroTol+sin(X*X));
        yrad := (constZeroTol+sin(Y*Y));

        T := X;
        X := Y*Y-(X)+(1-Y*sin(X*X));
        Y := T*T-(Y)+(1-T*sin(Y*Y));

        if (X) < xrad then
          XVal := (1 - (X)/(1e-15+xrad));

        if (Y) < yrad then
          YVal := (1 - (Y)/(1e-15+yrad));

        LSum_x := LSum_x+(1/(1+abs(XVal)));
        LSum_y := LSum_y+(1/(1+abs(YVal)));

        inc(i);
      end;

      MainForm.v.xtot := LSum_x/MainForm.v.Maxit;
      MainForm.v.ytot := LSum_y/MainForm.v.Maxit;

    end;

22:  begin
      i := 1;
      LSum_x := 0;
      LSum_y := 0;

      while i < MainForm.v.Iteration do
      begin
        xtmp := Abs(xArray[i]);
        ytmp := Abs(yArray[i]);

        xtmp := xtmp*MainForm.v.dFactor4;
        ytmp := ytmp*MainForm.v.dFactor4;

        length_z := xtmp*ytmp;

        xtmp := ytmp*xtmp-ytmp-sin(ytmp*ytmp);
        ytmp := ytmp*xtmp-xtmp-sin(xtmp*xtmp);

        if (xtmp) < MainForm.v.dStrands then
          LSum_x := sqrt(LSum_x + (1 - (xtmp)/MainForm.v.dStrands));

        if (ytmp) < MainForm.v.dStrands then
          LSum_y := sqrt(LSum_y + (1 - (ytmp)/MainForm.v.dStrands));

        inc(i);
      end;

      MainForm.v.xtot := LSum_x/MainForm.v.Maxit;
      MainForm.v.ytot := LSum_y/MainForm.v.Maxit;

    end;

23:  begin
      gaussian2;
      //MainForm.v.xtot := xmax;
      //MainForm.v.ytot := ymax;
      //MainForm.v.xtot := xtmp;
      //MainForm.v.ytot := ytmp;
      MainForm.v.xtot := (xmin + xmax)*0.5;;
      MainForm.v.ytot := (ymin + ymax)*0.5;

      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);

    end;

24:  begin
      gaussian3;
      //MainForm.v.xtot := xmax;
      //MainForm.v.ytot := ymax;

      //MainForm.v.xtot := xtmp;
      //MainForm.v.ytot := ytmp;

      MainForm.v.xtot := (xmin + xmax)*0.5;;
      MainForm.v.ytot := (ymin + ymax)*0.5;

      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);

     end;

25:  begin
      transform1;

      average;
      standard_deviation;
      min_max;
      dimension4;

      MainForm.v.xtot := dFDx;
      MainForm.v.ytot := dFDy;

      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);

     end;

26:  begin
      //transform1;
      dimension5;

      (*
      if dFDx > 1 then
        MainForm.v.xtot := 1 - 1/dFDx
      else
        MainForm.v.xtot := 1 - dFDx;

      if dFDy > 1 then
        MainForm.v.ytot := 1 - 1/dFDy
      else
        MainForm.v.ytot := 1 - dFDy;

      if MainForm.v.xtot = 0 then
         MainForm.v.xtot := constZeroTol;

      if MainForm.v.ytot = 0 then
         MainForm.v.ytot := constZeroTol;

      *)

      MainForm.v.xtot := dFDx;
      MainForm.v.ytot := dFDy;

      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);

     end;

27:  begin
      dimension6;
      MainForm.v.xtot := xrad;
      MainForm.v.ytot := yrad;

      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);

     end;

28:  begin
      il := 1/ln(2);                   // il = 1/log(@power)
      //lp := ln(ln(1e24)/2);          // lp = log(log(@bailout)/2.0)
      sum := 0;
      xsum := 0;
      ysum := 0;
      if MainForm.v.Iteration > 0 then
      begin
        i := 1;
        while i < MainForm.v.Iteration do
        begin
          xtmp := xArray[i];
          ytmp := yArray[i];

          if (abs(xtmp) < constZeroTol) then
            xtmp := constZeroTol;

          if (abs(ytmp) < constZeroTol) then
            ytmp := constZeroTol;

          length_zx := (xtmp*xtmp/(ytmp*ytmp));
          length_zy := 1/length_zx;

          if (abs(length_zx - 1) <= 1e-15) then
            length_zx := 1e-15 + 1.0;

          if (abs(length_zy - 1) <= 1e-15) then
            length_zy := 1e-15 + 1.0;

          //numerator_x := il*ln(abs(ln(length_zx)));
          numerator_x := ln(abs(ln(length_zx)));
          LTmp_x := (abs(ln(abs(numerator_x))));
          xsum := xsum - (1/(1 + LTmp_x));

          //numerator_y := il*ln(abs(ln(length_zy)));
          numerator_y := ln(abs(ln(length_zy)));
          LTmp_y := (abs(ln(abs(numerator_y))));
          ysum := ysum - (1/(1 + LTmp_y));

          inc(i);
        end;

        if (MainForm.v.Iteration > 1) then
        begin
          xsum := (xsum/(MainForm.v.iteration+ln(2+MainForm.v.iteration)));
          MainForm.v.xtot := xsum;

          ysum := (ysum/(MainForm.v.iteration+ln(2+MainForm.v.iteration)));
          MainForm.v.ytot := ysum;
        end;

        MainForm.v.xtot := abs(MainForm.v.xtot);
        MainForm.v.ytot := abs(MainForm.v.ytot);

      end;
    end;

29:  begin
      dimension8;
      MainForm.v.xtot := dFDx;
      MainForm.v.ytot := dFDy;

      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);

     end;

30:  begin
      dimension9;
      MainForm.v.xtot := 2*dFDx;
      MainForm.v.ytot := 2*dFDy;
      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);
     end;

31:  begin
      dimension10;
      MainForm.v.xtot := dFDx;
      MainForm.v.ytot := dFDy;
      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);
     end;

32:  begin
      dimension11;
      MainForm.v.xtot := dFDx;
      MainForm.v.ytot := dFDy;
      MainForm.v.xtot := abs(MainForm.v.xtot);
      MainForm.v.ytot := abs(MainForm.v.ytot);
     end;

  end;
end;

procedure min_max;
begin
  if MainForm.v.iteration >= 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      if xArray[i] > xmax then
        xmax := xArray[i];

      if xArray[i] < xmin then
        xmin := xArray[i];

      if yArray[i] > ymax then
        ymax := yArray[i];

      if yArray[i] < ymin then
        ymin := yArray[i];

      inc(i);
    end;
    xrange := xmax-xmin;
    yrange := ymax-ymin;
  end;
end;

procedure standard_deviation;
begin
  // calculate the standard deviation
  if MainForm.v.iteration >= 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := (xArray[i] - xmean);
      ydata := (yArray[i] - ymean);

      x_std := (x_std + xdata*xdata);
      y_std := (y_std + ydata*ydata);

      inc(i);
    end;

    x_std := sqrt(2*x_std/(MainForm.v.iteration+constZeroTol)+constZeroTol);
    y_std := sqrt(2*y_std/(MainForm.v.iteration+constZeroTol)+constZeroTol);

    x_std := x_std;
    y_std := y_std;

  end;
end;

procedure average;
begin
  if MainForm.v.iteration >= 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
      begin
      xmean := xmean + xArray[i];
      ymean := ymean + yArray[i];
      inc(i);
    end;
    xmean := ((xmean / (MainForm.v.Iteration+constZeroTol)));
    ymean := ((ymean / (MainForm.v.Iteration+constZeroTol)));
  end;
end;

procedure dimension11;
begin
  if MainForm.v.Iteration > 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xtmp := (abs(xArray[i]));
      ytmp := (abs(yArray[i]));

      xtmp := xtmp*MainForm.v.dFactor4;
      ytmp := ytmp*MainForm.v.dFactor4;

      r := (xtmp*xtmp+ytmp*ytmp);
      xtmp := abs(r/(1+xtmp*xtmp));
      ytmp := abs(r/(1+ytmp*ytmp));

      if xtmp > xmax then
        xmax := xtmp;

      if ytmp > ymax then
        ymax := ytmp;

      if xtmp < xmin then
        xmin := xtmp;

      if ytmp < ymin then
        ymin := ytmp;

      xmean := abs(xtmp)/(i+1+1e-15);
      ymean := abs(ytmp)/(i+1+1e-15);

      xsum := xtmp + xmean;
      ysum := ytmp + ymean;

      x_std := x_std + xsum*xsum;
      y_std := y_std + ysum*ysum;

      inc(i);
    end;

    x_std := sqrt(2*x_std/(MainForm.v.iteration+constZeroTol)+constZeroTol);
    y_std := sqrt(2*y_std/(MainForm.v.iteration+constZeroTol)+constZeroTol);

    denominator := abs(ln((MainForm.v.iteration+constZeroTol))+constZeroTol);

    xrange := xmax - xmin;
    yrange := ymax - ymin;

    if (xrange < constZeroTol) then
      xrange := constZeroTol;

    if (yrange < constZeroTol) then
      yrange := constZeroTol;

    if (Abs(xrange) > constZeroTol) then
      dFDx := abs((abs(ln(abs(Abs(xrange))/((MainForm.v.iteration+constZeroTol)*x_std)))))/denominator;

    if (Abs(yrange) > constZeroTol) then
      dFDy := abs((abs(ln(abs(Abs(yrange))/((MainForm.v.iteration+constZeroTol)*y_std)))))/denominator;

  end;
end;

procedure dimension10;
begin
  if MainForm.v.Iteration > 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := abs(xArray[i]);
      ydata := abs(yArray[i]);

      xdata := xdata*xdata;
      ydata := ydata*ydata;

      if xdata > xmax then
        xmax := xdata;

      if ydata > ymax then
        ymax := ydata;

      if (xdata < xmin) and (Abs(xdata) > constZeroTol) then
        xmin := xdata;

      if (ydata < ymin) and (Abs(ydata) > constZeroTol) then
        ymin := ydata;

      xmean := abs(xtmp)/(i+1e-15);
      ymean := abs(ytmp)/(i+1e-15);

      xsum := xdata + xmean;
      ysum := ydata + ymean;

      x_std := x_std + xsum*xsum;
      y_std := y_std + ysum*ysum;

      inc(i);
    end;

    x_std := sqrt(2*x_std/(MainForm.v.iteration+constZeroTol)+constZeroTol);
    y_std := sqrt(2*y_std/(MainForm.v.iteration+constZeroTol)+constZeroTol);

    denominator := abs(ln((MainForm.v.iteration+constZeroTol))+constZeroTol);

    temp := 1/(MainForm.v.maxit);

    xrange := xmax - xmin;
    yrange := ymax - ymin;

    if (xrange < constZeroTol) then
      xrange := constZeroTol;

    if (yrange < constZeroTol) then
      yrange := constZeroTol;

    dFDx := temp*ln((abs(xsum+1e-15))/xrange);
    dFDy := temp*ln((abs(ysum+1e-15))/yrange);

  end;
end;

procedure dimension9;
begin
  if MainForm.v.Iteration > 0 then
  begin
    i := 0;
    dFDx := 0;
    dFDy := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := abs(xArray[i]);
      ydata := abs(yArray[i]);

      if xdata > constOverFlow then
        xdata := constOverFlow;

      if ydata > constOverFlow then
        ydata := constOverFlow;

      if xdata < constZeroTol then
        xdata := constZeroTol;

      if ydata < constZeroTol then
        ydata := constZeroTol;

      xtmp := Round(xdata);
      ytmp := Round(ydata);

      xtmp := xtmp - xdata;
      ytmp := ytmp - ydata;

      r := (xtmp*xtmp+ytmp*ytmp);
      //xtmp := sin(r*xdata*xdata);
      //ytmp := sin(r*ydata*ydata);

      xtmp := (r+xdata*xdata);
      ytmp := (r+ydata*ydata);

      if xtmp > xmax then
        xmax := xtmp;

      if ytmp > ymax then
        ymax := ytmp;

      if (xtmp < xmin) and (Abs(xtmp) > constZeroTol) then
        xmin := xtmp;

      if (ytmp < ymin) and (Abs(ytmp) > constZeroTol) then
        ymin := ytmp;

      xmean := xtmp/(i+1e-15);
      ymean := ytmp/(i+1e-15);

      xsum := xtmp + xmean;
      xsumsum := xsumsum + xmean*xmean;

      ysum := ytmp + ymean;
      ysumsum := ysumsum + ymean*ymean;

      if abs(xsum) > constDMax then
        xsum := constDMax;

      if abs(ysum) > constDMax then
        ysum := constDMax;

      x_std := x_std + xsum*xsum;
      y_std := y_std + ysum*ysum;

      inc(i);
    end;

    x_std := sqrt(2*x_std/(MainForm.v.iteration+constZeroTol)+constZeroTol);
    y_std := sqrt(2*y_std/(MainForm.v.iteration+constZeroTol)+constZeroTol);

    denominator := abs(ln((MainForm.v.iteration+constZeroTol))+constZeroTol);

    temp := 1/(MainForm.v.maxit);

    xrange := xmax - xmin;
    yrange := ymax - ymin;

    if (xrange < constZeroTol) then
      xrange := constZeroTol;

    if (yrange < constZeroTol) then
      yrange := constZeroTol;

    dFDx := temp*ln((abs(xsum+1e-15))/xrange);
    dFDy := temp*ln((abs(ysum+1e-15))/yrange);

  end;
end;

procedure dimension8;
begin

  if MainForm.v.Iteration > 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := Abs(xArray[i]);
      ydata := Abs(yArray[i]);

      if xdata > constOverFlow then
        xdata := constOverFlow;

      if ydata > constOverFlow then
        ydata := constOverFlow;

      if xdata < constZeroTol then
        xdata := constZeroTol;

      if ydata < constZeroTol then
        ydata := constZeroTol;

      xtmp := sin(xdata*xdata);
      ytmp := sin(ydata*ydata);

      r := sqrt(xtmp*xtmp+ytmp*ytmp);
      xtmp := (r+sin(xdata));
      ytmp := (r+sin(ydata));

      if xtmp > xmax then
        xmax := xtmp;

      if ytmp > ymax then
        ymax := ytmp;

      if (xtmp < xmin) and (Abs(xtmp) > constZeroTol) then
        xmin := xtmp;

      if (ytmp < ymin) and (Abs(ytmp) > constZeroTol) then
        ymin := ytmp;

      xmean := xtmp/(i+1e-15);
      ymean := ytmp/(i+1e-15);

      xsum := xtmp + xmean;
      ysum := ytmp + ymean;

      x_std := x_std + xsum*xsum;
      y_std := y_std + ysum*ysum;

      inc(i);
    end;

    x_std := sqrt(2*x_std/(MainForm.v.iteration+constZeroTol)+constZeroTol);
    y_std := sqrt(2*y_std/(MainForm.v.iteration+constZeroTol)+constZeroTol);

    denominator := abs(ln((MainForm.v.iteration+constZeroTol))+constZeroTol);

    xmax := Abs(xmax);
    xmin := Abs(xmin);
    if xmax - xmin < constZeroTol then
      xmax := xmax +constZeroTol;

    ymax := Abs(ymax);
    ymin := Abs(ymin);
    if ymax - ymin < constZeroTol then
      ymax := ymax +constZeroTol;

    xrange := xmax - xmin;
    yrange := ymax - ymin;

    if (Abs(xrange) > constZeroTol) then
      dFDx := abs((abs(ln(abs(Abs(xrange))/((MainForm.v.iteration+constZeroTol)*x_std)))))/denominator;

    if (Abs(yrange) > constZeroTol) then
      dFDy := abs((abs(ln(abs(Abs(yrange))/((MainForm.v.iteration+constZeroTol)*y_std)))))/denominator;

  end;
end;

procedure fractal_dimension;
begin
  dFDx := 0;
  dFDy := 0;

  if MainForm.v.Iteration > 0 then
  begin
    denominator := abs(ln((MainForm.v.iteration+constZeroTol))+constZeroTol);

    if (Abs(xmax-xmin) > constZeroTol) then
      dFDx := abs(ln(abs(1/(abs(Abs(xrange))/((MainForm.v.iteration+constZeroTol)*x_std)))))/denominator;

    if (Abs(ymax-ymin) > constZeroTol) then
      dFDy := abs(ln(abs(1/(abs(Abs(yrange))/((MainForm.v.iteration+constZeroTol)*y_std)))))/denominator;

  end;
end;

procedure dimension4;
begin
  dFDx := 0;
  dFDy := 0;

  if MainForm.v.Iteration > 0 then
  begin
    denominator := abs(ln((MainForm.v.iteration+constZeroTol))+constZeroTol);

    if (Abs(xmax-xmin) > constZeroTol) then
      dFDx := abs((abs(ln(abs(Abs(xrange))/((MainForm.v.iteration+constZeroTol)*x_std)))))/denominator;

    if (Abs(ymax-ymin) > constZeroTol) then
      dFDy := abs((abs(ln(abs(Abs(yrange))/((MainForm.v.iteration+constZeroTol)*y_std)))))/denominator;

  end;
end;

procedure gaussian2;
begin
  if MainForm.v.Iteration > 1 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := Abs(xArray[i]);
      ydata := Abs(yArray[i]);

      xdata := xdata*MainForm.v.dFactor4;
      ydata := ydata*MainForm.v.dFactor4;

      if xdata > constOverFlow then
        xdata := constOverFlow;

      if ydata > constOverFlow then
        ydata := constOverFlow;

      if xdata < constZeroTol then
        xdata := constZeroTol;

      if ydata < constZeroTol then
        ydata := constZeroTol;

      xtmp := Round(xdata);
      ytmp := Round(ydata);

      xtmp := xtmp - xdata;
      ytmp := ytmp - ydata;

      r := (xtmp*xtmp+ytmp*ytmp);
      xtmp := r+xtmp*xtmp;
      ytmp := r+ytmp*ytmp;

      if xtmp > xmax then
        xmax := xtmp;

      if ytmp > ymax then
        ymax := ytmp;

      if (xtmp < xmin) and (xtmp > constZeroTol) then
        xmin := xtmp;

      if (ytmp < ymin) and (ytmp > constZeroTol) then
        ymin := ytmp;

      inc(i);
    end;
  end;
end;

procedure gaussian3;
begin
  if MainForm.v.Iteration > 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := Abs(xArray[i]);
      ydata := Abs(yArray[i]);

      xdata := xdata*MainForm.v.dFactor4;
      ydata := ydata*MainForm.v.dFactor4;

      if xdata > constOverFlow then
        xdata := constOverFlow;

      if ydata > constOverFlow then
        ydata := constOverFlow;

      if xdata < constZeroTol then
        xdata := constZeroTol;

      if ydata < constZeroTol then
        ydata := constZeroTol;

      xtmp := Round(xdata);
      ytmp := Round(ydata);

      xtmp := abs(xtmp - xdata);
      ytmp := abs(ytmp - ydata);

      r := sqrt(xtmp*xtmp+ytmp*ytmp);

      //xtmp := r+xtmp*xtmp;
      //ytmp := r+ytmp*ytmp;

      xtmp := r;
      ytmp := r;

      if xtmp > xmax then
        xmax := xtmp;

      if ytmp > ymax then
        ymax := ytmp;

      if (xtmp < xmin) and (xtmp > constZeroTol) then
        xmin := xtmp;

      if (ytmp < ymin) and (ytmp > constZeroTol) then
        ymin := ytmp;

      inc(i);
    end;
  end;
end;

procedure dimension5;
begin
  if MainForm.v.Iteration > 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := Abs(xArray[i]);
      ydata := Abs(yArray[i]);

      xdata := xdata*MainForm.v.dFactor4;
      ydata := ydata*MainForm.v.dFactor4;

      if xdata > constOverFlow then
        xdata := constOverFlow;

      if ydata > constOverFlow then
        ydata := constOverFlow;

      if xdata < constZeroTol then
        xdata := constZeroTol;

      if ydata < constZeroTol then
        ydata := constZeroTol;

      xtmp := Round(xdata);
      ytmp := Round(ydata);

      xtmp := xtmp - xdata;
      ytmp := ytmp - ydata;

      r := (xtmp*xtmp+ytmp*ytmp);
      //r := xdata*xdata+ydata*ydata;
      xtmp := r*xdata;
      ytmp := r*ydata;

      if xmin > xtmp then
      begin
        xmin := xtmp;
        MainForm.v.xtot := MainForm.v.iteration;
      end;

      if xmax < xtmp then
      begin
        xmax := xtmp;
      end;

      if ymin > ytmp then
      begin
        ymin := ytmp;
        MainForm.v.ytot := MainForm.v.iteration;
      end;

      if ymax < ytmp then
      begin
        ymax := ytmp;
      end;

        xsum := xsum + xtmp;
        xsumsum := xsumsum + xtmp*xtmp;

        ysum := ysum + ytmp;
        ysumsum := ysumsum + ytmp*ytmp;

      inc(i);
    end;

    temp := 1/(MainForm.v.maxit);
    if (xmax - xmin) > constZeroTol then
    begin
      dFDx := temp*ln((2*xsumsum-(xsum*sum)/Abs(xmax-xmin)));
    end;

    if (ymax - ymin) > constZeroTol then
    begin
      dFDy := temp*ln((2*ysumsum-(ysum*sum)/Abs(ymax-ymin)));
    end;
  end;
end;

procedure dimension6;
begin
  if MainForm.v.Iteration > 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := xArray[i];
      ydata := yArray[i];

      xdata := xdata*MainForm.v.dFactor4;
      ydata := ydata*MainForm.v.dFactor4;

      if abs(xdata) > constOverFlow then
        xdata := constOverFlow;

      if abs(ydata) > constOverFlow then
        ydata := constOverFlow;

      if abs(xdata) < constZeroTol then
        xdata := constZeroTol;

      if abs(ydata) < constZeroTol then
        ydata := constZeroTol;

      //xdata := xdata*xdata - ydata;
      //ydata := ydata*ydata - xdata;

      xdata := xdata - sin(ydata);
      ydata := ydata - sin(xdata);

      xtmp := Round(xdata);
      ytmp := Round(ydata);

      xtmp := (1e-15+xtmp - xdata);
      ytmp := (1e-15+ytmp - ydata);

      r := xtmp*xtmp+ytmp*ytmp;
      //xtmp := r;

      xtmp := r+xtmp*xtmp;
      ytmp := r+ytmp*ytmp;

      if xmin > xtmp then
      begin
        xmin := xtmp;
      end;

      if xmax < xtmp then
      begin
        xmax := xtmp;
        xrad := xtmp;
        xsum  := xtmp;
     end;

      if ymin > ytmp then
      begin
        ymin := ytmp;
      end;

      if ymax < ytmp then
      begin
        ymax := ytmp;
        yrad := ytmp;
        ysum  := ytmp;
      end;

      inc(i);
    end;

    xrad := xsum;
    yrad := ysum;

  end;
end;

procedure transform1;
begin
  i := 0;
  while i < MainForm.v.Iteration do
  begin
    xtmp := xArray[i];
    ytmp := yArray[i];

    xtmp := xtmp*MainForm.v.dFactor4;
    ytmp := ytmp*MainForm.v.dFactor4;

    if xtmp > constOverFlow then
      xtmp := constOverFlow;

    if Abs(xtmp) < constZeroTol then
      xtmp := constZeroTol;

    if ytmp > constOverFlow then
      ytmp := constOverFlow;

    if Abs(ytmp) < constZeroTol then
      ytmp := constZeroTol;

    xArray[i] := sin(xtmp);
    yArray[i] := sin(ytmp);

    //xArray[i] := xtmp/ytmp;
    //yArray[i] := ytmp/xtmp;

    inc(i);
  end;
end;

procedure gaussian;
begin
  if MainForm.v.iteration > 1 then
  begin
    i := 1;
    while i < MainForm.v.iteration do
    begin
      xdata := (xArray[i]);
      ydata := (yArray[i]);

      xdata := xdata*MainForm.v.dFactor4;
      ydata := ydata*MainForm.v.dFactor4;

      if abs(xdata) > constOverFlow then
        xdata := constOverFlow;

      if abs(ydata) > constOverFlow then
        ydata := constOverFlow;

      temp := xdata*xdata+ydata*ydata;
      if temp < constZeroTol then
        temp := constZeroTol;

      //xtmp := xtmp-sin(temp);
      //ytmp := ytmp-sin(temp);

      xtmp := temp;
      ytmp := temp;

      xdata := Round(xtmp);
      ydata := Round(ytmp);

      xmean := xmean + (Abs(xtmp - xdata))/(i+ln(i+2+1e-15));
      ymean := ymean + (Abs(ytmp - ydata))/(i+ln(i+2+1e-15));

      inc(i);
    end;
    xmean := (xmean / (ln(MainForm.v.iteration+1+1e-15)));
    ymean := (ymean / (ln(MainForm.v.iteration+1+1e-15)));

  end;
end;

procedure submenu1;
begin
  if MainForm.v.iteration > 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := Abs(xArray[i]);
      ydata := Abs(yArray[i]);

      xdata := xdata*MainForm.v.dFactor4;
      ydata := ydata*MainForm.v.dFactor4;

      if xdata > constOverFlow then
        xdata := constOverFlow;

      if ydata > constOverFlow then
        ydata := constOverFlow;

      xtmp := sin(xdata*(xdata+constZeroTol));
      ytmp := sin(ydata*(ydata+constZeroTol));

      if xtmp < constZeroTol then
        xtmp := constZeroTol;

      if ytmp < constZeroTol then
        ytmp := constZeroTol;

      xtmp := sin(Abs(xtmp));
      ytmp := sin(Abs(ytmp));

      xdata := Round(xtmp);
      ydata := Round(ytmp);

      xmean := xmean + (Abs(xtmp - xdata))/(i+ln(i+2+1e-15));
      ymean := ymean + (Abs(ytmp - ydata))/(i+ln(i+2+1e-15));

      inc(i);
    end;

    xmean := 3*(xmean / (ln(MainForm.v.iteration+1+1e-15)));
    ymean := 3*(ymean / (ln(MainForm.v.iteration+1+1e-15)));

  end;
end;

procedure submenu2;
begin
  if MainForm.v.iteration > 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := Abs(xArray[i]);
      ydata := Abs(yArray[i]);

      xdata := xdata*MainForm.v.dFactor4;
      ydata := ydata*MainForm.v.dFactor4;

      if xdata > constOverFlow then
        xdata := constOverFlow;

      if ydata > constOverFlow then
        ydata := constOverFlow;

      xtmp := abs(sin((xdata+constZeroTol)));
      ytmp := abs(sin((ydata+constZeroTol)));

      if xtmp < constZeroTol then
        xtmp := constZeroTol;

      if ytmp < constZeroTol then
        ytmp := constZeroTol;

      xtmp := sin(Abs(xtmp+(ytmp+constZeroTol)));
      ytmp := sin(Abs(ytmp+(xtmp+constZeroTol)));

      xdata := Round(xtmp);
      ydata := Round(ytmp);

      xmean := xmean + (Abs(xtmp - xdata))/(i+ln(i+2+1e-15));
      ymean := ymean + (Abs(ytmp - ydata))/(i+ln(i+2+1e-15));

      inc(i);
    end;

    xmean := 3*(xmean / (ln(MainForm.v.iteration+1+1e-15)));
    ymean := 3*(ymean / (ln(MainForm.v.iteration+1+1e-15)));

  end;
end;

procedure submenu3;
begin
  if MainForm.v.iteration > 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := Abs(xArray[i]);
      ydata := Abs(yArray[i]);

      xdata := xdata*MainForm.v.dFactor4;
      ydata := ydata*MainForm.v.dFactor4;

      if xdata > constOverFlow then
        xdata := constOverFlow;

      if ydata > constOverFlow then
        ydata := constOverFlow;

      xtmp := sin((xdata*ydata+constZeroTol));
      ytmp := cos((xdata*ydata+constZeroTol));

      if xtmp < constZeroTol then
        xtmp := constZeroTol;

      if ytmp < constZeroTol then
        ytmp := constZeroTol;

      xtmp := sin(Abs(xtmp));
      ytmp := sin(Abs(ytmp));

      xdata := Round(xtmp);
      ydata := Round(ytmp);

      xmean := xmean + (Abs(xtmp - xdata))/(i+ln(i+2+1e-15));
      ymean := ymean + (Abs(ytmp - ydata))/(i+ln(i+2+1e-15));

      inc(i);
    end;
    xmean := 3*(xmean / (ln(MainForm.v.iteration+1+1e-15)));
    ymean := 3*(ymean / (ln(MainForm.v.iteration+1+1e-15)));
  end;
end;

procedure submenu4;
begin
  ctemp := MainForm.v.cx*MainForm.v.cx+MainForm.v.cy*MainForm.v.cy;
  if MainForm.v.iteration > 0 then
  begin
    i := 0;
    while i < MainForm.v.Iteration do
    begin
      xdata := Abs(xArray[i]);
      ydata := Abs(yArray[i]);

      xdata := xdata*MainForm.v.dFactor4;
      ydata := ydata*MainForm.v.dFactor4;

      if xdata > constOverFlow then
        xdata := constOverFlow;

      if ydata > constOverFlow then
        ydata := constOverFlow;

      xtmp := Abs(Abs(MainForm.v.cx)-(xdata+constZeroTol));
      ytmp := Abs(Abs(MainForm.v.cy)-(ydata+constZeroTol));

      xtmp := (xtmp+constZeroTol);
      ytmp := (ytmp+constZeroTol);

      if xtmp < constZeroTol then
        xtmp := constZeroTol;

      if ytmp < constZeroTol then
        ytmp := constZeroTol;

      xmean := xmean + (Abs(xtmp - xdata))/(i+ln(i+2+1e-15));
      ymean := ymean + (Abs(ytmp - ydata))/(i+ln(i+2+1e-15));

      inc(i);
    end;

    xmean := 3*(xmean / (ln(MainForm.v.iteration+1+1e-15)));
    ymean := 3*(ymean / (ln(MainForm.v.iteration+1+1e-15)));

  end;
end;

end.