92

I found this plot on Wikipedia:

domain-colored plot of sine function

Domain coloring of $\sin(z)$ over $(-\pi,\pi)$ on $x$ and $y$ axes. Brightness indicates absolute magnitude, saturation represents imaginary and real magnitude.

Despite following the link and reading the page nothing I have tried is giving me the result shown. How should this be done?

J. M.'s missing motivation
  • 124,525
  • 11
  • 401
  • 574
Mr.Wizard
  • 271,378
  • 34
  • 587
  • 1,371

5 Answers5

109

Building on Heike's ColorFunction, I came up with this:

enter image description here

The white bits are the trickiest - you need to make sure the brightness is high where the saturation is low, otherwise the black lines appear on top of the white ones.

The code is below. The functions defined are:

  • complexGrid[max,n] simply generates an $n\times n$ grid of complex numbers ranging from $-max$ to $+max$ in both axes.

  • complexHSB[Z] takes an array $Z$ of complex numbers and returns an array of $\{h,s,b\}$ values. I've tweaked the colour functions slightly. The initial $\{h,s,b\}$ values are calculated using Heike's formulas, except I don't square $s$. The brightness is then adjusted so that it is high when the saturation is low. The formula is almost the same as $b2=\max (1-s,b)$ but written in a way that makes it Listable.

  • domainImage[func,max,n] calls the previous two functions to create an image. func is the function to be plotted. The image is generated at twice the desired size and then resized back down to provide a degree of antialiasing.

  • domainPlot[func,max,n] is the end user function which embeds the
    image in a graphics frame.


complexGrid = Compile[{{max, _Real}, {n, _Integer}}, Block[{r},
    r = Range[-max, max, 2 max/(n - 1)];
    Outer[Plus, -I r, r]]];

complexHSB = Compile[{{Z, _Complex, 2}}, Block[{h, s, b, b2},
    h = Arg[Z]/(2 Pi);
    s = Abs[Sin[2 Pi Abs[Z]]];
    b = Sqrt[Sqrt[Abs[Sin[2 Pi Im[Z]] Sin[2 Pi Re[Z]]]]];
    b2 = 0.5 ((1 - s) + b + Sqrt[(1 - s - b)^2 + 0.01]);
    Transpose[{h, Sqrt[s], b2}, {3, 1, 2}]]];

domainImage[func_, max_, n_] := ImageResize[ColorConvert[
    Image[complexHSB@func@complexGrid[max, 2 n], ColorSpace -> "HSB"],
    "RGB"], n, Resampling -> "Gaussian"];

