Essays, Posts & Presentations

Finding Zeros with Domain Coloring and Winding Numbers

Garrett Dubofsky

Author

Garrett Dubofsky

Title

Finding Zeros with Domain Coloring and Winding Numbers

Description

A colorful process for determining the zeros of any continuous 2D transformation.

Finding Zeros with Domain Coloring and Winding Numbers

Garrett Dubofsky

Western Washington University

Given a 2-dimensional transformation, we shall start by plotting its domain coloring. Specifically, every point on the plane takes the hue of where that point’s transformation lands on a color wheel centered at the origin. The winding number of any boundary on this domain is representative of how many loops of the color spectrum are completed when moving clockwise along the boundary. If the winding number is a nonzero integer, the boundary can be repeatedly halved as it converges onto a zero of the transformation. We seek to create this process visually and mathematically using the Wolfram Language.

2D Transformations and Domain Coloring

Determining the Winding Number

Decreasing the Boundary

Concluding Remarks

Keywords

Acknowledgements

References

2D Transformations and Domain Coloring

The Goal

In many areas of math, it is important to determine the points at which an equation outputs a zero. These zeros can represent extrema, intersections, optimal values, and much more. There are many methods available to locate the zeros of simplistic equations but it becomes more difficult as the equations increase in complexity. In this project, we seek to discover one more method of finding zeros for a particular type of equation – a continuous 2-dimensional transformation – using a process involving domain coloring and winding numbers.

Visualizing Four Dimensions

It is a challenge to visualize 2-dimensional transformations in our 3-dimensional world. This task might sound intuitive but these 2D transformations require two inputs and two outputs, resulting in a total of 4 components. In order to visualize our search for zeros, we must overcome this obstacle.

Let us define a 2D transformation to use as an example for our computations. The Q function below is designed to be continuous and differentiable while remaining difficult to algebraically manipulate. Notice that the zeros of this type of transformation are points that map to the coordinates

{0,0}

. For the sake of innovation, I did not determine the zeros of Q before starting this project.

In[]:=

Q[{x_,y_}]={

-5

-2x*y+3+Cos[y],

+y*(x-3)+2y-3-Sin[y]};

There exist a variety of potential solutions to visualizing the above transformation, some of which work better than others. A common representation of 2D transformations is a vector graph. As seen below, this graph uses an arbitrary number of vectors located at different points to give the direction that Q sends each point. Each vector’s magnitude corresponds to the euclidean distance the given point travels. To be consistent, we will focus our attention on Q around the domain

[-6,6]×[-6,6]

.
(Plot the vector graph of Q with vector scaling enabled.)

In[]:=

VectorPlot[Q[{x,y}],{x,-6,6},{y,-6,6},ImageSize350,AspectRatio1,VectorColorFunctionNone,VectorStyleBlack,VectorScalingAutomatic,VectorPoints->16,PlotLabelStyle["Vector Graph of Q",Black]]

Out[]=

While this graph clearly explains the directions that Q induces, we cannot determine where each point travels. This is an insufficient amount of information for our goal. To compensate, we can create an animation of where Q sends a collection of points on the plane. The x and y axes have been highlighted to show which coordinates are reasonably close to a zero.
(Plot several animated parameterized points without frame ticks, which would otherwise decenter the image.)

In[]:=

Animate[Graphics[Table[Point[{x,y}+t(Q[{x,y}]-{x,y})],{x,-6,6,0.6},{y,-6,6,0.6}],FrameTrue,PlotRangeAll,AspectRatio1,FrameTicksNone,Axes->True,AxesStyleDirective[Opacity[0.75],Gray,Dashed],PlotLabelStyle["Animated Transformation of Q",Black]],{{t,0,"Transformation time "},0,1},AnimationRate.2,AnimationRunningFalse,AnimationDirection->ForwardBackward]

Out[]=

The chaotic nature of this animation does not help narrow our search radius. While we can determine approximately where a point is taken, it is difficult to be precise around the origin (where multiple points coincide). The above results teach us that an ideal visualization would need to show both a transformed point’s direction and its distance from the origin. Therefore, we will introduce a new form of visualization: domain coloring.

Domain Coloring

The process of domain coloring is unorthodox, a little unintuitive, but ultimately beautiful. To begin, we consider placing a hue wheel centered at the origin and covering the entire



plane. Points closer to the origin shall be shaded darker to exaggerate their importance. This shall be accomplished by creating the function hueMap, which uses a normalized form of ArcTan to assign each point’s angle around the origin with a value between 0 and 1. Another function hueMapVal shall only provide the numerical hue value, which will be useful in later computations. By applying the Hue function to each output, the desired hue wheel is created.
(The If function ensures ArcTan[0,0] is not evaluated. We specify 1 saturation and change the hue value based on the distance from the origin.)

In[]:=

hueMap[{x_,y_}]:=Ifx!=0||y!=0,Hue

ArcTan[x,y]+

2π

,1,0.5Norm[{x,y}],Black;

(Another If function distinguishes between the positive and negative x-values, restricting the range between 0 and 1 and preserves the color orientation seen above.)

In[]:=

hueMapVal[{x_,y_}]:=Ifx!=0||y!=0,Ifx≥0,

ArcTan[x,y]+

2π

ArcTan

+

2π

+.5,0;

(We use ColorFunction to apply the hueMap to the plane.)

In[]:=

Grid[{{ParametricPlot[{x,y},{x,-6,6},{y,-6,6},ColorFunctionFunction[{x,y},hueMap[{x,y}]],ColorFunctionScalingFalse,FrameTrue,AxesFalse,ImageSize350,PlotLabelStyle["Hue Wheel at the Origin (with shading)",Black],PlotPoints200],ParametricPlot[{x,y},{x,-6,6},{y,-6,6},ColorFunctionFunction[{x,y},Hue[hueMapVal[{x,y}]]],ColorFunctionScalingFalse,FrameTrue,AxesFalse,ImageSize350,PlotLabelStyle["Hue Wheel at the Origin (without shading)",Black],PlotPoints200]}}]

