Welcome to Fractal Forums

// from JIT Mandelcup of Pupukuusikko :D
#include <math.h>

void __attribute__((fastcall)) Formula(
             double* x,      // [eax]
             double* y,      // [edx]
             double* arg,    // [ebp+8], points to TIteration3Dext.C1
             void* dummy     // so we end w/ ret 8 as delphi expects
             );

           
__attribute__((packed)) struct TIteration3Dext {
double Cw,Rold,RStopD,x,y,z,w; // use with neg indizes before C1. w is also used for 3d ifs and abox analytic DE
double Px, Py, Pz;             //      actual position, never change these!  Can be used as input.
double Cx, Cy, Cz;             //+24   these are the constants for adding. Pxyz or the julia seed.  Cw is @-56 (begin of struct)
void*  PVar;                   //+48   the actual formulas adress of constants and vars, constants at offset 0 increasing, user vars below -8 
float  SmoothItD;              //+52
double Rout;                   //+56   the square of the vector length, can be used as input
int    ItResultI;              //+64
int    maxIt;                  //+68
float  RStop;                  //+72
int    nHybrid[6];             //+76   all formulas iteration counts or single weights in interpol-hybrid
void*  fHPVar[6];              //+100  all formulas pointer to constants+vars, PVars-8=0.5, use PVar for the actual formula
void*  fHybrid[6];             //+124  the formulas adresses
int    CalcSIT;                //+148
int    DoJulia;                //+152
float  LNRStop;                //+156
int    DEoption;               //+160
float  fHln[6];                //+164  for SmoothIts
int    iRepeatFrom;            //+188
double OTrap;                  //+192
double VaryScale;              //+200  to use in vary by its
int    bFirstIt;               //+208  used also as iteration count, is set to 0 on it-start
int    bTmp;                   //+212  tmpBuf, free of use.
double Dfree1;                 //+216
double Dfree2;                 //+224
double Deriv1;                 //+232  for 4D first deriv or as full derivs
double Deriv2;                 //+240
double Deriv3;                 //+248
float  SMatrix4[4][4];         //+256  for 4d rotation, used like most other values only by the programs iteration loop procedure
};

__CRT_INLINE double __cdecl Hypot(double x,double y) { return sqrt(x*x+y*y); }
__CRT_INLINE double __cdecl length(double x,double y,double z) { return sqrt(x*x+y*y+z*z); }
__CRT_INLINE double __cdecl recip (long double _x)
{
  double res;
  __asm__ ("fld1
	"
	    "fdivp": "=t" (res) : "0" (_x));
  return res;
}

// fastcall is not quite delphi fastcall.
// first two args are ok, third is in ecx in delphi, on stack here.
void __attribute__((fastcall)) Formula(
             double* x,      // [eax]
             double* y,      // [edx]
             double* arg,    // [ebp+8], points to TIteration3Dext.C1
             void* dummy     // so we end w/ ret 8 as delphi expects
             ) {
  // Compute ptr to proper start of TIteration3Dext struct.
  struct TIteration3Dext* cfg = (struct TIteration3Dext*)(arg-7);
  int* UserV = (int*)(cfg->PVar-12); // beginning of int vars
  double* UserV2 = (double*)(cfg->PVar-12-16); // beginning of dbl vars
  // no cns used!
  double temp, z_y, z_x, r;
  double mag1, mag2;
  double Scale = UserV2[0];
  double Shape = UserV2[-1]; // AAARRGGHH neg indices but it workssss :P
  
  // stereographic projection, modified a bit
  r = length(cfg->x,cfg->y,cfg->z);
  if (UserV[0]!=0) {
     cfg->x =cfg->x/r; cfg->y =cfg->y/r;
  }
  cfg->z = (cfg->z/r)+Shape;
  
  if (UserV[-1]!=0) { 
     cfg->z =cfg->z*cfg->z; Shape = Shape*Shape;
  }
  
  if (UserV[-2]!=0) { 
    temp = cfg->z;
  } else {
    temp = recip(cfg->z);
  }
  
  z_x = cfg->x*temp; z_y = cfg->y*temp;
 
  //complex multiplication
  temp = z_x*z_x-z_y*z_y;
  z_y = z_x*z_y; z_y += z_y;
  z_x = temp;
  
  temp = Scale*(1.+(Shape*Shape));
  //z_x = z_x*temp; 
  //z_y = z_y*temp;

  
  //p = vec3(2.*z.x,2.*z.y,dot(z,z)-1)/dot(z.z+1);
  
  // inverse stereographic 
  mag1 = z_x*z_x+z_y*z_y;
  mag2 = recip(mag1+1.); mag1 = mag1-1.;
              
  cfg->x = temp*mag2;
  cfg->y = z_y*cfg->x;
  cfg->x *= z_x;
  cfg->z = mag1*mag2;
  
  temp = r*r;//+Offset3;
  
  cfg->x = cfg->x * temp + cfg->Cx;
  cfg->y = cfg->y * temp + cfg->Cy;
  cfg->z = cfg->z * temp + cfg->Cz;
}