Board logo

标题: 疑难分形征解 [打印本页]

作者: 柳烟    时间: 2012-5-3 18:13     标题: 疑难分形征解

选自叶瑞松论文:《复动力系统生成的拟3D图象》
该论文地址:http://115.com/file/e76qbelp#
用几何画板如何按:z=z^2-0.5*z+0.563,w=z=z^2-0.5*v+0.563,v=z,w=z.我刚才用几何画板白干一场。
我试着用UF做,结果好象正确。UF代码是:
叶瑞松J {
;
; Generic Julia set.
;
init:
  z = #pixel
loop:
  z = z^@power -0.5*z+@seed
  w=  z^@power -0.5*v+@seed
  v=z
  z=w
  
bailout:
  |z| <= @bailout
default:
  title = "叶瑞松范例J"
  helpfile = "Uf*.chm"
  helptopic = "Html\formulas\standard\julia.html"
$IFDEF VER50
  rating = recommended
$ENDIF
  param seed
    caption = "Julia seed"
    default = (0.563, 0)
    hint = "Use this parameter to create many different Julia sets. A good \
            way to set this parameter is with the Switch, Eyedropper, or \
            Explore features."
  endparam
  param power
    caption = "Power"
    default = (2,0)
    hint = "This parameter sets the exponent for the Julia formula. \
            Increasing the real part to 3, 4, and so on, will multiply the \
            symmetry of the Julia figure. Non-integer real values and non-zero \
            imaginary values will create distorted Julia fractals. Use (2, 0) \
            for the standard Julia set."
  endparam
  param bailout
    caption = "Bailout value"
    default = 128.0
    min = 1.0
$IFDEF VER40
    exponential = true
$ENDIF
    hint = "This parameter defines how soon an orbit bails out while \
            iterating. Larger values give smoother outlines; values around 4 \
            give more interesting shapes around the set. Values less than 4 \
            will distort the fractal."
  endparam
switch:
  type = "J233"
  power = power
  bailout = bailout
}
作者: 柳烟    时间: 2012-5-3 18:39

我迷惑的是:给了四个式子,究竟如何迭代?那里的v简直让人迷糊。
如何用几何画板实现这个分形?不一定非要论文中的3D先用GSP搞出普通的复分形来再说。
我把我UF中扫出的图片发在这:
未命名.jpg

图片附件: 未命名.jpg (2012-5-3 22:48, 93.35 KB) / 下载次数 2103
http://inrm3d.cn/attachment.php?aid=17179&k=3a3979950180eba64ba5bb8926a15e66&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-5-3 20:44

就是这个分形了(好象与UF中的phoenix极象,但一直没办法找到范例中的中心位置):
http://sprott.physics.wisc.edu/carlson/
作者: 柳烟    时间: 2012-5-3 21:50

我将v的初始值定为0,并让其参与迭代,与UF中的图有百分之九十九相当。些微差别不知何因。
叶氏范例探究.gsp (15.75 KB)
未命名.jpg

附件: 叶氏范例探究.gsp (2012-5-3 21:50, 15.75 KB) / 下载次数 3472
http://inrm3d.cn/attachment.php?aid=17182&k=dc0bc1ac10fd9d33ca5a531765b5cba6&t=1732387156&sid=xR1Q1q

图片附件: 未命名.jpg (2012-5-3 22:46, 58.74 KB) / 下载次数 2140
http://inrm3d.cn/attachment.php?aid=17183&k=4906569c120b5eef66d4183c46477db0&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-5-3 21:54

有点象但细看又觉得不象:

图片附件: phoenix.JPG (2012-5-3 21:54, 86.16 KB) / 下载次数 2004
http://inrm3d.cn/attachment.php?aid=17184&k=5d24130d32ad041ef5b5e7f5d46e8ed5&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-5-3 21:57

如果按叶氏论文中的方法着色会怎样呢?
作者: 柳烟    时间: 2012-5-3 22:04

改变初始值v将其定位在z=pixel处,局部放大一块,好看!
未命名.jpg

图片附件: 未命名.jpg (2012-5-3 22:47, 54.32 KB) / 下载次数 2080
http://inrm3d.cn/attachment.php?aid=17185&k=5d4d60887d69de97c921ca2c4d9a4425&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-5-3 22:36

6# 榕坚
调色板作色不只限于陷阱,按那论文的意思。
作者: 柳烟    时间: 2012-5-3 23:07

胡乱用调色板,设了两种颜色,照自已的意思整一幅,尝试尝试:
未命名.jpg

图片附件: 未命名.jpg (2012-5-6 18:45, 49.11 KB) / 下载次数 1861
http://inrm3d.cn/attachment.php?aid=17189&k=5fb810c7095f4c9f4a15d0cbd1468cd5&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-5-3 23:38

我晃眼看了一下论文中调色板一节中的(2),感觉用画板作问题不是特别大,明天再继续探索。论文中的k,好象相当于画板中的et,不知我这理解是否有理。大家一同想办法,众人捧柴火焰高呀。
作者: 榕坚    时间: 2012-5-4 07:45

好象有点象了,其实就是很早前的间隔et着色方法,但是没有颜色过渡(3D味道):

图片附件: 3D Phoenix Spirals.JPG (2012-5-4 07:45, 73.46 KB) / 下载次数 1969
http://inrm3d.cn/attachment.php?aid=17192&k=8643647e407efecb658ba5c91fc49b6f&t=1732387156&sid=xR1Q1q



图片附件: 3D Phoenix Spirals-5.JPG (2012-5-4 10:28, 105.06 KB) / 下载次数 2023
http://inrm3d.cn/attachment.php?aid=17197&k=34fab58667c36b93dfcec24c37006466&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-5-5 10:49

要使扫出的图有3D味,如何使这些粗尾巴的中间线特亮,然后向两边颜色暗淡下去呢?按前楼连接网址,找到的这图,看起来美:
cmphx03.gif

图片附件: cmphx03.gif (2012-5-5 10:49, 125.33 KB) / 下载次数 2058
http://inrm3d.cn/attachment.php?aid=17222&k=55893673025164727d62ea59a0a66bb9&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-5-5 17:06

12# 柳烟


中间线高亮倒也不是那么难办到,只要调整一下参数:

图片附件: 3D Phoenix Spirals-14.JPG (2012-5-5 17:06, 105.36 KB) / 下载次数 1977
http://inrm3d.cn/attachment.php?aid=17223&k=897037cfac60d9b65f64dea77a919bd2&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-5-5 17:42

目前,如何将背景色设置成想要的颜色,我弄了好一阵子,没办到。请问,这是如何做的?
作者: changxde    时间: 2012-5-5 18:51

13# 榕坚
榕老师这个定位在哪儿。我用叶文中的数据在UF中试试好像不对。
作者: 榕坚    时间: 2012-5-5 19:43

15# changxde


BackgroundLocation {
location:
  center=-0.0252539550257283315/1.1829770596358987 magn=14542012
}
作者: 榕坚    时间: 2012-5-5 20:23

