DubeauPhi1 {
; by Etienne Saint-Amant
; 2009/06
; Based on the article of Francois Dubeau that appeared on the Journal of
; Computational and Appied Mathematics 224 in 2009, pp. 66-76.
; doi:10.1016/j.cam.2008.04.014
global:
;$define DEBUG
init:
; factorial
int func factorial(int i)
int total
total = 1
while i >= 1
total = total * i
i = i - 1
endwhile
return total
endfunc
; num function
float func num(float delta, int i)
float total
float vardelta
total = 1
vardelta = delta
while i >= 1
total = total * vardelta
vardelta = vardelta - 1
i = i - 1
endwhile
return total
endfunc
; binomial delta function
float func binomialdelta(float delta, int i)
float result
if i == 0
result = 1
else
result = num(delta,i)/factorial(i)
endif
return result
endfunc
z = pixel
loop:
complex sum
int i
sum = 0
i = 0
while i <= @p - 1
print(1/cabs(@n),":",i,":",binomialdelta(1/cabs(@n),i))
sum = sum + binomialdelta(1/cabs(@n),i) * (@r/z^@n - 1)^i
i = i + 1
endwhile
zold = z
z = z * sum
bailout:
|z - zold| >= @bailout
default:
title = "Dubeau Phi 1"
periodicity = 0
$IFDEF VER50
rating = recommended
$ENDIF
maxiter = 100
param n
caption = "Exponent"
default = (3,0)
hint = "Specifies the exponent of the equation that is solved by \
Dubeau's method. Use real numbers (set the imaginary part \
to zero) to obtain the correct Dubeau's Method"
endparam
param r
caption = "Root"
default = (1,0)
hint = "Specifies the root of the equation that is solved. Use larger \
numbers for slower convergence."
endparam
int param p
caption = "Order of convergence"
default = 2
hint = "Must be a natural number greater than 2. 2 will generate a \
standard Newton's method."
min = 2
endparam
param bailout
caption = "Bailout value"
default = 0.00001
min = 0
$IFDEF VER40
exponential = true
$ENDIF
hint = "This parameter defines how soon a convergent orbit bails out while \
iterating. Smaller values give more precise results but usually \
require more iterations."
endparam
}
DubeauPhi1(1).gsp (21.07 KB)
|