Simulating waves and macroscopic phonons

Wave phenomena are fundamental for many branches of physics, not only for mechanics. Thus, students at university level should become familiar with the underlying theory, and especially with solutions of the wave equation. This paper refers to a multimedia program for PCs, tablets or smartphones, and introduces and discusses several animated illustrations. The introduced program can illustrate the power of the wave equation. Basic wave phenomena can be computed with numerical methods. The finite difference method which is used can be lead back to fundamental assumptions of physics. Yet, this numerical method is capable of calculating a great variety of wave phenomena. Several examples can demonstrate the power of this numerical procedure. Phenomena like the propagation of circular waves diverging from a point source, Huygens’s elementary waves, superposition of waves, reflexion, refraction and transmission at the interface between different media, scattering, interference phenomena, as well as dispersion and the Doppler effect can be studied in animated visualisations. This opens up a variety of teaching and learning strategies, e.g. the single concept principle and the inverted classroom concept, and the modelling of partial differential equations (the wave equation as well as the heat equation) can be supported. The conceptual approach is described here; an empirical study about learning will be a next step.

(Some figures may appear in colour only in the online journal) Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

Introduction
Wave phenomena appear in many contexts throughout physics. Thus, our attention should be concentrated not only on the knowledge of specific phenomena but also on causalities, fundamental principles, and methods for calculations based on theoretical approaches.
There are already many websites dealing with waves. A useful interactive tutorial including simulations and additional material is provided by Forinash and Wolfgang Christian [1]. Open source physics offers a significant number of simulations [2]. Though already older, the report by Bruce Mason [3] still provides an extensive selection of websites about waves with additional explanations. Meanwhile, most resources are updated and therefore still a good choice. Falstad provides some educational applets to visualize various concepts in math, physics, and engineering, mostly already converted to JavaScript, so that they can run in common browsers [4]. Vibrations and waves are also a topic for 'the physics classroom', providing additional textual information [5]. Also, Phet interactive simulations [6] offer nice opportunities to get experiences with basic phenomena based on ripple tank experiments, sound and also light. While in these simulations the underlying theoretical model cannot be varied by a user, Net-Logo provides a programmable modeling environment with an example that simulates wave motion in a membrane with the four edges of the membrane fixed to a frame [7].
The intention of this introduction is not to give an overview of existing web resources or to reproduce similar movies but to introduce a special approach and identify specialities of the presented program. Above all, an important intention is to provide an anchor to connect phenomena to numerical calculations based on the wave equation. To keep it familiar, the provided illustrations have an appearance similar to water waves. It is a tool for teachers to illustrate important aspects in lectures, and for learners it provides opportunities to examine the behaviour of waves under different conditions.
To pave the way to the applied numerical procedure the following fundamental assumption is feasible. Wave oscillations and their propagating in space can be studied based on the idea that the time dependent spread of local changes of a physical quantity is proportional to spatial differences of this quantity (beside further coupling parameter). To make this principle more vivid the applications to mechanics are helpful. Particularly for transversal waves this means that the vertical acceleration of a particle is proportional to the differences in vertical positions between neighbouring particles.
To keep it more tangible, the speed of mechanical waves is determined by the inertial and elastic properties of the medium. This assumption can be easily justified by Hooke's law and Newton's second law (see below). The variables of the wave equation, a second-order linear partial differential equation (see equation (1)), have very different meanings for acoustic waves or electromagnetic waves. Nevertheless, the pure numerical procedures to solve the secondorder linear partial differential equation are very similar. Thus, it makes sense to have a closer look at suitable numerical procedures.
There are several, also quite sophisticated methods to calculate the behaviour of waves. The 'finite volume method' is a suitable way to calculate the behavior of waves reaching shallow water [8]. A 'finite element method' is useful for complex situations, and Muliyati et al [9] used a mesh with quadratic triangular cells for calculating the behavior of ocean waves in coastal areas.
Finite element and finite volume methods have been developed to great generality and are more precise than finite differences. Nevertheless, simplicity, both from a mathematical and a programming perspective, makes the latter method quite interesting. At least for a basic numerical analysis finite-difference methods (FDM) are still suitable.
These methods use discretization for solving differential equations by approximating them with difference equations. Hence, basic ideas of the FDM can be easily associated to fundamental concepts of mechanics, making the derivation particularly relevant also for an introductory university course. In addition, this approach comprises the same ideas being used for calculating phonons in solid state physics. Especially the so called 'numerical dispersion' can be linked to dispersion relations of phonons.
Basic wave phenomena can be computed with numerical methods. The FDM which is used, can be lead back to fundamental assumptions of physics. Yet, this numerical method is capable of calculating a great variety of wave phenomena.
Based on the finite-difference method the mightiness of the basic wave equation (see equation (1)) for computing general wave phenomena can be visualised and studied for individually adjustable configurations. Special topics are: However, more sophisticated and complex behaviour, like those of Rayleigh or Lamb waves depending on the thickness of the material, or the coupling of magnetic and electric field components in electromagnetic waves go beyond what can be described by the classical two-dimensional wave equation.

