{"id":235,"date":"2012-08-06T12:32:42","date_gmt":"2012-08-06T17:32:42","guid":{"rendered":"https:\/\/blogs.mathworks.com\/cleve\/?p=235"},"modified":"2016-12-05T13:50:43","modified_gmt":"2016-12-05T18:50:43","slug":"can-one-hear-the-shape-of-a-drum-part-1-eigenvalues","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/cleve\/2012\/08\/06\/can-one-hear-the-shape-of-a-drum-part-1-eigenvalues\/","title":{"rendered":"Can One Hear the Shape of a Drum? Part 1, Eigenvalues"},"content":{"rendered":"&nbsp;\r\n<div class=\"content\"><!--introduction-->The title of this multi-part posting is also the title of a 1966 article by Marc Kac in the American Mathematical Monthly <a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">[1]<\/a>. This first part is about isospectrality.\r\n\r\n<!--\/introduction-->\r\n<h3>Contents<\/h3>\r\n<div>\r\n<ul>\r\n \t<li><a href=\"#1c0cb71f-6127-48af-b43d-d47e97ee9f04\">Isospectrality<\/a><\/li>\r\n \t<li><a href=\"#8116bc65-58f9-4361-b652-c4fd84229ec5\">Vertices<\/a><\/li>\r\n \t<li><a href=\"#f95b2356-f327-4e01-8624-d11253941088\">Finite difference grid<\/a><\/li>\r\n \t<li><a href=\"#13d23acf-93ea-4116-b982-03d6b9ddc6c4\">Finite difference Laplacian<\/a><\/li>\r\n \t<li><a href=\"#ea6b1e68-d22c-4b6a-b049-d460f734502e\">Compare eigenvalues<\/a><\/li>\r\n \t<li><a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">How about a proof?<\/a><\/li>\r\n \t<li><a href=\"#efc546e6-cc7b-4279-bedc-683e8056a689\">References<\/a><\/li>\r\n<\/ul>\r\n<\/div>\r\n<h4>Isospectrality<a name=\"1c0cb71f-6127-48af-b43d-d47e97ee9f04\"><\/a><\/h4>\r\nKac's article is not actually about a drum, which is three-dimensional, but rather about the two-dimensional drum head, which is more like a tambourine or membrane. The vibrations are modeled by the partial differential equation\r\n\r\n$$ \\Delta u + \\lambda u = 0 $$\r\n\r\nwhere\r\n\r\n$$ \\Delta u(x,y) = \\frac{\\partial^{2} u}{\\partial x^{2}} + \\frac{\\partial^{2} u}{\\partial y^{2}} $$\r\n\r\nThe boundary conditions are the key. Requiring $u(x,y) = 0$ on the boundary of a region in the plane corresponds to holding the membrane fixed on that boundary. The values of $\\lambda$ that allow nonzero solutions, the <i>eigenvalues<\/i>, are the squares of the frequencies of vibration, and the corresponding functions $u(x,y)$, the <i>eigenfunctions<\/i>, are the modes of vibration.\r\n\r\nThe MathWorks logo comes from this partial differential equation, on an L-shaped domain <a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">[2]<\/a>, <a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">[3]<\/a>, but this article is not about our logo.\r\n\r\nA region determines its eigenvalues. Kac asked about the opposite implication. If one specifies all of the eigenvalues, does that determine the region?\r\n\r\nIn 1991, Gordon, Webb and Wolpert showed that the answer to Kac's question is \"no\". They demonstrated a pair of regions that had different shapes but exactly the same infinite set of eigenvalues <a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">[4]<\/a>. In fact, they produced several different pairs of such regions. The regions are known as \"isospectral drums\". Wikipedia has a good background article on Kac's problem <a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">[5]<\/a>.\r\n\r\nI am interested in finite difference methods for membrane eigenvalue problems. I want to show that the finite difference operators on these regions have the same sets of eigenvalues, so they are also isospectral.\r\n\r\nI was introduced to isospectral drums by Toby Driscoll, a professor at the University of Delaware. A summary of Toby's work is available at his Web site <a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">[6]<\/a>. A reprint of his 1997 paper in SIAM Review is also available from his Web site <a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">[7]<\/a>. Toby developed methods, not involving finite differences, for computing the eigenvalues very accurately.\r\n\r\nThe isospectral drums are not convex. They have reentrant $270^\\circ$ corners. These corners lead to singularities in most of the eigenfunctions -- the gradients are unbounded. This affects the accuracy and the rate of convergence of finite difference methods. It is possible that for convex regions the answer to Kac's question is \"yes\".\r\n<h4>Vertices<a name=\"8116bc65-58f9-4361-b652-c4fd84229ec5\"><\/a><\/h4>\r\nI will look at the simplest known isospectral pair. The two regions are specified by the <i>xy<\/i>-coordinates of their vertices.\r\n<pre class=\"codeinput\">   drum1 = [0 0 2 2 3 2 1 1 0\r\n            0 1 3 2 2 1 1 0 0];\r\n   drum2 = [1 0 0 2 2 3 2 1 1\r\n            0 1 2 2 3 2 1 1 0];\r\n   vertices = {drum1,drum2};\r\n<\/pre>\r\nLet's first plot the regions.\r\n<pre class=\"codeinput\">   clf\r\n   shg\r\n   set(gcf,<span class=\"string\">'color'<\/span>,<span class=\"string\">'white'<\/span>)\r\n   <span class=\"keyword\">for<\/span> d = 1:2\r\n      <span class=\"comment\">% Plot the region.<\/span>\r\n      vs = vertices{d};\r\n      subplot(2,2,d)\r\n      plot(vs(1,:),vs(2,:),<span class=\"string\">'k.-'<\/span>);\r\n      axis([-0.1 3.1 -0.1 3.1])\r\n      axis <span class=\"string\">square<\/span>\r\n      title(sprintf(<span class=\"string\">'drum%d'<\/span>,d));\r\n   <span class=\"keyword\">end<\/span>\r\n<\/pre>\r\n<img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/drums_1_01.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/>\r\n<h4>Finite difference grid<a name=\"f95b2356-f327-4e01-8624-d11253941088\"><\/a><\/h4>\r\nI want to investigate simple finite difference methods for this problem. The MATLAB function <tt>inpolygon<\/tt> determines the points of a rectangular grid that are in a specified region.\r\n<pre class=\"codeinput\">   <span class=\"comment\">% Generate a coarse finite difference grid.<\/span>\r\n   ngrid = 5;\r\n   h = 1\/ngrid;\r\n   [x,y] = meshgrid(0:h:3);\r\n\r\n   <span class=\"comment\">% Loop over the two regions.<\/span>\r\n   <span class=\"keyword\">for<\/span> d = 1:2\r\n\r\n      <span class=\"comment\">% Determine points inside and on the boundary.<\/span>\r\n\r\n      vs = vertices{d};\r\n      [in,on] = inpolygon(x,y,vs(1,:),vs(2,:));\r\n      in = xor(in,on);\r\n\r\n      <span class=\"comment\">% Plot the region and the grid.<\/span>\r\n      subplot(2,2,d)\r\n      plot(vs(1,:),vs(2,:),<span class=\"string\">'k-'<\/span>,x(in),y(in),<span class=\"string\">'b.'<\/span>,x(on),y(on),<span class=\"string\">'k.'<\/span>)\r\n      axis([-0.1 3.1 -0.1 3.1])\r\n      axis <span class=\"string\">square<\/span>\r\n      title(sprintf(<span class=\"string\">'drum%d'<\/span>,d));\r\n\r\n      grid{d} = double(in);\r\n   <span class=\"keyword\">end<\/span>\r\n<\/pre>\r\n<img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/drums_1_02.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/>\r\n<h4>Finite difference Laplacian<a name=\"13d23acf-93ea-4116-b982-03d6b9ddc6c4\"><\/a><\/h4>\r\nDefining the 5-point finite difference Laplacian involves numbering the points in the region. The <tt>delsq<\/tt> function generates a sparse matrix representation of the operator and a <tt>spy<\/tt> plot of the nonzeros in the matrix shows its the band structure.\r\n<pre class=\"codeinput\">   <span class=\"keyword\">for<\/span> d = 1:2\r\n\r\n      <span class=\"comment\">% Number the interior grid points.<\/span>\r\n      G = grid{d};\r\n      p = find(G);\r\n      G(p) = (1:length(p))';\r\n\r\n      <span class=\"comment\">% Display the numbering.<\/span>\r\n      fprintf(<span class=\"string\">'grid%d ='<\/span>,d);\r\n      minispy(flipud(G))\r\n\r\n      <span class=\"comment\">% Discrete Laplacian.<\/span>\r\n      A = delsq(G);\r\n\r\n      <span class=\"comment\">% Spy plot<\/span>\r\n      subplot(2,2,d)\r\n      markersize = 6;\r\n      spy(A,markersize)\r\n      title(sprintf(<span class=\"string\">'delsq(grid%d)'<\/span>,d));\r\n   <span class=\"keyword\">end<\/span>\r\n<\/pre>\r\n<pre class=\"codeoutput\">grid1 = \r\n  .    .    .    .    .    .    .    .    .    .    .    .    .    .    .    .\r\n  .    .    .    .    .    .    .    .    .    .    .    .    .    .    .    .\r\n  .    .    .    .    .    .    .    .    .   56    .    .    .    .    .    .\r\n  .    .    .    .    .    .    .    .   48   55    .    .    .    .    .    .\r\n  .    .    .    .    .    .    .   41   47   54    .    .    .    .    .    .\r\n  .    .    .    .    .    .   35   40   46   53    .    .    .    .    .    .\r\n  .    .    .    .    .   30   34   39   45   52   60   63   65   66    .    .\r\n  .    .    .    .   26   29   33   38   44   51   59   62   64    .    .    .\r\n  .    .    .   18   25   28   32   37   43   50   58   61    .    .    .    .\r\n  .    .   11   17   24   27   31   36   42   49   57    .    .    .    .    .\r\n  .    5   10   16   23    .    .    .    .    .    .    .    .    .    .    .\r\n  .    4    9   15   22    .    .    .    .    .    .    .    .    .    .    .\r\n  .    3    8   14   21    .    .    .    .    .    .    .    .    .    .    .\r\n  .    2    7   13   20    .    .    .    .    .    .    .    .    .    .    .\r\n  .    1    6   12   19    .    .    .    .    .    .    .    .    .    .    .\r\n  .    .    .    .    .    .    .    .    .    .    .    .    .    .    .    .\r\n \r\ngrid2 = \r\n  .    .    .    .    .    .    .    .    .    .    .    .    .    .    .    .\r\n  .    .    .    .    .    .    .    .    .    .    .    .    .    .    .    .\r\n  .    .    .    .    .    .    .    .    .    .    .   57    .    .    .    .\r\n  .    .    .    .    .    .    .    .    .    .    .   56   62    .    .    .\r\n  .    .    .    .    .    .    .    .    .    .    .   55   61   65    .    .\r\n  .    .    .    .    .    .    .    .    .    .    .   54   60   64   66    .\r\n  .    5   11   18   26   30   34   38   42   46   50   53   59   63    .    .\r\n  .    4   10   17   25   29   33   37   41   45   49   52   58    .    .    .\r\n  .    3    9   16   24   28   32   36   40   44   48   51    .    .    .    .\r\n  .    2    8   15   23   27   31   35   39   43   47    .    .    .    .    .\r\n  .    1    7   14   22    .    .    .    .    .    .    .    .    .    .    .\r\n  .    .    6   13   21    .    .    .    .    .    .    .    .    .    .    .\r\n  .    .    .   12   20    .    .    .    .    .    .    .    .    .    .    .\r\n  .    .    .    .   19    .    .    .    .    .    .    .    .    .    .    .\r\n  .    .    .    .    .    .    .    .    .    .    .    .    .    .    .    .\r\n  .    .    .    .    .    .    .    .    .    .    .    .    .    .    .    .\r\n \r\n<\/pre>\r\n<img decoding=\"async\" src=\"https:\/\/blogs.mathworks.com\/images\/cleve\/drums_1_03.png\" alt=\"\" hspace=\"5\" vspace=\"5\" \/>\r\n<h4>Compare eigenvalues<a name=\"ea6b1e68-d22c-4b6a-b049-d460f734502e\"><\/a><\/h4>\r\nThe Arnoldi method implemented in the <tt>eigs<\/tt> function readily computes the eigenvalues and eigenvectors. Here are the first twenty eigenvalues.\r\n<pre class=\"codeinput\">   <span class=\"comment\">% How many eigenvalues?<\/span>\r\n\r\n   eignos = 20;\r\n\r\n   <span class=\"comment\">% A finer grid.<\/span>\r\n   ngrid = 32;\r\n   h = 1\/ngrid;\r\n   [x,y] = meshgrid(0:h:3);\r\n\r\n   inpoints = (7*ngrid-2)*(ngrid-1)\/2;\r\n   lambda = zeros(eignos,2);\r\n   V = zeros(inpoints,eignos,2);\r\n\r\n   <span class=\"keyword\">for<\/span> d = 1:2\r\n      vs = vertices{d};\r\n      [in,on] = inpolygon(x,y,vs(1,:),vs(2,:));\r\n      in = xor(in,on);\r\n\r\n      <span class=\"comment\">% Number the interior grid points.<\/span>\r\n      G = double(in);\r\n      p = find(G);\r\n      G(p) = (1:length(p))';\r\n      grid{d} = G;\r\n\r\n      <span class=\"comment\">% The discrete Laplacian<\/span>\r\n      A = delsq(G)\/h^2;\r\n\r\n      <span class=\"comment\">% Sparse matrix eigenvalues and vectors.<\/span>\r\n      [V(:,:,d),E] = eigs(A,eignos,0);\r\n      lambda(:,d) = diag(E);\r\n   <span class=\"keyword\">end<\/span>\r\n\r\n   format <span class=\"string\">long<\/span>\r\n   lambda = flipud(lambda)\r\n<\/pre>\r\n<pre class=\"codeoutput\">lambda =\r\n\r\n  10.165879621248989  10.165879621248980\r\n  14.630600866993412  14.630600866993376\r\n  20.717633982094981  20.717633982094902\r\n  26.115126153750705  26.115126153750676\r\n  28.983478457829740  28.983478457829762\r\n  36.774063407607322  36.774063407607287\r\n  42.283017757114585  42.283017757114699\r\n  46.034233949715308  46.034233949715386\r\n  49.213425509524733  49.213425509524775\r\n  52.126973962396342  52.126973962396370\r\n  57.063486161172904  57.063486161172769\r\n  63.350675017756465  63.350675017756259\r\n  67.491111510445251  67.491111510445194\r\n  70.371453210957910  70.371453210957810\r\n  75.709992784621988  75.709992784621818\r\n  83.153242199788906  83.153242199788821\r\n  84.673734481953971  84.673734481953886\r\n  88.554340162610174  88.554340162610174\r\n  94.230337192953158  94.230337192953257\r\n  97.356922250794767  97.356922250794653\r\n\r\n<\/pre>\r\n<h4>How about a proof?<a name=\"10036699-ceff-48c0-9db1-e39171e470b4\"><\/a><\/h4>\r\nVarying the number of eigenvalues, <tt>eignos<\/tt>, and the grid size, <tt>ngrid<\/tt>, in this script provides convincing evidence that the finite difference Laplacians on the two domains are isospectral. But this is not a proof. For the continuous problem, Chapman <a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">[8]<\/a> outlines an approach where any eigenfunction on one of the domains can be constructed from triangular pieces of the corresponding eigenfunction on the other domain. It is necessary to prove that these pieces fit together smoothly and that the differential equation continues to be satisfied across the boundaries. For this proof Chapman refers to a paper by Berard <a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">[9]<\/a>. I will explore the discrete analog of these arguments in a later post.\r\n<h4>References<a name=\"efc546e6-cc7b-4279-bedc-683e8056a689\"><\/a><\/h4>\r\n<div>\r\n<ol>\r\n \t<li>Marc Kac, Can one hear the shape of a drum?, Amer. Math. Monthly 73 (1966), 1-23.<\/li>\r\n \t<li><a href=\"https:\/\/www.mathworks.com\/company\/newsletters\/articles\/the-mathworks-logo-is-an-eigenfunction-of-the-wave-equation.html\">Cleve Moler, The MathWorks logo is an eigenfunction of the wave equation (2003).<\/a><\/li>\r\n \t<li><a href=\"http:\/\/people.maths.ox.ac.uk\/trefethen\/publication\/PDF\/2006_116.pdf\">Lloyd N. Trefethen and Timo Betcke, Computed eigenmodes of planar regions (2005).<\/a><\/li>\r\n \t<li><a href=\"http:\/\/www.ams.org\/journals\/bull\/1992-27-01\/S0273-0979-1992-00289-6\">Carolyn Gordon, David Webb, Scott Wolpert, One cannot hear the shape of a drum, Bull. Amer. Math. Soc. 27 (1992), 134-138.<\/a><\/li>\r\n \t<li><a href=\"http:\/\/en.wikipedia.org\/wiki\/Hearing_the_shape_of_a_drum\">Wikipedia, Hearing the shape of a drum.<\/a><\/li>\r\n \t<li><a href=\"http:\/\/www.math.udel.edu\/~driscoll\/research\/drums.html\">Tobin Driscoll, Isospectral Drums.<\/a><\/li>\r\n \t<li>Tobin Driscoll, Eigenmodes of isospectral drums, SIAM Review 39 (1997), 1-17.<\/li>\r\n \t<li><a href=\"http:\/\/www.math.ucdavis.edu\/~saito\/courses\/ACHA.READ.F03\/drum-chapman.pdf\">S. J. Chapman, Drums that sound the same, Amer. Math. Monthly 102 (1995), 124-138.<\/a><\/li>\r\n \t<li><a href=\"http:\/\/www.springer.com\/mathematics\/journal\/208\">Pierre Berard, Transplantation et isospectralite, Math. Ann. 292 (1992), 547-559.<\/a><\/li>\r\n<\/ol>\r\n<\/div>\r\nIf you are interested in pursuing this topic, see the PDE chapter of <a href=\"https:\/\/www.mathworks.com\/content\/dam\/mathworks\/mathworks-dot-com\/moler\/pdes.pdf\"><i>Numerical Computing with MATLAB<\/i><\/a>, and check out <a href=\"https:\/\/www.mathworks.com\/content\/dam\/mathworks\/mathworks-dot-com\/moler\/ncm\/pdegui.m\">pdegui.m<\/a>.\r\n\r\n<script>\/\/ <![CDATA[\r\nfunction grabCode_b3b44631f1934e41818144b719024275() {\r\n        \/\/ Remember the title so we can use it in the new page\r\n        title = document.title;\r\n\r\n        \/\/ Break up these strings so that their presence\r\n        \/\/ in the Javascript doesn't mess up the search for\r\n        \/\/ the MATLAB code.\r\n        t1='b3b44631f1934e41818144b719024275 ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' b3b44631f1934e41818144b719024275';\r\n    \r\n        b=document.getElementsByTagName('body')[0];\r\n        i1=b.innerHTML.indexOf(t1)+t1.length;\r\n        i2=b.innerHTML.indexOf(t2);\r\n \r\n        code_string = b.innerHTML.substring(i1, i2);\r\n        code_string = code_string.replace(\/REPLACE_WITH_DASH_DASH\/g,'--');\r\n\r\n        \/\/ Use \/x3C\/g instead of the less-than character to avoid errors \r\n        \/\/ in the XML parser.\r\n        \/\/ Use '\\x26#60;' instead of '<' so that the XML parser\r\n        \/\/ doesn't go ahead and substitute the less-than character. \r\n        code_string = code_string.replace(\/\\x3C\/g, '\\x26#60;');\r\n\r\n        copyright = 'Copyright 2012 The MathWorks, Inc.';\r\n\r\n        w = window.open();\r\n        d = w.document;\r\n        d.write('\r\n\r\n<pre>\\n');\r\n        d.write(code_string);\r\n\r\n        \/\/ Add copyright line at the bottom if specified.\r\n        if (copyright.length > 0) {\r\n            d.writeln('');\r\n            d.writeln('%%');\r\n            if (copyright.length > 0) {\r\n                d.writeln('% _' + copyright + '_');\r\n            }\r\n        }\r\n\r\n        d.write('<\/pre>\r\n\r\n\r\n\\n');\r\n\r\n        d.title = title + ' (MATLAB code)';\r\n        d.close();\r\n    }\r\n\/\/ ]]><\/script>\r\n<p style=\"text-align: right; font-size: xx-small; font-weight: lighter; font-style: italic; color: gray;\">\r\n<a><span style=\"font-size: x-small; font-style: italic;\">Get\r\nthe MATLAB code<noscript>(requires JavaScript)<\/noscript><\/span><\/a>\r\n\r\nPublished with MATLAB\u00ae 7.14<\/p>\r\n<p class=\"footer\">\r\nPublished with MATLAB\u00ae 7.14<\/p>\r\n\r\n<\/div>\r\n<!--\r\nb3b44631f1934e41818144b719024275 ##### SOURCE BEGIN #####\r\n%% Can One Hear the Shape of a Drum? Part 1, Eigenvalues\r\n% The title of this multi-part posting is also the title of a\r\n% 1966 article by Marc Kac in the American Mathematical Monthly <#7 [1]>.\r\n% This first part is about isospectrality.\r\n\r\n%% Isospectrality\r\n% Kac's article is not actually about a drum, which is three-dimensional,\r\n% but rather about the two-dimensional drum head, which is more like a\r\n% tambourine or membrane.  The vibrations are modeled by the\r\n% partial differential equation\r\n%\r\n% $$ \\Delta u + \\lambda u = 0 $$\r\n%\r\n% where\r\n%\r\n% $$ \\Delta u(x,y) = \\frac{\\partial^{2} u}{\\partial x^{2}} + \\frac{\\partial^{2} u}{\\partial y^{2}} $$\r\n%\r\n% The boundary conditions are the key.  Requiring $u(x,y) = 0$ on the boundary\r\n% of a region in the plane corresponds to holding the membrane fixed on that\r\n% boundary.  The values of $\\lambda$ that allow nonzero solutions, the\r\n% _eigenvalues_, are the squares of the frequencies of vibration, and the\r\n% corresponding functions $u(x,y)$, the _eigenfunctions_, are the modes of\r\n% vibration.\r\n%\r\n% The MathWorks logo comes from this partial differential equation,\r\n% on an L-shaped domain <#7 [2]>, <#7 [3]>, but this article is not\r\n% about our logo.\r\n%\r\n% A region determines its eigenvalues.  Kac asked about the opposite\r\n% implication.  If one specifies all of the eigenvalues, does that determine\r\n% the region?\r\n%\r\n% In 1991, Gordon, Webb and Wolpert showed that the answer to Kac's question\r\n% is \"no\".  They demonstrated a pair of regions that had different shapes\r\n% but exactly the same infinite set of eigenvalues <#7 [4]>.\r\n% In fact, they produced several different pairs of such regions.\r\n% The regions are known as \"isospectral drums\".\r\n% Wikipedia has a good background article on Kac's problem <#7 [5]>.\r\n%\r\n% I am interested in finite difference methods for membrane eigenvalue problems.\r\n% I want to show that the finite difference operators on these regions\r\n% have the same sets of eigenvalues, so they are also isospectral.\r\n%\r\n% I was introduced to isospectral drums by Toby Driscoll,\r\n% a professor at the University of Delaware.\r\n% A summary of Toby's work is available at his Web site <#7 [6]>.\r\n% A reprint of his 1997 paper in SIAM Review is also available from his\r\n% Web site <#7 [7]>.  Toby developed methods, not involving finite differences,\r\n% for computing the eigenvalues very accurately.\r\n%\r\n% The isospectral drums are not convex.  They have reentrant $270^\\circ$\r\n% corners.  These corners lead to singularities in most of the\r\n% eigenfunctions REPLACE_WITH_DASH_DASH the gradients are unbounded.\r\n% This affects the accuracy and the rate of convergence of finite\r\n% difference methods.\r\n% It is possible that for convex regions the answer to Kac's\r\n% question is \"yes\".\r\n\r\n%% Vertices\r\n% I will look at the simplest known isospectral pair.\r\n% The two regions are specified by the _xy_-coordinates of their vertices.\r\n\r\ndrum1 = [0 0 2 2 3 2 1 1 0\r\n0 1 3 2 2 1 1 0 0];\r\ndrum2 = [1 0 0 2 2 3 2 1 1\r\n0 1 2 2 3 2 1 1 0];\r\nvertices = {drum1,drum2};\r\n\r\n%%\r\n% Let's first plot the regions.\r\n\r\nclf\r\nshg\r\nset(gcf,'color','white')\r\nfor d = 1:2\r\n% Plot the region.\r\nvs = vertices{d};\r\nsubplot(2,2,d)\r\nplot(vs(1,:),vs(2,:),'k.-');\r\naxis([-0.1 3.1 -0.1 3.1])\r\naxis square\r\ntitle(sprintf('drum%d',d));\r\nend\r\n\r\n%% Finite difference grid\r\n% I want to investigate simple finite difference methods for this problem.\r\n% The MATLAB function |inpolygon| determines the points of a rectangular\r\n% grid that are in a specified region.\r\n\r\n% Generate a coarse finite difference grid.\r\nngrid = 5;\r\nh = 1\/ngrid;\r\n[x,y] = meshgrid(0:h:3);\r\n\r\n% Loop over the two regions.\r\nfor d = 1:2\r\n\r\n% Determine points inside and on the boundary.\r\n\r\nvs = vertices{d};\r\n[in,on] = inpolygon(x,y,vs(1,:),vs(2,:));\r\nin = xor(in,on);\r\n\r\n% Plot the region and the grid.\r\nsubplot(2,2,d)\r\nplot(vs(1,:),vs(2,:),'k-',x(in),y(in),'b.',x(on),y(on),'k.')\r\naxis([-0.1 3.1 -0.1 3.1])\r\naxis square\r\ntitle(sprintf('drum%d',d));\r\n\r\ngrid{d} = double(in);\r\nend\r\n\r\n%% Finite difference Laplacian\r\n% Defining the 5-point finite difference Laplacian involves numbering the\r\n% points in the region.  The |delsq| function generates a sparse matrix\r\n% representation of the operator and a |spy| plot of the nonzeros in the\r\n% matrix shows its the band structure.\r\n\r\nfor d = 1:2\r\n\r\n% Number the interior grid points.\r\nG = grid{d};\r\np = find(G);\r\nG(p) = (1:length(p))';\r\n\r\n% Display the numbering.\r\nfprintf('grid%d =',d);\r\nminispy(flipud(G))\r\n\r\n% Discrete Laplacian.\r\nA = delsq(G);\r\n\r\n% Spy plot\r\nsubplot(2,2,d)\r\nmarkersize = 6;\r\nspy(A,markersize)\r\ntitle(sprintf('delsq(grid%d)',d));\r\nend\r\n\r\n%% Compare eigenvalues\r\n% The Arnoldi method implemented in the |eigs| function readily computes\r\n% the eigenvalues and eigenvectors.  Here are the first twenty eigenvalues.\r\n\r\n% How many eigenvalues?\r\n\r\neignos = 20;\r\n\r\n% A finer grid.\r\nngrid = 32;\r\nh = 1\/ngrid;\r\n[x,y] = meshgrid(0:h:3);\r\n\r\ninpoints = (7*ngrid-2)*(ngrid-1)\/2;\r\nlambda = zeros(eignos,2);\r\nV = zeros(inpoints,eignos,2);\r\n\r\nfor d = 1:2\r\nvs = vertices{d};\r\n[in,on] = inpolygon(x,y,vs(1,:),vs(2,:));\r\nin = xor(in,on);\r\n\r\n% Number the interior grid points.\r\nG = double(in);\r\np = find(G);\r\nG(p) = (1:length(p))';\r\ngrid{d} = G;\r\n\r\n% The discrete Laplacian\r\nA = delsq(G)\/h^2;\r\n\r\n% Sparse matrix eigenvalues and vectors.\r\n[V(:,:,d),E] = eigs(A,eignos,0);\r\nlambda(:,d) = diag(E);\r\nend\r\n\r\nformat long\r\nlambda = flipud(lambda)\r\n\r\n%% How about a proof?\r\n% Varying the number of eigenvalues, |eignos|, and the grid size, |ngrid|,\r\n% in this script provides convincing evidence that the finite difference\r\n% Laplacians on the two domains are isospectral.  But this is not a proof.\r\n% For the continuous problem, Chapman <#7 [8]> outlines an approach where\r\n% any eigenfunction on one of the domains can be constructed from triangular\r\n% pieces of the corresponding eigenfunction on the other domain.\r\n% It is necessary to prove that these pieces fit together smoothly and\r\n% that the differential equation continues to be satisfied across the\r\n% boundaries.  For this proof Chapman refers to a paper by Berard <#7 [9]>.\r\n% I will explore the discrete analog of these arguments in a later post.\r\n\r\n%% References\r\n% # Marc Kac, Can one hear the shape of a drum?, Amer. Math. Monthly 73 (1966), 1-23.\r\n% # <https:\/\/www.mathworks.com\/company\/newsletters\/articles\/the-mathworks-logo-is-an-eigenfunction-of-the-wave-equation.html % Cleve Moler, The MathWorks logo is an eigenfunction of the wave equation % (2003).>\r\n% # <http:\/\/people.maths.ox.ac.uk\/trefethen\/publication\/PDF\/2006_116.pdf % Lloyd N. Trefethen and Timo Betcke, Computed eigenmodes of planar regions % (2005).>\r\n% # <http:\/\/www.ams.org\/journals\/bull\/1992-27-01\/S0273-0979-1992-00289-6 % Carolyn Gordon, David Webb, Scott Wolpert, % One cannot hear the shape of a drum, Bull. Amer. Math. Soc. 27 (1992), % 134-138.>\r\n% # <http:\/\/en.wikipedia.org\/wiki\/Hearing_the_shape_of_a_drum % Wikipedia, Hearing the shape of a drum.>\r\n% # <http:\/\/www.math.udel.edu\/~driscoll\/research\/drums.html % Tobin Driscoll, Isospectral Drums.>\r\n% #\r\n% Tobin Driscoll, Eigenmodes of isospectral drums, SIAM Review 39 (1997),\r\n% 1-17.\r\n% # <http:\/\/www.math.ucdavis.edu\/~saito\/courses\/ACHA.READ.F03\/drum-chapman.pdf % S. J. Chapman, Drums that sound the same, Amer. Math. Monthly 102 (1995), % 124-138.>\r\n% # <http:\/\/www.springer.com\/mathematics\/journal\/208 % Pierre Berard, Transplantation et isospectralite, Math. Ann. 292 (1992), % 547-559.>\r\n\r\n%%\r\n% If you are interested in pursuing this topic, see the PDE chapter of\r\n% <https:\/\/www.mathworks.com\/moler.htmlpdes.pdf _Numerical Computing with MATLAB_>,\r\n% and check out <https:\/\/www.mathworks.com\/moler.htmlncm\/pdegui.m pdegui.m>.\r\n\r\n##### SOURCE END ##### b3b44631f1934e41818144b719024275\r\n-->","protected":false},"excerpt":{"rendered":"<!--introduction-->The title of this multi-part posting is also the title of a 1966 article by Marc Kac in the American Mathematical Monthly <a href=\"#10036699-ceff-48c0-9db1-e39171e470b4\">[1]<\/a>. This first part is about isospectrality.\r\n\r\n<!--\/introduction-->... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/cleve\/2012\/08\/06\/can-one-hear-the-shape-of-a-drum-part-1-eigenvalues\/\">read more >><\/a><\/p>","protected":false},"author":78,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[13],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/235"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/users\/78"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/comments?post=235"}],"version-history":[{"count":18,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/235\/revisions"}],"predecessor-version":[{"id":2173,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/posts\/235\/revisions\/2173"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/media?parent=235"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/categories?post=235"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/cleve\/wp-json\/wp\/v2\/tags?post=235"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}