domainPlot[func_: Identity, max_: Pi, n_: 500] :=
  Graphics[{}, Frame -> True, PlotRange -> max, RotateLabel -> False, 
   FrameLabel -> {"Re[z]", "Im[z]", 
      "Domain Colouring of " <> ToString@StandardForm@func@"z"},
   BaseStyle -> {FontFamily -> "Calibri", 12},
   Prolog -> Inset[domainImage[func, max, n], {0, 0}, {Center, Center}, 2` max]];

domainPlot[Sin, Pi]

Other examples follow:

It's informative to plot the untransformed complex plane to understand what the colours indicate:

domainPlot[]

enter image description here

A simple example:

domainPlot[Sqrt]

enter image description here

Plotting a pure function:

domainPlot[(# + 2 I)/(# - 1) &]

enter image description here

I think this one is very pretty:

domainPlot[Log]

enter image description here

Simon Woods
  • 84,945
  • 8
  • 175
  • 324
  • Thanks for the accept. It's still not quite right - I have a feeling the original might have gone through more sophisicated processing than just computing HSB values from function samples. – Simon Woods Jun 27 '12 at 09:04
  • You're welcome. Re: processing, quite possibly, however I wasn't looking for an exact replication. Yours is a complete answer as far as I'm concerned. – Mr.Wizard Jun 27 '12 at 09:15
  • Lovely images! Trying to imagine where to use these for a proper Texture[]. – Yves Klett Jun 27 '12 at 14:12
  • 3
    +1, brilliant! In addition to just being lovely images, I love how the poles and branch cuts are clearly visible. – rcollyer Jun 27 '12 at 14:56
  • 2
    Incidentally, I posted this on Facebook. The response I got: "Holy crap, mathematica keeps getting more and more awesome." – rcollyer Jun 27 '12 at 15:53
  • Really nice work on the update. I'm sure this will earn a "Good Answer" badge soon. – Mr.Wizard Jun 27 '12 at 19:55
  • Simon, I'm having some trouble using domainPlot: I keep getting Coordinate 6.283185307179586 should be a pair of numbers, or a Scaled or Offset form. Since I'm seeing 2 Pi this looks like it's coming from Inset; has the syntax for Inset changed in version 8? – Mr.Wizard Jun 27 '12 at 20:09
  • Strange, it works if I use domainPlot[Log, N@Pi, 500] but fails if I use domainPlot[Log, Pi, 500] ; I'm not seeing why that should be. – Mr.Wizard Jun 27 '12 at 20:13
  • @Mr.Wizard, the documentation says that Inset was last modified in version 7, so perhaps there was a bug that was fixed. Can you confirm that using N[2 max] for the last argument of Inset remedies the problem in version 7? – Simon Woods Jun 27 '12 at 22:05
  • Simon, yes, using that, or as I prefer 2` max fixes it for version 7. Also for version 7 the option Resampling->"Gaussian" cannot be used but it looks fine without it. – Mr.Wizard Jun 28 '12 at 04:23
  • 2
    Thanks for that, I've updated the code to work in version 7. – Simon Woods Jun 28 '12 at 09:43
  • 2
    @rcollyer, that's great! Your friend has excellent taste. My friends would just think I was weird if I posted domain colouring plots on Facebook :-) – Simon Woods Jun 28 '12 at 09:44
  • The update is terrific! – rm -rf Jun 28 '12 at 14:22
  • I am wondering how could you work out such wonderful picture^-^ – withparadox2 Aug 09 '12 at 06:31
  • I would upvote this twice, if I could, just for its beauty. – Renan Sep 21 '12 at 17:31
  • 1
    On Mathematica V.9, the complexGrid[] function was wrong for me, the Outer[] was not giving the expected result (it puts z with Re[z]<0 and Im[z]<0 in the upper left quadrant). See e.g. the following code: Outer[Plus, I {-1, 0, 1}, {-1, 0, 1}] // MatrixForm I replaced it with Outer[Plus, r, I r // Reverse]\[Transpose] which works for me. –  Jan 31 '13 at 13:19
  • @guerom00, well spotted, thanks! A simpler correction is Outer[Plus, -I r, r] – Simon Woods Jan 31 '13 at 13:35
  • @SimonWoods by any chance, do you know what transform turns a unit circle into the real number line? – rcollyer Apr 16 '13 at 17:01
  • @rcollyer, I'm not sure what you mean. There must be many ways to map from the complex unit circle to the real number line, e.g. Re[z]/(1-Im[z]) but I don't get how it relates to this answer? – Simon Woods Apr 18 '13 at 09:23
  • @SimonWoods it doesn't have any relationship to the answer. I was grasping at straws for another problem, and I thought such a transformation might help. (Further investigation suggests not.) So, remembering this answer, I thought you might know. Sorry for any confusion. – rcollyer Apr 18 '13 at 11:45
  • Simon, I notice that presently the color orientation of the first graphic no longer matches that of the output of domainPlot[Sin] -- do you think this is something that should be corrected? – Mr.Wizard Mar 08 '16 at 12:49
  • @Mr.Wizard, I'll update it. The original contained two errors and I only fixed one. – Simon Woods Mar 08 '16 at 19:27
36

Not as pretty as the one in the original post, but it's getting in the right direction I think:

RegionPlot[True,
 {x, -Pi, Pi}, {y, -Pi, Pi},
 ColorFunction -> (Hue[Rescale[Arg[Sin[#1 + I #2]], {-Pi, Pi}],
     Sin[2 Pi Abs[Sin[#1 + I #2]]]^2,
     Abs@(Sin[Pi Re[Sin[#1 + I #2]]] Sin[Pi Im[Sin[#1 + I #2]]])^(1/
        4), 1] &),
 ColorFunctionScaling -> False, PlotPoints -> 200]

Mathematica graphics

It seems that the hue of the colour function is a function of Arg[Sin[z]], saturation is a function of Abs[Sin[z]] and the brightness is related to Re[Sin[z]] and Im[Sin[z]].

Heike
  • 35,858
  • 3
  • 108
  • 157
24

This is a good way :

DensityPlot[ Rescale[ Arg[Sin[-x - I y]], {-Pi, Pi}], {x, -Pi, Pi}, {y, -Pi, Pi}, 
             MeshFunctions -> Function @@@ {{{x, y, z}, Re[Sin[x + I y]]}, 
                                            {{x, y, z}, Im[Sin[x + I y]]},
                                            {{x, y, z}, Abs[Sin[x + I y]]}}, 
             MeshStyle -> {Directive[Opacity[0.8], Thickness[0.001]], 
                           Directive[Opacity[0.7], Thickness[0.001]], 
                           Directive[White, Opacity[0.3], Thickness[0.006]]}, 
             ColorFunction -> Hue, Mesh -> 50, Exclusions -> None, PlotPoints -> 100]

enter image description here

Another ways to tackle the problem, which apprears promising.

ContourPlot[ Evaluate @ {Table[Re @ Sin[x + I y] == 1/2 k, {k, -25, 25}], 
                         Table[Im @ Sin[x + I y] == 1/2 k, {k, -25, 25}]}, 
             {x, -Pi, Pi}, {y, -Pi, Pi}, PlotPoints -> 100, MaxRecursion -> 5]

enter image description here

and

RegionPlot[ Evaluate @ {Table[1/2 (k + 1) > Re @ Sin[x + I y] > 1/2 k, {k, -25, 25}],
                        Table[1/2 (k + 1) > Im @ Sin[x + I y] > 1/2 k, {k, -25, 25}]},
            {x, -Pi, Pi}, {y, -Pi, Pi}, PlotPoints -> 50, MaxRecursion -> 4, 
            ColorFunction -> Function[{x, y}, Hue[Re@Sin[x + I y]]]]

enter image description here

These plots seem to be good points for further playing around to get better solutions.

Artes
  • 57,212
  • 12
  • 157
  • 245
23

I already mentioned Bernd Thaller's package Graphics`ComplexPlot` in the comments; if one blends the ideas from Artes's and Heike's answers, and then use the function $ComplexToColorMap[] from Thaller's package (I won't include it here; again, see the package for that), we get this:

domain-colored plot

Needs["Graphics`ComplexPlot`"] (* Thaller's package; get it yourself *)

f1 = RegionPlot[True, {x, -Pi, Pi}, {y, -Pi, Pi}, 
 ColorFunction -> ($ComplexToColorMap[Abs[Sin[#1 + I #2]], 
     Arg[Sin[#1 + I #2]], {Pi, 1/10, 1, 1/10, 1}] &), 
 ColorFunctionScaling -> False, PlotPoints -> 200];

f2 = ContourPlot[
 Evaluate@{Table[Re@Sin[x + I y] == 1/2 k, {k, -25, 25}], 
   Table[Im@Sin[x + I y] == 1/2 k, {k, -25, 25}]}, {x, -Pi, 
  Pi}, {y, -Pi, Pi}, PlotPoints -> 100, ContourStyle -> Gray];

f3 = ContourPlot[
 Evaluate@Table[Abs@Sin[x + I y] == 1/2 k, {k, -25, 25}], {x, -Pi, 
  Pi}, {y, -Pi, Pi}, PlotPoints -> 100, ContourStyle -> White, 
 MaxRecursion -> 5];

Show[f1, f2, f3]

The $ComplexToColorMap[] function could probably be optimized a fair bit for new Mathematica, but I won't get into that for now. One might also consider tweaking the Opacity[] of the contour lines for the absolute value as well, but I'll leave that as an experiment for the reader.


Another thing you can try:

RegionPlot[True, {x, -Pi, Pi}, {y, -Pi, Pi}, 
 ColorFunction -> ($ComplexToColorMap[Abs[Sin[#1 + I #2]], 
     Arg[Sin[#1 + I #2]], {Pi, 1/50, 1, 1/50, 1}] &), 
 ColorFunctionScaling -> False, Mesh -> 51, 
 MeshFunctions -> {Re[Sin[#1 + I #2]] &, Im[Sin[#1 + I #2]] &}, 
 MeshStyle -> Gray, PlotPoints -> 95]
J. M.'s missing motivation
  • 124,525
  • 11
  • 401
  • 574
17

With Mathematica 12.0, there's now a ComplexPlot function that replaces user made solutions. As with other Plot functions, it allows us to specify a ColorFunction option to manipulate how to color the plot. This particular coloring is implemented natively in the "CyclicReImLogAbs" option.

So the modern equivalent is

ComplexPlot[Sin[z], {z, -Pi - Pi I, Pi + Pi I}, 
 ColorFunction -> "CyclicReImLogAbs", Frame -> False]

Plot Result

Mr.Wizard
  • 271,378
  • 34
  • 587
  • 1,371
lombardo2
  • 71
  • 2
  • 6