News: Support the forums via Fractalforums.com Online Store
 Welcome, Guest. Please login or register. August 02, 2014, 04:42:41 AM 1 Hour 1 Day 1 Week 1 Month Forever Login with username, password and session length ••• ••• ••• ••• ••• ••• ••• ••• ••• ••• ••• •••

 Pages: [1]   Go Down
 Author Topic: n dimensional extension of Mandelbrot set, based on the SVD  (Read 1057 times) Description: SVD being the singular value decomposition 0 Members and 1 Guest are viewing this topic.
kram1032
Fractal Senior

Posts: 1533

 « on: March 06, 2011, 01:44:31 AM »

So yesterday I tried out something I wondered for quite a while...
I took the Singular Value Decompoition of an aribitary 2D vector M={{x},{y}}, squared the results in a nearly trivial manner and multiplied it all together.
What I discovered, was exactly what I hoped for:

M²=U²Sigma²transpose(V)²,
where U and V in general are rotation-style matrices and Sigma is the Matrix of Singular values. It behaves a lot like a radius.
In this special case, V = 1 and Sigma = {sqrt(a^2+b^2),0,0}.

squaring a square matrix, which is the case for both U and V, is always defined. Sigma will always be as close to a diagonal matrix as possible, which also is defined: Simply square all the numbers.

Well, doing so for said 2D vector above, you actually get {{x²-y²},{2xy}}, which perfectly corresponds to complex number multiplication.

Now, what would stop anybody to do the very same thing for a 3D vector?
Nothing, luckily.
So I have a more or less natural extension of the Mset to 3D or to nD in general

What I got for 3D after svd, squaring and simplifying, is this:

x->x²-(y*z*(x+sqrt(x²+y²+z²)))/sqrt(x²+z²)
y->x*y+z*sqrt(x²+z²)
z->x*z+(y*(z²+x*sqrt(x²+y²+z²)))/sqrt(x²+z²)

now, just add the constants a,b,c and you should, in theory, get an Mset as close as possible to a "true" 3D version.
At least as far as vector stuff goes.

The Singular Value Decomposition is so to speak the "ideal" decomposition of any Matrix where all the information can be found: rotation+translation, scaling, back-rotation+translation
Using it should make it possible for a wide range of problems to be "propely" extended, as far as this even exists.

I already experimented a bit with this SV3D MSet and did an xy-slice buddhabrot.
Evident from that is, that this is NOT a extension in the sense of simply adding a dimension to the Mset: the slice has the same basic body vs. head structure but besides that doesn't look anywhere close to the basic Mset usually seen.
It's rather a port to 3D, based on directly perform a "double angle" operation in a 3D-vector sense.
How exactly this angular conversion works in this case, I'm not sure about.

I'll soon post the xy plane and after that start to render the xz plane. What I can see so far is, that it is symmetric to the y-axis....

I'm also curious what happens when looking at the full 3D object... How much whipped cream vs how much "actual detail"...
By now I'm kinda certain that the "true" nD Mset with n>2 will always include whipped cream. I guess, it's partially due to the whole thing being a kind of rotation. So I wouldn't be surprised if the above transform results in that. Still I'd love to see the results or how it works together with other variants. - The various hybrids look awesome and adding another "template" to play with and cross-section with other sets could give new interesting patterns...

xy plane (edited and updated)

xz plane
 « Last Edit: March 11, 2011, 08:59:02 PM by kram1032 » Logged
M Benesi
Fractal Bachius

Posts: 560

 « Reply #1 on: March 06, 2011, 02:08:55 AM »

It's... assuming I got the formula right (code follows after image)... a little.. solid.

Code:
r1=sqrt(sqr(sx)+sqr(sz));
r2=sqrt(sqr(sx)+sqr(sy)+sqr(sz));

nx=sqr(sx)-(sy*sz*(sx+r2))*r1^-1;    //then add in pixel values ....
ny=sx*sy+sz*r1;                             // tried it with r1^-1 in case the multiplication was a mistake... was swirly...
nz=sx*sz+(sy*(sqr(sz)+sx*r2))*r1^-1;
 Logged

kram1032
Fractal Senior

Posts: 1533

 « Reply #2 on: March 06, 2011, 02:39:43 AM »

yup, you got them right...
the multiplication wasn't an error.
Ok, I expected whipped cream... though, that's quite a lot...