Methods
where w represents the wave function, which is for mechanical waves the displacement at a position specified by the coordinates x and y; u is the speed of the wave. For fundamental considerations and basic examples, a square grid with mesh size a can be used as discretization of a two-dimensional plane. In the simplest case, only the next four neighbours of a grid point (x 0 , y 0 ) (see figure 1) are taken into account to solve the wave equation numerically.
Then, the second partial derivation of a function w(x, y, t) with respect to x at the location (x 0 , y 0 ) is to be replaced by An analogous expression applies to the partial derivatives along y. The time derivatives are also replaced by a difference operator, whereby the values must come from different 'time planes'. (See also [10] or [11] for further details.) A more precise treatment can also take into account the second neighbours with appropriate weighting. For example, the computer program presented works with a nine-point operator. This reduces a numerical anisotropy (see e.g. [12]).
The finite-difference method particularly offers an easy-to-understand numerical solution, which can also be derived from elementary physical laws. No continuum is assumed for the calculations but rather a discrete mass distribution. Then, Newton's second law can be applied to each mass element. The accelerating forces result from the different displacements of adjacent masses. For that, a linear law of force is assumed, which should already be a well-known procedure from Hooke's law. It is crucial to assume a linear law of forces depending on the local displacements. The basic assumptions are best illustrated with the help of a linear chain. From one grid point to the other, forces act which are proportional to the deflection differences. For transverse waves, it is assumed that the deflection differences from grid point to grid point remain small relative to the grid distance. Then, the restoring forces in the vertical direction are approximately proportional to the deflection differences. The following sketch (figure 2) serves as an orientation aid.
The vertical components of the coupling forces F l and F r are calculated by the following formulas (for small deflections and absolute values only): with a constant tension in the rope F l ≈ F r , and thus also the constant f is equivalent in both equations. The horizontal components are always directed in opposite direction to each other. They compensate each other in the first approximation.
Thus, the acceleration of an element n can be calculated in the following way: These assumptions, which are elementary from a physical point of view, directly prompt a use of finite differences.
It is thus possible to use the deflection values from the current and previous timestep to find the deflection in the next time step: or, rearranged for the new deflection: with f · Δt 2 = f * . Thus, the displacements w n new for the time t + Δt can be calculated numerically from the current values and the values before them. The numerical formulas open up a simple model that can be applied to all kinds of waves, giving them the rank of fundamental ideas.
However, the considerations above are not only suitable for making the numerical approach plausible. They can also help to uncover peculiarities of this method. It is of particular note that a grid structure is utilized. Thus, nearly the same procedures are used as for elementary solid state physics. In particular, this results in waves with dispersion, just as in calculations of phonons, and dispersion can even become visible (see 3.1.9). Furthermore, there is also a smallest wavelength that can be realized on this system. Nevertheless, with a grid of 100 × 100 points transversal waves (gravity waves) can be simulated quite well. In spite of new technologies, the computing power for such calculations is still limited. Thus, especially for smartphones also a smaller area of 80 × 80 or 60 × 60 points can be selected.

