I played a bit with the original quaternion adding iteration rotation, box fold and sphere fold.
The question I have is regarding the accuracy of DE formula in this case.
If the 'original' quaternion looked as follows:
float DE(vec3 pos) {
vec4 p = vec4(pos, 0.0);
vec4 dp = vec4(1.0, 0.0,0.0,0.0);
for (int i = 0; i < Iterations; i++) {
dp = 2.0* vec4(p.x*dp.x-dot(p.yzw, dp.yzw), p.x*dp.yzw+dp.x*p.yzw+cross(p.yzw, dp.yzw));
p = vec4(p.x*p.x-dot(p.yzw, p.yzw), vec3(2.0*p.x*p.yzw)) + C;
float p2 = dot(p,p);
orbitTrap = min(orbitTrap, abs(vec4(p.xyz,p2)));
if (p2 > Threshold) break;
}
float r = length(p);
return 0.5 * r * log(r) / length(dp);
}
And roughly we wanted to tweak as follows:
float DE(vec3 pos) {
vec4 p = vec4(pos, 0.0);
vec4 dp = vec4(1.0, 0.0,0.0,0.0);
for (int i = 0; i < Iterations; i++) {
dp = 2.0* vec4(p.x*dp.x-dot(p.yzw, dp.yzw), p.x*dp.yzw+dp.x*p.yzw+cross(p.yzw, dp.yzw));
p = vec4(p.x*p.x-dot(p.yzw, p.yzw), vec3(2.0*p.x*p.yzw)) + C;
====> roatate(p);
====> boxfold(p);
====> spherefold(p);
float p2 = dot(p,p);
orbitTrap = min(orbitTrap, abs(vec4(p.xyz,p2)));
if (p2 > Threshold) break;
}
float r = length(p);
return 0.5 * r * log(r) / length(dp);
}
So we get something like
float DEQT(vec3 pos) {
vec4 p = vec4(pos, 0.0);
vec4 dp = vec4(1.0, 0.0,0.0,0.0);
for (int i = 0; i <QTIterations; i++){
// Quaternion iterations
dp = 2.0* vec4(p.x*dp.x-dot(p.yzw, dp.yzw), p.x*dp.yzw+dp.x*p.yzw+cross(p.yzw, dp.yzw));
p = vec4(p.x*p.x-dot(p.yzw, p.yzw), vec3(2.0*p.x*p.yzw)) + QTC;
//rotation
p.xyz *= rot;
// Box Folding
if (foldBoxIts > 0){
if ( i >= foldBoxStartIt && i < foldBoxStartIt + foldBoxIts){
// Box Folding (v2) - separate axes + offset
p.x = clamp(p.x, -foldingLimitX, foldingLimitX) * foldingValueX - p.x + Trans.x;
p.y = clamp(p.y, -foldingLimitY, foldingLimitY) * foldingValueY - p.y + Trans.y;
p.z = clamp(p.z, -foldingLimitZ, foldingLimitZ) * foldingValueZ - p.z + Trans.z;
}
}
// Sphere Folding
if(foldSphereIts > 0){
// Sphere folding & scale
if ( i >= foldSphereStartIt && i < foldSphereStartIt + foldSphereIts){
float r2 = dot(p.xyz, p.xyz);
p *= clamp(max(MinRad2/r2, MinRad2), 0.0, 1.0);
p = p*scale + vec4(Julia,0);
}
}
float p2 = dot(p,p);
orbitTrap = min(orbitTrap, abs(vec4(p.xyz,p2)));
if (p2 > QTThreshold) break;
}
//QT
float r = length(p);
return 0.5 * r * log(r) / length(dp);
}
So we got p, we got dp then we introduce a few other transformations on p while leaving dp untouched for the final DE calculation.
Is the DE formula still valid if we do so? I have a feeling this is not exactly accurate and that might be a reason for much noise I am getting out of these images.
My math skills are really below my aspirations, so wrapping my head around all this gives me a hard time.
Just started thinking about it after having read syntopia's blogposts on DE fractals.