返回列表 回复 发帖

复分形难题再度征解

搞了好几遍,尝试了各种情形,仍不能得出结果:
Carr1965 {
; Updated for UF2 by Erik Reckase, March 2000
           ; Modified Sylvie Gallet frm.
init:
  z=c=pixel
  int compt=0, int limit=round(real(p1)), float bailout=4
loop:
  IF (compt>=limit)
    c = (-.743380900000982,.131850030300002)
  ELSE
    c = z + c/2.125
  ENDIF
  z=z*z+c
  compt=compt+1
bailout:
  |real(z)|<=bailout
default:
  title = "Carr 1965"
  periodicity = 0
  maxiter = 500
  magn = .9
  center = (0,0)
  method = multipass
  param p1
    caption = "Iter Limit"
    default = (100,0)
    hint = "The real part of this parameter signals a change \
            in the value of c. The imaginary part is not used."
  endparam
}
UF中效果图如下:
Fractal2.jpg
昨天干了五六遍,失败,今天干了一上午,失败了三次,终于我已经找到解决方案了,正在扫描。
未命名2.jpg
改变P1的横坐标为2,得:
未命名1.jpg
这个结果很有诱惑力,可是没办法完成:
pgd_Twofold2_Julia {
  ; Two-fold rotational symmetry and complementary critical
  ; point fates make for some very intertwining (and
  ; entertaining) Julia sets.
  ; f(z) = (z^2 - 1/3)/cz^3
global:
  bool foundatt1 = false
  bool foundatt2 = true
  complex att1 = 1
  complex att2 = -1
  float b = 1/@bailout
  int j = 0
  ; Locate attractor 1 goes to, if any.
  WHILE (j < #maxit)
    att1 = (sqr(att1) - 1/3)/(@seed * sqr(att1) * att1)
    j = j + 1
  ENDWHILE
  ; Is it really an attractor? Or did it just wander off
  ; and die somewhere? Is the attractor self-complementary
  ; or is there a distinct mirror attractor?
  j = 0
  complex a = att1
  WHILE (j < #maxit)
    att1 = (sqr(att1) - 1/3)/(@seed * sqr(att1) * att1)
    IF (|att1 - a| < b)
      foundatt1 = true      ; it repeated, yay!
    ENDIF
    att2 = (sqr(att2) - 1/3)/(@seed * sqr(att2) * att2)
    IF (|att2 - a| < b)
      foundatt2 = false      ; other cp went to same.
    ENDIF
    j = j + 1
  ENDWHILE
init:
  z = #pixel
  bool notdone = true
  bool landedright = false
  IF ((@what == 1) && (!foundatt1))
    notdone = false
  ENDIF
  IF ((@what == 2) && (!foundatt2))
    notdone = false
  ENDIF
  i = 0
loop:
  z = (sqr(z) - 1/3)/(@seed * sqr(z) * z)
  IF (|z| > @bailout)
    notdone = false
    landedright = (@what == 0)
  ENDIF
  IF (foundatt1)
    IF (|z - att1| < b)
      notdone = false
      landedright = (@what == 1)
    ENDIF
  ENDIF
  IF (foundatt2)
    IF (|z - att2| < b)
      notdone = false
      landedright = (@what == 2)
    ENDIF
  ENDIF
  IF (i == #maxit)
    IF (@what != 3)
      notdone = false
      landedright = false
    ENDIF
  ENDIF
  i = i + 1
  IF (!notdone)
    IF (!landedright)
      z = (501,10)  ; Tell Smooth (Generalized) 2 to color this thing solid.
    ENDIF
  ENDIF
bailout:
  notdone
default:
  title = "Twofold2 Julia"
  periodicity = 0
  maxiter = 1000
  method = onepass
  magn = 0.5
  param seed
    caption = "Julia seed"
    default = (0.75, 0)
  endparam
  param what
    caption = "Draw what"
    enum = "Main basin" "Other attractor #1" "Other attractor #2" "Trapped points"
    default = 0
  endparam
  param bailout
    caption = "Bailout"
    hint = "Escape radius (inverse is used for finite attractors)"
    default = 100000.0
    min = 0
  endparam
switch:
  type = "pgd_Twofold2_Mandel"
  bailout = bailout
}

Fractal1.jpg (52.69 KB)

Fractal1.jpg

这个代码读起来,太费脑细胞了,UF中有几个分形代码,要么太过漫长,要么逻辑判断太多,作起来太费劲,只好望洋兴叹了。
这个应该就是分形边界的一种构造方法,把内、外部都挖空了,要通过多次迭代才能完成。
实在挡不住诱惑,虽然伤了很多脑细胞但算是有成果:

pgd_Twofold2_Julia.JPG (39.51 KB)

pgd_Twofold2_Julia.JPG

这个代码,我可读不懂了。榕兄作得不错。
7# 柳烟


在对Z迭代之前先进行两次迭代获得一些准备数值,其它的就是按照代码依样画葫芦了。另一个它的“M"集那真是没办法了(其中有数组变量不知怎么转为几何画板能完成的),分开结构非常漂亮:
pgd_Twofold2_Mandel {
  ; Two-fold rotational symmetry and complementary critical
  ; point fates make for some very intertwining (and
  ; entertaining) Julia sets.
  ; f(z) = (z^2 - 1/3)/cz^3
init:
  z = 1
  bool notdone = true
  bool landedright = false
  int i = 0
  int j = 0
  bool odd = false
  complex z_values[@zmax]
  float b = 1/@bailout
loop:
  z = (sqr(z) - 1/3)/(#pixel * sqr(z) * z)
  IF (|z| > @bailout)
    notdone = false
    landedright = (@what == 0)
  ENDIF
  z_values = z
  IF (odd)
    IF (|z - z_values[j]| < b)
      notdone = false
      complex w = -1
      j = 0
      bool sep = true
      WHILE ((j < #maxit) && sep)
        w = (sqr(w) - 1/3)/(#pixel * sqr(w) * w)
        IF (|w - z| < b)
          sep = false
        ENDIF
        j = j + 1
      ENDWHILE
      IF (sep)
        landedright = (@what == 2)
      ELSE
        landedright = (@what == 1)
      ENDIF
    ENDIF
    j = j + 1
  ENDIF
  odd = !odd
  i = i + 1
  IF (i == @zmax)
    IF (@what != 3)
      notdone = false
      landedright = false
    ENDIF
  ENDIF
  IF (!notdone)
    IF (!landedright)
      z = (501,10)  ; Tell Smooth (Generalized) 2 to color this thing solid.
    ENDIF
  ENDIF
bailout:
  notdone
default:
  title = "Twofold2 Mandelbrot"
  periodicity = 0
  maxiter = 1000
  method = onepass
  param what
    caption = "Draw what"
    enum = "Main basin" "One finite attractor" "Two attractors" "Trapped points"
    default = 0
  endparam
  param zmax
    caption = "Max iters to record"
    default = 1000
  endparam
  param bailout
    caption = "Bailout"
    hint = "Escape radius (inverse is used for finite attractors)"
    default = 100000.0
    min = 0
  endparam
switch:
  type = "pgd_Twofold2_Julia"
  seed = #pixel
  bailout = bailout
}

捕获.JPG (94.01 KB)

捕获.JPG

8# 榕坚
谢谢。佩服老兄的钻研精神与思维的敏捷。UF中,每当我看到这种代码,就如珠峰挡在前面。对这种逻辑判断太多的分形,我是望而生畏。榕兄有空时,能否制成视频,详细反映出制作过程。注意有空时。另还有一个分形,是谢尔宾斯基三角形,代码也十分漫长,一时忘了在何处,我找着后发在这里。
位于gnd.ufm中,
gnd-slope-sierpinski2 {
; Based on a formula by Ramiro Perez
; Slope version by Gilles Nadeau, 2007
init:
  z1 = #pixel
  z2 = #pixel + @offset
  z3 = #pixel + flip(@offset)
  int done = 2
  float modz = 0.0
  float e1 = 0.0         ; potentials
  float e2 = 0.0
  float e3 = 0.0
  float vx = 0.0         ; normal vector
  float vy = 0.0
  float vz = 0.0
  float vd = 0.0
  float iterexp1 = 0.0
  float iterexp2 = 0.0
  float iterexp3 = 0.0

loop:

  IF ((imag(z1)>=.575*real(z1)) && (-.575*real(z1)<=imag(z1)))
    z1 = 2*z1-1i
  ELSEIF (real(z1)<=0)
    z1 = 2*z1+(.8660254, 0.5)
  ELSEIF (real(z1)>0)
    z1 = 2*z1+(-.8660254, 0.5)
  ELSE
    z1 = 2*z1
  ENDIF

  IF ((imag(z2)>=.575*real(z2)) && (-.575*real(z2)<=imag(z2)))
    z2 = 2*z2-1i
  ELSEIF (real(z2)<=0)
    z2 = 2*z2+(.8660254, 0.5)
  ELSEIF (real(z2)>0)
    z2 = 2*z2+(-.8660254, 0.5)
  ELSE
    z2 = 2*z2
  ENDIF

  IF ((imag(z3)>=.575*real(z3)) && (-.575*real(z3)<=imag(z3)))
    z3 = 2*z3-1i
  ELSEIF (real(z3)<=0)
    z3 = 2*z3+(.8660254, 0.5)
  ELSEIF (real(z3)>0)
    z3 = 2*z3+(-.8660254, 0.5)
  ELSE
    z3 = 2*z3
  ENDIF
  
  modz = |z1|
  
  iterexp1 = iterexp1 + exp(-cabs(z1))
  iterexp2 = iterexp2 + exp(-cabs(z2))
  iterexp3 = iterexp3 + exp(-cabs(z3))
  done = done + 1         ; increment iteration counter

  IF (modz > @bailout ||       \
      @everyiter ||            \
      done == #maxit + 2)      ; done, or every iteration, or last
    ; determine continuous iteration (height) for each point
      e1 = iterexp1 * @zscale
      e2 = iterexp2 * @zscale
      e3 = iterexp3 * @zscale

    ; apply transfer function
    ; a function is not used because these are floats
    ; and not all functions apply to floats
    IF (@xfer == 1)             ; log
      e1 = log(e1)
      e2 = log(e2)
      e3 = log(e3)
    ELSEIF (@xfer == 2)         ; sqrt
      e1 = sqrt(e1)
      e2 = sqrt(e2)
      e3 = sqrt(e3)
    ELSEIF (@xfer == 3)         ; cuberoot
      e1 = (e1)^(1/3)
      e2 = (e2)^(1/3)
      e3 = (e3)^(1/3)
    ELSEIF (@xfer == 4)         ; exp
      e1 = exp(e1)
      e2 = exp(e2)
      e3 = exp(e3)
    ELSEIF (@xfer == 5)         ; sqr
      e1 = sqr(e1)
      e2 = sqr(e2)
      e3 = sqr(e3)
    ELSEIF (@xfer == 6)         ; cube
      e1 = (e1)^3
      e2 = (e2)^3
      e3 = (e3)^3
    ELSEIF (@xfer == 7)         ; sin
      e1 = sin(e1)
      e2 = sin(e2)
      e3 = sin(e3)
    ELSEIF (@xfer == 8)         ; cos
      e1 = cos(e1)
      e2 = cos(e2)
      e3 = cos(e3)
    ELSEIF (@xfer == 9)         ; tan
      e1 = tan(e1)
      e2 = tan(e2)
      e3 = tan(e3)
    ENDIF

    ; apply post-scale
    e1 = e1 * @zscale2
    e2 = e2 * @zscale2
    e3 = e3 * @zscale2

      vx = e2-e1
      vy = e3-e1
      vz = -@offset
      ; normalize vector
      vd = 1/sqrt(sqr(vx)+sqr(vy)+sqr(vz))
      vx = vx*vd
      vy = vy*vd
      vz = vz*vd
      z = vx + flip(vy)         ; fudge z from vector
  ELSE               ; didn't compute z this time
    z = z1            ; use primary iteration value to keep
                      ; periodicity working
  ENDIF
  IF (modz > @bailout)         ; we're done
    done = 0
  ENDIF

bailout:
  (done > 0)
default:
  title = "Slope Sierpinski Triangle II"
  center = (0.0, 0.0)
  magn = 1.1538
  maxiter = 149
  method = multipass
  periodicity = 0

  float param version
    default = 1.1
    visible = false
  endparam

  param bailout
    caption = "Bailout value"
    default = 127
    min = 1
  endparam

  param offset
    caption = "Orbit Separation"
    default = 0.01
    hint = "Defines how far apart the simultaneous orbits are.  Smaller \
            distances will produce more accurate results."
  endparam
  param xfer
    caption = "Height Transfer"
    default = 0
    enum = "linear" "log" "sqrt" "cuberoot" "exp" "sqr" "cube" \
    "sin" "cos" "tan"
    hint = "This function will be applied to the height value \
            before a slope is calculated."
  endparam
  param zscale
    caption = "Height Pre-Scale"
    default = 1.0
    hint = "Specifies the ratio between height and distance.  Higher \
            values will exaggerate differences between high and low. \
       In general, you will want to use smaller numbers here."
  endparam
  param zscale2
    caption = "Height Post-Scale"
    default = 0.025
    hint = "Specifies the ratio between height and distance; like \
            Height Pre-Scale, except that this value is applied after \
       the transfer function."
  endparam
  param everyiter
    caption = "Every Iteration"
    default = false
    hint = "If set, the surface normal will be computed at every \
            iteration.  If you are using a coloring algorithm which \
       processes every iteration, you will need this."
  endparam
}
Fractal2.jpg
我从面板上看了,好象那些函数开关项,一点作用都不起,榕兄再看看,是不是这样。
返回列表