Since you already keep in mind signs of a, c in my opinion f should just be 1 (
you multiply by b then divide by a modulus -> equals to sgn(b) if c=0 etc. ... )
, in my case 0.25 because I moves out some constants, so this should be the correct;
(sorry for the strange stuff all around but it just works into assembly!
)
double x,y,z;
//vec3 lz=z+vec3(1.), llz=z+vec3(-1.);
double DF = load1();
x=loadX(); y=loadY(); z=loadZ();
double lzx = x+DF; double lzy = y+DF; double lzz = z+DF;
double llzx = x-DF; double llzy = y-DF; double llzz = z-DF;
double DE = loadNthDoubleCns(0,0);
double d2;
double oC = DE;
double a, b, AexJ;
//double bdiva = (fabs(b)>loadNthDoubleCns(1,0))?b/a:load1();
double bdiva, distJosfix;
double f; // f * 0.25
int Box_Iterations = loadNthIntVar(0,6*szof1dbl);
int Final_Iterations = loadNthIntVar(1,6*szof1dbl);
int MiscOps;
int SpI = (loadNthIntVar(2,6*szof1dbl))&4;
double bsx = loadNthDoubleVar(2,0);
double bsz = loadNthDoubleVar(3,0);
double bsx2 = bsx*loadNthDoubleCns(5,0);
double bsz2 = bsz*loadNthDoubleCns(5,0);
int i;
if (Box_Iterations>100) Box_Iterations = 100;
if (Final_Iterations>20) Final_Iterations = 20;
//Box_Iterations = 30; Final_Iterations = 4;
a = loadNthDoubleVar(0,0); b = loadNthDoubleVar(1,0);
AexJ = loadNthDoubleVar(10,3*szof1int);
bdiva = (b/a);
WORD oldcw=fegr_floor(); // wrap() works ok
// f = sig0p(b)*sig0p(AexJ)*loadNthDoubleCns(3,0); // f * 0.25
f = loadNthDoubleCns(3,0); // = 0.25
distJosfix = sqrt(b*b+AexJ*AexJ);
if (distJosfix>loadNthDoubleCns(1,0)) distJosfix = recip(distJosfix);
// vars:
// 0 = dbl a; 1 = dbl b; 2 = dbl box_size_x; 3 = dbl box_size_y;
// 4 = dbl clampY; 5 = dbl clampDF;
// 6+off(3int) = iCx 7+off(3int) = iCx 8+off(3int) = iCx
// 0+off(6dbl) = int Box_Iterations;
// 1+off(6dbl) = int Final_Iterations;
MiscOps = (loadNthIntVar(2,6*szof1dbl))&4;
if (MiscOps) {
//z=z-InvCenter;
x = x - loadNthDoubleVar(6,3*szof1int);
y = y - loadNthDoubleVar(7,3*szof1int);
z = z - loadNthDoubleVar(8,3*szof1int);
double iR = loadNthDoubleVar(9,3*szof1int); iR *= iR;
d2 = length2(x,y,z);
iR = iR/d2;
x = x*iR + loadNthDoubleVar(6,3*szof1int);
y = y*iR + loadNthDoubleVar(7,3*szof1int);
z = z*iR + loadNthDoubleVar(8,3*szof1int);
}
MiscOps = (loadNthIntVar(2,6*szof1dbl))&1;
for (i = 0; i < Box_Iterations ; i++)
{
//if(z.y<0. || z.y>a) break;
x += bdiva*y;
/*
if (FourGen)
z = wrap(z, vec3(2. * box_size_x, a, 2. * box_size_z), vec3(- box_size_x, 0., - box_size_z));
else
z.xz = wrap(z.xz, vec2(2. * box_size_x, 2. * box_size_z), vec2(- box_size_x, - box_size_z));
*/
if (MiscOps) {
wrap3D(&x,&y,&z,
bsx2,a, bsz2,
-bsx,loadzero(),-bsz);
} else {
wrap2D(&x,&z,
bsx2,bsz2,
-bsx,-bsz);
}
x -= bdiva*y;
//If above the separation line, rotate by 180° about (-b/2, a/2)
// if (z.y >= a * (0.5 + f * 0.25 * sign(z.x + b * 0.5)* (1. - exp( - 3.2 * abs(z.x + b * 0.5)))))
{
// Jos' rotated separation line for Aexion mod! THANKS JOS!
double x00 = (b * (x + b * loadNthDoubleCns(2,0))
-AexJ*(z + AexJ * loadNthDoubleCns(2,0))) * distJosfix;
// old separation line:
// double x00 = x + b * loadNthDoubleCns(2,0);
x00 = loadNthDoubleCns(2,0) + f * sign(x00) *
(load1() - exp( loadNthDoubleCns(4,0) * fabs(x00) ) ) ;
if (y >= a *x00 )
{
//z = vec3(-b, a, 0.) - z;
// x = -b - x; y = a - y; z = -z ;original
x = -b - x; y = a - y; z = AexJ-z; // Pupukuusikko mod THANKS
//z.xy = vec2(-b, a) - z.xy; // commented out in Jos orig src
}
}
//Apply transformation a (Orig)
// TransA(&x,&y,&z, &DF, a, b);
//Apply transformation a + add Aex to z
TransA_Aex(&x,&y,&z, &DF, a, b, AexJ);
/*
//For colouring
double oCn = fabs(fabs(y)-fabs(x)-fabs(z));
oC = min(oC, oCn);
*/
//If the iterated points enters a 2-cycle , bail out.
if ( length2(x-llzx,y-llzy,z-llzz) < loadNthDoubleCns(1,0) )
{
break;
}
//Store prévious iterates
llzx = lzx; llzy = lzy; llzz = lzz;
lzx = x; lzy = y; lzz = z;
}
bsx = loadNthDoubleVar(11,3*szof1int);
saveOTColor( (loadNthDoubleCns(6,0)-bsx)*lzy / ( bsx*sqrt((double)i)+loadNthDoubleCns(6,0)) );
bsx = loadNthDoubleVar(4,0); // clampy
bsz = loadNthDoubleVar(5,0); // clampdf
MiscOps = (loadNthIntVar(2,6*szof1dbl))&2;
MiscOps = (MiscOps != 0);
//WIP: Push the iterated point left or right depending on the sign of KleinI
i = 0;
// sorry for using GOTO but it just does its job perfectly here
do {
if (MiscOps) donothing();
else y = min(y, a-y); // Show balls (Strongly reccomended! Really! :D )
DE=min(DE,min(y,bsx)/max(DF,bsz));
if (i>=Final_Iterations) break;
//Apply inv-transformation a (Orig)
// TransAInv(&x,&y,&z, &DF, a, b);
//Apply inv-transformation a + add Aex to z
TransAInv_Aex(&x,&y,&z, &DF, a, b, AexJ);
i++;
} while (1);
MiscOps = (loadNthIntVar(2,6*szof1dbl))&4;
if (MiscOps) { // fixde finale sphereinv
DE=DE*d2/(loadNthDoubleVar(9,3*szof1int)+sqrt(d2)*DE);
}
fesetround(oldcw); // do this right before returning
saveDE (DE);
/*
cns; loadNthDoubleCns(n,0);
0 = 1e10
1 = 1e-5
2 = 0.5
3 = 0.25
4 = -3.2
5 = 2.
---------
*/
To Aexion post codeeeeeeeeeeee