Well, I'll try to derive the power 8 formulas or, if possible, the power n ones (with n being a real)
Just a moment or three

Edit:
Ok, this seems to be too much for wolfram alpha to take in a simplified form...

what needs to be done should be something like this this:

MatrixPower[{{a/x, -c/y, -(a b)/(x y)}, {b/x, 0, y/x}, {c/x, a/y, -(b c)/(x y)}}, n] * {{x^n},{0},{0}}, assume x=sqrt(a^2+b^2+c^2), y=sqrt(a^2+c^2) // Full Simplify

if somebody could put that into mathematica (and possibly correct it as I myself can't be sure if that actually is the right thing to put in Mathematica, not owning a copy of the program myself), we'd have the general solutions.
By using that, we could go on defining the exponential or different stuff...

Directly using the singular value decomposition should also be possible inside mathematica. I heavily simplified the results of the svd to get to something that's actually possible to work with, but there is always the possibility that I assumed something that only allows positive reals or something... (I hope, I avoided that, though. Most simplifications where done by assuming x*conjugate(x)=x², e.g. the values you put in are real. So the formulas above are not valid but never where meant to be if you use complex x,y or z...)

should be something like:

M={{x},{y},{z}}

SingularValueDecomposition[M]

U^n*Sigma^n*transpose(V)^n // Full Simplify
:= M^n

however, I'm not sure wether mathematica will like the vector power of Sigma, as that's not square... in that case, you might have to directly calculate the x coordinate over the nth power, as the other two coordinates will be 0 here, anyway...

if we have the nth power of M, we can easily define the, let's for now call it "pseudo-exponential" by the power series sum(1/n!*M^n,n,0,infinity)
 « Last Edit: March 06, 2011, 03:07:54 AM by kram1032 » Logged
M Benesi
Fractal Bachius

Posts: 560

 « Reply #3 on: March 06, 2011, 04:36:59 AM »

Tried loading LAPACK in maxima, wouldn't work.... sheesh.  Gotta find a working version.  Suppose I could update maxima and search the add on packs... but that may not happen for a while.

BTW-  Maxima is a freeware program that is similar in functionality to Mathematica (although Mathematica has more features- it should as you pay out the 4 for it...).
 Logged

kram1032
Fractal Senior

Posts: 1533

 « Reply #4 on: March 06, 2011, 11:47:53 AM »

Maxima is nice, just tried it out

However, what we want to do here, an algebratic SVD for ANY vector, can't be done in Maxima, apparently...
I searched around a bit and it only supports numberic evaluation of the SVD, as it seems...

So I fear, we do need mathematica. OR I try to use the wolfram alpha results which I simplified partially by hand, partially by wolfram alpha itself... let's see if Maxima handles that
*testing*

Edit: Ok, Maxima gets just as far as wolfram alpha.
It can't generally evalute the nth power of the matrix
It could evalute the 8th power, however, the result is "too long to be viewed" :-/

So I could use it but only inside Maxima...

I tried simplifying it but the result was still too long...

Wolfram Alpha actually was able to compute the 8th power result but it's a mess of factors of indefinite length...

I'll try one more trick. If I let WA evaluate a power n generic 3x3 matrix I could just substitute the other values...

Edit: darn it. WA *can* do a generic 2x2 matrix to the nth, but it can't evaluate a generic 3x3 one. - it pretends to not even understand the input -.-
That means, we gotta have to try mathematica or a different symbolic CAS. Maxima is great but not that great
 « Last Edit: March 06, 2011, 12:26:50 PM by kram1032 » Logged
msltoe
Iterator

Posts: 159

 « Reply #5 on: March 06, 2011, 03:19:42 PM »

kram1032,
Empathy over here. I'm trying to take the inverse square root of a 3x3 matrix. While I wrote a C program, it's iterative and slowwwww.
For 8th order power, try M2=M*M;M3=M2*M2;M4=M3*M3;
-mike
 Logged
kram1032
Fractal Senior

Posts: 1533

 « Reply #6 on: March 06, 2011, 05:41:56 PM »

(Edit: just for reference, here the generic 3x3 ^-2 matrix http://www.wolframalpha.com/input/?i=simplify+{{a%2Cb%2Cc}%2C{d%2Cf%2Cg}%2C{h%2Cj%2Ck}}^%28-2%29 )

Ok, by doing so, I get the following and hope, I don't mess up anything... it's REALLY long... I'm almost certain, it could be simplified a lot, though.

$V=\begin{pmatrix} a\\b\\c \end{pmatrix}$
$x=sqrt{a^2+b^2+c^2}$
$y=sqrt{a^2+c^2}$

$V^8=\begin{pmatrix}{{{2ac^4x^3y^7-2a^3c^2x^3y^7+6a^2c^4x^2y^7-4a^4c^2x^2y^7-6a^5c^2xy^7+a^8y^7-a^3bcx^4y^6-3bc^5x^3y^6+15a^2bc^3x^3y^6-4a^4bcx^3y^6-9abc^5x^2y^6+38a^3bc^3x^2y^6-8a^5bcx^2y^6+25a^4bc^3xy^6-12a^6bcxy^6-7a^7bcy^6+3a^2b^2c^2x^4y^5-20ab^2c^4x^3y^5+26a^3b^2c^2x^3y^5-2a^5b^2x^3y^5+3b^2c^6x^2y^5-54a^2b^2c^4x^2y^5+64a^4b^2c^2x^2y^5-4a^6b^2x^2y^5-40a^3b^2c^4xy^5+50a^5b^2c^2xy^5-6a^7b^2xy^5+21a^6b^2c^2y^5-3ab^3c^3x^4y^4+3b^3c^5x^3y^4-34a^2b^3c^3x^3y^4+15a^4b^3cx^3y^4+22ab^3c^5x^2y^4-90a^3b^3c^3x^2y^4+38a^5b^3cx^2y^4+30a^2b^3c^5xy^4-80a^4b3c^3xy^4+25a^6b^3cxy^4-35a^5b^3c^3y^4+b^4c^4x^4y^3+6ab^4c^4x^3y^3-20a^3b^4c^2x^3y^3+2a^5b^4x^3y^3-2b^4c^6x^2y^3+38a^2b^4c^4x^2y^3-54a^4b^4c^2x^2y^3+6a^6b^4x^2y^3-10ab^4c^6xy^3+60a^3b^4c^4xy^3-40a^5b^4c^2xy^3+35a^4b^4c^4y^3+3a^2b^5c^3x^3y^2-3a^4b^5cx^3y^2-4ab^5c^5x^2y^2+22a^3b^5c^3x^2y^2-9a^5b^5cx^2y^2+b^5c^7xy^2-20a^2b^5c^5xy^2+30a^4b^5c^3xy^2-21a^3b^5c^5y^2-2a^2b^6c^4x^2y+3a^4b^6c^2x^2y+2ab^6c^6xy-10a^3b^6c^4xy+7a^2b^6c^6y+a^2b^7c^5x-ab^7c^7}\over{y^7}}\cr {{\left(cy+ab \right)\left(a^3x^3y^6+c^4x^2y^6-4a^2c^2x^2y^6+a^4x^2y^6-4a^3c^2xy^6+a^5xy^6+a^6y^6-6a^2bcx^3y^5+14abc^3x^2y^5-14a^3bcx^2y^5+12a^2bc^3xy^5-14a^4bcxy^5-6a^5bcy^5+9ab^2c^2x^3y^4-3a^3b^2x^3y^4-6b^2c^4x^2y^4+40a^2b^2c^2x^2y^4-6a^4b^2x^2y^4-12ab^2c^4xy^4+38a^3b^2c^2xy^4-5a^5b^2xy^4+15a^4b^2c^2y^4-4b^3c^3x^3y^3+6a^2b^3cx^3y^3-24ab^3c^3x^2y^3+24a^3b^3cx^2y^3+4b^3c^5xy^3-40a^2b^3c^3xy^3+16a^4b^3cxy^3-20a^3b^3c^3y^3-3ab^4c^2x^3y^2+3b^4c^4x^2y^2-16a^2b^4c^2x^2y^2+3a^4b^4x^2y^2+17ab^4c^4xy^2-18a^3b^4c^2xy^2+15a^2b^4c^4y^2+2ab^5c^3x^2y-2a^3b^5cx^2y-2b^5c^5xy+8a^2b^5c^3xy-6ab^5c^5y-ab^6c^4x+b^6c^6\right)}\over{y^6}}\cr -{{3a^2c^3x^3y^7-a^4cx^3y^7-3ac^5x^2y^7+6a^3c^3x^2y^7-a^5cx^2y^7+5a^4c^3xy^7-a^6cxy^7-a^7cy^7-a^4bx^4y^6-7abc^4x^3y^6+14a^3bc^2x^3y^6-a^5bx^3y^6+3bc^6x^2y^6-30a^2bc^4x^2y^6+21a^4bc^2x^2y^6-a^6bx^2y^6-20a^3bc^4xy^6+16a^5bc^2xy^6-a^7bxy^6+7a^6bc^2y^6+a^3b^2cx^4y^5+3b^2c^5x^3y^5-30a^2b^2c^3x^3y^5+15a^4b^2cx^3y^5+30ab^2c^5x^2y^5-76a^3b^2c^3x^2y^5+19a^5b^2cx^2y^5+30a^2b^2c^5xy^5-55a^4b^2c^3xy^5+11a^6b^2cxy^5-21a^5b^2c^3y^5-3a^2b^3c^2x^4y^4+14ab^3c^4x^3y^4-34a^3b^3c^2x^3y^4+4a^5b^3x^3y^4-6b^3c^6x^2y^4+76a^2b^3c^4x^2y^4-64a^4b^3c^2x^2y^4+4a^6b^3x^2y^4-20ab^3c^6xy^4+80a^3b^3c^4xy^4-35a^5b^3c^2xy^4+35a^4b^3c^4y^4+ab^4c^3x^4y^3-b^4c^5x^3y^3+15a^2b^4c^3x^3y^3-12a^4b^4cx^3y^3-21ab^4c^5x^2y3+64a^3b^4c^3x^2y^3-15a^5b^4cx^2y^3+5b^4c^7xy^3-55a^2b^4c^5xy^3+50a^4b^4c^3xy^3-35a^3b^4c^5y^3-ab^5c^4x^3y^2+4a^3b^5c^2x^3y^2-a^5b^5x^3y^2+b^5c^6x^2y^2-19a^2b^5c^4x^2y^2+15a^4b^5c^2x^2y^2+16ab^5c^6xy^2-35a^3b^5c^4xy^2+21a^2b^5c^6y^2+ab^6c^5x^2y-4a^3b^6c^3x^2y-b^6c^7xy+11a^2b6c^5xy-7ab^6c^7y-ab^7c^6x+b^7c^8 }\over{y^7}}\cr }\end{pmatrix}$
 « Last Edit: March 06, 2011, 07:01:55 PM by kram1032 » Logged
M Benesi
Fractal Bachius

Posts: 560

 « Reply #7 on: March 07, 2011, 02:33:08 AM »

Ok... those long equations made me laugh.  On another note I remember getting those "too long to display" equations in maxima for some work I was doing in the past.  Gotta love it.

Another nickel to throw in:  I found that the simplification in Maxima and Mathematica (both, my father used his copy of Mathematica to try and simplify something for me a while back as Maxima couldn't do it) to sometimes miss simple tricks, although a lot of the time their simplification routines are totally awesome.  I believe it was simple trig tricks that Maxima missed, and I can't recall exactly what it was that both couldn't handle (I think, once again, it was trig stuff, substituting cosine for sqrt(1-sin^2) or something like that).
 Logged

KRAFTWERK
Fractal Supremo

Posts: 1145

Virtual Surreality

 « Reply #8 on: March 07, 2011, 01:34:47 PM »

This is interesting! Beautiful Buddhabrot, love that formula Mark, it covers all my three screens.
Now we need a formula that makes formulas....
 Logged

Fractalis Surrealis: http://mandelwerk.deviantart.com/gallery/28152444
And my 3.3 GIGA PIXEL zoomable mandelbulb: http://mandelwerk.com/
kram1032
Fractal Senior

Posts: 1533

 « Reply #9 on: March 11, 2011, 08:52:01 PM »

Either way, triple-squaring should work just fine, I guess...
this way of taking powers should be power-associative or how it's called (what I'm trying to say is, that x*x should be x² and x*x*x should be x²*x or x*x² or x³ and so on)

It's proably not associative, though. (x*y)!=(y*x), but nobody here should really have a problem with that

(updated first post with the xz plane)
 « Last Edit: March 11, 2011, 08:59:28 PM by kram1032 » Logged
 Pages: [1]   Go Down