|
Size: 2366
Comment: Started section with "Winding number of a plane curve"
|
← Revision 7 as of 2020-06-01 18:42:11 ⇥
Size: 2538
Comment:
|
| Deletions are marked like this. | Additions are marked like this. |
| Line 7: | Line 7: |
| by Pablo Angulo. Computes winding number as an integral, and also as a intersection number with a half line through the origin. | by Pablo Angulo. Computes winding number (with respect to the origin!) as an integral, and also as a intersection number with a half line through the origin. |
| Line 9: | Line 9: |
| {{{ | {{{#!sagecell |
| Line 13: | Line 13: |
| r2 = x^2 + y^2 | r2 = x**2 + y**2 |
| Line 16: | Line 16: |
| integrando = (x*yp -y*xp)/r2 | integrando = (x*yp -y*xp) / r2 |
| Line 18: | Line 18: |
| return round(i/(2*pi)) | return round(i / (2 * pi)) |
| Line 24: | Line 24: |
| delta= (b-a)/N | delta = (b - a) / N |
| Line 26: | Line 26: |
| for t in srange(a, b ,delta): | for t in srange(a, b, delta): |
| Line 28: | Line 28: |
| zeros.append(find_root(f, t-epsilon, t+delta+epsilon)) | zeros.append(find_root(f, t - epsilon, t + delta + epsilon)) |
| Line 32: | Line 32: |
| if not zeros: return zeros if abs(zeros[0] + 2*pi - zeros[-1])<epsilon: |
if not zeros: return zeros if abs(zeros[0] + 2*pi - zeros[-1]) < epsilon: |
| Line 37: | Line 38: |
| if abs(c - zeros_cleaned[-1])>epsilon: | if abs(c - zeros_cleaned[-1]) > epsilon: |
| Line 39: | Line 40: |
| if abs(zeros[0] + 2*pi - zeros[-1])<epsilon: | if abs(zeros[0] + 2*pi - zeros[-1]) < epsilon: |
| Line 46: | Line 47: |
| x = x.function(t); y = y.function(t); | x = x.function(t) y = y.function(t) |
| Line 48: | Line 50: |
| raise ValueError, "Curve is not closed!" | raise ValueError("Curve is not closed!") |
| Line 55: | Line 57: |
| html(r'$\int \frac{1}{x^2 + y^2}(xdy-ydx)=%d$'%winding_number_integral(x,y,a,b)) | pretty_print(html(r'$\int \frac{1}{x^2 + y^2}(xdy-ydx)=%d$'%winding_number_integral(x,y,a,b))) |
| Line 57: | Line 59: |
| zeros = all_the_zeros(x,a, b) | zeros = all_the_zeros(x, a, b) |
| Line 70: | Line 72: |
| print 'Winding number = (# of red points) - (# of green points): ', wn | print('Winding number = (number of red points) - (number of green points): {}'.format(wn)) |
| Line 73: | Line 75: |
| arrow((x(0),y(0)), (x(0) + xp1(0), y(0) + yp1(0))) ) | arrow((x(0),y(0)), (x(0) + xp1(0), y(0) + yp1(0))) + point2d([(0,0)], color = 'black', pointsize = 70)) |
Sage Interactions - Topology
goto interact main page
Winding number of a plane curve
by Pablo Angulo. Computes winding number (with respect to the origin!) as an integral, and also as a intersection number with a half line through the origin.