Out[]=

Returning to the earlier animation, we can color each point according to where its transformation lands on the hue wheel. In other words, the transformed point is “painted” by the hue wheel before returning to its starting location. Applying the hueMap function provides another layer of information and improves our ability to distinguish between points.
(Recreating the previous animation with the Style and hueMap functions.)

In[]:=

Animate[Graphics[Table[Style[Point[{x,y}+t(Q[{x,y}]-{x,y})],hueMap[Q[{x,y}]]],{x,-6,6,0.6},{y,-6,6,0.6}],FrameTrue,PlotRangeAll,AspectRatio1,FrameTicksNone,Axes->True,AxesStyleDirective[Opacity[0.75],Gray,Dashed],PlotLabelStyle["Animated Transformation of Q using hueMap",Black]],{{t,0,"Transformation time "},0,1},AnimationRate.1,AnimationRunningFalse,AnimationDirection->ForwardBackward]

Out[]=

This process can be improved further by increasing the density of points. With enough detail, the animation can be replaced by a static image that encodes both a transformed point’s direction and distance from the origin.
(Coloring each point in a region based on its Q-transformation.)

In[]:=

ParametricPlot[{x,y},{x,-6,6},{y,-6,6},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],ColorFunctionScalingFalse,FrameTrue,AxesFalse,ImageSize350,PlotLabelStyle["Domain Coloring of Q",Black],PlotPoints200]

Out[]=

The above graphic is called the domain coloring of Q. The shading makes is clear which points the Q transformation sends close to the origin. While this approaches our goal, domain coloring not only required evaluating hueMap over an entire area of the plane but it does not produce numerical estimates. Next, we will uncover a method that overcomes these obstacles.

Determining the Winding Number

Defining a Boundary

Instead of investigating an entire area, let us select just a closed boundary. Any arbitrary shape will suffice but few will make later processes significantly easier to compute. Therefore, we will limit our analysis to rectangular boundaries. Since a parametric plot used for the domain coloring, it is convenient to define a function rectangleParametric that parameterizes the perimeter of a rectangle defined by two corners.
(Parameterizing an arbitrary rectangle between t = 0 and t = 1.)

In[]:=

rectangleParametric[{{x1_,y1_},{x2_,y2_}},t_]:=Module{MaxX,MaxY,MinX,MinY,Xdiff,Ydiff,Perim},MaxX=Max[x1,x2];MaxY=Max[y1,y2];MinX=Min[x1,x2];MinY=Min[y1,y2];Xdiff=MaxX-MinX;Ydiff=MaxY-MinY;Perim=2Xdiff+2Ydiff;Piecewise1-

Perim

Xdiff

tMinX+

Perim

Xdiff

tMaxX,MaxY,0<=t<

Xdiff

Perim

,MaxX,1-

Perim

Ydiff

Xdiff

Perim

MaxY+

Perim

Ydiff

Xdiff

Perim

MinY,

Xdiff

Perim

<=t<

Xdiff+Ydiff

Perim

,

Perim

Xdiff

Xdiff+Ydiff

Perim

MinX+1-

Perim

Xdiff

Xdiff+Ydiff

Perim

MaxX,MinY,

Xdiff+Ydiff

Perim

≤t<

2Xdiff+Ydiff

Perim

,MinX,

Perim

Ydiff

2Xdiff+Ydiff

Perim

MaxY+1-

Perim

Ydiff

2Xdiff+Ydiff

Perim

MinY,

2Xdiff+Ydiff

Perim

≤t≤1;

This parameterized rectangle can be projected atop the above graphic to isolate the domain coloring of Q along its perimeter. This is a significant computational improvement; analyzing a perimeter instead of an area dramatically reduces the number of processed points.
(Plotting the domain coloring of Q along an arbitrary rectangle’s perimeter.)

In[]:=

Manipulate[ParametricPlot[rectangleParametric[{Pnt1,Pnt2},tt],{tt,0,1},FrameTrue,ImageSize350,AxesFalse,ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],ColorFunctionScalingFalse,PlotStyleThickness[0.01],PlotRange{{-6,6},{-6,6}},PlotStyle->Directive[CapForm["Round"],Thickness[0.01]],PlotLabelStyle["Domain Coloring of Q on Rectangle Boundary",Black]],{{Pnt1,{4,4}},{-6,-6},{6,6}},{{Pnt2,{-4,-4}},{-6,-6},{6,6}},ControlPlacementLeft]

Out[]=

What is Winding?

To understand the meaning behind a winding number, it is important to first describe what we are “winding”. Simply enough, we are winding around the hue wheel. This can be interpreted as the angular movement of an arrow that points towards its color on a surrounding hue wheel, as seen below. Moving in the “ROYGBIV” direction will increase the winding number, whereas moving in the opposite direction will decrease it.
(Interactive arrow and hue wheel interpretation of a winding number.)

In[]:=

hueDisk=ParametricPlot[{Cos[y],Sin[y]},{y,0,2π},ColorFunctionFunction[{x,y},Hue[hueMapVal[{x,y}]]],ColorFunctionScalingFalse,PlotStyleThickness[0.03]];ManipulateColumnShow[Graphics[{Arrowheads[0.1],Thickness[0.015],Style[Arrow[{{0,0},{Cos[θ],Sin[θ]}}],Hue[hueMapVal[{Cos[θ],Sin[θ]}]]]},ImageSize350,PlotRange{{-1.3,1.3},{-1.3,1.3}}],hueDisk],"Winding number: "<>ToStringN

2π

,{{θ,0},-2π,2π}

Out[]=

This arrow interpretation allows us to understand the winding number along our rectangular boundary. Since we restricted our analysis to continuous 2D transformations, we are assured that the arrow shall follow a continuous path. In the next section, we will introduce numeric values to accompany this visualization.
(Interactive rectangle boundary with cycling hue wheel arrow.)

In[]:=