14# 柳烟


判断函数啊。
作者: 柳烟    时间: 2012-5-5 22:02

终于解决了背景色问题:

图片附件: 未命名.jpg (2012-5-6 18:42, 77.07 KB) / 下载次数 1905
http://inrm3d.cn/attachment.php?aid=17225&k=5296e80f196a89d555c44fd65c1d4ad7&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-5-5 23:26

精典M集:
未命名.jpg

图片附件: 未命名.jpg (2012-5-6 18:40, 37.7 KB) / 下载次数 2128
http://inrm3d.cn/attachment.php?aid=17226&k=d7e963e3a0532d94ccaf0257cba6495b&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-5-6 16:15

没用调色板,按老办法扫一幅,自我感觉不错。
未命名.jpg

图片附件: 未命名.jpg (2012-5-6 18:41, 53.33 KB) / 下载次数 2319
http://inrm3d.cn/attachment.php?aid=17230&k=b5e7c806bde0af5e4be103a90e48601a&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-6-14 01:43

UF中的KochCurve代码,内含三个开关项,代码太长,解读困难,为降低难度,我将其代码精简为只整一个开关项,此精简后的代码在UF中运行后无误,我今晚试做了一下,结果不成功,扫出的图怪得很,大家帮忙看看,能否作成。榕坚兄有空时,看看,如何用GSP实现。
KochCurve {
init:
  z = #pixel
  x = real(z)
  y = imag(z)
  sq3 = sqrt(3)
  i = 0
loop:
  i = i + 1
  if i == 2
    arg = atan2(z)
      if (y + 1/sq3 > 0) && (sq3*x - y + 2/sq3 > 0) \
         && (sq3*x + y - 2/sq3 < 0)
        bail = true
      endif
      if (-y + 1/sq3 > 0) && (sq3*x + y + 2/sq3 > 0) \
         && (sq3*x - y - 2/sq3 < 0)
        bail = true
      endif

      if (arg > 5/6*pi) || (arg < -pi/2)
        z = z*exp(1i*4/3*pi)
      endif
      if (arg < pi/6) && (arg > -pi/2)
        z = z*exp(1i*2/3*pi)
      endif
      z = z - 1i*1/sq3
      
  elseif i > 2
    z = 3*z
    x = real(z)
    y = imag(z)
     if (y > 0) && (sq3*x - y + sq3 > 0) \
       && (sq3*x + y - sq3 < 0)
      bail2 = true
    endif
    z = z/3
    x = real(z)
    y = imag(z)
    if x < -1/3
      z = 3*z + 2
    elseif x > 1/3
      z = 3*z - 2
    else
      if x < 0
        z = z + 1/3
        z = z*exp(-1i*pi/3)
        z = 3*z - 1
      else
        z = z - 1/3
        z = z*exp(1i*pi/3)
        z = 3*z + 1
      endif
    endif
  endif
bailout:
  bail == false && bail2 == false
default:
  title = "Koch Curve1"
  helpfile = "sam-help/kochcurves.htm"
  helptopic = "kcurve"
  magn = 1.5
  center = (0.0002,0)
  maxiter = 50
}
作者: 柳烟    时间: 2012-6-14 01:48

这是UF中的效果:
未命名.jpg

图片附件: 未命名.jpg (2012-6-14 01:48, 16.69 KB) / 下载次数 1736
http://inrm3d.cn/attachment.php?aid=17693&k=6c375e0798c5a84caf43db3b78fe401b&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-6-14 21:00

这个分形的判断函数太多了,转来转去头都晕了。
作者: 柳烟    时间: 2012-6-14 21:31

23# 榕坚
问题丢在这,有兴趣时可弄弄。我昨晚算了一大屏数字,确如你所说,判断繁多弄昏了头,扫出的图十万八千里。又下面这个分形位于sam.ufm中的Sierpinski Triangle II,代码似比上楼的轻松了许多,整了一天,也整来不对劲。我也打算等脑壳清醒点,再来看看。
未命名.jpg

图片附件: 未命名.jpg (2012-6-14 21:31, 16.12 KB) / 下载次数 1742
http://inrm3d.cn/attachment.php?aid=17703&k=5ac858ee3109d2ac0ae5e4ecb783fe6e&t=1732387156&sid=xR1Q1q


作者: zwh2010    时间: 2012-6-14 21:41

这玩意好像很久很久以前弄过类似的,似乎原理差不多。
作者: 榕坚    时间: 2012-6-14 21:49

25# zwh2010


但柳老师想解释UF代码。
作者: 榕坚    时间: 2012-6-15 08:07

柳老师能否把reb系列中的indra's promise这个分形的每个开关项给分离出来,这个分形的代码太长了,因为他集成了很多个我们之前讨论的各种分形图:

图片附件: Fractal1.jpg (2012-6-15 08:07, 69.31 KB) / 下载次数 1557
http://inrm3d.cn/attachment.php?aid=17706&k=8693eed8111453eedff0bcc17febbd51&t=1732387156&sid=xR1Q1q



图片附件: Fractal2.jpg (2012-6-15 08:07, 64.12 KB) / 下载次数 1570
http://inrm3d.cn/attachment.php?aid=17707&k=778d91704eb2c30dda7b32defdb74b2b&t=1732387156&sid=xR1Q1q



图片附件: Fractal3.jpg (2012-6-15 08:07, 39.94 KB) / 下载次数 1412
http://inrm3d.cn/attachment.php?aid=17708&k=4c1d2e1d6e26656392578aa0d8afcbee&t=1732387156&sid=xR1Q1q



图片附件: Fractal4.jpg (2012-6-15 08:07, 26.64 KB) / 下载次数 1411
http://inrm3d.cn/attachment.php?aid=17709&k=2faf1b6d4d49b8881c693ca5e346a73a&t=1732387156&sid=xR1Q1q



图片附件: Fractal5.jpg (2012-6-15 08:07, 31.01 KB) / 下载次数 1430
http://inrm3d.cn/attachment.php?aid=17710&k=1895258fbf64141302e2eb9b0f2f9865&t=1732387156&sid=xR1Q1q



图片附件: Fractal6.jpg (2012-6-15 08:07, 47.84 KB) / 下载次数 1416
http://inrm3d.cn/attachment.php?aid=17711&k=543c4d3b596749040ff772acb22d6e0d&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-6-15 10:41

24# 柳烟


其中这段代码是多余的:
if m1*x + m1 < y || -m1*x + m1 < y \
       || y < 0
    endif
作者: 柳烟    时间: 2012-6-15 11:37

27# 榕坚
应该能分离,我今晚空了,办这个事。
作者: 柳烟    时间: 2012-6-15 13:35

28# 榕坚
按你的提示,将这段代码删掉后,UF中仍然显示此图,于是用GSP造,干了三四个钟点,不对,找原因,修改,往复了好几回,终于成了,扫出的图比UF中的图还要美!谢谢。
未命名.jpg
未命名.jpg
Sierpinski Triangle II.gsp (17.5 KB)

