Welcome to Fractal Forums

Ok, here's a first cut.  tabs=4.  Even though there may be some 
rudimentary C++, this is mainly C.  This code is cobbled together from 
various files.  Originally compiled for Win16, never ported to Win32.


I had a custom program that would send 'messages' to methods like this.  
FN_RUN is the one you're interested in.  To make inverse iteration work, 
and to avoid most of the ISS-like set being sparse, you have to prune 
iterations with a cumulative derivative that is too small.  That's what the 
jtree stuff is all about.


I'm not sure whether anything of this is useful.  Your best bet would be to 
get inverse iteration working to your satisfaction for the classic quadratic 
Julia, and then work from there.  A Glynn set is just a Julia set for an 
exponent other than 2.  Usually the exponent is fractional and between 1 and 2. 


parmA,B,C,.. correspond to the numbers the pictures.


Good luck,


Peter


#define GLOBAL extern


GLOBAL HWND hWndDynamics ;
GLOBAL HWND hWndParticles ;
GLOBAL HWND hWndMain ; //cheat


GLOBAL double parmA ;
GLOBAL double parmB ;
GLOBAL double parmC ;
GLOBAL double parmD ;
GLOBAL double parmE ;
GLOBAL double parmF ;
GLOBAL double parmG ;
GLOBAL double parmH ;
GLOBAL double temperature ;
GLOBAL DWORD lifetime ;
GLOBAL DWORD joltPeriod ;
//GLOBAL double frequency ;
GLOBAL double x ;
GLOBAL double y ;
GLOBAL BOOL bGo ;
GLOBAL BOOL bBounded ;
const int MAX_DEPTH=1000 ;


class CTree {
public:
double x[MAX_DEPTH] ;
double y[MAX_DEPTH] ;
double dp[MAX_DEPTH] ;
int choice[MAX_DEPTH] ;
int index ;
private:
} ; 


void glynn (int mode)
{
switch (mode){
case FN_RUN:
{
if (parmC==0.0) return ;
int branches = (int)ceil(fabs(parmC)) ;
// if cumulative derivative is larger than threshold, backtrack
if (jtree.dp[jtree.index] > cutoff || jtree.index > MAX_DEPTH-2){
for ( ; ; ){
jtree.choice[jtree.index] = 0 ;
jtree.index-- ;
if (jtree.index==0){
jtree.choice[0] = (jtree.choice[0]+1)%branches ;
break ;
}
if (jtree.choice[jtree.index] < branches-1){
jtree.choice[jtree.index]++ ;
break ;
}
}
}
double xt = jtree.x[jtree.index] - parmA ;
double yt = jtree.y[jtree.index] - parmB ;
double r = pow (::sqrt (xt*xt+yt*yt), 1.0/parmC) ;
double theta  = atan2 (yt,xt) ;
BOOL bNoRoot = FALSE ;
int branch = jtree.choice[jtree.index] ;
double limit = fabs(parmC)*PI ;
theta = theta + TWOPI*(branch - floor((theta+limit)/TWOPI)) ;
if (theta > -limit && theta < limit){
theta = theta/parmC ;
} else {
bNoRoot = TRUE ;
}


if (bNoRoot){
// don't need to assign x, y - just use existing values
jtree.index++ ;
assert(jtree.index<MAX_DEPTH) ;
jtree.dp[jtree.index] = 2.0*cutoff ; 
} else {
x = r*cos(theta) ;
y = r*sin(theta) ;
jtree.index++ ;
assert(jtree.index<MAX_DEPTH) ;
jtree.x[jtree.index] = x ;
jtree.y[jtree.index] = y ;
jtree.dp[jtree.index] = 
  fabs(parmC)*pow(::sqrt(x*x+y*y),parmC-1.0)*jtree.dp[jtree.index-1];
}
}
break ;
case FN_RESET:
{
//jtree.x[0] = -1 ; jtree.y[0] = y ; jtree.index = 0 ; jtree.dp[0]=1.0;
jtree.x[0] = -1 ; jtree.y[0] = 0 ; jtree.index = 0 ; jtree.dp[0]=1.0;
for (int i = 0; i<MAX_DEPTH; i++){jtree.choice[i]=0 ; }
}
break ;
case FN_UPDATEPARM:
  cutoff = pow (10.0, max(parmD, 0.0)) ;
dynamics (FN_RESET) ;
break ; 
case FN_INIT:
joltPeriod = lifetime = 100 ; temperature = .0 ;
parmA = -0.74543 ; parmB = 0.11301 ; parmC = 2.0 ; parmD=10.0 ;
dynamics (FN_UPDATEPARM) ;
x = y = 0.0 ;
break ;
}
}

z= (R^(@pow/2))*sin(@pow*ph)*exp(i*@pow*th)+c1
   zz=(R^(@pow/2))*cos(@pow*ph)+cz
     R=(|z|+zz*zz)
     th=atan2(z)
     ph=atan2(zz+i*cabs(z))
     if ph>=#pi/2,ph=#pi-ph,endif