Curva Y Copo De Nieve De Von Koch Con Python 3

por | 8 enero, 2014

Hola
El otro día vi un videotutorial sobre cómo dibujar la curva de Von Koch con Python, y también su copo de nieve. Al verlo me decidí a implementarlo en Python 3 utilizando las librerías gráficas para matemáticas Matplotlib.
La curva de Von Koch es un fractal que para construirlo partimos de un segmento, en el primer paso dividimos el segmento de partida en 3 partes iguales. La parte central la eliminamos y hacemos aparecer en dicho hueco un triángulo equilátero sin su base. Cada uno de los lados del triángulo equilátero mide la tercera parte del segmento de partida. Por lo tanto, en el primer paso tenemos una figura compuesta de 4 lados.
En el segundo paso lo que hacemos es volver a realizar lo que hemos hecho en el primer paso pero a cada uno de los 4 lados anteriores. Y esto lo podemos realizar de forma indefinida. De esa forma es como obtenemos un fractal, ya que cada porción es semejante al resto de porciones (autosemejanza).
El copo de nieve de Von Koch consiste, en grandes rasgos, en cerrar la curva de Von Koch. Es decir, obtener una curva cerrada, la cual se asemeja a un copo de nieve conforme vamos aumentando el número de pasos. Aquí una imagen sacada de la Wikipedia:
Von Koch curve

Von Koch lo que trataba era encontrar una curva continua que nunca fuese diferenciable, puedes encontrar más información sobre ello en el siguiente enlace:


Bien, nos vamos a fijar para construir el código en Python 3 en la primera iteración:
Von koch 1 etape
De izquierda a derecha vamos a nombrar los 4 puntos por:
(x_i,\: y_i),\: (x_1,\: y_1)\: (x_2,\: y_2),\: (x_3,\: y_3),\: (x_f,\: y_f)

En nuestro código debemos determinar los puntos 1, 2 y 3 respecto el punto inicial y el punto final. Los puntos 1 y 2 son los más sencillos de determinar. De forma general podemos considerar cada una de sus ordenadas como la fórmula que se utiliza para determinar un segmento [a,\: b]:
(1-t)a+bt\: ; \forall\: t\in [0,\: 1]

