• August 02, 2021, 08:23:13 AM

### Author Topic:  rational function  (Read 1619 times)

0 Members and 1 Guest are viewing this topic.

• Fractal Frankfurter
• Posts: 541
##### rational function
« on: June 17, 2021, 05:16:22 PM »
I want to draw   Julia set for rational function (image bo36.png)  by  Michael Becker
Quote
f(z)=1/((0,15+0,15i)z5+z3+(-3+3i)z), dargestellt auf [-3;3]x[-3;3].

It look like connected Julia set. To speed up computaions I want to find critical points and periodic points

Code: [Select]
kill(all);remvalue(all);display2d:false;f(z):= 1/(a*z^5+ z^3 + b*z);/* converts complex number z = x*y*%i to the list in a draw format:  [x,y] */d(z):=[float(realpart(z)), float(imagpart(z))]$ToPoints(myList):= points(map('d,myList))$/* give Draw List from one point*/ToPoint(z):=points([d(z)])$/* f(z) is used as a global function I do not know how to put it as a argument */GiveOrbit(z0):=block( [z,Orbit, iLength ], iLength : 100, z:z0, Orbit:[z0], for i:1 thru iLength step 1 do ( z:f(z), z: float(z), z:rectform(z), Orbit:endcons(z,Orbit)), Orbit:ToPoints(Orbit), return(Orbit) )$/* find fixed points  returns a list */GiveFixedPoints():= block(  [s,z],  s:solve(f(z)=z),  /* remove "z="  from list s */  s:map('rhs,s),  s:map('rectform,s),  s:map('float,s),  return(s))$/*compile(all);*/a: (0.15+0.15*%i);b: (-3+3*%i);dz: diff(f(z) ,z,1)$dz:ratsimp(float(ev(dz)));  /* dz = -((150*%i+150)*z^4+600*z^2+600*%i-600) /(9*%i*z^10+(60*%i+60)*z^8-160*z^6+(1200*%i-1200)*z^4-3600*%i*z^2)*/s: solve(dz=0)$s:map('rhs,s)$s:map('float,s)$s:map('rectform,s)$s:map('float,s);/* [0.389374737779177*%i-0.9400337727919567,        0.9400337727919567-0.389374737779177*%i,        (-1.816005797447402*%i)-0.7522142306508819,        1.816005797447402*%i+0.7522142306508819] *//* r1:solve(f(z)=z); */iLength:100;orbit1 :GiveOrbit(s[1])$orbit2 :GiveOrbit(s[2])$orbit3 :GiveOrbit(s[3])$orbit4 :GiveOrbit(s[4])$path:"~/Dokumenty/ijon/3/36/"\$  /*  if empty then file is in a home dir */load(draw); /* ( interface to gnuplot ) by Mario Rodriguez Riotorto http://www.telefonica.net/web2/biomates */draw2d(    title = concat("Critical orbit for f(z)=z^2 +", string(c)),    terminal  = png,    user_preamble = "set size square", /*    */    file_name = concat(path ,"8"),    dimensions = [1500,1500],    /* Since Maxima 5.23, pic_width and pic_height are deprecated. */    xrange = [-3,3],    yrange = [-3,3],    xlabel     = "z.re ",    ylabel     = "z.im",            point_type    = 7,     points_joined = false,    point_size    = 0.8,    key=" critical orbit ",    color             =red,    orbit1,    orbit2,    orbit3,    orbit4,       point_size    = 1.2,    key="",    ToPoint(s[1]),    ToPoint(s[2]),    ToPoint(s[3]),    key= "critical point",    color           = blue,    ToPoint(s[4])        /*        key= "fixed point",    color           = black,    alfa*/ );
My image of critical orbits  shows that critical orbits tend to z=0, but on the image it seems to be in the exterior
How to find periodic points ?

« Last Edit: June 17, 2021, 06:15:18 PM by Adam Majewski »

#### claude

• Posts: 1938
##### Re: rational function
« Reply #1 on: June 17, 2021, 05:27:02 PM »
critical orbits tend to z=0
Are you sure about that?  Possibly tending to fixed points near (but not equal to) 0?  I have not investigated, but it looks like there is a gap in the middle, in both images.