图片附件: 未命名.jpg (2012-6-15 14:01, 93.67 KB) / 下载次数 1813
http://inrm3d.cn/attachment.php?aid=17714&k=d082c2b1d40fc21ab11c6c353d0b15ed&t=1732387156&sid=xR1Q1q



图片附件: 未命名.jpg (2012-6-15 14:12, 102.87 KB) / 下载次数 1806
http://inrm3d.cn/attachment.php?aid=17715&k=698bde9b2531cd7e2abe189f67f580da&t=1732387156&sid=xR1Q1q



附件: Sierpinski Triangle II.gsp (2012-6-15 14:16, 17.5 KB) / 下载次数 3690
http://inrm3d.cn/attachment.php?aid=17716&k=726e747ca3653c089b8a89d35498211b&t=1732387156&sid=xR1Q1q
作者: 榕坚    时间: 2012-6-15 16:00

30# 柳烟


确实比之前的哪种做法好看多了,这个好象三角形内部所有点的et值都一样。
作者: 柳烟    时间: 2012-6-15 17:53

请给出27楼各图片所用的开关项及参数,这样好去除代码中的枝丫部分。
作者: 榕坚    时间: 2012-6-15 18:03

我看就先拿这个开刀吧,不过即使单独开关项我估计也是有难度的:

附件: Fractal1.rar (2012-6-15 18:03, 195.79 KB) / 下载次数 2010
http://inrm3d.cn/attachment.php?aid=17717&k=b93497ba24e25d07da942bc81e580144&t=1732387156&sid=xR1Q1q

图片附件: Fractal1.jpg (2012-6-15 18:03, 62.88 KB) / 下载次数 1603
http://inrm3d.cn/attachment.php?aid=17718&k=02d745bd0e9be0cd186662801f7df4c7&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-6-15 21:02

解压后,双击打开,图与上楼一样,但我按右边面板上的参数,打开重新打开 UF中的该范例文件,按你面板上的开关项,整了几次,怎么没有图形呢?
未命名.jpg
从你的文件来看,好象动了图片中的打钩部分与画圈部分,但我把reb系列中的文件调出按这弄,没图形,啥原因?

图片附件: 未命名.jpg (2012-6-15 21:20, 15.82 KB) / 下载次数 1540
http://inrm3d.cn/attachment.php?aid=17720&k=fa7ca6ce12222c3c1d7a450475826004&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-6-15 21:53

是这个参数,可能要等一会儿才会出图形,注意下方的进度条。或者检查一下中心与缩放倍数:(0,0),缩放1倍。
作者: 柳烟    时间: 2012-6-16 18:28

将前面的KochCurve曲线代码,删除了一段没用的代码,再浓缩成下面的代码:
KochCurve {
init:
  z = #pixel
  x = real(z)
  y = imag(z)
  sq3 = sqrt(3)
  bail =  false
  bail2 =false
  i = 1
loop:
  i = i + 1
  if i == 2
    arg = atan2(z)
      if (-y + 1/sq3 > 0) && (sq3*x + y + 2/sq3 > 0) \
         && (sq3*x - y - 2/sq3 < 0)
        bail = true
      endif

      if (arg > 5/6*pi) || (arg < -pi/2)
        z = z*exp(1i*4/3*pi)
      endif
      if (arg < pi/6) && (arg > -pi/2)
        z = z*exp(1i*2/3*pi)
      endif
      z = z - 1i*1/sq3

  elseif i > 2
    z = 3*z
    x = real(z)
    y = imag(z)
     if (y > 0) && (sq3*x - y + sq3 > 0) \
       && (sq3*x + y - sq3 < 0)
      bail2 = true
    endif
    z = z/3
    x = real(z)

    if x < -1/3
      z = 3*z + 2
    elseif x > 1/3
      z = 3*z - 2
    else
if x < 0
        z = z + 1/3
        z = z*exp(-1i*pi/3)
        z = 3*z - 1
      elseif x> 0
        z = z - 1/3
        z = z*exp(1i*pi/3)
        z = 3*z + 1
      endif
    endif
  endif
bailout:
  bail == false && bail2 == false
default:
  title = "Koch Curve1"
  helpfile = "sam-help/kochcurves.htm"
  helptopic = "kcurve"
  magn = 1.5
  center = (0.0002,0)
  maxiter = 50
}
作了几次,每次都不一样,查起来费劲,每干一次,扫出的图均不一样,怪。不过,倒是取得了进展,扫出的图有些地方可看到局部的科赫曲线影子。不知是不是我解读上面的代码,有误。
未命名.jpg

图片附件: 未命名.jpg (2012-6-16 19:58, 12.62 KB) / 下载次数 1577
http://inrm3d.cn/attachment.php?aid=17723&k=c77584f594ed814cf1899fd154d696ee&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-6-16 18:46

if x < -1/3
      z = 3*z + 2
    elseif x > 1/3
      z = 3*z - 2
    else
if x < 0
        z = z + 1/3
        z = z*exp(-1i*pi/3)
        z = 3*z - 1
      elseif x> 0
        z = z - 1/3
        z = z*exp(1i*pi/3)
        z = 3*z + 1
这段代码,理解起来有些费劲,可能我的解读有偏差。
作者: 柳烟    时间: 2012-6-17 23:51

未命名.jpg
未命名.jpg
干了几天,终于有了成果,源文件整理好后发,因算式多而混乱。

图片附件: 未命名.jpg (2012-6-18 09:24, 20.15 KB) / 下载次数 1738
http://inrm3d.cn/attachment.php?aid=17730&k=545a53de7c98b6d385976adbd02f682b&t=1732387156&sid=xR1Q1q



图片附件: 未命名.jpg (2012-6-18 09:23, 11.68 KB) / 下载次数 1742
http://inrm3d.cn/attachment.php?aid=17731&k=cb1ace370001863ab30629c4c72b9e52&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-6-18 08:20

38# 柳烟


漂亮,功夫不负有心人哪。
作者: 柳烟    时间: 2012-6-18 09:04

39# 榕坚
谢过榕兄支持并鼓励,并帮助。
文件整理好了,发于此,希提出批评并改进意见,使我进步。
Koch Curve(1).gsp (29.35 KB)
未命名.jpg

附件: Koch Curve(1).gsp (2012-6-18 09:04, 29.35 KB) / 下载次数 2845
http://inrm3d.cn/attachment.php?aid=17732&k=4ea26bc4214d5e22f67135dfab3b1ae0&t=1732387156&sid=xR1Q1q

图片附件: 未命名.jpg (2012-6-18 09:22, 16.96 KB) / 下载次数 1752
http://inrm3d.cn/attachment.php?aid=17733&k=6ce850725620395ce49bc78dcbefe68a&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-6-18 14:40