Calculations at boundaries.
Of special interest is the behavior of waves at boundaries. Several types of interfaces can be realized: • Fixed end points can be implemented easily. After each calculation step, the values of the boundary points are reset to a fixed value (e.g. zero). • For simulating a free end, virtual points beyond the boundary are adjusted to the same value as the points on the boundary. Thus there is no backaction from outside the lattice. • Far more sophisticated is simulating an infinite continuation. No reflections are permitted, and all the wave energy has to dissipate through the boundary. To realize this, the field must be terminated by the wave impedance/specific impedance of the medium. (For the numerical calculation it is helpful to realize that the wave impedance or specific impedance indicates how large the force exerted on a mass point becomes in relation to its vertical velocity. This can be derived for a simple string by dividing the tension of the string by the velocity of the wave along the string-not transverse to it, respectively by the square root of the tension divided by the mass density of the rope.) Thus, a virtual point beyond the boundary has to exert a force which is proportional to the vertical speed of the end point and the wave impedance.  figure 3).

Closer inspection of complex scenarios.
Animations can highlight different features and provide valuable insights especially into time-dependent processes. This is helpful for complex processes, e.g. for the partial transmission and reflection of waves at the boundary between different media (see figure 11). The changing directions of the wave front as well as different speeds become visible. Another example is dispersion, and an illustration that the speed of waves varies depending on their wavelengths (see figure 13).
In an experimental setup waves cannot be stopped for analysing special phases. Features like 'stop', 'backward', 'forward', 'save' belong to the benefits of the introduced program. The possibility to stop the animation or even to step backwards can enhance detailed observations. Particular phases can be studied in appropriate speed and as often as needed. Also the possibility to save pictures is helpful for comparing different states and the results of various starting conditions. Basic terms and phenomena become clear in scenarios which are not abstract but easily observed in visualizations. This is a first step, and also a guideline for further analogical thinking. The next step would be to proceed to more abstract waves, e.g. electromagnetic waves (see also [13]). Certainly, simulations shall not substitute real experiments, but they can extend possibilities of experimental studies (anytime and anywhere, even with smartphones). A setup for complex configurations, e.g. to study interference patterns, can be realized quickly and easily, allowing to accomplish quasi-experimental activities with less effort than for a real experimental setup. Special details emerge clearly on the screen, and can be reproduced multiple times. Results of calculations can be traced step by step. Specific questions can be pursued picture by picture, and repeated animations uncover the time-dependent behaviour.
In this way, dealing with inaccuracy of measurements or calculations is also practised and can appear similar to experiences in real experiments.