Animate[Manipulate[{xR,yR}=rectangleParametric[{Pnt1,Pnt2},t];Show[ParametricPlot[rectangleParametric[{Pnt1,Pnt2},t],{t,0,1},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],PlotStyleThickness[0.01],ColorFunctionScalingFalse,PlotRange{{-8,8},{-8,8}},AspectRatio1,AxesFalse,ImageSize350,FrameTrue,PlotStyle->Directive[CapForm["Round"],Thickness[0.01]],PlotLabelStyle["Domain Coloring of Rectangle Boundary",Black]],Graphics[{Arrowheads[0.08],Thickness[0.015],Style[Arrow[{{xR,yR},{xR,yR}+3Normalize[Q[{xR,yR}]]}],hueMap[Q[{xR,yR}]]]}]],{{Pnt1,{4,4}},{-6,-6},{6,6}},{{Pnt2,{-1,-4}},{-6,-6},{6,6}},ControlPlacementLeft],{{t,0,"Distance around perimeter: "},0,1},AnimationRunningFalse,AnimationRate0.075]

Out[]=

Finding the Winding Number

Now that we have established the concept of winding, we need to calculate its value. Recall that the previously defined function hueMapVal can be used to plot the hue value of the boundary’s color. By coloring the plot accordingly, the correspondence between color and value becomes clear:
(A generalized function that plots the value of hueMap along a rectangular boundary.)

In[]:=

hueMapCurve[Functn_,{x1_,y1_},{x2_,y2_},t_]:=hueMapVal@Functn[rectangleParametric[{{x1,y1},{x2,y2}},t]];

(Plotting hueMapCurve alongside its respective rectangular boundary.)

In[]:=

Animate[Manipulate[Grid[{{ParametricPlot[rectangleParametric[{Pnt1,Pnt2},tt],{tt,0,1},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],PlotStyle->Directive[CapForm["Round"],Thickness[0.01]],ColorFunctionScalingFalse,PlotRange{{-8,8},{-8,8}},AspectRatio1,AxesFalse,ImageSize350,FrameTrue,Epilog{PointSize[0.02],Point[rectangleParametric[{Pnt1,Pnt2},t]]},PlotLabelStyle["Domain Coloring of Rectangle Boundary",Black]],Plot[hueMapCurve[Q,Pnt1,Pnt2,tt],{tt,0,1},ColorFunctionFunction[{x,y},Hue[y]],ColorFunctionScalingFalse,AspectRatio->1,PlotStyleThickness[0.0075],ImageSize->350,PlotRange{{0,1},{0,1}},Epilog{PointSize[0.02],Point[{t,hueMapCurve[Q,Pnt1,Pnt2,t]}]},PlotLabelStyle["Values of Hue along Rectangle Boundary",Black],AxesLabel{"t","Hue"},LabelStyleBlack,PlotPoints60,Exclusions->hueMapCurve[Q,Pnt1,Pnt2,tt]==.98]}}],{{Pnt1,{4,4}},{-6,-6},{6,6}},{{Pnt2,{-4,-4}},{-6,-6},{6,6}},ControlPlacement{Left,Left,Top}],{{t,0,"Distance around perimeter: "},0,1},AnimationRunningFalse,AnimationRate0.1]

Out[]=

Asymptotic discontinuities occur because the hue value is restricted between 0 and 1. However, both endpoints correspond to red, so moving directly between 0 and 1 does not constitute as a discontinuity when moving around the hue wheel. In other words, the hue value may be discontinuous but the color order is not. Since we start and end at the same location, the number of complete hue wheel rotations (the winding number) must be an integer.

The winding number is equivalent to the net change of the hue value across the rectangle’s perimeter. Since the curve generated by hueMapCurve is typically discontinuous and non-differentiable, we rely on a discretization to find the winding number. That is, the t-axis is partitioned into a large number of equal-length intervals and hueMapCurve is evaluated at end endpoints of each interval. The change in the function across the interval is found simply by subtracting the value of hueMapCurve at the smaller t-value from the larger t-value. The summation of these differences produces the net change of the hue value across the

[0,1]

domain. The below function applies this process to a rectangular boundary for a specified total number of intervals.
(Find the winding number when (100*perc)% around a rectangular boundary, also removing outliers)

In[]:=