UF中科赫雪花的第二个开关项抽出后的代码如下:
KochCurve {
init:
  z = #pixel
  x = real(z)
  y = imag(z)
  sq3 = sqrt(3)
  bool bail2 = false
  bool bail = false
  i = 0
loop:
  i = i + 1
  if i == 2
    arg = atan2(z)
      if |x| > 1 || x/sq3 + y - 2*sq3/3 > 0 \
         || x/sq3 + y + 2*sq3/3 < 0 || x/sq3 \
         - y + 2*sq3/3 < 0 || x/sq3 - y - \
         2*sq3/3 > 0
        bail = true
      endif
      if (abs(arg) < pi/6)
        z = z*exp(-1i*pi/2)
      elseif (arg > pi/6) && (arg < pi/2)
        z = z*exp(-1i*5*pi/6)
      elseif (arg > pi/2) && (arg < 5*pi/6)
        z = z*exp(1i*5*pi/6)
      elseif (abs(arg) > 5*pi/6)
        z = z*exp(1i*pi/2)
      elseif (arg < -pi/6) && (arg > -pi/2)
        z = z*exp(-1i*pi/6)
      elseif (arg < -pi/2) && (arg > -5*pi/6)
        z = z*exp(1i*pi/6)
      endif
      z = z + 1i
      z = sq3*z

  elseif i > 2
    z = 3*z
    x = real(z)
    y = imag(z)
    if (y > 0) && (sq3*x - y + sq3 > 0) \
       && (sq3*x + y - sq3 < 0)
      bail2 = true
    endif
    z = z/3
    x = real(z)
    y = imag(z)
    if x < -1/3
      z = 3*z + 2
    elseif x > 1/3
      z = 3*z - 2
    elseif (x < 0)&&(x>-1/3)

        z = z + 1/3
        z = z*exp(-1i*pi/3)
        z = 3*z - 1
    elseif (x > 0)&&(x<1/3)
        z = z - 1/3
        z = z*exp(1i*pi/3)
        z = 3*z + 1
    endif
   endif

bailout:
  bail == false && bail2 == false

default:
  title = "Koch Curve"
  helpfile = "sam-help/kochcurves.htm"
  helptopic = "kcurve"
  magn = 1.5
  center = (0.0002,0)
  maxiter = 50

}
第一个开关项的科氏雪花作出后,这第二个开关项,是顺理成章的事,估计难度不是很大了。
作者: xiaongxp    时间: 2012-6-18 15:04

能将UF程序解读并用画板实现,是要有智慧和毅力的,非常钦佩柳老师的能力和钻劲。
作者: 榕坚    时间: 2012-6-18 15:47

38#的第二图不是第二个开关项吗?UF好象是三个开关项。
作者: 柳烟    时间: 2012-6-18 15:50

43# 榕坚
UF中共有三个开关项,已发表的几张科赫雪花图片均是第一个开关项,第二个开关项代码41楼,GSP版的,还未生产出。
作者: 榕坚    时间: 2012-6-18 16:00

有一个想法:如何象M集那样提取出它的边界?
作者: 柳烟    时间: 2012-6-18 21:20

45# 榕坚
扫描法搞边界,怕有些困难,MJ集边界成熟算法dem法,涉及导数,这科赫雪花好象找不到一个复函数,这导数没法办。
作者: 榕坚    时间: 2012-6-18 22:00

46# 柳烟


应该还是可以办到的,不一定要DEM法。
作者: 柳烟    时间: 2012-6-18 22:22

这第二个开关项,整出的不对劲,检查了一遍,没发现错误,算式多,眼睛都看花了,干脆重新整一遍,结果对头了。
作者: 榕坚    时间: 2012-6-19 07:38

Koch线,线条不是很清晰:

图片附件: K.JPG (2012-6-19 07:38, 31.72 KB) / 下载次数 1972
http://inrm3d.cn/attachment.php?aid=17738&k=9a6532b3dcac532eeb760b827d788ed7&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-6-19 08:08

49# 榕坚
作得不错,这是用的什么办法弄出的边界呢
作者: 柳烟    时间: 2012-6-19 08:10

UF中的效果图是这样的:
Fractal2.jpg
用GSP扫一幅:
未命名.jpg
再扫一幅别样味道的:
未命名2.jpg
Koch Curve(二).gsp (35.88 KB)

图片附件: 未命名.jpg (2012-6-19 09:03, 14.64 KB) / 下载次数 2153
http://inrm3d.cn/attachment.php?aid=17739&k=85013661ddf195f0ec276a38496150b4&t=1732387156&sid=xR1Q1q



图片附件: Fractal2.jpg (2012-6-19 09:03, 22.02 KB) / 下载次数 2138
http://inrm3d.cn/attachment.php?aid=17740&k=dfa74f27dfaecaeae166ef705d34e973&t=1732387156&sid=xR1Q1q



图片附件: 未命名2.jpg (2012-6-19 09:03, 30.54 KB) / 下载次数 2148
http://inrm3d.cn/attachment.php?aid=17741&k=4953df380558e2f386f551ab6565c9a0&t=1732387156&sid=xR1Q1q



附件: Koch Curve(二).gsp (2012-6-19 08:49, 35.88 KB) / 下载次数 3492
http://inrm3d.cn/attachment.php?aid=17742&k=b20b851c0c284b706b5c8a7c6f57baca&t=1732387156&sid=xR1Q1q
作者: 榕坚    时间: 2012-6-19 08:42

50# 柳烟


最早做等势线的办法。
作者: 柳烟    时间: 2012-6-19 08:47

52# 榕坚
谢谢!
作者: 柳烟    时间: 2012-6-19 10:34

