So far I have a code to handle 16-component numbers (sedenions).

I'm making a fractal marcher, using C++. Is there any way to avoid 16 nested for loops, or is there some C++ magic?

#include <cmath>
#include <iostream>
#include <vector>
#include <complex>
using namespace std;


class sedenion
{
public:

	sedenion(void)
	{
		vertex_data.resize(vertex_length, 0);
	}

	size_t vertex_length = 16;
	vector<float> vertex_data;
};


sedenion pow_sedenion(sedenion &in, float exponent)
{
	const float beta = exponent;

	const float fabs_beta = fabsf(beta);

	float all_self_dot = 0;

	for (size_t i = 0; i < in.vertex_length; i++)
		all_self_dot += (in.vertex_data[i] * in.vertex_data[i]);

	float imag_self_dot = 0;

	for (size_t i = 1; i < in.vertex_length; i++)
		imag_self_dot += (in.vertex_data[i] * in.vertex_data[i]);

	sedenion out;

	if (all_self_dot == 0)
	{
		for (size_t i = 0; i < out.vertex_length; i++)
			out.vertex_data[i] = 0;

		return out;
	}

	const float all_len = sqrtf(all_self_dot);
	const float imag_len = sqrtf(imag_self_dot);
	const float self_dot_beta = powf(all_self_dot, fabs_beta / 2.0f);

	out.vertex_data[0] = self_dot_beta * std::cos(fabs_beta * std::acos(in.vertex_data[0] / all_len));

	for (size_t i = 1; i < out.vertex_length; i++)
		out.vertex_data[i] = in.vertex_data[i] * self_dot_beta * sin(fabs_beta * acos(in.vertex_data[0] / all_len)) / imag_len;

	return out;

	//if (beta < 0)
	//	inverse(&out, 0, &out);
}


complex<float> pow_complex(const complex<float>& in, const float beta)
{
	float fabs_beta = fabsf(beta);

	float self_dot = in.real() * in.real() + in.imag() * in.imag();

	if (self_dot == 0)
		return complex<float>(0, 0);

	float len = std::sqrtf(self_dot);
	float self_dot_beta = std::powf(self_dot, fabs_beta / 2.0f);

	complex<float> out = complex<float>(
		self_dot_beta * std::cos(fabs_beta * std::acos(in.real() / len)),
		in.imag() * self_dot_beta * std::sin(fabs_beta * std::acos(in.real() / len)) / sqrtf(in.imag() * in.imag()));

	//if (beta < 0)
	//	out = conj(out) / powf(abs(out), 2.0f);

	return out;

}


int main(void)
{
	sedenion s;

	s.vertex_data[0] = 0.2f;
	s.vertex_data[1] = 0.5f;
	s = pow_sedenion(s, 2.0f);

	cout << s.vertex_data[0] << " " << s.vertex_data[1] << endl;


	complex<float> c(0.2f, 0.5f);
	c = pow_complex(c, 2.0);

	cout << c.real() << " " << c.imag() << endl;

	return 0;
}

taby said:
Is there any way to avoid 16 nested for loops, or is there some C++ magic?

Recently i saw an interesting library which gives auto simd vectorization for such arrays of numbers. But can't find it again, and don't remember how its code looked like. :(
It was not something trivial like the counter less loop below , but some new language feature to execute some code for each element. Maybe somebody else knows what i could mean… : )

for (auto/*sorry*/ n : in.vertex_data) all_self_dot += n*n;

Though, that's shorter than your classic for loops, so maybe preferred.

However, some performance worries on my side:

class sedenion
{
public:

	sedenion(void)
	{
		vertex_data.resize(vertex_length, 0);
	}

	size_t vertex_length = 16;
	vector<float> vertex_data;
};

This is really bad imo. By using a vector, you have no control where the data is, so it will scatter in memory randomly. There is a need to allocate the class and the data, and each access to the data is chasing a pointer.
If your size of 16 is constant, you should use something like this instead:

struct sedenion
{
	static conxtexpr size_t vertex_length = 16; // Make compile time constant. Though, we coudl just use vertex_data.size() instead.
	std::array<float, vertex_length> vertex_data = { 0 }; // Use array over vector. So if we had a vector of sedenions, It's nice sequential blocks of 16 numbers in memory.
};

Even if you need different lengths than 16, you could make the length a template argument.

You said often you like vectors, but it seems you tend to overuse them.
Using vectors when it's (small) size is known at compile time is a bad mistake. Imagine a math lib which uses std::vector for vec3, vec4, matrices… everything. Does not make anything easier, but really slow.

JoeJ said:
but some new language feature to execute some code for each element.

Kept searching, but did not find.
Ofc. we have this alternative, using a function object:

float all_self_dot = 0;
std::for_each(vertex_data.begin(), vertex_data.end(), [&all_self_dot](float n){ all_self_dot += n*n; });

But personally i shy away because i expect some overhead cost?
(Syntax may be wrong too because i use it rarely, mostly if i want to share a complex iteration mechanism for multiple functions.)

Once again you come to the rescue lol

I reduced the dimensions from 16 to 5, to make things manageable.

so now my question is: do you have a favourite way to reduce 5D data down to 3D?

taby said:
do you have a favourite way to reduce 5D data down to 3D?

Yes! Like this (4D to 3D):

But idk how it works :D

For 3D to 2D, we could slice 3D stuff with a 3D plane to get a 2D intersection. But for higher dimension? hmmm…

Thanks for the advice!

I ended up doing this:

vertex_3 v3;

v3.x = v4.x / v4.w;

v3.y = v4.y / v4.w;

etc

Use 3D cameras.

Check it out: it's a 5D point cloud converted to a 3D point cloud, then converted to a mesh using Marching Cubes.

It looks like the quaternion fractal, but with an extra parameter (5 in total).

P.S. I can now make the Mandelbrot set for any dimension n, where n is a positive integer. It's not just n = {2, 4, 8, 16} anymore. I consider this a breakthrough, if it's not already well known.