findWR[Functn_,{{x1_,y1_},{x2_,y2_}},perc_,num_:2500]:=Total[Select[Table[hueMapVal[Functn[rectangleParametric[{{x1,y1},{x2,y2}},n]]]-hueMapVal[Functn[rectangleParametric[{{x1,y1},{x2,y2}},n-

-1

num

]]],{n,

-1

num

,perc,

-1

num

}],Abs[#]<0.5&]];

Since the winding number is an integer, we should only consider the rounded value of findWR. The parameter that corresponds to the number of intervals, num, shall determine the value’s accuracy.

We can combine all previous visualizations into a singular graphic that determines the winding number of Q on any rectangular boundary.
(The previous animation in tandem with the hue wheel and findWR function. The number of intervals is low to improve latency at the cost of accuracy.)

In[]:=

AnimateManipulateGridParametricPlot[rectangleParametric[{Pnt1,Pnt2},tt],{tt,0,1},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],PlotStyleThickness[0.01],ColorFunctionScalingFalse,PlotRange{{-8,8},{-8,8}},AspectRatio1,AxesFalse,ImageSize350,FrameTrue,PlotStyle->Directive[CapForm["Round"],Thickness[0.01]],Epilog{PointSize[0.02],Point[rectangleParametric[{Pnt1,Pnt2},t]]},PlotLabelStyle["Domain Coloring of Rectangle Boundary",Black]],Plot[hueMapCurve[Q,Pnt1,Pnt2,tt],{tt,0,1},ColorFunctionFunction[{x,y},Hue[y]],ColorFunctionScalingFalse,AspectRatio->1,PlotStyleThickness[0.0075],ImageSize->350,PlotRange{{0,1},{0,1}},Epilog{PointSize[0.02],Point[{t,hueMapCurve[Q,Pnt1,Pnt2,t]}]},PlotLabelStyle["Values of Hue along Rectangle Boundary",Black],AxesLabel{"t","Hue"},LabelStyleBlack,PlotPoints40,Exclusions->hueMapCurve[Q,Pnt1,Pnt2,tt]==.98],ColumnShowGraphicsArrowheads[0.1],Thickness[0.015],StyleArrow{0,0},Cos2πhueMapCurve[Q,Pnt1,Pnt2,t]-

,Sin2πhueMapCurve[Q,Pnt1,Pnt2,t]-

,HuehueMapValCos2πhueMapCurve[Q,Pnt1,Pnt2,t]-

,Sin2πhueMapCurve[Q,Pnt1,Pnt2,t]-

,ImageSize250,PlotRange{{-1.3,1.3},{-1.3,1.3}},PlotLabelStyle["Hue Wheel Locator",Black],hueDisk,Framed[Column[{"Distance around perimeter: "<>ToString[Round[t*100,0.01]]<>"%","Current Winding Number: "<>ToString[Round[findWR[Q,{Pnt1,Pnt2},t,350],0.01]],"Total Winding Number: "<>ToString[Round[findWR[Q,{Pnt1,Pnt2},1,350]]]}]],{{Pnt1,{4,4}},{-6,-6},{6,6}},{{Pnt2,{-4,-4}},{-6,-6},{6,6}},ControlPlacement{Left,Left,Top},{{t,0,"Distance around perimeter: "},0,1},AnimationRunningFalse,AnimationRate0.035

Out[]=

Decreasing the Boundary

Properties of the Winding Number

Now we can elaborate why the winding number is an important tool for finding zeros of 2D transformations. Notice that zeros are the only points with the property that if we zoom infinitely close to the point, we are guaranteed that the surrounding domain coloring shall pass through the entire hue spectrum. This property directly follows from the requirement that the 2D transformation is continuous, which we discovered also implies that the domain coloring is color-continuous (no jumps in hue between two points). Since the domain coloring of the plane places the origin in the center of the hue wheel, any point that maps to the origin must also retain this quality.

Halving a boundary with a nonzero winding number is guaranteed to produce at least one smaller boundary with a nonzero winding number. This is a result from the fact that winding numbers are integers along with the addition properties discussed later. Since a nonzero winding number implies that the boundary passes through the entire hue spectrum, we can repeatedly apply this process until we obtain an infinitesimally small boundary that still passes through all colors. As we stated earlier, this infinitely small boundary (which can be thought of as a point) must be a zero. Therefore, any closed boundary with a nonzero winding number must contain a zero. We specify a nonzero winding number because it is possible to create a boundary that passes through the entire spectrum but, when split in two, neither of its halves does the same. On the other hand, a boundary with a zero winding number may actually contain a zero. This scenario is captured with a large boundary around the domain coloring of the Q below:
(Determining the winding number of the 11x11 square centered at the origin.)

In[]:=

Grid[{{ParametricPlot[rectangleParametric[{{-5.5,-5.5},{5.5,5.5}},tt],{tt,0,1},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],PlotStyleThickness[0.01],ColorFunctionScalingFalse,PlotRange{{-8,8},{-8,8}},AspectRatio1,AxesFalse,ImageSize350,FrameTrue,PlotStyle->Directive[CapForm["Round"],Thickness[0.01]],PlotLabelStyle["Domain Coloring of Q on Large Boundary",Black]],Framed["Total Winding Number: "<>ToString[Round[findWR[Q,{{-5.5,-5.5},{5.5,5.5}},1,2500]]]]}}]

Out[]=

Total Winding Number: 0

Although this large boundary has a winding number of zero, we will discover shortly that the region it encompasses surrounds at least one zero. This phenomenon is created by two regions with nonzero winding numbers, one positive and one negative, cancelling the other’s value. The process of winding number addition can be depicted below, where it is important to remember that the hueMapVal of every boundary is evaluated clockwise. The red and blue colors only identify the winding orientation and do not represent the domain coloring.
(Animation of the winding number addition of two rectangular boundaries using a table of points.)

In[]:=

Animate[Show[{ParametricPlot[rectangleParametric[{{-4,-3},{0,3}},tt],{tt,0,1},PlotStyle{Thickness[0.01],Red,Opacity[0.2]},PlotRange{{-5,5},{-4,4}},AspectRatio1,AxesFalse,ImageSize350,FrameTrue,PlotStyle->Directive[CapForm["Round"],Thickness[0.01]],PlotLabelStyle["Visualizing the Addition of Boundaries",Black],Epilog{PointSize[0.02],Red,Table[Point[rectangleParametric[{{-4,-3},{0,3}},Mod[t+n,1]]],{n,0,1,0.05}],Blue,Table[Point[rectangleParametric[{{0,3},{4,-3}},Mod[t+n,1]]],{n,0,1,0.05}]}],ParametricPlot[rectangleParametric[{{0,3},{4,-3}},tt],{tt,0,1},PlotStyle{Thickness[0.01],Blue,Opacity[0.2]}]}],{t,0,1},AnimationRate0.05]

Out[]=

In the above animation, it is clear that the winding of the shared edge will cancel when the boundaries and combined. Therefore, the winding number of the large boundary is equivalent to the sum of the winding number of the two smaller boundaries.
(Computing the winding numbers of the red, blue, and red-blue boundaries from above.)

In[]:=

Grid[{{Column[{ParametricPlot[rectangleParametric[{{-4,-3},{0,3}},tt],{tt,0,1},PlotStyleThickness[0.01],ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],ColorFunctionScalingFalse,PlotRange{{-5,5},{-4,4}},AspectRatio1,AxesFalse,ImageSize300,FrameTrue,PlotStyle->Directive[CapForm["Round"],Thickness[0.01]],PlotLabelStyle["Red Boundary",Black]],Framed["Winding Number: "<>ToString[Round[findWR[Q,{{-4,-3},{0,3}},1,2500]]],FrameStyleRed]},Center],Column[{ParametricPlot[rectangleParametric[{{0,3},{4,-3}},tt],{tt,0,1},PlotStyleThickness[0.01],ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],ColorFunctionScalingFalse,PlotRange{{-5,5},{-4,4}},AspectRatio1,AxesFalse,ImageSize300,FrameTrue,PlotStyle->Directive[CapForm["Round"],Thickness[0.01]],PlotLabelStyle["Blue Boundary",Black]],Framed["Winding Number: "<>ToString[Round[findWR[Q,{{0,-3},{4,3}},1,2500]]],FrameStyleBlue]},Center],Column[{ParametricPlot[rectangleParametric[{{-4,-3},{4,3}},tt],{tt,0,1},PlotStyleThickness[0.01],ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],ColorFunctionScalingFalse,PlotRange{{-5,5},{-4,4}},AspectRatio1,AxesFalse,ImageSize300,FrameTrue,PlotStyle->Directive[CapForm["Round"],Thickness[0.01]],PlotLabelStyle["Red-Blue Combined Boundary",Black]],Framed["Winding Number: "<>ToString[Round[findWR[Q,{{-4,-3},{4,3}},1,2500]]],FrameStylePurple]},Center]}}]