Para el punto 1 consideramos t=\frac{1}{3}, que es la distancia que hay del punto inicial al punto 1, así:
(x_1,\: y_1)\: = ((1-\frac{1}{3} )x_i+\frac{1}{3}x_f,\: (1-\frac{1}{3} )y_i+\frac{1}{3}y_f)\:
Para el punto 3 consideramos t=\frac{2}{3}, que es la distancia que hay del punto inicial al punto 3, así:
(x_1,\: y_1)\: = ((1-\frac{2}{3} )x_i+\frac{2}{3}x_f,\: (1-\frac{2}{3} )y_i+\frac{2}{3}y_f)\:
Te dejo los cálculos algebraicos de simplificación, son muy sencillos
Nos queda determinar el punto 2, en el vídeo que vi la fórmula que pone es correcta pero no he utilizado esa porque no se obtenerla, no tengo porqué saberlo todo . Pensando un poco lo que he hecho es considerar que el punto 2 pertenece a la circunferencia de centro el punto (x_1,\: y_1) y radio la distancia entre (x_1,\: y_1) y el punto (x_3,\: y_3). Es decir;
r = \sqrt{(x_3-x_1)^2+(y_3-y_1)^2}
Pero esto lo vamos a considerar en coordenadas polares, por lo que necesitamos determinar el ángulo.
\alpha = atan(\dfrac{y_3-y_1}{x_3-x_1})
En el código Python utilizaremos hypot para calcular el radio y atan2 para calcular el ángulo. Para iterar debemos aumentar en nuestro código el ángulo:
\alpha = \alpha +\frac{\pi}{3}
Los pi tercios es el ángulo que forma el triángulo equilátero de la primera iteración, al ser equilátero cada uno de sus ángulos mide 60 grados que en radianes son los pi tercios.
Recapitulando, nuestro punto 2 en coordenadas polares es:
(x_2,\: y_2)=(x_1+r\cdot cos(\alpha ),\: y_1+r\cdot sin(\alpha ))
Con esto ya podemos dibujar la primera iteración de la curva. Al ser la figura un fractal podemos calcular las siguientes iteraciones de la curva de forma recursiva. Te dejo el código abajo.

  1. #/usr/bin/env python3.3.2
  2. # -*- coding: utf-8 -*-
  3.  
  4.  
  5. from math import sin, cos, pi, atan2
  6. from pylab import *
  7.  
  8.  
  9. def copoVonKoch(lado, n):
  10.     x_vertice1 = 0
  11.     y_vertice1 = 0
  12.  
  13.     x_vertice2 = lado * cos(2 * pi / 3)
  14.     y_vertice2 = lado * sin(2 * pi / 3)
  15.  
  16.     x_vertice3 = lado * cos(pi / 3)
  17.     y_vertice3 = lado * sin(pi / 3)
  18.  
  19.     curvaVonKoch(x_vertice1, y_vertice1, x_vertice2, y_vertice2, n)
  20.     curvaVonKoch(x_vertice2, y_vertice2, x_vertice3, y_vertice3, n)
  21.     curvaVonKoch(x_vertice3, y_vertice3, x_vertice1, y_vertice1, n)
  22.  
  23.  
  24. def curvaVonKoch(xi, yi, xf, yf, n):
  25.     if n == 0:
  26.         plot([xi, xf], [yi, yf], lw=1.0, color='b')
  27.     elif n > 0:
  28.         x1 = xi + (xf - xi) / 3.0
  29.         y1 = yi + (yf - yi) / 3.0
  30.  
  31.         x3 = xf - (xf - xi) / 3.0
  32.         y3 = yf - (yf - yi) / 3.0
  33.  
  34.         radio = hypot(x3 - x1, y3 - y1)
  35.         alpha = atan2((y3 - y1), (x3 - x1))
  36.         alpha += pi / 3.0
  37.         x2 = x1 + radio * cos(alpha)
  38.         y2 = y1 + radio * sin(alpha)
  39.  
  40.         curvaVonKoch(xi, yi, x1, y1, n - 1)
  41.         curvaVonKoch(x1, y1, x2, y2, n - 1)
  42.         curvaVonKoch(x2, y2, x3, y3, n - 1)
  43.         curvaVonKoch(x3, y3, xf, yf, n - 1)
  44.  
  45.  
  46. def dibuja(lado, n):
  47.     axes().set_xlim(0, lado)
  48.     axes().set_ylim(-2, lado / 2.0)
  49.     axes().set_aspect(1.0)
  50.     title('Curva De Von Koch')
  51.     xlim(0, lado)
  52.     curvaVonKoch(0, 0, lado, 0, n)
  53.     show()
  54.     axes().set_xlim(-lado, lado)
  55.     axes().set_ylim(-2, lado + 0.5 * lado)
  56.     axes().set_aspect(1.0)
  57.     title('Copo De Nieve De Koch')
  58.     xlim(-lado, lado)
  59.     copoVonKoch(lado, n)
  60.     show()
  61.  
  62.  
  63. def principal():
  64.     print('CURVA DE VON KOCH Y SU COPO DE NIEVE.\n')
  65.     print('-------------------------------------')
  66.     n = int(raw_input('¿Hasta qué iteración quieres dibujar? '))
  67.     lado = float(raw_input('Indica la longitud del lado inicial. '))
  68.     dibuja(lado, n)
  69.  
  70.  
  71. if __name__ == "__main__":
  72.     principal()

Al parecer WP me esta haciendo la puñeta con los signos de desigualdad, dejo un enlace de descarga del código:


El vídeo en que vi la explicación son 2, pongo los enlaces:

Abierto a cualquier mejora del código
Saludos

2 pensamientos en “Curva Y Copo De Nieve De Von Koch Con Python 3

  1. Pingback: Bitacoras.com

  2. Pingback: Curva Y Copo De Nieve De Von Koch Con Python 3 ...

Los comentarios están cerrados.