科赫雪花开关项3的代码如下:
KochCurve {
init:
  z = #pixel
  zz=0
  x = real(z)
  y = imag(z)
  sq3 = sqrt(3)
  bool bail2 = false
  bool bail = false
  i = 2
loop:
  i = i + 1
  if i == 2
    arg = atan2(z)
      if (-y + 1/sq3 > 0) && (sq3*x + y + 2/sq3 > 0) \
         && (sq3*x - y - 2/sq3 < 0)
        bail = true
      endif
      if (arg > 5/6*pi) || (arg < -pi/2)
        z = z*exp(1i*4/3*pi)
      endif
      if (arg < pi/6) && (arg > -pi/2)
        z = z*exp(1i*2/3*pi)
      endif
      z = z - 1i*1/sq3

        zz = z
        z = #pixel
      if |x| > 1 || x/sq3 + y - 2*sq3/3 > 0 \
         || x/sq3 + y + 2*sq3/3 < 0 || x/sq3 \
         - y + 2*sq3/3 < 0 || x/sq3 - y - \
         2*sq3/3 > 0
        bail = true
      endif
      if (abs(arg) < pi/6)
        z = z*exp(-1i*pi/2)
      elseif (arg > pi/6) && (arg < pi/2)
        z = z*exp(-1i*5*pi/6)
      elseif (arg > pi/2) && (arg < 5*pi/6)
        z = z*exp(1i*5*pi/6)
      elseif (abs(arg) > 5*pi/6)
        z = z*exp(1i*pi/2)
      elseif (arg < -pi/6) && (arg > -pi/2)
        z = z*exp(-1i*pi/6)
      elseif (arg < -pi/2) && (arg > -5*pi/6)
        z = z*exp(1i*pi/6)
      endif
      z = z + 1i
      z = sq3*z
  elseif i > 2
      oldz = z
      z = zz
      zz = oldz
    z = 3*z
    x = real(z)
    y = imag(z)
    if (y > 0) && (sq3*x - y + sq3 > 0) \
       && (sq3*x + y - sq3 < 0)
      bail2 = true
    endif
    z = z/3
    x = real(z)
    y = imag(z)
    if x < -1/3
      z = 3*z + 2
    elseif x > 1/3
      z = 3*z - 2
    else
      if x < 0
        z = z + 1/3
        z = z*exp(-1i*pi/3)
        z = 3*z - 1
      else
        z = z - 1/3
        z = z*exp(1i*pi/3)
        z = 3*z + 1
      endif
    endif
  endif
bailout:
  bail == false && bail2 == false
default:
  title = "Koch Curve"
  helpfile = "sam-help/kochcurves.htm"
  helptopic = "kcurve"
  magn = 1.5
  center = (0.0002,0)
  maxiter = 50
}
好象是综合前两个开关项,此开关项算起来更繁,更具挑战性,更要有十足耐心与高精度的细致才行。
作者: 柳烟    时间: 2012-6-19 23:05

科赫雪花开关项3,计算量特大,我硬着头皮弄完后,扫出的图不对劲,数据太多了,查都没法查起,等脑子清醒占再弄,也许要放弃了。此分形,我算是服了。
作者: 柳烟    时间: 2012-6-20 23:54

此谢儿宾斯基三角形,原来我与榕老师研究过,榕老师最先造出。原代码太长,且繁杂不堪。今对原代码大势减肥,浓缩成下面代码,再用GSP作之,十分容易,好象效果图没啥两样,大家可对照对照,看看有无问题。
gnd-slope-sierpinski2 {
init:
  z1 = #pixel
  int done = 2
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
  modz = |z1|
  done = done + 1

  IF (modz > @bailout)
    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

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

}
作者: 柳烟    时间: 2012-6-20 23:55

未命名.jpg
Slope Sierpinski Triangle II.gsp (17.53 KB)

图片附件: 未命名.jpg (2012-6-20 23:55, 87.54 KB) / 下载次数 2316
http://inrm3d.cn/attachment.php?aid=17755&k=484bf79da6a7e6e9fa43406c47540009&t=1732387156&sid=xR1Q1q



附件: Slope Sierpinski Triangle II.gsp (2012-6-20 23:58, 17.53 KB) / 下载次数 3274
http://inrm3d.cn/attachment.php?aid=17756&k=2df3615b925b8b2f0956a8f887d23858&t=1732387156&sid=xR1Q1q
作者: 柳烟    时间: 2012-6-21 17:21