Out[]=

Winding Number: 1

Winding Number: 2

Generating New Coordinates

Based on what we discovered in the previous section, the first step of our analysis requires finding a boundary with a nonzero winding number. Next, we need to halve the boundary and determine which half still has the desired nonzero winding number. This process is repeated until we obtain the desired boundary area, which results in a small region that must contain a zero. At this stage, we can approximate the exact location of the zero by using the center of the rectangle.

However, a rectangular boundary shape means that halving the boundary repeatedly using only one axis will result in a very elongated strip. This is undesirable because we will only be able to approximate one coordinate of the zero while the other remains unknown. Thus, our halving process will alternate between the axes. In effect, we need to create new functions: newXCoords and newYCoords that halve the given boundary and test the resulting winding numbers, along with quarterArea that combines both newXCoords and newYCoords to return a rectangle that is only a quarter of the original boundary area. For instance, we can visualize newXCoords as slicing the given region halfway along the y-axis and selecting which half we should continue investigating. Once a sufficiently small region has been found, we need to determine its midpoint using the appropriately named rectangleMidpoint.
(The described functions are defined below.)

In[]:=

newXCoords[Functn_,{{x1_,y1_},{x2_,y2_}},num_:2500]:=If0!=RoundfindWRFunctn,{x1,y1},

x1+x2

,y2,1,num,{x1,y1},

x1+x2

,y2,

x1+x2

,y1,{x2,y2};newYCoords[Functn_,{{x1_,y1_},{x2_,y2_}},num_:2500]:=If0!=RoundfindWRFunctn,{x1,y1},x2,

y1+y2

,1,num,{x1,y1},x2,

y1+y2

,x1,

y1+y2

,{x2,y2};

In[]:=

quarterArea[Functn_,{{x1_,y1_},{x2_,y2_}},num_:2500]:=newXCoords[Functn,newYCoords[Functn,{{x1,y1},{x2,y2}},num],num];

In[]:=

rectangleMidpoint[{{x1_,y1_},{x2_,y2_}}]:=

x1+x2

y1+y2

;

Although these functions initially evaluate using 2500 intervals, the examples in the next section only required 500 intervals for the same level of accuracy.

To apply quarterArea, we must first generate a rectangular boundary with nonzero winding number. We will select the square region defined by corners

{-4,-4}

and

{4,4}

.
(Confirming this region satisfies our winding number requirement.)

In[]:=

Framed["Winding number on initial boundary: "<>ToString[Round[findWR[Q,{{-4,-4},{4,4}},1,500]]]]

Out[]=

Winding number on initial boundary: 1

Recall that rounding is used because findWR gives an approximation of an integer-valued winding number. Next, as seen below, both newXCoords and newYCoords produce boundaries that still have nonzero winding numbers:
(Testing the outputs of the new functions.)

In[]:=

Grid[{{Framed[Column[{"New y-cut region: "<>ToString[newXCoords[{{-4,-4},{4,4}}]],"New winding number: "<>ToString[Round[findWR[Q,newXCoords[Q,{{-4,-4},{4,4}}],1,500]]]}]],Framed[Column[{"New x-cut region: "<>ToString[newYCoords[{{-4,-4},{4,4}}]],"New winding number: "<>ToString[Round[findWR[Q,newYCoords[Q,{{-4,-4},{4,4}}],1,500]]]}]]}}]

Out[]=

New y-cut region: newXCoords[{{-4, -4}, {4, 4}}]

New winding number: 1

New x-cut region: newYCoords[{{-4, -4}, {4, 4}}]

New winding number: 1

Zeroing In on Zeros

Nested application of quarterArea quickly narrows our boundary around a zero. To determine when to stop evaluating, we shall define an arbitrary minimum required boundary area of

-10

. Only 50 iterations were necessary for the below examples:
(A list of rectangular boundary corners found by nesting quarterArea until the area fell below the cutoff value.)

In[]:=