EDIT fixed points and period 2 points:
Code: [Select]
(%i10) ToPoints(solve(z = f(z), z));rat: replaced 0.15 by 3/20 = 0.15rat: replaced 0.15 by 3/20 = 0.15(%o10) points([[z = - 0.1728968618954043, 0.0 = 0.4517976427751508], [z = 0.1728968618954043, 0.0 = - 0.4517976427751508], [z = - 1.542012716025067, 0.0 = 0.5826223205347855], [z = 1.542012716025067, 0.0 = - 0.5826223205347855], [z = - 1.052616733765589, 0.0 = - 2.511067897204381], [z = 1.052616733765589, 0.0 = 2.511067897204381]])(%i18) map('allroots, solve(f(f(z))=z));rat: replaced 0.15 by 3/20 = 0.15rat: replaced 0.15 by 3/20 = 0.15rat: replaced 0.15 by 3/20 = 0.15rat: replaced 0.15 by 3/20 = 0.15rat: replaced 0.15 by 3/20 = 0.15rat: replaced 0.15 by 3/20 = 0.15rat: replaced 0.15 by 3/20 = 0.15rat: replaced 0.15 by 3/20 = 0.15(%o18) [[z = 0.0], [z = 0.1004507524894826 - 0.04160806403172667 %i, z = 0.04160806403172667 %i - 0.1004507524894826, z = (- 0.1972126334382917 %i) - 0.4417257204225227, z = 0.1728968618954057 - 0.4517976427751519 %i, z = 0.1972126334382916 %i + 0.4417257204225231, z = 0.4517976427751519 %i - 0.172896861895406, z = (- 2.008427151852299 %i) - 0.8319177653355924, z = 1.486804283840914 - 0.6158544989534238 %i, z = 0.5826223205446455 %i - 1.542012716044921, z = 0.6783914544446792 %i - 1.502343841878277, z = 1.545643025048946 - 0.6402263035701801 %i, z = 1.542012716021183 - 0.5826223205577208 %i, z = 0.6158544989143758 %i - 1.48680428382238, z = 0.6402263036121009 %i - 1.545643025068914, z = 1.502343841903455 - 0.678391454434491 %i, z = 2.511067896520314 %i + 1.052616752123668, z = 2.008427151854269 %i + 0.83191776533551, z = (- 2.519905568763936 %i) - 1.031280709377238, z = (- 2.511067897025448 %i) - 1.052616732010429, z = 2.519905568548988 %i + 1.03128068923336, z = 2.496853258644555 %i + 1.034230487296836, z = (- 2.532331263153173 %i) - 1.048925954202523, z = (- 2.496853270435978 %i) - 1.034230487368577, z = 2.532331275662722 %i + 1.048925954304978]]
« Last Edit: June 17, 2021, 05:41:13 PM by claude, Reason: points »

• Fractal Frankfurter
• Posts: 541
##### Re: rational function
« Reply #2 on: June 17, 2021, 05:37:10 PM »
Are you sure about that?  Possibly tending to fixed points near (but not equal to) 0?  I have not investigated, but it looks like there is a gap in the middle, in both images.

Ok. It looks that z=0 will be inside ( ?center) of very small component which is smth like limit of white components.

So period 1 point will be not helpfull for me.

I ask about it because usually testing for attracting periodic points makes the code faster.

#### claude

• Posts: 1938
##### Re: rational function
« Reply #3 on: June 17, 2021, 05:42:42 PM »
see edit, period 2 points include 0 and +/-(0.04160806403172667 %i - 0.1004507524894826)

• Fractal Frankfurter
• Posts: 541
##### Re: rational function
« Reply #4 on: June 17, 2021, 07:29:01 PM »
Code: [Select]
/*https://web.archive.org/web/20161024194536/http://www.ijon.de/mathe/julia/some_julia_sets_3.htmla: (0,15+0,15i)b: (-3+3i)f: 1/(a*z^5+ b*z^3 + z); dz = -(5*a*z^4+3*b*z^2+1)/(a*z^5+z^3+b*z)^2  Adam Majewski  adammaj1 aaattt o2 dot pl  // o like oxygen not 0 like zero         Structure of a program or how to analyze the program       ============== Image X ========================    DrawImageOf -> DrawPointOf -> ComputeColorOf ( FunctionTypeT FunctionType , complex double z) -> ComputeColor      check only last function  which computes color of one pixel for given Function Type         ==========================================    ---------------------------------  indent d.c   default is gnu style   -------------------  c console progam     export  OMP_DISPLAY_ENV="TRUE"   gcc d.c -lm -Wall -march=native -fopenmp  time ./a.out > b.txt  gcc d.c -lm -Wall -march=native -fopenmp  time ./a.out  time ./a.out >i.txt  time ./a.out >e.txt              convert -limit memory 1000mb -limit disk 1gb dd30010000_20_3_0.90.pgm -resize 2000x2000 10.png      */#include <stdio.h>#include <stdlib.h> // malloc#include <string.h> // strcat#include <math.h> // M_PI; needs -lm also#include <complex.h>#include <omp.h> // OpenMP#include <limits.h> // Maximum value for an unsigned long long int// https://sourceforge.net/p/predef/wiki/Standards/#if defined(__STDC__)#define PREDEF_STANDARD_C_1989#if defined(__STDC_VERSION__)#if (__STDC_VERSION__ >= 199409L)#define PREDEF_STANDARD_C_1994#endif#if (__STDC_VERSION__ >= 199901L)#define PREDEF_STANDARD_C_1999#endif#endif#endif/* --------------------------------- global variables and consts ------------------------------------------------------------ *///FunctionTypetypedef enum  {Fatou = 0, IntLSM =1 , ExtLSM = 2, LSM = 3, DEM = 4, Unknown = 5 , BD = 6, MBD = 7 , SAC = 8, DLD = 9, ND = 10 , NP= 11, POT = 12 , Blend = 13 } FunctionTypeT; // FunctionTypeT FunctionType;// virtual 2D array and integer ( screen) coordinate// Indexes of array starts from 0 not 1 //unsigned int ix, iy; // varstatic unsigned int ixMin = 0; // Indexes of array starts from 0 not 1static unsigned int ixMax; //static unsigned int iWidth; // horizontal dimension of arraystatic unsigned int iyMin = 0; // Indexes of array starts from 0 not 1static unsigned int iyMax; //static unsigned int iHeight = 1000; //  // The size of array has to be a positive constant integer static unsigned long long int iSize; // = iWidth*iHeight; // memmory 1D array unsigned char *data;unsigned char *edge;//unsigned char *edge2;// unsigned int i; // var = index of 1D array//static unsigned int iMin = 0; // Indexes of array starts from 0 not 1static unsigned int iMax; // = i2Dsize-1  = // The size of array has to be a positive constant integer // unsigned int i1Dsize ; // = i2Dsize  = (iMax -iMin + 1) =  ;  1D array with the same size as 2D array// see SetPlanedouble radius = 3.0; complex double center = 0.0 ;double  DisplayAspectRatio  = 1.0; // https://en.wikipedia.org/wiki/Aspect_ratio_(image)// dx = dy compare setup : iWidth = iHeight;double ZxMin; //= -1.3; //-0.05;double ZxMax;// = 1.3; //0.75;double ZyMin;// = -1.3; //-0.1;double ZyMax;// = 1.3; //0.7;double PixelWidth; // =(ZxMax-ZxMin)/ixMax;double PixelHeight; // =(ZyMax-ZyMin)/iyMax;// demdouble BoundaryWidth ; //= 1.0*iWidth/2000.0  ; //  measured in pixels ( when iWidth = 2000) double distanceMax ; //= BoundaryWidth*PixelWidth;double ratio; /*  ER = pow(10,ERe);  AR = pow(10,-ARe);*///int ARe ; // increase ARe until black ( unknown) points disapear //int ERe ;double ER;double ER2; //= 1e60;double AR; // bigger values do not worksdouble AR2;double AR_max;//double AR12;int IterMax = 100000;int IterMax_LSM = 1000;int IterMax_DEM = 100000;/* colors = shades of gray from 0 to 255    unsigned char colorArray[2][2]={{255,231},    {123,99}};   color = 245;  exterior */unsigned char iColorOfExterior = 245;unsigned char iColorOfInterior = 99;//unsigned char iColorOfInterior2 = 183;unsigned char iColorOfBoundary = 0;unsigned char iColorOfUnknown = 5;// pixel countersunsigned long long int uUnknown = 0;unsigned long long int uInterior = 0;unsigned long long int uExterior = 0;// critical points// critical points//complex double zcr = -0.35*I; // only one critical point //complex double zc2 = -2.2351741790771484375e-08+9.4296410679817199707e-09*I;/*  0.389374737779177*%i-0.9400337727919567,        0.9400337727919567-0.389374737779177*%i,        (-1.816005797447402*%i)-0.7522142306508819,        1.816005797447402*%i+0.7522142306508819]*/const complex double zp2 = 0.389374737779177*I-0.9400337727919567; //-0.04160806403172667*I - 0.1004507524894826;/* a: (0.15+0.15*%i);b: (-3+3*%i);f: 1/(a*z^5+ b*z^3 + z); dz = -(5*a*z^4+3*b*z^2+1)/(a*z^5+b*z^3+z)^2*/const complex double a =  0.15 + 0.15*I;const complex double b = -3.0  + 3.0*I;//const int period = 1;// periodic points = attractors//complex double z1 =  0.0 ; //fixed point (period  1) =  attracting cycle/* ------------------------------------------ functions -------------------------------------------------------------*/// complex functioncomplex double f(const complex double z0) { double complex z = z0; z = 1/(a*z*z*z*z*z + z*z*z + b*z); return  z; } // complex double dfz(const complex double z0) { // dz= -(5*a*z^4+3*z^2+b) / (a*z^5+z^3+b*z)^2 double complex z = z0; complex double z2= z*z; complex double numerator = -5.0*a*z2*z2 - 3.0*z2 - b ; complex double denom = a*z2*z2*z + z2*z + b*z; denom = denom*denom; return  numerator/denom; } //------------------complex numbers -----------------------------------------------------// from screen to world coordinate ; linear mapping// uses global consdoubleGiveZx (int ix){  return (ZxMin + ix * PixelWidth);}// uses globaal consdoubleGiveZy (int iy){  return (ZyMax - iy * PixelHeight);} // reverse y axiscomplex doubleGiveZ (int ix, int iy){  double Zx = GiveZx (ix);  double Zy = GiveZy (iy);  return Zx + Zy * I;}double cabs2(complex double z){  return creal(z)*creal(z)+cimag(z)*cimag(z);}// ****************** DYNAMICS = trap tests ( target sets) ****************************/* -----------  array functions = drawing -------------- *//* gives position of 2D point (ix,iy) in 1D array  ; uses also global variable iWidth */unsigned int Give_i (unsigned int ix, unsigned int iy){  return ix + iy * iWidth;}// unsigned char ComputeColorOfFatou (complex double z){   double r2;  int i; // number of iteration  for (i = 0; i < IterMax; ++i)    {         r2 =cabs2(z);       if (r2 > ER2) // esaping = exterior {   uExterior += 1;   return iColorOfExterior; } if (cabs2(zp2-z) < AR2) // {     return iColorOfInterior; }       z = f(z); // complex iteration f(z)=z^3 + c    }    return iColorOfUnknown;}// f(z)=1+z−3z2−3.75z3+1.5z4+2.25z5unsigned char ComputeColorOfLSM (complex double z){   double r2;  int i; // number of iteration  for (i = 0; i < IterMax_LSM; ++i)    {       // complex iteration f(z)=z^3 + c            r2 =cabs2(z);       if (r2 > ER2) // esaping = exterior {   uExterior += 1;   return 255- ((i*15) % 255); }        if (cabs2(zp2-z) < AR2)  // interior {     return iColorOfInterior; } z = f(z);     }  return iColorOfUnknown;}// ***************************************************************************************************************************// ************************** DEM/J*****************************************// ****************************************************************************************************************************unsigned char ComputeColorOfDEMJ(complex double z){  // https://en.wikibooks.org/wiki/Fractals/Iterations_in_the_complex_plane/Julia_set#DEM.2FJ    int nMax = IterMax_DEM;  complex double dz = 1.0; //  is first derivative with respect to z.  double distance;  double cabsz;   int n;  for (n=0; n < nMax; n++){ //forward iteration    cabsz = cabs(z);    if (cabsz > 1e60 || cabs(dz)> 1e60) break; // big values     if (cabs2(zp2-z) < AR2)  return iColorOfInterior; // falls into finite attractor = interior      dz = dfz(z)*dz;     z = f(z) ; /* forward iteration : complex cubic polynomial */   }      distance = 2.0 * cabsz* log(cabsz)/ cabs(dz);  if (distance <distanceMax) return iColorOfBoundary; // distanceMax = BoundaryWidth*PixelWidth;  // else    return iColorOfExterior; }/* ==================================================================================================   ============================= Draw functions ===============================================================   =====================================================================================================*/ unsigned char ComputeColor(FunctionTypeT FunctionType, complex double z){  unsigned char iColor;   switch(FunctionType){    case Fatou :{iColor = ComputeColorOfFatou(z); break;}     // case IntLSM :{iColor = ComputeColorOfIntLSM(z); break;}   // case ExtLSM :{iColor = ComputeColorOfExtLSM(z); break;}    case LSM :{iColor = ComputeColorOfLSM(z); break;}      case DEM : {iColor = ComputeColorOfDEMJ(z); break;} /*   case Unknown : {iColor = ComputeColorOfUnknown(z); break;}   case BD : {iColor = ComputeColorOfBD(z); break;}   case MBD : {iColor = ComputeColorOfMBD(z); break;}   case SAC : {iColor = ComputeColorOfSAC(z); break;}    case DLD : {iColor = ComputeColorOfDLD(z); break;}   case ND : {iColor = ComputeColorOfND(z); break;}   case NP : {iColor = ComputeColorOfNP(z); break;}   case POT : {iColor = ComputeColorOfPOT(z); break;}   case Blend : {iColor = ComputeColorOfBlend(z); break;}  */   default: {}   }   return iColor;}// plots raster point (ix,iy) int DrawPoint ( unsigned char A[], FunctionTypeT FunctionType, int ix, int iy){  int i; /* index of 1D array */  unsigned char iColor;  complex double z;  i = Give_i (ix, iy); /* compute index of 1D array from indices of 2D array */    z = GiveZ(ix,iy);    iColor = ComputeColor(FunctionType, z);  A[i] = iColor ; //     return 0;}int DrawImage ( unsigned char A[], FunctionTypeT FunctionType){  unsigned int ix, iy; // pixel coordinate   fprintf (stdout, "compute image \n");  // for all pixels of image #pragma omp parallel for schedule(dynamic) private(ix,iy) shared(A, ixMax , iyMax, uUnknown, uInterior, uExterior)  for (iy = iyMin; iy <= iyMax; ++iy)    {      fprintf (stdout, " %d from %d \r", iy, iyMax); //info       for (ix = ixMin; ix <= ixMax; ++ix) DrawPoint(A, FunctionType, ix, iy); //      }  fprintf (stdout, "\n"); //info   return 0;}int PlotPoint(complex double z, unsigned char A[]){   unsigned int ix = (creal(z)-ZxMin)/PixelWidth;  unsigned int iy = (ZyMax - cimag(z))/PixelHeight;  unsigned int i = Give_i(ix,iy); /* index of _data array */   A[i]= 0; //255-A[i]; // Mark point with inveres color   return 0; }// ***********************************************************************************************// ********************** edge detection usung Sobel filter ***************************************// ***************************************************************************************************// from Source to Destinationint ComputeBoundaries(unsigned char S[], unsigned char D[]){   unsigned int iX,iY; /* indices of 2D virtual array (image) = integer coordinate */  unsigned int i; /* index of 1D array  */  /* sobel filter */  unsigned char G, Gh, Gv;   // boundaries are in D  array ( global var )   // clear D array  memset(D, iColorOfExterior, iSize*sizeof(*D)); // for heap-allocated arrays, where N is the number of elements = FillArrayWithColor(D , iColorOfExterior);   // printf(" find boundaries in S array using  Sobel filter\n");   #pragma omp parallel for schedule(dynamic) private(i,iY,iX,Gv,Gh,G) shared(iyMax,ixMax)  for(iY=1;iY<iyMax-1;++iY){     for(iX=1;iX<ixMax-1;++iX){       Gv= S[Give_i(iX-1,iY+1)] + 2*S[Give_i(iX,iY+1)] + S[Give_i(iX-1,iY+1)] - S[Give_i(iX-1,iY-1)] - 2*S[Give_i(iX-1,iY)] - S[Give_i(iX+1,iY-1)];      Gh= S[Give_i(iX+1,iY+1)] + 2*S[Give_i(iX+1,iY)] + S[Give_i(iX-1,iY-1)] - S[Give_i(iX+1,iY-1)] - 2*S[Give_i(iX-1,iY)] - S[Give_i(iX-1,iY-1)];      G = sqrt(Gh*Gh + Gv*Gv);      i= Give_i(iX,iY); /* compute index of 1D array from indices of 2D array */      if (G==0) {D[i]=255;} /* background */      else {D[i]=0;}  /* boundary */    }  }       return 0;}// copy from Source to Destinationint CopyBoundaries(unsigned char S[],  unsigned char D[]){   unsigned int iX,iY; /* indices of 2D virtual array (image) = integer coordinate */  unsigned int i; /* index of 1D array  */    //printf("copy boundaries from S array to D array \n");  for(iY=1;iY<iyMax-1;++iY)    for(iX=1;iX<ixMax-1;++iX)      {i= Give_i(iX,iY); if (S[i]==0) D[i]=0;}     return 0;}// *******************************************************************************************// ********************************** save A array to pgm file ****************************// *********************************************************************************************int SaveArray2PGMFile (unsigned char A[],  char * n, char *comment){  FILE *fp;  const unsigned int MaxColorComponentValue = 255; /* color component is coded from 0 to 255 ;  it is 8 bit color file */  char name[100]; /* name of file */  snprintf (name, sizeof name, "%s", n ); /*  */  char *filename = strcat (name, ".pgm");  char long_comment[200];  sprintf (long_comment, "fc(z)=z^2+z*c  %s", comment);  // save image array to the pgm file   fp = fopen (filename, "wb"); // create new file,give it a name and open it in binary mode   fprintf (fp, "P5\n # %s\n %u %u\n %u\n", long_comment, iWidth, iHeight, MaxColorComponentValue); // write header to the file  size_t rSize = fwrite (A, sizeof(A[0]), iSize, fp); // write whole array with image data bytes to the file in one step   fclose (fp);  // info   if ( rSize == iSize)   {  printf ("File %s saved ", filename);  if (long_comment == NULL || strlen (long_comment) == 0)    printf ("\n");  else { printf (". Comment = %s \n", long_comment); }  }  else {printf("wrote %zu elements out of %llu requested\n", rSize,  iSize);}  return 0;}int PrintCInfo (){  printf ("gcc version: %d.%d.%d\n", __GNUC__, __GNUC_MINOR__, __GNUC_PATCHLEVEL__); // https://stackoverflow.com/questions/20389193/how-do-i-check-my-gcc-c-compiler-version-for-my-eclipse  // OpenMP version is displayed in the console : export  OMP_DISPLAY_ENV="TRUE"  printf ("__STDC__ = %d\n", __STDC__);  printf ("__STDC_VERSION__ = %ld\n", __STDC_VERSION__);  printf ("c dialect = ");  switch (__STDC_VERSION__)    { // the format YYYYMM     case 199409L:      printf ("C94\n");      break;    case 199901L:      printf ("C99\n");      break;    case 201112L:      printf ("C11\n");      break;    case 201710L:      printf ("C18\n");      break;      //default : /* Optional */    }  return 0;}intPrintProgramInfo (){  // display info messages  printf ("Numerical approximation of Julia set for F(z) =   z*c) \n");  //printf ("parameter C = ( %.16f ; %.16f ) \n", creal (C), cimag (C));    printf ("Image Width = %f in world coordinate\n", ZxMax - ZxMin);  printf ("PixelWidth = %.16f \n", PixelWidth);  //printf ("AR = %.16f = %f *PixelWidth = %f %% of ImageWidth \n", AR, AR / PixelWidth, AR / ZxMax - ZxMin);  // image corners in world coordinate  // center and radius  // center and zoom  // GradientRepetition  printf ("Maximal number of iterations = iterMax = %d \n", IterMax);  printf ("ratio of image  = %f ; it should be 1.000 ...\n", ratio);  //  return 0;}int SetPlane(complex double center, double radius, double a_ratio){  ZxMin = creal(center) - radius*a_ratio;   ZxMax = creal(center) + radius*a_ratio; //0.75;  ZyMin = cimag(center) - radius; // inv  ZyMax = cimag(center) + radius; //0.7;  return 0;}// Check Orientation of z-plane image : mark first quadrant of complex plane // it should be in the upper right position// uses global var :  ...int CheckZPlaneOrientation(unsigned char A[] ){  double Zx, Zy; //  Z= Zx+ZY*i; unsigned i; /* index of 1D array */ unsigned int ix, iy; // pixel coordinate fprintf(stderr, "compute image CheckOrientation\n");  // for all pixels of image #pragma omp parallel for schedule(dynamic) private(ix,iy, i, Zx, Zy) shared(A, ixMax , iyMax) for (iy = iyMin; iy <= iyMax; ++iy){    //fprintf (stderr, " %d from %d \r", iy, iyMax); //info     for (ix = ixMin; ix <= ixMax; ++ix){    // from screen to world coordinate     Zy = GiveZy(iy);    Zx = GiveZx(ix);   i = Give_i(ix, iy); /* compute index of 1D array from indices of 2D array */   if (Zx>0 && Zy>0) A[i]=255-A[i];   // check the orientation of Z-plane by marking first quadrant */    }    }        return 0;}// *****************************************************************************//;;;;;;;;;;;;;;;;;;;;;;  setup ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;// **************************************************************************************int setup (){  fprintf (stderr, "setup start\n");  /* 2D array ranges */  iWidth = iHeight* DisplayAspectRatio ;  iSize = iWidth * iHeight; // size = number of points in array   // iy  iyMax = iHeight - 1; // Indexes of array starts from 0 not 1 so the highest elements of an array is = array_name[size-1].  //ix  ixMax = iWidth - 1;  /* 1D array ranges */  // i1Dsize = i2Dsize; // 1D array with the same size as 2D array  iMax = iSize - 1; // Indexes of array starts from 0 not 1 so the highest elements of an array is = array_name[size-1].    SetPlane( center, radius,  DisplayAspectRatio );   /* Pixel sizes */  PixelWidth = (ZxMax - ZxMin) / ixMax; //  ixMax = (iWidth-1)  step between pixels in world coordinate   PixelHeight = (ZyMax - ZyMin) / iyMax;  ratio = ((ZxMax - ZxMin) / (ZyMax - ZyMin)) / ((double) iWidth / (double) iHeight); // it should be 1.000 ...  ER = 200.0; //   ER2 = ER*ER;  //AR_max = 5*PixelWidth*iWidth/2000.0 ; // adjust first number   // GiveTunedAR(const int i_Max, const complex double zcr, const double c, const double zp){  AR = 0.003; // GiveTunedAR(200, C, AR_max);  AR2 = AR * AR;  //AR12 = AR/2.0;        BoundaryWidth = 1.0*iWidth/2000.0  ; //  measured in pixels ( when iWidth = 2000)   distanceMax = BoundaryWidth*PixelWidth;  /* create dynamic 1D arrays for colors ( shades of gray ) */  data = malloc (iSize * sizeof (unsigned char));  edge = malloc (iSize * sizeof (unsigned char));  if (data == NULL || edge == NULL)    {      fprintf (stderr, " Could not allocate memory");      return 1;    }   fprintf (stderr, " end of setup \n");  return 0;} // ;;;;;;;;;;;;;;;;;;;;;;;;; end of the setup ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;int end (){  fprintf (stderr, " allways free memory (deallocate )  to avoid memory leaks \n"); // https://en.wikipedia.org/wiki/C_dynamic_memory_allocation  free (data);  free(edge);  PrintProgramInfo ();  PrintCInfo ();  return 0;}// ********************************************************************************************************************/* -----------------------------------------  main   -------------------------------------------------------------*/// ********************************************************************************************************************int main (){  setup ();   DrawImage (data, Fatou); // first   SaveArray2PGMFile (data,  "Fatou_1" , "Fatou ");    DrawImage (data, LSM); // first   SaveArray2PGMFile (data,  "LSM_1" , "LSM ");        DrawImage (data, DEM); // first   SaveArray2PGMFile (data,  "DEM_1" , "DEM ");      end ();  return 0;}
« Last Edit: June 19, 2021, 11:38:38 AM by Adam Majewski, Reason: b »

• Fractal Frankfurter
• Posts: 541
##### Re: rational function
« Reply #5 on: June 17, 2021, 08:20:56 PM »

• Fractal Frankfurter
• Posts: 541
##### Re: rational function
« Reply #6 on: June 17, 2021, 09:44:20 PM »
Wolfram alfa

Wolfram Alpha

solve{ 1/((0.15*I+0.15)*z^5 +z^3 + (3*I-3)*z)==z}

z≈ ± (1.0526 + 2.5111 i)
z≈ ± (0.17290 - 0.45180 i)
z≈ ± (1.5420 - 0.5826 i)

So there are 6 fixed points. Stability index = abs(f'(z))
[1.028219433424695, 1.028219433424695, 31.94911214885403,  31.94911214885403, 235.8426231893345, 235.8426231893345]

So none of fixed point is attracting.

Period 2 points ( numerically):
period 2 points :  [ 0.0, // 1+ 0.1004507524894826 - 0.04160806403172667 %i,  //2- 0.1004507524894826 + 0.04160806403172667 %i , //2//- 0.1972126334382917 %i - 0.4417257204225227, + 0.1972126334382916 %i + 0.4417257204225231, //- 0.172896861895406  + 0.4517976427751519 %i , // period 1+ 0.1728968618954057 - 0.4517976427751519 %i, // period 1//- 1.542012716044921 + 0.5826223205446455 %i , // period 1+ 1.542012716021183 - 0.5826223205577208 %i, // period 1//- 1.48680428382238  + 0.6158544989143758 %i , + 1.486804283840914 - 0.6158544989534238 %i, //- 1.545643025068914 + 0.6402263036121009 %i , + 1.545643025048946 - 0.6402263035701801 %i, //- 1.502343841878277 + 0.6783914544446792 %i , + 1.502343841903455 - 0.678391454434491 %i, //- 2.008427151852299 %i - 0.8319177653355924, + 2.008427151854269 %i + 0.83191776533551, //+ 2.511067896520314 %i + 1.052616752123668, // period 1 - 2.511067897025448 %i - 1.052616732010429, // period 1 //- 2.519905568763936 %i - 1.031280709377238, + 2.519905568548988 %i + 1.03128068923336, //  2.496853258644555 %i + 1.034230487296836, - 2.496853270435978 %i - 1.034230487368577, //- 2.532331263153173 %i - 1.048925954202523,  2.532331275662722 %i + 1.048925954304978]

Rational Julia Sets by  Mark McClure

use ( comma separated lists of coefficients in the order of ascending degree from 0 to degree of numerator/denominator ) :
1
0,  -3 +3*i,0,1, 0,0.15 + 0.15*i

or xaos user fomule : 1/((0.15*I+0.15)*z^5+(3*I-3)*z+z^3)
It seem the image is upside down in Xaos ?

« Last Edit: June 19, 2021, 03:39:56 PM by Adam Majewski, Reason: not fails »

#### pauldelbrot

• 3f
• Posts: 2994
##### Re: rational function
« Reply #7 on: June 19, 2021, 03:20:49 PM »
The morphology of the Julia shapes makes it clear that the finite attractors are 2-cycles, so it's not surprising you found no attracting finite fixed points.

#### xenodreambuie

• Fractal Friar
• Posts: 149
##### Re: rational function
« Reply #8 on: June 20, 2021, 12:40:17 AM »
I did this one in Jux some years ago, with a small change to the z^5 coefficient to 0.1591+0.1591i for more circular shapes. It appears to actually be period 4, with odd cycles and even cycles in separate clusters (not affected by the coefficient change).

#### gerson

• Fractal Freak
• Posts: 676
##### Re: rational function
« Reply #9 on: June 20, 2021, 02:42:59 AM »
@ Adam Majewski very interesting equation
@ xenodreambuie liked a lot the image

• Fractal Frankfurter
• Posts: 541
##### Re: rational function
« Reply #10 on: June 20, 2021, 08:59:33 AM »
The morphology of the Julia shapes makes it clear that the finite attractors are 2-cycles, so it's not surprising you found no attracting finite fixed points.
Probably {0, infinity} is superattractive period 2 cycle.
Are ther any other attractive cycles ?
Is this Julia set connected ?

Here is my naive numerical anlysis of periodic points ( c code below). Comments are welcome
Code: [Select]
k = 0i = 0 z = +0.000000 +0.000000   z = +inf -nan cabs(zn) > 3.0 = escaping to infinity or period 2 cycle = {z0, inf} or repelling period 2 k = 1i = 0 z = +0.100451 -0.041608  i = 1 z = -0.831918 -2.008427 cabs(z) = 2.173906 d= 2.176623i = 2 z = +0.100451 -0.041608 cabs(z) = 0.108727 d= 2.176623period 2 cycle probably attracting :  fabs(dp-dn) < 0.000001 = distance is not changing or slowly changing k = 2i = 0 z = -0.100451 +0.041608  i = 1 z = +0.831918 +2.008427 cabs(z) = 2.173906 d= 2.176623i = 2 z = -0.100451 +0.041608 cabs(z) = 0.108727 d= 2.176623period 2 cycle probably attracting :  fabs(dp-dn) < 0.000001 = distance is not changing or slowly changing k = 3i = 0 z = -0.441726 -0.197213  i = 1 z = +0.441726 +0.197213 cabs(z) = 0.483750 d= 0.967501i = 2 z = -0.441726 -0.197213 cabs(z) = 0.483750 d= 0.967501period 2 cycle probably attracting :  fabs(dp-dn) < 0.000001 = distance is not changing or slowly changing k = 4i = 0 z = +0.441726 +0.197213  i = 1 z = -0.441726 -0.197213 cabs(z) = 0.483750 d= 0.967501i = 2 z = +0.441726 +0.197213 cabs(z) = 0.483750 d= 0.967501period 2 cycle probably attracting :  fabs(dp-dn) < 0.000001 = distance is not changing or slowly changing k = 5i = 0 z = -0.172897 +0.451798  i = 1 z = -0.172897 +0.451798 cabs(z) = 0.483750 d= 0.000000i = 2 z = -0.172897 +0.451798 cabs(z) = 0.483750 d= 0.000000 z = -0.172897 +0.451798   dp < 0.000001 = probably period 1 k = 6i = 0 z = +0.172897 -0.451798  i = 1 z = +0.172897 -0.451798 cabs(z) = 0.483750 d= 0.000000i = 2 z = +0.172897 -0.451798 cabs(z) = 0.483750 d= 0.000000 z = +0.172897 -0.451798   dp < 0.000001 = probably period 1 k = 7i = 0 z = -1.542013 +0.582622  i = 1 z = -1.542013 +0.582622 cabs(z) = 1.648409 d= 0.000000i = 2 z = -1.542013 +0.582622 cabs(z) = 1.648409 d= 0.000000 z = -1.542013 +0.582623   dp < 0.000001 = probably period 1 k = 8i = 0 z = +1.542013 -0.582622  i = 1 z = +1.542013 -0.582622 cabs(z) = 1.648409 d= 0.000000i = 2 z = +1.542013 -0.582622 cabs(z) = 1.648409 d= 0.000000 z = +1.542012 -0.582623   dp < 0.000001 = probably period 1 k = 9i = 0 z = -1.486804 +0.615854  i = 1 z = +1.034230 +2.496853 cabs(z) = 2.702575 d= 3.145437i = 2 z = -1.486804 +0.615854 cabs(z) = 1.609305 d= 3.145437i = 3 z = +1.034215 +2.496862 cabs(z) = 2.702577 d= 3.145430i = 4 z = -1.487475 +0.614659 cabs(z) = 1.609468 d= 3.146682i = 5 z = +0.938842 +2.545302 cabs(z) = 2.712930 d= 3.100709i = 6 z = -0.171898 -0.241483 cabs(z) = 0.296417 d= 2.999986i = 7 z = +0.771033 -0.122995 cabs(z) = 0.780781 d= 0.950346i = 8 z = -0.177912 -0.300455 cabs(z) = 0.349179 d= 0.965395i = 9 z = +0.637916 -0.158362 cabs(z) = 0.657279 d= 0.828110unknown  or if distance is constant probably attracting period 2 cycle  or repelling period 2 cycle if distance is changing  k = 10i = 0 z = +1.486804 -0.615854  i = 1 z = -1.034230 -2.496853 cabs(z) = 2.702575 d= 3.145437i = 2 z = +1.486805 -0.615855 cabs(z) = 1.609306 d= 3.145437i = 3 z = -1.034234 -2.496871 cabs(z) = 2.702592 d= 3.145450i = 4 z = +1.488172 -0.616124 cabs(z) = 1.610673 d= 3.146385i = 5 z = -1.055143 -2.607916 cabs(z) = 2.813281 d= 3.230432i = 6 z = -0.239382 +0.188797 cabs(z) = 0.304873 d= 2.913257i = 7 z = +0.101878 +0.781313 cabs(z) = 0.787927 d= 0.683765i = 8 z = -0.203446 +0.173630 cabs(z) = 0.267466 d= 0.680075i = 9 z = +0.079584 +0.889911 cabs(z) = 0.893462 d= 0.770171unknown  or if distance is constant probably attracting period 2 cycle  or repelling period 2 cycle if distance is changing  k = 11i = 0 z = -1.545643 +0.640226  i = 1 z = -1.048926 -2.532331 cabs(z) = 2.740976 d= 3.211207i = 2 z = -1.545643 +0.640227 cabs(z) = 1.672992 d= 3.211207i = 3 z = -1.048899 -2.532346 cabs(z) = 2.740979 d= 3.211226i = 4 z = -1.544285 +0.642759 cabs(z) = 1.672709 d= 3.213518i = 5 z = -0.799562 -2.636953 cabs(z) = 2.755507 d= 3.363201i = 6 z = +0.064716 +0.091997 cabs(z) = 0.112479 d= 2.862541i = 7 z = -2.059373 +0.355826 cabs(z) = 2.089888 d= 2.140411i = 8 z = -0.085504 +0.016463 cabs(z) = 0.087074 d= 2.002830i = 9 z = +1.518575 +2.246207 cabs(z) = 2.711368 d= 2.746785unknown  or if distance is constant probably attracting period 2 cycle  or repelling period 2 cycle if distance is changing  k = 12i = 0 z = +1.545643 -0.640226  i = 1 z = +1.048926 +2.532331 cabs(z) = 2.740976 d= 3.211207i = 2 z = +1.545643 -0.640226 cabs(z) = 1.672992 d= 3.211207i = 3 z = +1.048933 +2.532362 cabs(z) = 2.741007 d= 3.211236i = 4 z = +1.542793 -0.639610 cabs(z) = 1.670123 d= 3.210188 z = +1.107056 +2.828193 cabs(zn) > 3.0 = escaping to infinity or period 2 cycle = {z0, inf} or repelling period 2 k = 13i = 0 z = -1.502344 +0.678391  i = 1 z = +1.502344 -0.678391 cabs(z) = 1.648409 d= 3.296818i = 2 z = -1.502344 +0.678391 cabs(z) = 1.648409 d= 3.296818period 2 cycle probably attracting :  fabs(dp-dn) < 0.000001 = distance is not changing or slowly changing k = 14i = 0 z = +1.502344 -0.678391  i = 1 z = -1.502344 +0.678391 cabs(z) = 1.648409 d= 3.296818i = 2 z = +1.502344 -0.678391 cabs(z) = 1.648409 d= 3.296818period 2 cycle probably attracting :  fabs(dp-dn) < 0.000001 = distance is not changing or slowly changing k = 15i = 0 z = -0.831918 -2.008427  i = 1 z = +0.100451 -0.041608 cabs(z) = 0.108727 d= 2.176623i = 2 z = -0.831918 -2.008427 cabs(z) = 2.173906 d= 2.176623period 2 cycle probably attracting :  fabs(dp-dn) < 0.000001 = distance is not changing or slowly changing k = 16i = 0 z = +0.831918 +2.008427  i = 1 z = -0.100451 +0.041608 cabs(z) = 0.108727 d= 2.176623i = 2 z = +0.831918 +2.008427 cabs(z) = 2.173906 d= 2.176623period 2 cycle probably attracting :  fabs(dp-dn) < 0.000001 = distance is not changing or slowly changing k = 17i = 0 z = +1.052617 +2.511068  i = 1 z = +1.052616 +2.511064 cabs(z) = 2.722764 d= 0.000004i = 2 z = +1.051601 +2.511174 cabs(z) = 2.722473 d= 0.001022i = 3 z = +1.082008 +2.772741 cabs(z) = 2.976379 d= 0.263328i = 4 z = +0.070973 -0.059033 cabs(z) = 0.092315 d= 3.006848i = 5 z = -0.236804 -2.546474 cabs(z) = 2.557461 d= 2.506409i = 6 z = +0.049184 -0.002415 cabs(z) = 0.049243 d= 2.560084 z = -3.214658 -3.549163 cabs(zn) > 3.0 = escaping to infinity or period 2 cycle = {z0, inf} or repelling period 2 k = 18i = 0 z = -1.052617 -2.511068  i = 1 z = -1.052617 -2.511068 cabs(z) = 2.722768 d= 0.000000i = 2 z = -1.052715 -2.511071 cabs(z) = 2.722809 d= 0.000098i = 3 z = -1.052554 -2.488103 cabs(z) = 2.701579 d= 0.022969i = 4 z = +0.462528 -0.999277 cabs(z) = 1.101130 d= 2.124165i = 5 z = +0.014954 -0.195738 cabs(z) = 0.196309 d= 0.919781i = 6 z = +0.771090 -0.911883 cabs(z) = 1.194198 d= 1.041444i = 7 z = -0.070214 -0.203606 cabs(z) = 0.215373 d= 1.099749i = 8 z = +0.972784 -0.475440 cabs(z) = 1.082752 d= 1.077839i = 9 z = -0.130633 -0.302244 cabs(z) = 0.329266 d= 1.116926unknown  or if distance is constant probably attracting period 2 cycle  or repelling period 2 cycle if distance is changing  k = 19i = 0 z = -1.031281 -2.519906  i = 1 z = +1.031281 +2.519906 cabs(z) = 2.722768 d= 5.445536i = 2 z = -1.031188 -2.519888 cabs(z) = 2.722717 d= 5.445485i = 3 z = +1.035915 +2.498359 cabs(z) = 2.704610 d= 5.427312i = 4 z = -1.610454 +0.770400 cabs(z) = 1.785239 d= 3.160555i = 5 z = +0.134719 -0.407770 cabs(z) = 0.429448 d= 2.105639i = 6 z = +0.223892 -0.495799 cabs(z) = 0.544008 d= 0.125304i = 7 z = +0.124084 -0.414032 cabs(z) = 0.432226 d= 0.129026i = 8 z = +0.235618 -0.484984 cabs(z) = 0.539190 d= 0.132190i = 9 z = +0.114235 -0.422698 cabs(z) = 0.437862 d= 0.136431unknown  or if distance is constant probably attracting period 2 cycle  or repelling period 2 cycle if distance is changing  k = 20i = 0 z = +1.031281 +2.519906  i = 1 z = -1.031281 -2.519901 cabs(z) = 2.722764 d= 5.445531i = 2 z = +1.032305 +2.519975 cabs(z) = 2.723221 d= 5.445984i = 3 z = -0.993050 -2.779869 cabs(z) = 2.951919 d= 5.673660i = 4 z = -0.040063 +0.082568 cabs(z) = 0.091774 d= 3.016908i = 5 z = -0.836507 +2.428865 cabs(z) = 2.568877 d= 2.477788i = 6 z = -0.004093 +0.026235 cabs(z) = 0.026553 d= 2.542743 z = -5.232646 +7.169638 cabs(zn) > 3.0 = escaping to infinity or period 2 cycle = {z0, inf} or repelling period 2 k = 21i = 0 z = +1.034230 +2.496853  i = 1 z = -1.486803 +0.615854 cabs(z) = 1.609305 d= 3.145436i = 2 z = +1.034230 +2.496787 cabs(z) = 2.702513 d= 3.145396i = 3 z = -1.481624 +0.615777 cabs(z) = 1.604491 d= 3.141293i = 4 z = +1.011372 +2.146284 cabs(z) = 2.372637 d= 2.925317i = 5 z = -0.106695 +0.068182 cabs(z) = 0.126620 d= 2.359785i = 6 z = +0.404193 +1.823882 cabs(z) = 1.868133 d= 1.828522i = 7 z = -0.085208 +0.044710 cabs(z) = 0.096226 d= 1.845256i = 8 z = +0.731458 +2.343232 cabs(z) = 2.454744 d= 2.439293i = 9 z = -0.120305 -0.002622 cabs(z) = 0.120334 d= 2.495703unknown  or if distance is constant probably attracting period 2 cycle  or repelling period 2 cycle if distance is changing  k = 22i = 0 z = -1.034230 -2.496853  i = 1 z = +1.486804 -0.615854 cabs(z) = 1.609305 d= 3.145437i = 2 z = -1.034230 -2.496860 cabs(z) = 2.702580 d= 3.145440i = 3 z = +1.487300 -0.615818 cabs(z) = 1.609749 d= 3.145860i = 4 z = -1.031051 -2.535870 cabs(z) = 2.737463 d= 3.166811i = 5 z = -0.369923 +1.279341 cabs(z) = 1.331749 d= 3.872070i = 6 z = -0.032201 +0.134449 cabs(z) = 0.138252 d= 1.193664i = 7 z = -0.884125 +1.454884 cabs(z) = 1.702458 d= 1.571408i = 8 z = +0.039601 +0.087176 cabs(z) = 0.095749 d= 1.650423i = 9 z = -2.299897 +0.862691 cabs(z) = 2.456372 d= 2.464685unknown  or if distance is constant probably attracting period 2 cycle  or repelling period 2 cycle if distance is changing  k = 23i = 0 z = -1.048926 -2.532331  i = 1 z = -1.545643 +0.640226 cabs(z) = 1.672992 d= 3.211207i = 2 z = -1.048927 -2.532322 cabs(z) = 2.740968 d= 3.211198i = 3 z = -1.546480 +0.640104 cabs(z) = 1.673719 d= 3.211207i = 4 z = -1.058996 -2.455627 cabs(z) = 2.674243 d= 3.133879i = 5 z = +0.278573 -0.455770 cabs(z) = 0.534162 d= 2.405934i = 6 z = +0.075496 -0.441610 cabs(z) = 0.448017 d= 0.203570i = 7 z = +0.280187 -0.431527 cabs(z) = 0.514510 d= 0.204939i = 8 z = +0.068681 -0.461364 cabs(z) = 0.466448 d= 0.213600i = 9 z = +0.276134 -0.407608 cabs(z) = 0.492335 d= 0.214304unknown  or if distance is constant probably attracting period 2 cycle  or repelling period 2 cycle if distance is changing  k = 24i = 0 z = +1.048926 +2.532331  i = 1 z = +1.545642 -0.640226 cabs(z) = 1.672991 d= 3.211207i = 2 z = +1.048928 +2.532431 cabs(z) = 2.741068 d= 3.211305i = 3 z = +1.536481 -0.640006 cabs(z) = 1.664446 d= 3.209682 z = +0.900158 +3.704194 cabs(zn) > 3.0 = escaping to infinity or period 2 cycle = {z0, inf} or repelling period 2
here is c src code
Code: [Select]
#include <stdio.h>#include <stdlib.h> // malloc#include <math.h> // M_PI; needs -lm also#include <complex.h>/* a: (0.15+0.15**I);b: (-3+3**I);define(f(z), 1/(a*z^5+ z^3 + b*z));period 2 points : */double complex zp2[25] = { 0.0, // period 2 = found by analysis of the function//+ 0.1004507524894826 - 0.04160806403172667 *I,  //2-0.1004507524894826 + 0.04160806403172667 *I , //2//- 0.1972126334382917 *I - 0.4417257204225227, //2+ 0.1972126334382916 *I + 0.4417257204225231, //- 0.172896861895406  + 0.4517976427751519 *I , // period 1+ 0.1728968618954057 - 0.4517976427751519 *I, // period 1//- 1.542012716044921 + 0.5826223205446455 *I , // period 1+ 1.542012716021183 - 0.5826223205577208 *I, // period 1//- 1.48680428382238  + 0.6158544989143758 *I , // probably repelling period 2 + 1.486804283840914 - 0.6158544989534238 *I, // probably repelling period 2 //- 1.545643025068914 + 0.6402263036121009 *I ,  // probably repelling period 2 + 1.545643025048946 - 0.6402263035701801 *I,  // probably repelling period 2 //- 1.502343841878277 + 0.6783914544446792 *I , // period 2+ 1.502343841903455 - 0.678391454434491 *I, //- 2.008427151852299 *I - 0.8319177653355924, //period 2+ 2.008427151854269 *I + 0.83191776533551, //+ 2.511067896520314 *I + 1.052616752123668, // period 1 - 2.511067897025448 *I - 1.052616732010429, // period 1 //- 2.519905568763936 *I - 1.031280709377238, + 2.519905568548988 *I + 1.03128068923336, //  2.496853258644555 *I + 1.034230487296836, // ??- 2.496853270435978 *I - 1.034230487368577, //- 2.532331263153173 *I - 1.048925954202523,  2.532331275662722 *I + 1.048925954304978};complex double a = 0.15+0.15*I;complex double b = -3.0+3.0*I;complex double f(complex double z){ complex double z3 = z*z*z; complex double z5 = z3*z*z; return  1.0/(a*z5+ z3 + b*z); } void TestZ(complex double z0){ complex double zm = z0; complex double zn; double dp; double dn; int i=0; printf("i = %d \t z = %+f %+f  \n",  i, creal(zm), cimag(zm) ); for (i = 1; i < 10; ++i)  {        zn = f(zm);    if (cabs(zn) > 3.0){printf(" z = %+f %+f \t cabs(zn) > 3.0 = escaping to infinity or period 2 cycle = {z0, inf} or repelling period 2 \n\n", creal (zn), cimag(zn)); return;}    dn = cabs(zm - zn);        if (i>2) {        if ( dn > 6.0 ){ printf(" dp>6.0\n"); break;}    if ( dn < 0.000001 ) { printf(" z = %+f %+f \t  dp < 0.000001 = probably period 1 \n\n", creal (zn), cimag(zn)); return;}    if ( fabs(dp-dn) < 0.000001 ) { printf("period 2 cycle probably attracting :  fabs(dp-dn) < 0.000001 = distance is not changing or slowly changing \n\n"); return;}    }         //    printf("i = %d \t z = %+f %+f \t cabs(z) = %f \t d= %f\n", i, creal(zn), cimag(zn), cabs(zn), dn );    zm= zn;    dp = dn;      } printf("unknown \n or if distance is constant probably attracting period 2 cycle \n or repelling period 2 cycle if distance is changing  \n\n");}int main (void){ int k; int kMax = sizeof(zp2)/sizeof(zp2[0]); for (k = 0; k < kMax; ++k) {printf("k = %d\n\n",k); TestZ(zp2[k]);} return 0;}
« Last Edit: June 20, 2021, 04:06:17 PM by Adam Majewski, Reason: code »

• Fractal Frankfurter
• Posts: 541
##### Re: rational function
« Reply #11 on: June 20, 2021, 09:13:35 PM »
the best image

#### xenodreambuie

• Fractal Friar
• Posts: 149
##### Re: rational function
« Reply #12 on: June 20, 2021, 11:35:17 PM »
Jux detects this as two basins each with period 4. Here they are coloured with warm colours in one and cool in the other. That seems to be at odds with the analysis of period 2 points by everyone else. Any ideas as to what might cause this disparity?

#### Chris_M_Thomasson

• Fractal Feline
• Posts: 154
• Programmer by trade; love fractals and math! :^)
##### Re: rational function
« Reply #13 on: June 21, 2021, 01:10:30 AM »
Wolfram alfa

Wolfram Alpha

[...]

or xaos user fomule : 1/((0.15*I+0.15)*z^5+(3*I-3)*z+z^3)
It seem the image is upside down in Xaos ?

Back when I used Xaos, it always erroneously inverted the y axis where +y pointed down, and -y pointed up. This is wrong. Btw, I found some interesting bugs in Xaos:

https://element90.wordpress.com/2013/10/09/another-unexpected-difference/
It's a Fractal Life, 247... ;^)

• Fractal Frankfurter
• Posts: 541
##### Re: rational function
« Reply #14 on: June 21, 2021, 02:50:50 PM »
Critical points are centers of biggest components ( gray dots on a_Fatou_cr.png)

Periodic points ( without zere and infinity) are on the boundary of componnets ( gray black  dots on a_Fatou-zp2.png are hard to seee, sorry) so this is parabolic case.

Zero and infinity is a period 2 cycle in the exterior of Julia set, so here there is no  gap in the middle. See that components also end on 2 periodic points near zero.

« Last Edit: June 22, 2021, 03:08:15 PM by Adam Majewski, Reason: spell »

### Similar Topics

###### Definitions of degree of rational function
3 Replies
339 Views
December 06, 2020, 07:30:43 AM
###### How to convert a function to a recursive function

Started by Deliberate Dendrite on Noob's Corner

2 Replies
481 Views
December 13, 2018, 05:32:30 PM
by Deliberate Dendrite
###### Rounding function

Started by mclarekin on Color Snippets

0 Replies
241 Views
March 26, 2019, 01:46:13 AM
by mclarekin
###### Lambert W function

Started by superheal on Share a fractal

3 Replies
808 Views
December 06, 2018, 05:38:25 AM
by FractalDave
###### Programming a recursive function

Started by Fraktalist on Fractal Humor

0 Replies
1105 Views
March 21, 2018, 11:48:01 AM
by Fraktalist