跳转到内容

分形/复平面上的迭代/三角不等式

来自维基教科书,开放的书籍,开放的世界
TIA

用于为集合(曼德布罗特或朱利亚)外部着色的算法

算法族:分支平均着色算法(“因为计算了在截断轨道中所有点上评估的函数的平均值”)[1]

  • 三角不等式平均着色算法 = TIA = Trippie。“线性插值,专门为曼德布罗特集合设计。在映射到零的某些点处不连续”
  • 曲率平均着色算法(由 Damien M. Jones 于 1999 年为 Ultra Fractal 开发,在映射到零的点处连续)
  • 纤维着色方法(来自 Jux 程序):“首先,我更喜欢带权重参数 (0 - 0.999) 的运行平均值,它避免了任何可能的精度问题,并且允许你调整细节。其次,我使用 sin 或 cos 而不是简单的角度,并添加一个频率参数。以及一些其他变体” xenodreambuie
  • Kerry Mitchell 1999.02.05
  • Damien M. Jones 2000.04.02

“三角不等式基本上是一种计算角度的廉价方法。TIA 对所有迭代的角进行平均以获得平滑的结果。因此,需要一些初始化以及每次迭代一些计算才能完成求和。它必须在添加下一个之前存储上一次迭代的求和,以便你可以在迭代带之间插值以获得连续函数(与通常对连续电势进行插值的方式相同)。这种插值是在迭代退出后完成的。” xenodreambuie [2]

“创建平滑着色的... 从分形向外延伸的巨大火焰状图案”(来自 Frederik Slijkerman 的 Ultra Fractal)[3]

  • 逃逸半径: “使用非常大的退出值 (1E20 = 10^20) 以获得良好的结果。” Ultra Fractal
  • 迭代最大 = 2500

Ultrafractal