examplePoints=NestWhileList[quarterArea[Q,#,50]&,{{-4,-4},{4,4}},Area[Rectangle[First@#,Last@#]]>

-10

&];Framed[Column[{"List of boundary corners: ",N[examplePoints]}]]

Out[]=

List of boundary corners:

{{{-4.,-4.},{4.,4.}},{{-4.,-4.},{0.,0.}},{{-2.,-2.},{0.,0.}},{{-1.,-1.},{0.,0.}},{{-1.,-1.},{-0.5,-0.5}},{{-0.75,-1.},{-0.5,-0.75}},{{-0.75,-0.875},{-0.625,-0.75}},{{-0.75,-0.875},{-0.6875,-0.8125}},{{-0.75,-0.875},{-0.71875,-0.84375}},{{-0.75,-0.875},{-0.734375,-0.859375}},{{-0.742188,-0.867188},{-0.734375,-0.859375}},{{-0.738281,-0.863281},{-0.734375,-0.859375}},{{-0.738281,-0.863281},{-0.736328,-0.861328}},{{-0.738281,-0.862305},{-0.737305,-0.861328}},{{-0.737793,-0.862305},{-0.737305,-0.861816}},{{-0.737793,-0.862061},{-0.737549,-0.861816}},{{-0.737793,-0.862061},{-0.737671,-0.861938}},{{-0.737732,-0.862},{-0.737671,-0.861938}},{{-0.737701,-0.862},{-0.737671,-0.861969}},{{-0.737686,-0.862},{-0.737671,-0.861984}},{{-0.737679,-0.861992},{-0.737671,-0.861984}}}

The above computation gives an estimation of a zero at the following coordinates:
(Find the midpoint of the last region using rectangleMidpoint.)

In[]:=

Framed["Zero estimate: "<>ToString[N[rectangleMidpoint[Last@examplePoints]]]]

Out[]=

Zero estimate: {-0.737675, -0.861988}

We can determine the accuracy of this estimation by evaluating the output of Q and its euclidean norm at this coordinate. Ideally, the norm value of the zero estimate is close to 0.
(Plugging in the estimated zero point into the Q function, plus its norm.)

In[]:=

Column[{Framed["Output of Q at zero estimate: "<>ToString[N[Q[rectangleMidpoint[Last@examplePoints]]]]],Framed["Norm of Q at zero estimate: "<>ToString[N[Norm[Q[rectangleMidpoint[Last@examplePoints]]]]]]}]

Out[]=

Output of Q at zero estimate: {-0.0000165564, 0.0000164511}

Norm of Q at zero estimate: 0.0000233399

This is very accurate!
We visually depict the utility of quarterArea by plotting the rectangular search boundaries it outputs as the function is repeatedly applied. As the name implies, we can see that each boundary is a quarter of the previous iteration. The norm of the output of Q at each rectangle’s midpoint is also included to check that the accuracy increases.
(Graphing the outputs of quarterArea based on the iteration number. Minimum required boundary area is

-10

In[]:=

domainColor=ParametricPlot[{tx,ty},{tx,-5,5},{ty,-5,5},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],ColorFunctionScalingFalse,PlotStyleOpacity[0.05],PlotPoints200];Animate[Column[{Show[Graphics[Style[Point[rectangleMidpoint[examplePoints[[n]]]],Black,PointSize[0.01]],PlotRange{{-5,5},{-5,5}},AspectRatio1,AxesFalse,ImageSize350,FrameTrue,PlotLabelStyle["Process of quarterArea Boundary Reduction",Black]],domainColor,ParametricPlot[rectangleParametric[examplePoints[[n]],tt],{tt,0,1},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],PlotStyleThickness[0.005],ColorFunctionScalingFalse,PlotStyle->Directive[CapForm["Round"],Thickness[0.01]]]],Framed[Column[{"Number of iterations: "<>ToString[n-1],"Midpoint (in black): "<>ToString[N[rectangleMidpoint[examplePoints[[n]]]]],"Norm of Q at midpoint: "<>ToString[N[Norm[Q[rectangleMidpoint[examplePoints[[n]]]]]]]}]]},Center],{{n,1,"Process timeline:"},1,Length[examplePoints],1},AnimationRate1,AnimationRunningFalse]

Out[]=

The process of shrinking boundaries can be summarized in a single image:
(Graphing all quarterArea boundaries.)

In[]:=