Modelling the wave equation.
Linking to theoretical considerations with respect to the wave equation helps to develop a deeper understanding. For example, the superposition principle is to be connected to the linearity and isotropy of the wave equation, and can be reconstructed quite well by the numerical nine-point operator (see equation (10)).
The described program is written in html5, and is an open source application. Especially for modelling and the option for changing the underlying mathematical model, the theoretically relevant function has been placed at the top of the source code and particularly marked. Thus, parameters can be changed easily in order to study their influence. For programming a changed notation compared to the usually applied mathematical notation was necessary, and the variables were renamed according to figure 4.   w n,m new = 2w n,m − w n,m old + f * m · (w m+1,n+1 + 4 · w m+1,n + w m+1,n−1 + 4 · w m,n+1 − 20 · w m,n + 4 · w m,n−1 + w m−1,n+1 + 4 · w m−1,n + w m−1,n−1 .

(10)
An example with only a few, but nevertheless quite interesting changes is shown in figure 6. Thus, even the heat equation can be studied.
Considerations concerning interactions between neighbouring points/particles or boundaries initiate interesting discussions, e.g. about their influence on reflection and transmission at a boundary between two different media. In that case of course, the simulations cannot stand alone and should be integrated into guiding lecture notes and exercises.
However, we found that also the opposite approach can be quite stimulating and motivating for theoretical considerations: the fact that basic wave phenomena come out so clearly just by numerically solving the wave-equation was quite impressive for students as well as for teachers. (Yet we do not have empirical data about that.) Inspecting the progression of waves in combination with a derivation of the numerical formula based on fundamental laws of physics (see above) can motivate and lead to further insights.

Results and applications
The user interface of the program provides options • to select various sources for excitation, like drops, up to ten point sources, moving sources, sinusoidal oscillating sources generating linear wavefronts; • to define different structures for the membrane like barriers, slits, scattering centres; • to change the stiffness of the medium, also generating a transition between different media; • to select different boundaries; as well as • to change the resolution of the mesh and the perspective.
The following examples show the results of different settings. All the numerical calculations are executed step by step-no analytical solution for the wave equation nor a simple mathematical wave function is used. Via the program menu only the initial values and the boundary conditions are set. All the calculations are based on a discretization of the wave equation according to the finite difference method described in section 2.1. (Nevertheless, there is also an option to change the underlying model by changing four lines in the program, like it is described in section 2.2.5, and thus for example also the propagation of thermal waves for different boundary conditions can be studied.) Beside the program, a clickable mind map and a selection of animated gifs can be downloaded from: https://www.en.didaktik.physik.uni-muenchen.de/multimedia/waves.

Simple circular waves
A simple starting configuration can be chosen producing a circular wave, similar to a water wave caused by a drop falling into water. Several sources can be started by positioning excitation pulses via mouse click (see figure 7).

Huygens' principle and elementary waves
According to Huygens's principle every point on a wave front is the source of a spherical wave. Thus, a small slit (significantly smaller than the wavelength) on a two-dimensional plane can be used to create circular waves (see figure 8).

The superposition of waves
A fascinating characteristic shared by all waves is that an overlapping of waves in a region produces quite complex disturbances. Afterwards, travelling waves move on as if nothing had happened, and the waves maintain their 'integrity'. The local deflection of two or more overlapping waves is the algebraic sum of their contributions at that point. Figure 9 shows the superposition of two waves and a picture of their further progress in time.

Interference and double-slit experiment
When waves emitted by phase-coupled sources overlap, interference structures can be observed. (Details depend on the sources and the structure of the medium.) Twenty coherent point-like sources can be placed at different positions. Furthermore, different and also periodic boundaries can be set up. Hence, a huge variety of interference pattern can be studied, like interference effects at barriers, scattering objects. A few examples are shown in figures 10 and 11.

Scattering
The scattering of waves hitting an obstacle can be observed. Figure 12 shows an example where a circular wave (which could be created by a drop of water) encounters an obstacle.

Resonances of a vibrating membrane and Chladni figures
Standing waves (eigenresonances of a two-dimensional surface) can be seen, being the solution of the two-dimensional wave-equation with the given boundary conditions. In figure 13 there is an oscillator in the center of a plate with open/free boundaries. Chladni figures can be realized by selecting a fixed center and periodic swinging edges.

Waves passing through different media (reflection, transmission and refraction)
When waves pass through the boundary between two media, where waves have different phase velocities, then transmission, reflection and refraction can be observed. Snell's law and Fresnel's formula provide deeper insights into the quantitative relations. Additionally, the introduced program can offer visual insights into the varying behaviour. Two examples are shown. In figure 14(a), circular wave passes the boundary between two different media. The flow of energy will change. Figure 14 shows that the wave pulse is partly reflected (due to the different tension and 'inertia per area' of the two media).
In figure 15 two periodical wave fronts propagate through two media. In this example the periodic waves naturally have the same frequency in both media, however a different speed and thus a different wavelength (see figure 15). The stronger the linkage between the masses and the smaller the density/inertia of the medium, the greater the wavelength.

Wave packets and dispersion
In section 2.1.1 the numerical dispersion was already mentioned, showing similarities to dispersion relations in solid state physics. Related to this, the pictures of figure 16, taken from an animation, need some explanations. The starting configuration is a uniform displacement in a square region (first picture). According to Fourier this can be regarded as a superposition of a broad spectrum of waves with different wavelengths.
The second picture shows the situation after 94 calculation steps. It becomes clear that the propagation speed of waves depends on their wavelengths. At the end-just like in solid state physics-those with a wavelength of twice the lattice constant remain at the initial position. (Their propagation speed is zero.)

Doppler effect
A point source moving through a homogeneous medium creates waves with different wavelengths on either side of the source. The distance between two phase-fronts differ depending   on the direction of the moving source. In figure 17 shorter wavelengths can be seen in the direction of the moving source and longer wavelengths behind it. Figure 18 shows situations where the source is moving with a speed higher than the propagations speed of waves. The left picture refers to a periodic oscillator creating a shock-wave cone (Mach's cone).
Yet, a ship's keel is not a periodical oscillator but creates a bow wave starting with a constant amplitude. This is simulated in figure 18 (right picture). However, this is a special, simplified model. The situation is much more complex in reality, as there are extensive wave sources, and water surface waves show other kinds of dispersion.

Energy conservation and reversibility
A feature with interesting consequences is the reverse button implemented in the simulation program. The possibility to change the arrow of time goes in line with the conservation of energy. It is possible to reverse the numerical calculation-except for all the configurations where no energy conservation persists. This concerns boundaries with infinite continuation or the different kinds of periodic excitations.

Options for working
(a) Already prepared animated gifs can be downloaded, simply to study phenomena. A clickable mind map (see figure 3) might be useful for offering structure or guidance. (b) More detailed studies are possible with a special recorder, providing the option to speed up, slow down, stop or even step back. (c) Work with the full html5-program (which also can be used offline), to create your own starting configurations, or even modulate the underlying numerical procedure (as described in section 2.2.5).

Discussion
On the one hand, similarities to water waves can make visual displays more familiar. On the other hand, it must be clear that there are also fundamental differences. Surface water waves can be attributed to surface tension, gravitational and inertial forces. Wind or seismic excitations create different forms. Furthermore, also the familiar gravity-dominated surface waves are a mixture of longitudinal and transverse waves. (This is described in many textbooks, and an illustration can be seen in [14].) Hence, beside a numerical inaccuracy also the classical wave equation is only an approximation. The same is true for an infinite continuation or a complete absorption respectively for the corresponding calculations at boundaries. Only an approximation for the wave impedance or specific impedance can be realized. Actually in this program a combined calculation with a specific impedance and a kind of further calculations in a virtual area beyond the borderline is used. Nevertheless, corners or angles in the borderline cause problems. E.g. calculations for a total absorber with slits (single or double slit) lead to visible errors and unrealistic disturbances after several thousand calculation steps.
The single concept principle and mapping with linked pictures might be helpful for inverted classroom strategies. However, additional guidance and further explanations are required. Also, studying concepts in complete isolation is not realistic. On closer inspection more aspects will always become relevant.
Modelling and changing the program is limited. Only changing the basic function for the wave equation makes sense, and this might also offer a motivation to think about numerical calculations. Other changes, e.g. the calculations at boundaries, are sophisticated or require deeper knowledge of JavaScript and html-programming.
The main focus is on the physics. Hence, only a few mouse clicks are necessary to change initial situations and frame conditions. It makes it quick and easy to examine and investigate the behaviour of the simulated waves. This can reduce the need for some sophisticated experimental activities but should not replace basic ones. A comparison with real experiments is always necessary to match theory and reality. It should always be clear that only a simplified numerical calculation is provided. Pointing to the fact that a differential equation is substituted by a difference equation should make this clear.

Conclusions
The animations discussed above can be tools to vividly illustrate various wave phenomena. The described derivation of a numerical calculation which is connected to basic laws from physics in a straightforward way can help to gain deeper insights. Nevertheless, as already mentioned, critical thinking is always required.