pgd_Twofold2_Mandel {
  z = 1
  bool notdone = true
  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

  ENDIF
  z_values[i] = z
  IF (odd)
    IF (|z - z_values[j]| < b)
      notdone = false
      complex w = -1
      j = 0

      WHILE ((j < #maxit) )
        w = (sqr(w) - 1/3)/(#pixel * sqr(w) * w)

        j = j + 1
      ENDWHILE

    ENDIF
    j = j + 1
  ENDIF
  odd = !odd
  i = i + 1
bailout:
  notdone
default:
  title = "Twofold2 Mandelbrot"
  periodicity = 0
  maxiter = 1000
  method = onepass

  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
}
这是一个分形的浓缩代码,通过UF验证,此代码与原长代码的该开关项,图形差别甚微(不知究竟差别在那,大家可帮我验证验证),因本人将原代码删掉太多,难免要走样。此分形代码中:z_values[@zmax]
是什么意思?这句搞懂后,分形易造。
作者: 榕坚    时间: 2012-6-21 20:06

58# 柳烟


这个是UF中的数组变量,要在几何画板中实现应该比较难。就看能否用其它的方法替代,就象判断语句用sgn替代那样。而且这还是双重循环的,难。
下面红色部分你的代码中有错:
loop:
  z = (sqr(z) - 1/3)/(#pixel * sqr(z) * z)
  IF (|z| > @bailout)
    notdone = false

  ENDIF
z_values = z
IF (odd)
    IF (|z - z_values[j]| < b)
      notdone = false
      complex w = -1
      j = 0
应该是是z_values[I] = z,那个I为小写。
作者: 柳烟    时间: 2012-6-21 20:31

59# 榕坚
对的,因为加了颜色后,没有小i了,现已更正。谢谢,既然如此,摆在这,看有没有不世出的高人弄出。
作者: 柳烟    时间: 2012-6-23 09:48

UF中的另一特效:
作色程序:
TangentCircles {
; Mark Townsend, 29 Apr 1999
;
; ****** NOTE *******
; Instead of using this coloring method use
; the Tangent Balls mode of the Paul Carlson's
; Orbit Traps method in pwc.ucl.
; *******************
;
; Paul W. Carlson's Tangent Circles oribit trap
; method.
;
init:
  int iter = 0
  float x = 0
  float y = 0
  float Xabs = 0
  float Yabs = 0
  float Dsgd0 = 0
  float Dsgd1 = 0
  float Dsgd2 = 0
  float Circle = 0
  float ZtoPsqd = 0
  float Rc = @r
  float Phi = #pi * 0.125
  float Rm = Rc / sin(Phi)
  float RcSqd = Rc^2
  float Py = Rm * sin(2 * Phi)
  float Px = Rm * cos(2 * Phi)
  bool Trapped = false
loop:
  iter = iter + 1
  x = real(#z)
  y = imag(#z)
  if abs((cabs(#z) - Rm)) < Rc && iter > @skip && !Trapped
    Xabs = abs(x)
    Yabs = abs(Y)
    Dsgd0 = Xabs^2 + (Yabs - Rm)^2
    Dsgd1 = (Xabs - Px)^2 + (Yabs - Py)^2
    Dsgd2 = (Xabs - Rm)^2 + Yabs^2
    if Dsgd0 < RcSqd
      Trapped = true
      ZtoPsqd= Dsgd0
      if y > 0
        Circle = 0
      else
        Circle = 4
      endif
    elseif Dsgd1 < RcSqd
      Trapped = true
      ZtoPsqd= Dsgd1
      if y >  0&& x > 0
        Circle = 1
      elseif y < 0 && x > 0
        Circle = 3
      elseif y < 0 && x < 0
        Circle = 5
      else
        Circle = 7
      endif
    elseif Dsgd2 < RcSqd
      trapped = true
      ZtoPsqd = Dsgd2
      if x > 0
        Circle = 2
      else
        Circle = 6
      endif
   endif
endif
final:
  if !trapped
    #solid = true
  else
    Circle = (Circle + @off) % 8
    float Ratio = sqrt(ZtoPsqd/Rcsqd)
    float ColorIndex = 29 * Ratio + Circle * 30
    #index = (ColorIndex + 1) % 256 /256
  endif
default:
  title = "Carlson Tangent Circles"
  param r
    caption = "Circle radius"
    default = 0.2
  endparam
  param skip
    caption = "Iters to skip"
    default = 1
    hint = "Iterations to skip."
  endparam
  param off
    caption = "Hue cycle"
    default = 0
    min = 0
    max = 7
    hint = "This rotates the coloring order of the balls."
  endparam
}
感觉到程序不是太难,从UF的效果来看,类似作过的圆陷阱,我正在用画板整,感觉到可能要失败。大家想想办法。
Fractal1.jpg

图片附件: Fractal1.jpg (2012-6-23 09:48, 137.06 KB) / 下载次数 1872
http://inrm3d.cn/attachment.php?aid=17765&k=d84ceb9f7a8585bbfab1ce3c067424d3&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-6-23 15:13

61# 柳烟


UF代码中的#solid 是什么意思一直摸不透。
作者: 柳烟    时间: 2012-6-23 18:39

62# 榕坚
此句好象是调除M集肚内与众多圆陷阱外的颜色用的,我试了,此值若为true,则为黑,若为false,则为默认的天蓝色。我将
if !trapped
    #solid = true
  else
    Circle = (Circle + @off) % 8
    float Ratio = sqrt(ZtoPsqd/Rcsqd)
    float ColorIndex = 29 * Ratio + Circle * 30
    #index = (ColorIndex + 1) % 256 /256
中if后的!号去掉,并将#solid = true
  else
拿掉,变成这样:
if trapped

        Circle = (Circle + @off) % 8
    float Ratio = sqrt(ZtoPsqd/Rcsqd)
    float ColorIndex = 29 * Ratio + Circle * 30
    #index = (ColorIndex + 1) % 256 /256
在UF中试,结果影响不大。另外,将#solid改为任一变量,如a等等,不影响。
作者: 榕坚    时间: 2012-6-23 20:48

63# 柳烟


好象是如果trapped=1,按#index着色。否则调用UF的内部着色机制着色。也就是说落入陷井的点按#index着色。是这样理解的吗?
作者: 柳烟    时间: 2012-6-23 22:27

64# 榕坚
应该是这样理解的。我刚才乱整一通,整出了几个圈圈,有些相象。可能我的迭代有问题,好象是将着色程序中的x,y,也就是将M集的#z,可能是z0,替代x,y,让其与M集一同迭代。这我还没试。刚才乱整出的圈圈预言,这能用画板实现。
未命名.JPG
另外,从UF中看,那些圆圈被黑M集压住,应该是用dxy老师的M集外部q法用#index作色。

图片附件: 未命名.JPG (2012-6-23 22:27, 14.68 KB) / 下载次数 1674
http://inrm3d.cn/attachment.php?aid=17777&k=8cbedbcc9e556e0efbe316250401e247&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-6-23 23:47

好象我已经做出来了,正在扫图:
未命名.JPG
扫了一个黑白的,白天再整一个彩色的。

图片附件: 未命名.JPG (2012-6-24 00:51, 70.78 KB) / 下载次数 1925
http://inrm3d.cn/attachment.php?aid=17778&k=3f39444d38d542f1685b2d8f9c76057f&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-6-24 06:49

调一彩色的。
未命名.jpg
未命名.jpg
Carlson Tangent Circles下M集.gsp (30.13 KB)

图片附件: 未命名.jpg (2012-6-24 06:49, 128.02 KB) / 下载次数 1815
http://inrm3d.cn/attachment.php?aid=17779&k=39f6b8ff4d28547cd738fe9402627f96&t=1732387156&sid=xR1Q1q



图片附件: 未命名.jpg (2012-6-24 06:58, 149.99 KB) / 下载次数 1801
http://inrm3d.cn/attachment.php?aid=17780&k=fd0b47fdbe4f567547a5d5e7d508b998&t=1732387156&sid=xR1Q1q



附件: Carlson Tangent Circles下M集.gsp (2012-6-25 20:46, 30.13 KB) / 下载次数 3330
http://inrm3d.cn/attachment.php?aid=17792&k=76fd98744abbe06c88847053d269a1b6&t=1732387156&sid=xR1Q1q
作者: 榕坚    时间: 2012-6-24 07:30

其实就是八圆陷井了,不同的是分别标注着色。

图片附件: TangentCircles-7.JPG (2012-6-24 10:28, 78.02 KB) / 下载次数 1455
http://inrm3d.cn/attachment.php?aid=17781&k=0da30140314722dde72d8a94d525b771&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-6-24 10:29

67# 柳烟


如何让每个圆的中心透亮使其更有立方体感?
作者: 柳烟    时间: 2012-6-24 14:02

69# 榕坚
可能有些麻烦,好象与原来造陷阱的算法不一样,不知有何更好的办法。
作者: 柳烟    时间: 2012-6-25 13:01

bezier-curve {
init:
  float x0=real(@z0)
  float y0=imag(@z0)
  float x1=real(@z1)
  float y1=imag(@z1)
  float x2=real(@z2)
  float y2=imag(@z2)
  float x3=real(@z3)
  float y3=imag(@z3)
  float cx=3*(x1-x0)
  float bx=3*(x2-x1)-cx
  float ax=x3-x0-cx-bx
  float cy=3*(y1-y0)
  float by=3*(y2-y1)-cy
  float ay=y3-y0-cy-by

  float r=0.0
  float x=0.0
  float y=0.0
  float u=0.0
  float v=0.0
  float t=0.0
  float rmin=1.0e20
loop:
  u=real(#z)
  v=imag(#z)
   t=0.0
  while(t<=1.0)
    x=((ax*t+bx)*t+cx)*t+x0
    y=((ay*t+by)*t+cy)*t+y0
    r=(x-u)*(x-u)+(y-v)*(y-v)
    if(r<rmin)
      rmin=r
    endif
    t=t+@dt
  endwhile
final:
    #index=rmin^@nexp
default:
  title="Bezier Curve"
  helpfile="lkm-help\lkm-bezier.html"
  param z0
    caption="1st anchor point"
    default=(1.0,0.0)
    hint="Curve starts at this point."
  endparam
  param z1
    caption="1st control point"
    default=(1.0,1.0)
    hint="Influences the shape of the curve."
  endparam
  param z2
    caption="2nd control point"
    default=(0.0,0.0)
    hint="Influences the shape of the curve."
  endparam
  param z3
    caption="2nd anchor point"
    default=(0.0,1.0)
    hint="Curve ends at this point."
  endparam
  param dt
    caption="step size"
    default=0.1
    hint="Decrease for smoother line, increase \
      to see dots.  Should be between 0 & 1."
    min=0.0
    max=1.0
  endparam
  param nexp
    caption="power"
    default=0.1
    min=0.0
    hint="Decrease to make thinner lines. Use \
      with 'minimum distance' coloring."
  endparam

}
此代码,经UF中验证,没问题。此特效一直没整成功,今天去除另三个不很美的开关项,打算重新整这特效,力争整出原样。我干了一遍,结果不对。打算换一个方法来弄。
作者: 柳烟    时间: 2012-6-25 16:28

越做越不象样,不知这些值:
float r=0.0
  float x=0.0
  float y=0.0
  float u=0.0
  float v=0.0
  float t=0.0
参不参予M集的迭代?
作者: 柳烟    时间: 2012-6-25 21:48

未命名.jpg
未命名.jpg

图片附件: 未命名.jpg (2012-6-27 07:50, 30.41 KB) / 下载次数 1498
http://inrm3d.cn/attachment.php?aid=17794&k=73cc3e7ff72cc2bac744c3ab158fd75a&t=1732387156&sid=xR1Q1q



图片附件: 未命名.jpg (2012-6-27 07:49, 29.18 KB) / 下载次数 1444
http://inrm3d.cn/attachment.php?aid=17795&k=088e55e194fd956bac0d4a46c1a398aa&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-6-27 00:18

未命名.jpg
将上一楼的第一个图片中的八切圆往下放一点,好看。
未命名.jpg

图片附件: 未命名.jpg (2012-6-27 07:50, 30.81 KB) / 下载次数 1483
http://inrm3d.cn/attachment.php?aid=17805&k=f262b660878555a13eace21b30b71b6b&t=1732387156&sid=xR1Q1q



图片附件: 未命名.jpg (2012-6-27 07:51, 30.28 KB) / 下载次数 1481
http://inrm3d.cn/attachment.php?aid=17806&k=d2830238eddbc5627ef3d3a65a4a4f23&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-6-27 07:43

74# 柳烟


如果把圆陷井改为M集陷井就可能在谢氏三角形内部做M集了。
作者: 柳烟    时间: 2012-6-27 08:34

75# 榕坚
M集作陷阱,一直找不到方向。老兄给个实例。
作者: 榕坚    时间: 2012-6-27 08:42

76# 柳烟


常老师的贴子:http://www.inrm3d.cn/viewthread. ... age%3D2&page=24
UF中也有范例:julia trap:

图片附件: Fractal1.jpg (2012-6-28 07:56, 32.98 KB) / 下载次数 1216
http://inrm3d.cn/attachment.php?aid=17827&k=42f06031c875ba9d8199c2dc96288210&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-6-28 19:51

77# 榕坚
谢谢,常老师的帖子,与圆陷阱等,作法相类,能看懂。感觉到与M集作陷阱,好象要进入另一片天。我今天试了试,想弄一个J集的M陷阱,有点胡拼海整,结果纯色一片。问题是,如何将J集迭代的终点x 坐标与y坐标是M集上的点或者不是M集上的点,?我在这点上老是迷不穿。我今天将这x 坐标与y坐标代入P=sgn(0.5-spn(x^2+y^2-baiout)),若p=1,说明这终点落入M集内吗?好象不对,好象是落入以原点为心的圆的陷阱内,所以越整越不知从何着手,好象这种陷阱不及图形陷阱好处理。
作者: 柳烟    时间: 2012-6-28 20:07

77# 榕坚
此楼特效不错,开关项多,代码长,我看能否将各开关项抽出来,再用画板搞搞。
怎么我改变开关项,及参数,不能得到你那帖图呢
作者: 榕坚    时间: 2012-6-28 22:03

79# 柳烟


中心选(1,1)点,陷井大小选0.5,旋转角为-30
作者: 柳烟    时间: 2012-6-28 23:45

此分形漂亮!我在UF中取了两张,好看。
qqq.jpg
未命名.jpg
我先将这开关项的代码抽出来,看能否用画板演绎。

图片附件: qqq.jpg (2012-6-29 00:39, 41.04 KB) / 下载次数 1587
http://inrm3d.cn/attachment.php?aid=17828&k=dd1cb5cd3a6c3662fa40360e0552ff17&t=1732387156&sid=xR1Q1q



图片附件: 未命名.jpg (2012-6-29 00:39, 98.45 KB) / 下载次数 1758
http://inrm3d.cn/attachment.php?aid=17830&k=8b82459ed259758c76cefeb07c85bb0c&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-6-29 00:42

juliatrap(BOTH) {
; By Samuel Monnier, 1999
init:
  bool first = true
  count = 0
  float stdz = 0.0
  start = @seed

loop:
  trz = exp(flip(-pi/180*@rot))*(#z-@center)/@size
  x = 1/sqrt(@r)*real(trz)
  y = sqrt(@r)*imag(trz)
  if @freq != 0
    x = 2*sin(2*pi*@freq*x)
    y = 2*sin(2*pi*@freq*y)
  endif
  trz = x + flip(y)
  zz = trz
  float i = 0
  iter = @niter
  float dz = 0.0

    start = zz
    zz = @seed


  while i < @niter
    i = i + 1
    zz = zz^2 + start
    if |zz| > 1e20
      iter = i
      i = @niter
    endif

  endwhile
  float logp = 1/log(2)
  float logb = log(log(1e20))

    dz = real(iter + logp*logb - logp*log(log(cabs(zz))))


    if first
      stdz = dz
      first = false
    elseif dz > stdz
      stdz = dz
    endif

final:
  #index = .1*stdz
default:
  title = "Julia Trap"
  helpfile = "sam-help/juliatrap.htm"

  param seed
    caption = "Seed (for 'custom')"
    default = (0,0)
  endparam
  param center
    caption = "Center"
    default = (1,1)
  endparam
  param rot
    caption = "Rotation"
    default = -30.0
  endparam
  param size
    caption = "Size"
    default = 0.38750
  endparam
  param r
    caption = "Ratio Width/Heigh"
    default = 1.0
  endparam
  param @niter
    caption = "Julia Iterations"
    default = 100
  endparam
  param freq
    caption = "Trap Frequency"
    default = 0.0
  endparam
}
更新了代码。
作者: 柳烟    时间: 2012-6-29 07:36

里面有个嵌套循环,可能有占麻烦。
作者: 榕坚    时间: 2012-6-29 16:54

83# 柳烟


不一定要按UF做,一般陷井构造步骤:
一:确定陷井区域,之前的表达式表示的区域比较简单。M集陷井的话需要先确定中心,旋转角,这些坐标变换公式都好办。然后按一般M集作法先进行迭代(这就是UF中的嵌套循环),完了后得到终点Zn的同时也得到落入M集内部点Zn满足的条件,这时就与表达式区域一样了。
二:做主体分形,此时可以按上一步的条件确定迭代点是否落入M集区域内进入判断并标记该点作色,对于落入与没有落入区域的点分别赋于不同着色参数。
三:对坐标平面上的点按以上着色参数进行着色。
作者: 榕坚    时间: 2012-6-29 19:52

这个也很好玩的:

图片附件: Fractal1.jpg (2012-6-29 19:52, 70.08 KB) / 下载次数 1470
http://inrm3d.cn/attachment.php?aid=17837&k=9e778dae815767d950dc955fd9a4a4f8&t=1732387156&sid=xR1Q1q


作者: xiaongxp    时间: 2012-6-29 20:02

85# 榕坚
好看,眼前一亮!
作者: 柳烟    时间: 2012-6-29 20:20

84# 榕坚
就是你说的第二条,老是迷不穿。
作者: 榕坚    时间: 2012-6-30 09:06

30# 柳烟
柳老师好象偷懒得更利害了,把当中的开关项:Transversal twist也给删了吧?它是一个沿x轴方向的切变变换,我怎么一直不能实现。有空请试试看。
作者: 柳烟    时间: 2012-6-30 12:39

88# 榕坚
呵呵,当时觉得这种变换的图不怎么样,还是那个等边三角形的谢氏图案美,所以删掉了,诚如你所说,人是有点懒了,不过,这功能也易,添一个判断进去。问好。
顺便问问:M集陷阱第二条,当一个复分形迭代终点为zn',如何建立与M集的迭代终点zn建立判断,以知晓其落入M集上?
作者: 柳烟    时间: 2012-6-30 12:51

我昨天研究了一天嵌套循环,好象用GSP能够实现嵌套循环时满足某条件,嵌套循环终止并退出的问题,我调看了迭代表后得知,可惜作以前一个分形,扫出的图不对劲,昨天一天没找到原因。这个分形是:
http://www.inrm3d.cn/viewthread.php?tid=3225&extra=page%3D2
一二楼。
此分形一直没作出,成为一块心病,此分形意义重大。
作者: 榕坚    时间: 2012-6-30 14:30

89# 柳烟


Zn'做为Zn的初象进行M集的迭代啊。
作者: 榕坚    时间: 2012-6-30 21:53

89# 柳烟


就是这个判断被折磨了半天哪,真见鬼了。看似很容易的事情却怎么也办不到这才气人。就那么一行的代码。
作者: 柳烟    时间: 2012-7-1 08:45

92# 榕坚
仅将前面的那个谢氏画板文件中的z的初值横坐标x[z]编辑计算:x[pixel] - @trt*y[pixel] ,其余的都原封不动,立马得到这效果。x[pixel]、y[pixel]  为pixel的横标,纵标。
未命名.jpg

图片附件: 未命名.jpg (2012-7-1 09:21, 102.47 KB) / 下载次数 2441
http://inrm3d.cn/attachment.php?aid=17842&k=bea3c22e82182e84b9a303d45cf9ac01&t=1732387156&sid=xR1Q1q


作者: 榕坚    时间: 2012-7-1 10:36

93# 柳烟


真是一语点醒梦中人哪,有的时候头脑总是跟着代码走,且误读了代码。其实它就是一个变换而已。真是水到渠成,之前就一直想把M集弄到空白三角形中但一直没有实现,因此,学习分形切不可急于求成哪。

图片附件: SierpinskiTriangleII-4'.JPG (2012-7-1 10:36, 45.69 KB) / 下载次数 1930
http://inrm3d.cn/attachment.php?aid=17843&k=dfd4465458ec7d36aaf95f700a2299c7&t=1732387156&sid=xR1Q1q


作者: xiaongxp    时间: 2012-7-1 10:53

两位老师真是锲而不舍,快到光辉顶点了,为你们喝彩。
作者: 柳烟    时间: 2012-7-1 21:02

if first
      stdz = dz
      first = false
    elseif dz > stdz
      stdz = dz
    endif

这段代码不知何意,但是这几句代码又重要。做出的图相去甚远,连陷阱的影子都没见着,思维越来越理不清头绪了。榕兄此分形作出来没有?
作者: 柳烟    时间: 2012-7-1 21:05

95# xiaongxp
问好向老师,太费脑筋了,有时忙乎老半天,结果得到虚无,烦。
作者: 榕坚    时间: 2012-7-2 09:11

96# 柳烟


之前弄过一个按代码做的,但逃逸区非常杂不是很满意:
http://www.inrm3d.cn/viewthread. ... age%3D1&page=55
作者: 柳烟    时间: 2012-7-2 12:58

修改了那段费解代码,按下面代码作,在UF中对照了一下,好象与UF中的效果一致。
juliatrap1(BOTH) {
; By Samuel Monnier, 1999
init:
  float dz = 0.0
  float stdz = 0.0
  start = @seed

loop:
  trz = exp(flip(-pi/180*@rot))*(#z-@center)/@size
  x = 1/sqrt(@r)*real(trz)
  y = sqrt(@r)*imag(trz)
  if @freq != 0
    x = 2*sin(2*pi*@freq*x)
    y = 2*sin(2*pi*@freq*y)
  endif

  start = x + flip(y)

  iter = @niter
  zz = @seed

  float i = 0
  while i < @niter
    i = i + 1
    zz = zz^2 + start
    if |zz| > 1e20
      iter = i
      i = @niter
    endif

  endwhile
  float logp = 1/log(2)
  float logb = log(log(1e20))

    dz = real(iter + logp*logb - logp*log(log(cabs(zz))))

    if dz >stdz
      stdz = dz
    endif

final:
  #index = .1*stdz
default:
  title = "Copy of Julia Trap"
  helpfile = "sam-help/juliatrap.htm"

  param seed
    caption = "Seed (for 'custom')"
    default = (0,0)
  endparam
  param center
    caption = "Center"
    default = (1,1)
  endparam
  param rot
    caption = "Rotation"
    default = -30.0
  endparam
  param size
    caption = "Size"
    default = 0.38750
  endparam
  param r
    caption = "Ratio Width/Heigh"
    default = 1.0
  endparam
  param @niter
    caption = "Julia Iterations"
    default = 100
  endparam
  param freq
    caption = "Trap Frequency"
    default = 0.0
  endparam
}
我按此代码作出图后,逃逸区,参数值与面板上的值不完全一致,有问题。这变种图是:
未命名.jpg

图片附件: 未命名.jpg (2012-7-2 12:58, 131.09 KB) / 下载次数 2392
http://inrm3d.cn/attachment.php?aid=17848&k=0428269fc7c62983000dc6423fa723ac&t=1732387156&sid=xR1Q1q


作者: 柳烟    时间: 2012-7-2 18:08

之所以与UF有差别,原因在于那循环,用画板如何实现,是个迷,也许终难实现。我是按逃逸时间算法造那循环的M集陷阱,所以等势圈过于花梢。




欢迎光临 inRm3D: 画板论坛 (http://inrm3d.cn/) Powered by Discuz! 7.0.0