[编辑 | 编辑源代码]
/* 
http://formulas.ultrafractal.com/

comment {
  dmj-pub.ucl 2.1
  Coloring methods for Ultra Fractal 2
  by Damien M. Jones
  April 2, 2000

  For more information about this formula collection,
  please visit its home page:
  
      http://www.fractalus.com/ultrafractal/dmj-pub-uf.htm

  Many of these coloring algorithms are general-purpose tools
  which work with a variety of fractals. A few are more
  specialized. Huge portions of this compilation build upon
  the techniques explored by others.
}
*/
dmj-Triangle {
;
; This coloring method implements a variation of Kerry
; Mitchell's Triangle Inequality Average coloring
; method.  Because the smoothing used in this formula
; is based on the dmj-Smooth formula, which only works
; for z^n+c and derivates, the smoothing here will only
; work for those types as well.
;
init:
  float sum = 0.0
  float sum2 = 0.0
  float ac = cabs(#pixel)
  float il = 1/log(@power)
  float lp = log(log(@bailout)/2.0)
  float az2 = 0.0
  float lowbound = 0.0
  float f = 0.0
  BOOL first = true
  float ipower = 1/@apower

loop:
  sum2 = sum
  IF (!first)
    az2 = cabs(#z - #pixel)
    lowbound = abs(az2 - ac)
    IF (@aflavor == 0)
      sum = sum + ((cabs(#z) - lowbound) / (az2+ac - lowbound))^@apower
    ELSEIF (@aflavor == 1)
      sum = sum + 1-(1-(cabs(#z) - lowbound) / (az2+ac - lowbound))^ipower
    ENDIF
  ELSE
    first = false
  ENDIF

final:
  sum = sum / (#numiter)
  sum2 = sum2 / (#numiter-1)
  f = il*lp - il*log(log(cabs(#z)))
  #index = sum2 + (sum-sum2) * (f+1)
  
default:
  title = "Triangle Inequality Average"
  helpfile = "dmj-pub\dmj-pub-uf-tia.htm"

  param apower
    caption = "Average Exponent"
    default = 1.0
    hint = "This skews the values averaged by raising them to \
            this power. Use 1.0 for the classic coloring."
  endparam
  param aflavor
    caption = "Average Flavor"
    default = 0
    enum = "normal" "reversed"
    hint = "Controls whether values are reversed before being \
            raised to a power. Has no effect if Average Exponent \
	    is 1.0."
  endparam
  param power
    caption = "Exponent"
    default = 2.0
    hint = "This should be set to match the exponent of the \
            formula you are using.  For Mandelbrot, this is 2."
  endparam
  param bailout
    caption = "Bailout"
    default = 1e20
    min = 1
    hint = "This should be set to match the bailout value in \
            the Formula tab.  Use a very high bailout!"
  endparam
}
// http://www.kerrymitchellart.com/tutorials/formulas2/uf2-2.htm
// http://formulas.ultrafractal.com/cgi/formuladb?view;file=lkm.ufm;type=.txt
comment { ; copyright Kerry Mitchell 07feb99

Triangle Inequality

The triangle inequality method is based on a simple characteristic of
complex numbers:  the magnitude of the sum of two complex numbers, |a+b|,
is strictly limited to a range determined by a and b:

|a+b| >= ||a| - |b||, and
|a+b| <= |a| + |b|,

where |z| is the square root of the sum of the squares of the
components, not the sum of the squares, as in Ultra Fractal.  The
extremes of this inequality are easily seen with a few examples.
If a=1 and b=2, then:

|a| = 1, |b| = 2;
||a| - |b|| = |1-2| = |-1| = 1;
|a| + |b| = 1+2 = 3;
|a+b| = |3| = 3;
1 <= 3 <= 3.

The upper bound occurs when both addends have the same polar angle.
The geometric interpretation of this is that the complex numbers add
up, and the length of the sum is simply the sum of the individual
lengths.

The lower bound occurs when the polar angles of the complex numbers
differ by 180 degrees; the two numbers are diametrically opposed.
Then, the length of the sum is the difference of the lengths.  For
example, if a=3i and b=-5i, then:

|a| = 3, |b| = 5;
||a| - |b|| = |3-5| = |-2| = 2;
|a| + |b| = 3 + 5 = 8;
|a+b| = |-2i| = 2;
2 <= 2 <= 8.

In general, the length of the sum is somewhere inbetween, and can be
thought of in terms of a triangle, which is how the inequality gets
its name.  If |a| is the length of one side of a triangle, and |b|
is the length of the second side, then |a+b| is the length of the
third side, and lies somewhere within the range shown above.

Back to fractals.  The two numbers of interest are z^n and c.  Given
z (the previous iterate) and c (the Mandelbrot or Julia parameter),
the range for the magnitude of the new iterate can then be determined.
With this range, the magnitude of the new iterate can be rescaled to
a fraction between 0 and 1 inclusive:

min = ||z_old| - |c||, max = |z_old| + |c|,
z_new = z_old * z_old + c
fraction = (|z_new| - min) / (max - min).

During the iterating, these fractions are average together.  After the
last iteration, this average fraction is stored in the real part of z,
and the actual fraction for the last iteration is stored in the imaginary
part of the z.  If the orbit diverges (bails out), then renormalization
techniques are used to color the outside pixels smoothly without iteration
bands.  This works best with large bailout values.  The large bailout
needs to be reduced with larger exponents, to prevent color gaps due to
overflow errors.

For best results, use the "basic" coloring method, with "real(z)" to show
the average triangulation fraction, or "imag(z)" to show the actual
fraction for the last iteration.

}

triangle-general-julia { ; Kerry Mitchell 05feb99
;
; z^n+c Julia set, colors by
; triangle inequality
;
init:
  zc=#pixel
  c=@julparam
  float rc=cabs(c)
  float rzc=0.0
  zn=(0.0,0.0)
  float rnz=0.0
  int iter=0
  bool done=false
  float logn=log(@nexp)
  float llbail=log(log(@bailout))+log(0.5)
  float count=0.0
  float count1=0.0
  float count2=0.0
  float countx=0.0
  float county=0.0
  float fac1=0.0
  float fac2=0.0
  float min=0.0
  float max=0.0
  float k=0.0
  float kfrac=0.0
loop:
  iter=iter+1
  if(@nexp==2.0)
    zn=sqr(zc)
    rnz=|zc|
  else
    zn=zc^@nexp
    rnz=cabs(zn)
  endif
  zc=zn+c
  min=cabs(rnz-rc)
  max=rnz+rc
  rzc=cabs(zc)
  county=(rzc-min)/(max-min)
  count=count+county
  countx=count/iter
  z=countx+flip(county)
  if(|zc|>@bailout)
    done=true
    count1=count/iter
    fac1=log(log(rzc))-llbail
    iter=iter+1
    if(@nexp==2.0)
      zn=sqr(zc)
      rnz=|zc|
    else
      zn=zc^@nexp
      rnz=cabs(zn)
    endif
    min=cabs(rnz-rc)
    max=rnz+rc
    zc=zn+c
    rzc=cabs(zc)
    county=(rzc-min)/(max-min)
    count=count+county
    count2=count/iter
    fac2=log(log(rzc))-llbail-logn
    k=exp((fac1+fac2)/2.0)
    kfrac=(k-1.0)/(@nexp-1.0)
    countx=kfrac*count1+(1.0-kfrac)*count2
    z=countx+flip(county)
  endif
bailout:
  done==false
default:
  title="Triangle General Julia"
  maxiter=100
  periodicity=0
  center=(0.0,0.0)
  magn=1.0
  angle=0
  param julparam
    caption="Julia parameter"
    default=(1.0,0.0)
  endparam
  param bailout
    caption="bailout value"
    default=1000000.0
    min=0.0
  endparam
  param nexp
    caption="exponent"
    default=2.0
    hint="z exponent, > 1.0"
    min=1.0
  endparam
switch:
  type="triangle-general-mandelbrot"
  bailout=bailout
  nexp=nexp
}

triangle-general-mandelbrot { ; Kerry Mitchell 05feb99
;
; z^n+c Mandelbrot set, colors by
; triangle inequality
;
init:
  c=#pixel
  zc=@manparam+c
  float rc=cabs(c)
  float rzc=0.0
  zn=(0.0,0.0)
  float rnz=0.0
  int iter=0
  bool done=false
  float logn=log(@nexp)
  float llbail=log(log(@bailout))+log(0.5)
  float count=0.5
  float count1=0.0
  float count2=0.0
  float countx=0.0
  float county=0.0
  float fac1=0.0
  float fac2=0.0
  float min=0.0
  float max=0.0
  float k=0.0
  float kfrac=0.0
loop:
  iter=iter+1
  if(@nexp==2.0)
    zn=sqr(zc)
    rnz=|zc|
  else
    zn=zc^@nexp
    rnz=cabs(zn)
  endif
  zc=zn+c
  min=cabs(rnz-rc)
  max=rnz+rc
  rzc=cabs(zc)
  county=(rzc-min)/(max-min)
  count=count+county
  countx=count/iter
  z=countx+flip(county)
  if(|zc|>@bailout)
    done=true
    count1=count/iter
    fac1=log(log(rzc))-llbail
    iter=iter+1
    if(@nexp==2.0)
      zn=sqr(zc)
      rnz=|zc|
    else
      zn=zc^@nexp
      rnz=cabs(zn)
    endif
    min=cabs(rnz-rc)
    max=rnz+rc
    zc=zn+c
    rzc=cabs(zc)
    county=(rzc-min)/(max-min)
    count=count+county
    count2=count/iter
    fac2=log(log(rzc))-llbail-logn
    k=exp((fac1+fac2)/2.0)
    kfrac=(k-1.0)/(@nexp-1.0)
    countx=kfrac*count1+(1.0-kfrac)*count2
    z=countx+flip(county)
  endif
bailout:
  done==false
default:
  title="Triangle General Mandelbrot"
  maxiter=100
  periodicity=0
  center=(0.0,0.0)
  magn=1.0
  angle=0
  param manparam
    caption="Mandelbrot start"
    default=(0.0,0.0)
    hint="use (0,0) for standard Mandelbrot"
  endparam
  param bailout
    caption="bailout value"
    default=1000000.0
    min=0.0
  endparam
  param nexp
    caption="exponent"
    default=2.0
    hint="z exponent, > 1.0"
    min=1.0
  endparam
switch:
  type="triangle-general-julia"
  julparam=pixel
  bailout=bailout
  nexp=nexp
}

来自 dmj.ucl

comment {
  dmj-pub.ucl 2.1
  Coloring methods for Ultra Fractal 2
  by Damien M. Jones
  April 2, 2000

  For more information about this formula collection,
  please visit its home page:
  
      http://www.fractalus.com/ultrafractal/dmj-pub-uf.htm

  Many of these coloring algorithms are general-purpose tools
  which work with a variety of fractals. A few are more
  specialized. Huge portions of this compilation build upon
  the techniques explored by others.
}

dmj-Curvature {
;
; This formula averages curvature over the course
; of all iterations by averaging the difference in
; absolute value of angles between orbit steps. The
; technique was suggested by Prof. Javier Barrallo.
;
; Note that the results are often very similar to
; Kerry Mitchell's Triangle Inequality Average method,
; but are slightly more angular.
;
init:
  complex zold = (0,0)
  complex zold2 = (0,0)
  float a = 0.0
  float a2 = 0.0
  int i = 0

loop:
  a2 = a
  IF (i >= 2 && @aflavor != 2)			; zold and zold2 are valid
    a = a + abs(atan2((#z-zold)/(zold-zold2)))	; update average
  ENDIF
  i = i + 1					; count the iteration
  zold2 = zold					; save the orbit values
  zold = #z

final:
  IF (@aflavor == 0)
    float il = 1/log(@power)
    float lp = log(log(@bailout)/2.0)
    float f = il*lp - il*log(log(cabs(#z)))
    a = a / i
    a2 = a2 / (i-1)
    #index = (a2 + (a-a2) * (f+1)) / #pi
  
  ELSEIF (@aflavor == 1)
    #index = (a/i) / #pi + 1
    
  ELSE
    a = atan2((#z-zold)/(zold-zold2))		; update average
    #index =  a    / #pi + 1
  ENDIF

混沌之影 = VOC

[编辑 | 编辑源代码]
#version 400
// http://www.fractalforums.com/programming/triangle-inequality-average-coloring/
// code by Softology
uniform vec2 resolution;
uniform vec3 palette[256];
uniform double xmin;
uniform double xmax;
uniform double ymin;
uniform double ymax;
uniform double bailout;
uniform int maxiters;

double bailout_squared=double(bailout*bailout);
double magnitude,r1,r2,g1,g2,b1,b2,tweenval;
double realiters;
vec4 finalcol,col;
int superx,supery;
double stepx=(xmax-xmin)/resolution.x;
double stepy=(ymax-ymin)/resolution.y;
int colval,colval1,colval2;
dvec2 z,c;
//triangle inequality average coloring
double sum,sum2,ac,il,lp,az2,lowbound,f,index,tr,ti;
int mandelbrotPower;
double rval,gval,bval,rval1,gval1,bval1,rval2,gval2,bval2;

void main(void)
{
	sum=0;
	sum2=0;
	ac=0;
	il=0;
	lp=0;
	mandelbrotPower=2;

	finalcol=vec4(0,0,0,0);
			c.x = xmin+gl_FragCoord.x/resolution.x*(xmax-xmin);
			c.y = ymin+gl_FragCoord.y/resolution.y*(ymax-ymin);
			int i;
			z = dvec2(0.0,0.0);

				//triangle inequality average coloring
				sum = 0;
				sum2 = 0;
				ac = sqrt(c.x * c.x + c.y * c.y);
				il = 1.0 / log(mandelbrotPower);
				lp = log(float(log(float(bailout)) / mandelbrotPower));
				az2 = 0.0;
				lowbound = 0.0;
				f = 0.0;
				index = 0.0;

			for(i=0; i<maxiters; i++) 
			{
				//START OF FRACTAL FORMULA
				double x = (z.x * z.x - z.y * z.y) + c.x;
				double y = (z.y * z.x + z.x * z.y) + c.y;
				//END OF FRACTAL FORMULA
				
				magnitude=(x * x + y * y);
				if(magnitude>bailout_squared) break;
				z.x = x;
				z.y = y;

					//tia
					sum2=sum;
					if ((i!=0)&&(i!=maxiters-1)) {
						tr=z.x-c.x;
						ti=z.y-c.y;
						az2=sqrt(tr * tr + ti * ti);
						lowbound=abs(az2 - ac);
						sum+=((sqrt(z.x * z.x + z.y * z.y)-lowbound)/(az2+ac-lowbound));
					}

			}

			if (i==maxiters) {
				col=vec4(0.0,0.0,0.0,1.0);
			} else {
				//triangle inequality average
				sum=sum/i;
				sum2=sum2/(i-1.0);
				f=il*lp - il*log(log(float(length(z))));
				index=sum2+(sum-sum2)*(f+1.0);
				realiters=255*index;
				colval1= int(mod(realiters,255));
				colval2=int(mod((colval1+1),255));
				tweenval=fract(realiters);
				if (colval1<0) { colval1=colval1+255; }
				if (colval2<0) { colval2=colval2+255; }
				rval1 =palette[colval1].r;
				gval1 =palette[colval1].g;
				bval1 =palette[colval1].b;
				rval2 =palette[colval2].r;
				gval2 =palette[colval2].g;
				bval2 =palette[colval2].b;
				rval =rval1 +((rval2 - rval1)*tweenval);
				gval =gval1 +((gval2 - gval1)*tweenval);
				bval =bval1 +((bval2 - bval1)*tweenval);
				col=vec4(rval,gval,bval,1.0);
			}
	gl_FragColor = vec4(col);
}

参考资料

[编辑 | 编辑源代码]

参考资料

[编辑 | 编辑源代码]
华夏公益教科书