Show[ParametricPlot[{tx,ty},{tx,-5,5},{ty,-5,5},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],ColorFunctionScalingFalse,PlotStyleOpacity[0.05],PlotPoints200,PlotRange{{-5,5},{-5,5}},AspectRatio1,AxesFalse,ImageSize350,FrameTrue,PlotLabelStyle["Finalized quarterArea Boundary Reduction",Black]],ParametricPlot[rectangleParametric[#,tt]&/@examplePoints,{tt,0,1},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],PlotStyleThickness[0.005],ColorFunctionScalingFalse,PlotStyle->Directive[CapForm["Round"],Thickness[0.01]]]]

Out[]=

Remember that darker hues correlate to points that Q sends near the origin. Thus, darker shading corresponds to smaller norm values.
As desired, the norm value decreases quickly. To investigate how quickly, we can create the following plots:
(Plots of the norm value against the number of iterations, one with log scaling.)

In[]:=

Grid[{{ListLinePlot[Norm[Q[rectangleMidpoint[#]]]&/@examplePoints,ImageSize450,PlotMarkersAutomatic,PlotStyleBlack,PlotRangeAll,PlotLabelStyle["Norm of Q at Each Midpoint",Black],AxesLabel{"Iteration #","Norm of Q"},DataRange{0,Length[examplePoints]-2}],ListLogPlot[Norm[Q[rectangleMidpoint[#]]]&/@examplePoints,ImageSize450,PlotStyleBlack,PlotRangeAll,DataRange{0,Length[examplePoints]-2},PlotLabelStyle["Norm of Q at Each Midpoint (logrithmic scale)",Black],AxesLabel{"Iteration #","Log of Norm of Q"}]}}]

Out[]=

Notice that the logarithmically scaled plot is linear, implying that the estimation’s error decreases exponentially!
For a different initial boundary, such as one with corners

{0,0}

and

{4,4}

, we see similar results:
(Repeating the quarterArea investigation with a new initial boundary.)

In[]:=

Framed["Winding number of initial boundary: "<>ToString[Round[findWR[Q,{{0,0},{4,4}},1,1500]]]]

Out[]=

Winding number of initial boundary: 1

In[]:=

domainColorSmall=ParametricPlot[{tx,ty},{tx,-1,5},{ty,-1,5},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],ColorFunctionScalingFalse,PlotStyleOpacity[0.05],PlotPoints200];examplePointsTwo=NestWhileList[quarterArea[Q,#,50]&,{{0,0},{4,4}},Area[Rectangle[First@#,Last@#]]>

-10

&];Animate[Module[{point},point=examplePointsTwo[[n]];Grid[{{Show[Graphics[Style[Point[rectangleMidpoint[point]],Black,PointSize[0.01]],PlotRange{{-1,5},{-1,5}},AspectRatio1,AxesFalse,ImageSize350,FrameTrue,PlotLabelStyle["Process of quarterArea Boundary Reduction",Black]],domainColorSmall,ParametricPlot[rectangleParametric[point,tt],{tt,0,1},ColorFunctionFunction[{x,y},hueMap[Q[{x,y}]]],PlotStyleThickness[0.005],ColorFunctionScalingFalse,PlotStyle->Directive[CapForm["Round"],Thickness[0.01]]]],ListLogPlot[Norm[Q[rectangleMidpoint[#]]]&/@examplePointsTwo,ImageSize410,PlotStyleLightGray,PlotRangeAll,AspectRatio1,DataRange{0,Length[examplePoints]-2},PlotLabelStyle["Norm of Q at Each Midpoint (logrithmic scale)",Black],AxesLabel{"Iteration #","Log of Norm of Q"},EpilogStyle[Point[{n-1,Log[Norm[Q[rectangleMidpoint[point]]]]}],Black,PointSize[0.015]]],Framed[Column[{"Number of iterations: "<>ToString[n-1],"Estimation (in black): "<>ToString[N[rectangleMidpoint[point]]],"Output of Q at estimation: "<>ToString[DecimalForm[N[Q[rectangleMidpoint[point]]],3]],"Norm of Q at estimation: "<>ToString[N[Norm[Q[rectangleMidpoint[point]]]]]}]]}}]],{{n,1,"Process timeline:"},1,Length[examplePoints]-1,1},AnimationRate1,AnimationRunningFalse]

Out[]=

It should be noted that multiple (potentially infinite) initial boundaries can converge to the same zero:
(Confirming winding number and finding zero estimate on new boundary corners

{0,0}

and

{3,3}

In[]:=

Framed[Column[{"Initial boundary corners: {0,0} and {3,3}",StringJoin["Winding number of initial boundary: ",ToString[Round[findWR[Q,{{0,0},{3,3}},1,1500]]]],"Zero estimate: "<>ToString[N[rectangleMidpoint[Last@NestWhileList[quarterArea[Q,#,50]&,{{0,0},{3,3}},Area[Rectangle[First@#,Last@#]]>

-10

&]]]],"Norm of Q at zero estimate: "<>ToString[Norm[Q[{0.8728952407836914`,2.0376405715942383`}]]]}]]

Out[]=

Initial boundary corners: {0,0} and {3,3}

Winding number of initial boundary: 1

Zero estimate: {0.872895, 2.03764}

Norm of Q at zero estimate: 0.0000180174

One function can be created to culminate the entire process of domain coloring, winding number computation, and boundary reduction:
(Defining a function that outputs the zero estimate for a given initial boundary.)

In[]:=

zeroEstimate[Functn_,{point1_,point2_},num_:2500,sizeLimit_:

-10

]:=rectangleMidpoint[NestWhile[quarterArea[Functn,#,num]&,{point1,point2},Area[Rectangle[First@#,Last@#]]>sizeLimit&]]

This function generates the same zero estimations as found above:
(Confirming that zeroEstimate achieves the same results as the above processes.)

In[]:=

Framed[Column[{"Zero estimate of {-4,-4} to {4,4} boundary: "<>ToString[N[zeroEstimate[Q,{{-4,-4},{4,4}},50]]],"Zero estimate of {0,0} to {4,4} boundary: "<>ToString[N[zeroEstimate[Q,{{0,0},{4,4}},50]]]}]]

Out[]=

Zero estimate of {-4,-4} to {4,4} boundary: {-0.737675, -0.861988}

Zero estimate of {0,0} to {4,4} boundary: {0.87289, 2.03764}

Concluding Remarks

Future Endeavours

Although continuousness is required, the given transformation does not need be differentiable. For instance, we can define the following transformation that is not differentiable when

x=0

In[]:=

A[{x_,y_}]={3Sin[Abs[x]],2y};

Let us select a boundary with a nonzero winding number and find the respective domain coloring:
(Plotting the domain coloring of A.)

In[]:=

Column[{Show[ParametricPlot[{x,y},{x,-5,5},{y,-5,5},ColorFunctionFunction[{x,y},hueMap[A[{x,y}]]],ColorFunctionScalingFalse,FrameTrue,AxesFalse,ImageSize350,AspectRatio1,PlotLabelStyle["Domain Coloring of A",Black],PlotPoints200],ParametricPlot[rectangleParametric[{{-2,-3},{4,4}},t],{t,0,1},PlotStyle{Dashed,Gray}]],Framed["Winding number of gray boundary: "<>ToString[Round[findWR[A,{{-2,-3},{4,1}},1,2500]]]]}]

Out[]=

Winding number of gray boundary: 1

Even though the boundary passes through two non-differentiable points, it is still possible to find a zero in this region:
(Zero estimation for the boundary on A.)

In[]:=

Module[{point},point=zeroEstimate[A,{{-2,-3},{4,4}},500];Framed[Column[{StringForm["Zero estimate: ``",N[point]],StringForm["Norm of A at zero estimate: ``",Round[Norm[A[point]],

-8

]]}]]]

Out[]=

Zero estimate: {3.14159,-1.43051×

-6

}

Norm of A at zero estimate:

151

50000000

By revisiting the definition of A, note that a zero exists at

{π,0}

. For future exploration, one could investigate a 2-dimensional transformation that is not differentiable at any point. Such a transformation could be created by modifying the Weierstrass Function.

Other areas of exploration include finding winding numbers of non-rectangular boundaries. We chose to investigate rectangular boundaries for their convenience, but any closed curve has a winding number. Future endeavors would involve redefining findWR for a curve that is formed from any number of points.

Additionally, recall that the newXCoords, newYCoords, and quarterArea functions require the initial boundary to have a nonzero winding number. These functions also only produce the coordinates of one boundary, even though both boundaries might have nonzero winding numbers. Addressing these two issues would improve the reliability of finding all zeros in a given region.

Potential Issues

The methodology that led to estimateZeros is not a perfect process. For instance, lots of evaluations are required for finding a winding number depending on the level of accuracy. The nature of winding number addition also leads to instances where an initial region has a winding number of zero but encloses multiple zeros of the transformation. The function FindRoot is faster in many situations, but its reliance on Newton’s method may fail; for initial coordinates

{0,0}

, FindRoot cannot find the roots of Q whereas the winding number process succeeds.

Not all initial boundaries will produce desirable results. For instance, findWR determines that a rectangle on the domain coloring of A with corners

{3,0}

and

{4,1}

has a noninteger winding number:
(Winding number of a specific boundary on A.)

In[]:=

Framed[StringForm["Winding number: ``",findWR[A,{{3,0},{4,1}},1,2500]]]

Out[]=

Winding number: 0.5

How is this possible? Plotting the domain coloring around this boundary makes the issue clear: the boundary passes through a zero of the transformation at

{π,0}

.
(Domain coloring of A around the specific boundary. The zero is highlighted in black.)

In[]:=

Show[ParametricPlot[{x,y},{x,2.9,4.1},{y,-0.1,1.1},ColorFunctionFunction[{x,y},Hue[hueMapVal[A[{x,y}]]]],ColorFunctionScalingFalse,FrameTrue,AxesFalse,ImageSize350,AspectRatio1,PlotLabelStyle["Domain Coloring of A",Black],PlotPoints200,Epilog->Style[Point[{π,0}],Black,PointSize[0.015]]],ParametricPlot[rectangleParametric[{{3,0},{4,1}},t],{t,0,1},PlotStyle{Dashed,Gray}]]

Out[]=

Therefore, a boundary that passes through a zero does not produce a sensible winding number. Although the objective is to find zeros, one must be careful that the search boundaries surround the zero instead of intersecting it. For this reason, any boundary with a non-integer winding number should be thoroughly investigated.

Conclusion

Out of the many methods available to find zeros of 2D transformations, none are as visually stunning as domain coloring and winding numbers. By taking an unorthodox approach to zero estimation, we have revealed new ways to visualize 2D transformations, uncover complexities of relatively simplistic concepts, and estimate zeros with an exponentially-decaying error. This process brings more color to the world of math and highlights how the same result can be approached in a myriad of ways. Many more avenues of exploration remain for this intriguing and beautiful method, and it is very exciting to see what awaits.

Keywords

◼

2D transformation:

A mapping from



→



of the form

{x,y}→{f(x,y),g(x,y)}

for some functions f and g.For example, we primarily use

Q(x,y)=

-5

-2xy+Cos[y]+3,

+y(x-3)+2y-Sin[y]-3

above.

◼

Zero:

Most occurrences of a “zero” are referring to a point in 2-dimensional space that maps to

{0,0}

when the given transformation is applied. Many processes create a “zero estimation”, which is a point that approximates the actual value of this type of “zero”. Sometimes referred to as a “root”.

◼

Color wheel:

Otherwise referred to as a “hue wheel”, the color wheel is depicted as the full spectrum of colors arranged in ROYGBIV order around the perimeter of a circle. In this document, we aligned the colors such that red is in the southern region and the remaining colors are arranged in counter-clockwise order. See Section 1 for a visualization.

◼

Domain coloring:

Every point on the



plane takes the hue of where that point’s transformation lands on a color wheel centered at the origin. To paraphrase Grant Sanderson, we can also interpret domain coloring as the process of applying the 2D transformation, ‘painting’ the plane with the color wheel, and reversing the transformation.

◼

Boundary:

In general, a “boundary” refers to the outline of a 2-dimensional closed curve with 1-dimensional width. Our exploration focuses on boundaries of rectangular shape.

◼

Winding number:

In this context, a winding number of a given boundary represents the number of how many loops of the color spectrum are completed when moving clockwise along the boundary. The direction moved along the boundary has no negative effect on the uses of the winding number. Likewise, positive winding numbers correspond to moving in the “ROYGBIV” spectrum order.

Acknowledgments

In no particular order, I would like to acknowledge those who inspired, encouraged, and assisted me throughout the process of building this document:

◼

My mentor, Clayton Shonkwiler, give very useful feedback, suggested innovative visualizations, and provided crucial code assistance. His expertise in Mathematica visualizations was a huge inspiration and crucial for the many pitfalls I encountered.

◼

Although Paul Abbott started my journey in the Wolfram Summer School program, his involvement did not stop there. His dedication to my success drove me to achieve something beyond what I thought I could accomplish. His supply of tips, tools, and resources kept me both excited about my project and the world of math as a whole.

◼

I appreciate Erin Cherry’s patience, tenacity, and commitment to making the WSS run smoothly. She kept everyone on track and made sure to remind me of important deadlines. Many important events would have been missed without her guidance.

◼

It was a pleasure to meet Stephen Wolfram and hear his feedback on my project topic. Dr. Wolfram inspired me during my educational journey and it was amazing to participate in his program.

◼

Branko Ćurgus made me the mathematician I am today. He taught me many fundamental math concepts, one of which being the Wolfram Language. Without his guidance, I would not have had the opportunity to participate in the WSS.

◼

I would like to thank my girlfriend Samantha Kirlin and my parents for listening, assisting, and supporting me throughout the duration of this project.

References

◼

Grant Sanderson’s (3Blue1Brown) YouTube video: https://youtu.be/b7FxPsqfkOY

Cite this as: Garrett Dubofsky, "Finding Zeros with Domain Coloring and Winding Numbers" from the Notebook Archive (2021), https://notebookarchive.org/2021-07-60wgb9i

Finding Zeros with Domain Coloring and Winding Numbers

Table of Contents

2D Transformations and Domain Coloring

The Goal

Visualizing Four Dimensions

Domain Coloring

Determining the Winding Number

Defining a Boundary

What is Winding?

Finding the Winding Number

Decreasing the Boundary

Properties of the Winding Number

Generating New Coordinates

Zeroing In on Zeros

Concluding Remarks

Future Endeavours

Potential Issues

Conclusion

Keywords

Acknowledgments

References