Maxima is mostly an algebraic system , however using various
simplifications
numerical calculations are also possible ; even those
that have a tendency to cause
numerical overflow. Yes numerical
overflow is possible in maxima unless various
precautions are taken.
This is the maxima equivalent to the fractint formula d2jbMandelbrot
that I posted
recently .
You may not need to compile() the function ,
I tend to do this as I'm
always assuming that a speed increase will
result.
The file iterdiff2c.wxm
/* [wxMaxima batch file version 1] [
DO NOT EDIT BY HAND! ]*/
/* [ Created with wxMaxima version 11.08.0 ]
*/
/* [wxMaxima: comment start ]
.
iterdiff2c.wxm
(c) Copyright 2017
, sciwise(a)ihug.co.nz
Methodology for expressing differential equations
in terms of an
iterated function.
In particular a first order
differential equation and the derivative
thereof .
Note that f(z,c)
might be any iterated function .
.
[wxMaxima: comment end ] */
/*
[wxMaxima: comment start ]
'
Differential equation we're examining .
diff(f(z),z) + 3*f(z)*z = sin(z)
for an iterated function z = f(z) ,
g(z) =f(f(z))
diff(g(z),z) + 3*g(z)*f(z) = sin(f(z))
.
[wxMaxima:
comment end ] */
/* [wxMaxima: input start ]
*/
kill(all)$
eq1:diff(g(z,c),z) + 3*g(z,c)*f(z,c) -
sin(f(z,c))$
expand(eq1);
expand(diff(eq1,z)) ;
/* [wxMaxima: input end
] */
/* [wxMaxima: comment start ]
For the iterated Mandelbrot fractal
function .
[wxMaxima: comment end ] */
/* [wxMaxima: input start ]
*/
f(z,c):=z*z+c$
g(z,c):=f(f(z,c),c)$
f(z,c);
g(z,c);
/* [wxMaxima:
input end ] */
/* [wxMaxima: comment start ]
[wxMaxima: comment end ]
*/
/* [wxMaxima: input start ] */
eq1:diff(g(z,c),z) + 3*g(z,c)*f(z,c)
- sin(f(z,c))$
ed:expand(eq1)$
edp:expand(diff(eq1,z))$
/* [wxMaxima:
input end ] */
/* [wxMaxima: comment start ]
[wxMaxima: comment end ]
*/
/* [wxMaxima: input start ] */
ed;
edp;
/* [wxMaxima: input end ]
*/
/* [wxMaxima: comment start ]
This may well simplify the
construction of fractint formulas for this
type of equation .
[wxMaxima: comment end ] */
/* [wxMaxima: comment start ]
z=Pixel
c=Pixel
:
ed =
-sin(z^2+c)+3*z^6+9*c*z^4+4*z^3+9*c^2*z^2+3*c*z^2+4*c*z+3*c^3+3*c^2
edp
= -2*z*cos(z^2+c)+18*z^5+36*c*z^3+12*z^2+18*c^2*z+6*c*z+4*c
[wxMaxima:
comment end ] */
/* [wxMaxima: input start ] */
kill(all)$
/*
[wxMaxima: input end ] */
/* [wxMaxima: input start ] */
ratprint :
false$
/* [wxMaxima: input end ] */
/* [wxMaxima: comment start
]
Iterated function , Differential and Newton's method.
[wxMaxima:
comment end ] */
/* [wxMaxima: input start ]
*/
iterdiffn(c,n):=block([z,ed,edp,mag,j,i],
mag:0.001,
j:0,
z:c,
ed:1,
for i:0 while( i < n) and ( abs(ed) > mag ) do(
z:expand(z),
z:radcan(z),
z:realpart(z) + %i*imagpart(z),
z:ev(z,numer),
ed :
-sin(z^2+c)+3*z^6+9*c*z^4+4*z^3+9*c^2*z^2+3*c*z^2+4*c*z+3*c^3+3*c^2,
ed:expand(ed),
ed:radcan(ed),
ed:realpart(ed) + %i*imagpart(ed),
ed:ev(ed,numer),
edp :
-2*z*cos(z^2+c)+18*z^5+36*c*z^3+12*z^2+18*c^2*z+6*c*z+4*c,
edp:expand(edp),
edp:radcan(edp),
edp:realpart(edp) +
%i*imagpart(edp),
edp:ev(edp,numer),
z : z - (ed/edp),
j : j+1
),
realpart(ed) + %i*imagpart(ed)
)$
/* [wxMaxima: input end ] */
/*
[wxMaxima: comment start ]
[wxMaxima: comment end ] */
/* [wxMaxima:
input start ] */
compile(iterdiffn)$
/* [wxMaxima: input end ] */
/*
[wxMaxima: comment start ]
This returns the value of the differential
equation for the intial value of c.
If this is close to zero then the
last value of z might be a complex root of the
Differential Equation .
[wxMaxima: comment end ] */
/* [wxMaxima: input start ] */
n:6$
c:0.01
+ %i*0.1$
iterdiffn(c,n);
/* [wxMaxima: input end ] */
/* [wxMaxima:
comment start ]
[wxMaxima: comment end ] */
/* Maxima can't
load/batch files which end with a comment! */
"Created with
wxMaxima"$
Also the file iterdiff2c.mac
/*
.
iterdiff2c.mac
(c)
Copyright 2017 , sciwise(a)ihug.co.nz
Methodology for expressing
differential equations in terms of an
iterated function.
In particular
a first order differential equation and the derivative
thereof .
Note
that f(z,c) might be any iterated function .
.
*/
/*
'
Differential
equation we're examining .
diff(f(z),z) + 3*f(z)*z = sin(z)
for an
iterated function z = f(z) , g(z) =f(f(z))
diff(g(z),z) + 3*g(z)*f(z)
= sin(f(z))
.
*/
kill(all)$
eq1:diff(g(z,c),z) + 3*g(z,c)*f(z,c) -
sin(f(z,c))$
expand(eq1);
expand(diff(eq1,z)) ;
/*
For the iterated
Mandelbrot fractal function
.
*/
f(z,c):=z*z+c$
g(z,c):=f(f(z,c),c)$
f(z,c);
g(z,c);
/*
*/
eq1:diff(g(z,c),z)
+ 3*g(z,c)*f(z,c) -
sin(f(z,c))$
ed:expand(eq1)$
edp:expand(diff(eq1,z))$
/*
*/
ed;
edp;
/*
This may well simplify the construction of fractint formulas for
this
type of equation .
*/
/*
z=Pixel
c=Pixel
:
ed =
-sin(z^2+c)+3*z^6+9*c*z^4+4*z^3+9*c^2*z^2+3*c*z^2+4*c*z+3*c^3+3*c^2
edp
=
-2*z*cos(z^2+c)+18*z^5+36*c*z^3+12*z^2+18*c^2*z+6*c*z+4*c
*/
kill(all)$
ratprint
: false$
/*
Iterated function , Differential and Newton's
method.
*/
iterdiffn(c,n):=block([z,ed,edp,mag,j,i],
mag:0.001,
j:0,
z:c,
ed:1,
for i:0 while( i < n) and ( abs(ed) > mag ) do(
z:expand(z),
z:radcan(z),
z:realpart(z) + %i*imagpart(z),
z:ev(z,numer),
ed :
-sin(z^2+c)+3*z^6+9*c*z^4+4*z^3+9*c^2*z^2+3*c*z^2+4*c*z+3*c^3+3*c^2,
ed:expand(ed),
ed:radcan(ed),
ed:realpart(ed) + %i*imagpart(ed),
ed:ev(ed,numer),
edp :
-2*z*cos(z^2+c)+18*z^5+36*c*z^3+12*z^2+18*c^2*z+6*c*z+4*c,
edp:expand(edp),
edp:radcan(edp),
edp:realpart(edp) +
%i*imagpart(edp),
edp:ev(edp,numer),
z : z - (ed/edp),
j : j+1
),
realpart(ed) + %i*imagpart(ed)
)$
/*
*/
compile(iterdiffn)$
/*
This
returns the value of the differential equation for the intial value of
c.
If this is close to zero then the last value of z might be a
complex root of the
Differential Equation .
*/
n:6$
c:0.01 +
%i*0.1$
iterdiffn(c,n);
/*
*/
