David Todd and Dr. Blake Barker, Department of Mathematics
The goal of this project was to create code for a generic numerical solver utilizing a finite difference scheme. Prior to the creation of this code, those wishing to study the evolution of traveling waves in STABLAB, a MATLAB based program specifically designed for studying the stability of traveling waves, have had to seek out external numerical solvers for the individual systems that they were studying, or develop them themselves. By developing this code, we have helped to save these individuals valuable time and ensure that they are utilizing an efficient numerical method for studying the evolution of their system. In addition, we have improved the convenience and value of STABLAB, by providing another easily accessible resource for users of this mathematics software package.
The first step of the project was to determine which finite difference scheme to utilize in the code. We tested the explicit, implicit, and Crank-Nicolson schemes against several known systems, including Burgers Equation and the Isentropic Gas Equations.
Each finite difference method was sampled individually to create a numerical solver in Python that can take various initial boundary conditions and run a time evolution for each example system. To start, we initialized a matrix for each function that satisfied the boundary conditions whose indices represented each time and space step. Then we created the Jacobian Matrix for these functions given the finite difference equations we established above. Using Newton’s Method, we solved iteratively for the behavior of the functions at each time step. The end result was a grid representing numerical approximations of the value of the function at each point.
Once we had tested several different sample systems with various boundary conditions, we determined that the Crank-Nicolson finite difference scheme was more accurate and efficient than either the implicit or explicit methods. From this point, we sought to generalize our code so that the process we employed could be repeated for any given system. Once we developed generalized code in Python that prints out the needed MATLAB code to run time evolution studies on the viscous conservation systems we were interested in studying, we uploaded the code to the internet for public access.
The first several months of research were spent exploring the theory of finite difference and the options that could be used in the final code that we would produce. After much testing and computation, we determined that the Crank- Nicolson finite difference scheme, in which the average of the implicit and explicit numerical finite difference approximations for partial derivatives is used, is the most computationally efficient and accurate method for our purposes.
The most tangible results of the research project was the development of new code accessible to anyone who wishes to use it on GITHUB under the repository nonlinear-waves/stablab_matlab. In addition to the code that was produced, this research also provided several opportunities for professional development. I had the opportunity to present preliminary research findings at the BYU College of Physical and Mathematical Sciences Research Conference earlier this year in front of many colleagues and local figures in mathematics. I also am seeking to publish my research and have submitted my research findings for publication in collaboration with fellow undergraduate Jalen Morgan to the Minnesota Journal of Undergraduate Mathematics.
In the completion of this project, we encountered several difficult challenges that required creative thinking and learning new material and processes. For example, in testing the finite difference schemes that we would use, we experimented with several different sample systems. The first system, Burgers Equation, was a simple partial differential equation of a single variable that was not computationally difficult. However, the next system that we worked with, the Isentropic Gas equations, was a system of two variables in which each variable represented its own function of time and space we needed to solve for. This meant that there were twice as many variables to deal with in the code and that the matrix used to solve for each time step required four times as many indices. This was a consideration we had to keep when generalizing the code, as we needed it to work for systems with any given number of variables.
In addition, as some of the systems we would be working with were nonlinear and we chose to use the Crank- Nicolson finite difference method, it became impossible to develop our numerical solver without the use of Newton’s Method. Some of the challenges we faced in adapting our code to work with a Newton solver included calculating the Jacobian Matrix of a system with finite difference components and the necessity of solving a matrix equation by taking the inverse matrix in Python.
The completion of this research project and the development of generalized time evolution code will simplify the process of studying perturbations of traveling waves in STABLAB. Researchers will be able to focus on understanding the dynamics of the system instead of spending time and effort on the numerical methods. In addition, it allows people who have limited programming experience the opportunity to study these systems.
There are opportunities for future research in further developing this coding package. While the Crank- Nicolson finite difference method proved to be the best option overall for any general system of the form studied, further development of the code could allow the user to choose which finite difference scheme to utilize for specific systems that might benefit better from another method. It could also be diversified to allow for time evolution of additional types of systems, such as reaction diffusion equations.