01_Introduction.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "![title](figs/title.png)\n",
    "\n",
    "<h1><center>Module 01: Introduction to Geophysical Computing</center></h1>\n",
    "\n",
    "\n",
    "For the modern professional geophysicists, the use of computers and algorithms is ubiquitous in all aspects of the daily routine.  Whether using a \"canned\" commercial geophysical software application to or developing new code to address a particular task, employing/implementing some sort of computational algorithm is almost always at the heart of any geophysical task.\n",
    "\n",
    "## The challenge of canned software\n",
    "\n",
    "One of the key challenges geophysicsts face when using commercial geophysical software is that it is often difficult to determine what is going on \"under the hood\".  In many instances it may seem that geophysical results are found in the following way:\n",
    "\n",
    "$$ \n",
    "\\boxed{\\begin{array}{c}\\mathrm{Field}\\\\\\mathrm{Geophysical}\\\\\\mathrm{Data}\\end{array}}\n",
    "\\Longrightarrow \n",
    "\\boxed{\\begin{array}{c}\\mathrm{Magic\\,Black}\\\\\\mathrm{Box\\,of\\, Geophysical}\\\\\\mathrm{Data\\,Processing}\\end{array}}\n",
    "\\Longrightarrow \n",
    "\\boxed{\\begin{array}{c}\\mathrm{Top-notch}\\\\\\mathrm{Geophysical}\\\\\\mathrm{Result}\\end{array}}\n",
    "$$\n",
    "\n",
    "Often there are numerous parameters that go into the blackbox, various combinations of which can greatly affect the output results! It can be extremely frustrating for users - especially when the documentation and usage examples cannot be found, and there is no open-source code to help figure it out.\n",
    "\n",
    "## Python, Scientific Libraries\n",
    "\n",
    "One of the major goals of this course is to strip away some of the **magic black box** nature of how students are commonly introduced to geophysical computing, and to provide a better conceptual idea of what is going on when software is applied to tackle geophysical problems.  That said, it is not our goal to \"reinvent the wheel\" on absolutely every algorithm because this would take far too long and one would never finish finding solutions to the problem at hand.  \n",
    "\n",
    "The approach we are taking in this course tries to find a middle way.  In particular, the notes will focus on the **Python** language and leverage its open-source, diverse and well-validated computational toolkits including **Numpy** and **Scipy** as well as the strong plotting library **Matplotlib**.  \n",
    "\n",
    "However, to fully appreciate the theory behind some of these tools, it is important to code up some algorithms by oneself to help deepen the understanding of just what the algorithm is doing. Moreover, there will be scenarios where one has to develop a full computer program for a specialized task for which no community-based solution exists.  Having a deeper understanding of the geophysical algoritms will be greatly helpful in these scenarios. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## C and Parallel Computing (GPGN536)\n",
    "\n",
    "In many cases, Python may not be the optimal language in which to develop and thus there will be the need to include algorithms directly embedded into the code. In particular, when handling large-scale 3D geophysical computing problems, it is highly desirable to develop **parallel computing** solutions that allow you to leverage modern computational hardware to solve these problem in a fraction of the time that **serial computing** solutions will take. Having a deeper understanding of the geophysical algoritms will be greatly helpful in these scenarios. \n",
    "\n",
    "In this course we will begin to explore these more advanced concepts using the student's choice of [\"low-level\" programming language](https://en.wikipedia.org/wiki/High-level_programming) C and [OpenMP](https://www.openmp.org/). Overall, it is our hope that this will assist students in designing better geophysical computing solutions for real-world three-dimensional (3-D) problems that will facilitate achievement of the student's geophysical research and data processing goals."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Course Aims\n",
    "\n",
    "The course is intended to provide senior-level undergraduate and first-year graduate students with material that aimed at improving the algorithmic, programming and computing skills, while simultaneously enhancing skills in handling geophysical data sets.  In some ways, these skills are those that are not really taught in any one particular geophysics course; however, they are often the \"glue\" that students need to help stick together disparate components learned in specific courses into a workable solution.\n",
    "\n",
    "If we were to broadly classify the material presented in this course, it could be broken down into the following three skill sets:\n",
    "\n",
    "* Computing Skills: \n",
    "    * **Scientific computing using Python (GPGN436) and C (GPGN536)**\n",
    "    * **Working in a Linux Environment (GPGN536)**\n",
    "    * **Introduction to Parallelism with OpenMP (GPGN536)**\n",
    "    \n",
    "    \n",
    "* Numerical Solutions of Geophysical Partial Differential Equations (PDEs): \n",
    "    * **Elliptical PDEs (2D Laplace and Poisson Equations)**\n",
    "    * **Parabolic PDEs (2D Diffusion-Advection Flow)**\n",
    "    * **Hyperbolic PDEs (2D Acoustic wave equation)**\n",
    "    \n",
    "    \n",
    "* Applied Numerical Algorithms and Data Handling\n",
    "    * **Applied Linear Algebra for Geophysical Problems**\n",
    "    * **Differentiation and Discretization**\n",
    "    * **Regression**\n",
    "    * **Introduction to Machine Learning**\n",
    "    \n",
    "We briefly summarize each of these sections below."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Module 02: Numerical Discretization / Differentiation\n",
    "\n",
    "We begin our exploration of numerical methods used to solve partial differential equations (PDEs) such as the multi-dimensional advection+convection equation:\n",
    "\n",
    "$$\\frac{\\partial u}{\\partial t} = \\kappa \\nabla^2 u -\\nabla \\cdot ({\\bf v} u)  + F\\tag{1}$$\n",
    "\n",
    "where ${\\bf v}$ is velocity, $u$ is a scalar field, $\\kappa$ is the diffusivity field, and $F$ represents sources or sinks of quantity $u$. Examples in nature include the movement of, say, a pollutant introduced into a slowing flowing stream.\n",
    "\n",
    "To solve these equations numerically in a computer, we need to solve PDEs on a discrete solution grid that necessarily requires the process of **discretizing the physical system** including the partial differential operators (e.g., $\\frac{\\partial^2 }{\\partial x^2}$) comprising the PDE.  This is commonly done using **finite-difference** approximations, which commonly lead to numerical **stencils** and very **efficient** solution algorithms.\n",
    "\n",
    "For those of you who have taken digital signal processing, you should be familiar with the concept of discretization. To remind you, the left panel below shows an image of the moon. The center panel shows this picture after applying a **convolutional Laplacian filter** for the purpose of emphasizing the discontinuous structure. The right panel shows the **2D finite-difference approximation** of the continuous Laplacian operators $\\nabla^2$.\n",
    "\n",
    "<table><tr>\n",
    "<td> <img src=\"figs/Moon_discretization.jpg\" alt=\"Drawing\" style=\"width: 500px;\"/> </td>\n",
    "<td> <img src=\"figs/Lapfilt.png\" alt=\"Drawing\" style=\"width: 200px;\"/> </td>\n",
    "</tr></table>\n",
    "\n",
    "**Figure 1. Moon image before (left) and after (right) applying a numerical filter representing the 2D Laplacian operator, which is a fundamental component of many of the PDEs studied in this course.**  \n",
    "\n",
    "However, we will go beyond this to explore how different discretizations affect the **stabilty** and **accuracy** of numerical solution of **time-dependent** PDEs.  In this module we will be taking a deeper look at the numerical approximations and discretization with the goal of developing the skills required to solve PDEs of interest like those in Modules 5, 7 and 8."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Module 03: Applied Linear Algebra for Geophysical Problems\n",
    "\n",
    "In your undergraduate prepatory career you have probably taken a semester-long course in Linear Algebra. It is equally likely that this course was pretty abstract and didn't really make a strong connection with geophysics.  However, applied linear algebra comes up all of the time in numerical geophysical computing - especially in the context of generating the solution to a system of linear equations. \n",
    "\n",
    "In many situations we will investigate this semester, the geophysical problem can be represeted by a straightforward matrix equation:\n",
    "\n",
    "$${\\bf Gm}= {\\bf d}$$\n",
    "\n",
    "where ${\\bf m}$ is some geophysical **model parameter** (e.g., acoustic wave speed, thermal conductivity, electric charge), ${\\bf d}$ is some sort of geophysical **data** (e.g., acoustic pressure, heat distibution, electric potential), and ${\\bf G}$ is a numerical representation of the **physics** and often **experimental geometry** that is used to forward model ${\\bf d}$ given ${\\bf m}$.  The following figure illustrates this concept using graphical depictions of ${\\bf G}$, ${\\bf m}$ and ${\\bf d}$.\n",
    "\n",
    "<img src=\"figs/Gm_d_Inversion.jpeg\" width=\"750\">\n",
    "\n",
    "**Figure 2.  Illustration of a modeling operator (left) acting on a model vector (center) to develop data (right).  You will see this paradigm constantly throughout this course!**\n",
    "\n",
    "A good portion of this course looks at how one can solve these types of forward modeling equations in an efficient and straightforward manner using different methods of **applied linear algebra**.  We will also examine scenarios where solutions cannot actually be generated due to various numerical challenges. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Module 04: Elliptical PDEs\n",
    "\n",
    "Having tackled issues pertaining to discretization of partial differential operators and thereby PDEs, we are now ready to move on to actually solving some PDEs of interest in geophysics: **elliptical PDEs**.  One of the most straightforward PDEs to solve is the 2D Poisson's equation for a potential surface $U(x,y,z)$:\n",
    "\n",
    "$$\\nabla^2 U = \\frac{\\partial^2 U}{\\partial x^2} + \\frac{\\partial^2 U}{\\partial y^2}  = f(x,y) \\tag{2}$$\n",
    "\n",
    "which reduces to Laplace's equation in instances where $f(x,y)=0$ throughout the solution domain.  You may recall that one must specify the **boundary conditions** prior to calcuating solutions to these types of equations. \n",
    "\n",
    "The example below shows the potential surface solution related to the distribution of either four point masses (gravitational potential) or the response to four positive point charges (electrical potential).\n",
    "\n",
    "<img src=\"figs/solution_mixt.png\" width=\"600\">\n",
    "\n",
    "**Figure 3. Numerical solution of Poisson's equation solution for four point charges within a domain held at zero potential.**\n",
    "\n",
    "In this module we will explore 2D and 3D numerical approaches that will be couched in the framework of - and be used to motivate - applied linear algebra (Module 03)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Module 05: Parabolic PDEs \n",
    "\n",
    "The next set of PDEs we will study fall into the family of **parabolic** equations.  In particular, we are looking to generate numerical solutions of 2D heat flow, 2D advection equations and 2D diffusion-advection equation, all of which can be written in the following form:\n",
    "\n",
    "$${\\frac {\\partial \\phi (x,y,t)}{\\partial t}}=\\nabla \\cdot {\\big [}D(x,y)\\ \\nabla \\phi (x,y,t){\\big ]} \\tag{3}$$\n",
    "\n",
    "where $\\phi (x,y,z,t)$ is, e.g., the distribution of heat through the 2D solution domain, and $D(x,y,z)$ is the heterogeneous thermal conductivity field.  Computing solutions requires setting both the **initial condition** (i.e., at time $t=0$s) and the **boundary conditions** at the edges of the computational domain.\n",
    "\n",
    "In Figure 5, the left panels show an example of the time evolution of the heat distribution $U(x,y,t)$ for a square computational domain where the four edges are held at $0^\\circ$C and all interior points initially start at $100^\\circ$C. The right panel shows a cross-section through the solution and is easier to see how the $U(x,y=0.5,t)$ evolves over the computed solution time.\n",
    "\n",
    "<img src=\"figs/2dheat.gif\" width=\"800\">\n",
    "\n",
    "**Figure 4. Illustrating of the numerical solution of the 2D heat equation.**\n",
    "\n",
    "In this module we will be developing numerical solutions to the 2D diffusion+advection equation, validating them against analytical solutions, and then using your validated code to solve some interesting geophysical problems for which **no analytical solutions exist**.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Module 06: Hyperbolic PDEs\n",
    "\n",
    "The final set of PDEs that we will study fall into the family of **hyperbolic** equations.  In particular, we are looking to generate numerical solutions of 2D acoustic wave equation:\n",
    "\n",
    "$$\\left[ \\frac{ \\partial^2 }{\\partial x^2}+\\frac{ \\partial^2 }{\\partial y^2} -\\frac{1}{c^2(x,y)} {\\frac {\\partial^2}{\\partial t^2}}\\right] U(x,y,t) = f(x,y,t) \\tag{4}$$\n",
    "\n",
    "where $U(x,y,t)$ is the displacement of the acoustic wave disturbance, and $c(x,y)$ is the heterogeneous acoustic wave speed, and $f(x,y,t)$ is the force source distribution (e.g., what is causing the acoustic wave disturbance).  Computing solutions requires setting the initial conditions and the boundary conditions at the edges of the computational domain.\n",
    "\n",
    "The example below shows a numerical solution to the 2D acoustic wave equation that is modeling the temporal and spatial evolution of a taut square drum head (i.e., clamped boundaries with zero displacement) and some initial amplitude and/or velocity distribution.\n",
    "\n",
    "<img src=\"figs/2dwave.gif\" width=\"400\">\n",
    "\n",
    "**Figure 5. Illustration of the numerical solution of the 2D acoustic wave equation for a taut membrane on a square frame.**\n",
    "\n",
    "In this module we will be developing numerical solutions to the 2D acoustic wave equation, validating them against analytical solutions to a PDE system like in the example above, and then using our validated code to solve some interesting geophysical wave propagation problems for which there are no analytical solution."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Module 07: Interpolation\n",
    "\n",
    "Geophysicists handle data that are acquired at regular and irregular (spatial and/or temporal intervals); however, most of the time the locations where data are acquired are insufficient in number or not in all of the desired locations required to make the corresponding geophysical or geological interpretation.   \n",
    "\n",
    "Let's look a the following example.  For simplicitly, let's say that we have 30+ elevation data points that were acquired over a 1 km$^2$ area.  An important question that will no doubt be asked is: *Based on this data, can we determine the elevation profile throughout the entire area?*\n",
    "\n",
    "<img src=\"figs/2D_Interpolation_intro.png\" width=\"500\">\n",
    "\n",
    "**Figure 6. Illustration of an underlying interpolation solution developed from the given points in red.**\n",
    "\n",
    "The good news is that yes we can create high-density maps using the acquired data points; however, the challenge is that there is no one unique map that we can create from the data that would be arguably correct.  To do this, we we will explore some of the fundamentals of **interpolation** and see a number of different ways in which plausible maps can be created in 1D and 2D.  Many of these algorithms even can be exported to higher dimensions for those cases where 3D, 4D and even 5D interpolation is required; however, these fall beyond the scope of this module.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# Module 08: Regression\n",
    "\n",
    "In previous modules you have learned how to stimulate the physics using numerical models. However, in reality, models can never fit the data perfectly, due to the various kinds of noise existing in the measurement. Finding the model that \"best\" fits the data is called a regression problem, which is the core part of many geophysicists' job responsibility. \n",
    "\n",
    "<img src=\"figs/regression_intro.png\" width=\"300\">\n",
    "\n",
    "**Figure 7: Fitting the production rate of an oil producer using a physical model (exponential decay).**\n",
    "\n",
    "Most regression problems can be presented like this:\n",
    "\n",
    "   * Find $m$ such that $||G(m)-d|| + \\lambda \\Gamma(m)$ is minimized. \n",
    "\n",
    "Where $G$ is the model function that links the model parameter $m$ with the observations $d$. $G(m)$ is the model-predict data, while $d$ is the actually measured one. Symbol $||\\dots||$ represents some norm of the difference between the predicted data and real data. $\\Gamma(m)$ is some regularization term, representing our a priori knowledge of the model. \n",
    "\n",
    "Two catalogs of regression problems will be discussed in this section: linear and non-linear regression. For the linear regression problem, both model function $G$ and regularization function $\\Gamma$ can be presented as a system of linear functions, and the problem usually has an analytical solution. On the other hand, the non-linear regression problem has to be solved numerically and several algorithms will be introduced. \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Module 09: Introduction to Machine Learning\n",
    "\n",
    "In this chapter, we will introduce some fundamentals of machine learning, as well as the associated python packages. \n",
    "\n",
    "Machine-learning models, or statistical-learning models, can be considered as highly adaptive empirical functions that project the input data X with output y. In mathematics, it can be presented as:\n",
    "\n",
    "   * Find $p$ such that $||G(X,p)-y||$ is minimized. \n",
    "\n",
    "where $p$ is the parameters to be \"trained\". You may find that this formula is quite similar to the regression problem from the last chapter: they are fundamentally the same, except machine-learning algorithms use statistical models instead of physical ones.\n",
    "\n",
    "We will learn several popular machine learning models and basic machine learning workflow: data processing, feature engineering, hyperparameter selection, and model validation.\n",
    "\n",
    "The following figure shows the same data fitted by a machine-learning based model, except in this case, no physical knowledge is required to model the data. \n",
    "\n",
    "<img src=\"figs/machine_learning_intro.png\" width=\"300\">\n",
    "\n",
    "**Figure 8: Fitting the production rate of an oil producer using a machine-learning model (multiple-layer neural network).**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.6"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": false,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {
    "height": "calc(100% - 180px)",
    "left": "10px",
    "top": "150px",
    "width": "241.455px"
   },
   "toc_section_display": true,
   "toc_window_display": true
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}