{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Finite difference modeling\n",
"\n",
"*Prof. dr.ir.T.N.Olsthoorn*\n",
"\n",
"*Heemstede, Sept. 2016, 24 May 2017*"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Approach\n",
"\n",
"In this chapter we set-up a 3D steady-state finite difference model from scratch. We do this by computing a numerical groundwater problem step by step, by hand, using finite difference, building up the pieces of the model, which we will assemble in the next chapter."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Setup of the model by specifying its dimensions"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The 3D steady state FDM will be based on a regular grid consisting of rows an columns and layers. The column widths and the row heigts are constant on a per column and per row basis, but the layer thickness can vary on a cell by cell basis. The grid of a full 3D model will thus be specified in general by a vector of x cell boundary coordinates, a vector y row boundary coordinates and a full 3D array of cell top and bottom coordinates.\n",
"\n",
"Notice that the arrays are interpreted as [z, y, x] or [layer row col]. This is a convenience in Python where when Phi is a 3D array of the shape of the grid [Nz, Ny, Nz] we have\n",
"\n",
"Phi[k].shape is [Ny, Nx], the entire layer number i.\n",
"Phi[k][j] = Nx, the entire row j of layer i.\n",
"Phi[k][j][i] = the head in cell [k, j, i] which is the same as Phi[k, j, i]"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"import numpy as np\n",
"\n",
"# specify a rectangular grid\n",
"x = np.arange(-1000., 1000., 25.)\n",
"y = np.arange(-1000., 1000., 25.) # backward, i.e. first row grid line has highest y\n",
"z = np.arange(-100., 0., 20.) # backward, i.e. from top to bottom"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"From these coordinates we obtain the number of cells along each axis and the cell sizes and"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"# as well as the number of cells along the three axes\n",
"Nx = len(x)-1\n",
"Ny = len(y)-1\n",
"Nz = len(z)-1\n",
"\n",
"sz = (Nz,Ny,Nx) # the shape of the model\n",
"Nod = np.prod(sz) # total number of cells in the model\n",
"\n",
"# from this we have the width of columns, rows and layers\n",
"dx = np.diff(x).reshape(1, 1, Nx)\n",
"dy = np.diff(y).reshape(1, Ny,1)\n",
"dz = np.abs(np.diff(z)).reshape(Nz, 1,1)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## IBOUND array - telling which cells are active and which have a prescribed head"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let's first specify which of the cells have their head prescrided and which cells are inactive. We have to tackle inactive cells early to make sure their conductance is made zero (in case there conductivities might be specified as non-zeros).\n",
"\n",
"We do that my means of a so-called boundary array IBOUND (MODFLOW terminology), which is an integer array of the shape of the model grid that tells which cells have a prescribed head, which cells are inactive (i.e. which cells does not take part of the computation, such as cells that represent impermeable rock) and for which cells the head should be computed.\n",
"\n",
"* IBOUND > 0, means heads will be computed\n",
"* IBOUND == 0, means cells are inactive\n",
"* IBOUND <0 , means heads prescribed\n",
"\n",
"In this particular example we specify that the vertical zx plane at the last row of the model will have prescribed heads equal to zero."
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"IBOUND = np.ones(sz)\n",
"IBOUND[:,-1,:] = -1 # last row of model heads are prescribed\n",
"IBOUND[:, 40:45, 20:70]=0 # these cells are inactive"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"This boundary array makes it easy telling which cells cells are active (head computed), inactive, and fixed-head."
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"active = (IBOUND>0).reshape(Nod) # active is now a vector of booleans of length Nod\n",
"inact = (IBOUND==0).reshape(Nod) # dito for inact\n",
"fxhd = (IBOUND<0).reshape(Nod) # dito for fxhd"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Cell conductancies: defining the ease of flow between adjacent cells"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The first thing to define based on the properties of the cells is the flow resistance of each cell in the 3 grid directions, `x`, `y` and `z`. For that we need the cell sizes from the coordinates and the hydraulic conductivities in the `x`, `y` and `z` direction. The latter are given as full 3D arrays `kx`, `ky`, `kz` whose shapes correspond to that of the model mesh."
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"k = 10.0 # m/d uniform conductivity\n",
"kx = k * np.ones(sz) # [L/T] 3D kx array\n",
"ky = k * np.ones(sz) # [L/T] 3D ky array with same values as kx \n",
"kz = k * np.ones(sz) # [L/T] 3D kz array with same values as kx"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The flow resistances for each cell is the head loss across opposite cell faces due to a unit flux through the cell along the axis perperndicular to them. These resistances are cell properties that can immediately be computed for the entire grid of the model. Because we always need the resistance between the cell center and its outer faces, we use the factor 0.5 (flow over half the lenght of the cell in each direction)"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"# half cell flow resistances\n",
"Rx = 0.5 * dx / (dy * dz) / kx # [T/L2], flow resistance half cell in x-direction\n",
"Ry = 0.5 * dy / (dz * dx) / ky # same in y-direction\n",
"Rz = 0.5 * dz / (dx * dy) / kz # same in z-direction"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Make inactive cells inactive by setting their resistance to np.Inf (infinite):"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"Rx = Rx.reshape(Nod,); Rx[inact] = np.Inf; Rx=Rx.reshape(sz)\n",
"Ry = Ry.reshape(Nod,); Ry[inact] = np.Inf; Ry=Ry.reshape(sz)\n",
"Rz = Rz.reshape(Nod,); Rz[inact] = np.Inf; Rz=Rz.reshape(sz)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"From this we compute the conductance between each pair of adjacent cells across their connecting cell face. The conductance is just the reverse of the resistance of the two connected half cells. This resistance is the sum of the resistances of the two connected half cells because these resistances are placed in series with respect to the flow."
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"# conductances between adjacent cells\n",
"Cx = 1 / (Rx[:, :, :-1] + Rx[:, :,1:]) # [L2/T] in x-direction\n",
"Cy = 1 / (Ry[:, :-1,:] + Ry[:, 1:,:]) # idem in y-direction\n",
"Cz = 1 / (Rz[:-1,:,:] + Rz[1:,:,:]) # idem in z-direction"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Setting up the system matrix - set of water balance equations"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The system matrix has size of (`Nod`, `Nod`) allowing a connection between each pair of cells. Of course only cells that share their cell face are connected in reality. In a 3D model this means that each cell is connected to its 6 neighbors instead of to all other cells in the model. This means that most of the matrix entries will be zero.\n",
"\n",
"To be able to indentify adjacent cells we generate cell numbers in an array that has the size of the model grid:"
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"NOD = np.arange(Nod).reshape(sz) # this is a full 3D array of node numbers (cell numbers)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"With this array it's easy to identify adjacent cells by their cell number. Thus we generate arrays with the cel numbers of right hand neigbor of the cells (east neighbor), the left hand neighbor (the west neigbor), the north neighbor, south neighbor, the top neighbor and the bottom neighbor as follows"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"IE = NOD[:, :, 1:] # numbers of the eastern neighbors of each cell\n",
"IW = NOD[:, :, :-1] # same western neighbors\n",
"IN = NOD[:, :-1,:] # same northern neighbors\n",
"IS = NOD[:, 1:,:] # southern neighbors\n",
"IT = NOD[:-1,:,:] # top neighbors\n",
"IB = NOD[1:,:,:] # bottom neighbors"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Notice that the shape of the `IE` and `IW` is the same as that of `Cx`, the size of `IN` and `IS` is the same as that of `Cy` and the size of `IT` and `IB` is the same as that of `Cx`.\n",
"\n",
"To put the conductances into the system matrix we need their row and column indices together with their value, so that we can say `a`[`j`,`i`] = `value`. Because we have the numbers of adjacent cells in the arrays `IE`, `IW` etc, we can immediately place all the system matrix coefficiencts at the place into a sparse matrix."
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {
"collapsed": true,
"scrolled": false
},
"outputs": [],
"source": [
"import scipy.sparse as sp\n",
"\n",
"R = lambda x : x.ravel() # define short hand for x.ravel()\n",
"\n",
"# notice the call signature:\n",
"# csc_matrix( (data, (row_index, col_index) ), (M,N)); This is a tuple within tuple.\n",
"A = sp.csc_matrix(( np.concatenate(( R(Cx), R(Cx), R(Cy), R(Cy), R(Cz), R(Cz)) ),\\\n",
" (np.concatenate(( R(IE), R(IW), R(IN), R(IS), R(IB), R(IT)) ),\\\n",
" np.concatenate(( R(IW), R(IE), R(IS), R(IN), R(IT), R(IB)) ),\\\n",
" )),(Nod,Nod))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We now have to define the diagonal elements of the system matrix `A`, i.e. the values `a`[`i`,`i`] for `i`=[0:`Nod`].\n",
"\n",
"These are just the negative sum of the row coefficients. Hence we sum `A` over the second axis (`axis`=1) to get them in a [`Nod`,1] sized vector. (Notice stat sparace matrix derived vectors keep their orientation, contrary to vectors obained from numpy arrays, which produce dimensionless vectors).\n",
"\n",
"Generate the diagonal values:"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"# to use the vector of diagonal values int a call of sp.diags() we need to have it aa a \n",
"# standard nondimensional numpy vector.\n",
"# To get this:\n",
"# - first turn the matrix obtained by A.sum(axis=1) into a np.array by np.array( .. )\n",
"# - then take the whole column to loose the array orientation (to get a dimensionless numpy vector)\n",
"adiag = np.array(-A.sum(axis=1))[:,0]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Then generate a diagonal array from these values, that we can add it to `A`, to complete `A`."
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"Adiag = sp.diags(adiag)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"More complex alternative: Generate diagonal array by calling the csr_matrix constructor:\n",
"\n",
"Adiag = sp.csr_matrix((adiag,(np.arange(Nod),np.arange(Nod))),(Nod,Nod))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Boundary conditions\n",
"\n",
"For this chapter we only use fixed flow and fixed head boundary conditions."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Fixed flows"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Fixed flow boundary conditions are specified by an 3D array of the size of the grid. Each values specifies the inflow for the corresponding cell (injections are positive). Cells without a specified flow are, in fact, cells where the specified flow is zero. Hence the fixed-flows array is a full 3D array with flow values that are zero where no flow enters or leaves the cells and have non-zero values elsewhere.\n",
"\n",
"For this example, we specify a single extraction of Q=-1200 m3/d in cell [30,25,2]:"
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"FQ = np.zeros(sz) # all flows zero. Note sz is the shape of the model grid\n",
"FQ[2, 30, 25] = -1200 # [m3/d] extraction in this cell"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The righ-hand size of the matrix equation to be solved, the vector RHS, contains the flows. So we can generate it by assignment of `FQ` and converting it to a numpy vector"
]
},
{
"cell_type": "code",
"execution_count": 20,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"RHS = FQ.reshape(Nod)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"See further down how we use `RHS` for only the active and non-fixed head rows.\n",
"\n",
"The next step is to add fixed head boundary conditions."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### Fixed heads"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Fixed heads are known heads. This implies that in the set of equations that represent the model, i.e\n",
"\n",
"$A \\times Phi = RHS$\n",
"\n",
"Some of the Phis are prescribed and should not be computed as defined by `IBOUND` and contained in the boolean vectors `active`, `fxhd` and `inact` specified and computed above."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Now that we know which cell have fixed heads, we can multiply out these heads with the corresponding columns of the system matrix, which yields a vector of constant values with dimension flow [m3/d] that can be added to the fixed flow vector in the `RHS` vector. The `RHS` vector is now the sum of the `FQ` and the contribution from the fixed heads.\n",
"\n",
"Notice that the fixed heads will be obtained from the given array `HI` of the initial heads, where the head in the cells where `IBOUND`>0 correspond with the fixed heads."
]
},
{
"cell_type": "code",
"execution_count": 21,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"HI = np.zeros(sz)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We reshape `FQ` and `HI` to a column vector to allow matrix multiplication"
]
},
{
"cell_type": "code",
"execution_count": 22,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"RHS = FQ.reshape(Nod,1) - A[:,fxhd].dot(HI.reshape(Nod,1)[fxhd])"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We have now the complete `RHS` of the matrix equation to solve:\n",
"\n",
"$A \\times Phi = RHS$"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Solving the matrix equation for the unknown heads"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"We use the sparse matrix solver in module `scipy.sparse.linalg` to compute the unknown heads."
]
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"from scipy.sparse.linalg import spsolve # import with from to use its short name"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Of course we only need the active rows and columns of `A` and the active rows from `RHS`.\n",
"\n",
"But first allocate a full-fledged vector of heads to store the result."
]
},
{
"cell_type": "code",
"execution_count": 24,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"Phi = HI.flatten() "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Then compute the unknown heads (i.e. the active cels only).\n",
"\n",
"Remark: If we want to select a submatrix from `A` defiend by a given vectors of row and column indices, we can do so in sequence: Rows (I) first, columns (J) next, like so:\n",
"\n",
"$A[I][:,J]$\n",
"\n",
"which we apply in the next line"
]
},
{
"cell_type": "code",
"execution_count": 25,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
" Phi[active] = spsolve( (A+Adiag)[active][:,active] ,RHS[active] )"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"At this point we solved the problem and now have the heads for all cells in the vector Phi.\n",
"\n",
"We didn't touch the rows and columns that are inactive. So the heads of these inactive cells whatever they are in `HI` are know still in `Phi`. Just to make sure we detect them and won't use them, set them to `NaN` (Not a Number)."
]
},
{
"cell_type": "code",
"execution_count": 26,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"Phi[inact] = np.NaN"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Finally we reshape the head vector to that of the model grid."
]
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"Phi=Phi.reshape(sz) # reshape vector Phi to 3D shape of the grid"
]
},
{
"cell_type": "code",
"execution_count": 28,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"array([[[ 1.77107246, 1.77132861, 1.77183761, ..., 1.56399507,\n",
" 1.56320139, 1.56280385],\n",
" [ 1.77081631, 1.77107575, 1.77159133, ..., 1.56360172,\n",
" 1.56280525, 1.56240631],\n",
" [ 1.77030071, 1.77056677, 1.77109553, ..., 1.56281367,\n",
" 1.56201158, 1.56160982],\n",
" ..., \n",
" [ 0.05529699, 0.05524583, 0.0551436 , ..., 0.03419641,\n",
" 0.03426788, 0.03430386],\n",
" [ 0.02763571, 0.02761016, 0.02755909, ..., 0.01708962,\n",
" 0.01712512, 0.01714299],\n",
" [ 0. , 0. , 0. , ..., 0. ,\n",
" 0. , 0. ]],\n",
"\n",
" [[ 1.77107246, 1.77132861, 1.77183761, ..., 1.56399507,\n",
" 1.56320139, 1.56280385],\n",
" [ 1.77081631, 1.77107575, 1.77159133, ..., 1.56360172,\n",
" 1.56280525, 1.56240631],\n",
" [ 1.77030071, 1.77056677, 1.77109553, ..., 1.56281367,\n",
" 1.56201158, 1.56160982],\n",
" ..., \n",
" [ 0.05529699, 0.05524583, 0.0551436 , ..., 0.03419641,\n",
" 0.03426788, 0.03430386],\n",
" [ 0.02763571, 0.02761016, 0.02755909, ..., 0.01708962,\n",
" 0.01712512, 0.01714299],\n",
" [ 0. , 0. , 0. , ..., 0. ,\n",
" 0. , 0. ]],\n",
"\n",
" [[ 1.77107246, 1.77132861, 1.77183761, ..., 1.56399507,\n",
" 1.56320139, 1.56280385],\n",
" [ 1.77081631, 1.77107575, 1.77159133, ..., 1.56360172,\n",
" 1.56280525, 1.56240631],\n",
" [ 1.77030071, 1.77056677, 1.77109553, ..., 1.56281367,\n",
" 1.56201158, 1.56160982],\n",
" ..., \n",
" [ 0.05529699, 0.05524583, 0.0551436 , ..., 0.03419641,\n",
" 0.03426788, 0.03430386],\n",
" [ 0.02763571, 0.02761016, 0.02755909, ..., 0.01708962,\n",
" 0.01712512, 0.01714299],\n",
" [ 0. , 0. , 0. , ..., 0. ,\n",
" 0. , 0. ]],\n",
"\n",
" [[ 1.77107246, 1.77132861, 1.77183761, ..., 1.56399507,\n",
" 1.56320139, 1.56280385],\n",
" [ 1.77081631, 1.77107575, 1.77159133, ..., 1.56360172,\n",
" 1.56280525, 1.56240631],\n",
" [ 1.77030071, 1.77056677, 1.77109553, ..., 1.56281367,\n",
" 1.56201158, 1.56160982],\n",
" ..., \n",
" [ 0.05529699, 0.05524583, 0.0551436 , ..., 0.03419641,\n",
" 0.03426788, 0.03430386],\n",
" [ 0.02763571, 0.02761016, 0.02755909, ..., 0.01708962,\n",
" 0.01712512, 0.01714299],\n",
" [ 0. , 0. , 0. , ..., 0. ,\n",
" 0. , 0. ]]])"
]
},
"execution_count": 28,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Phi # show Phi"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Plotting the heads as contours"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Import the required plotting module and setup the plot."
]
},
{
"cell_type": "code",
"execution_count": 29,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"%matplotlib notebook\n",
"import matplotlib.pyplot as plt # combines namespace of numpy and pyplot"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"For coordinates of the cells use their centers."
]
},
{
"cell_type": "code",
"execution_count": 30,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"xm = 0.5 * (x[:-1] + x[1:]) # [L] coordinates of column centers\n",
"ym = 0.5 * (y[:-1] + y[1:]) # [L] coordinates of row centers\n",
"layer = 2 # contours for this layer\n",
"nc = 50 # number of contours in total"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Plot the results using plt functions like in Matlab"
]
},
{
"cell_type": "code",
"execution_count": 31,
"metadata": {},
"outputs": [
{
"data": {
"application/javascript": [
"/* Put everything inside the global mpl namespace */\n",
"window.mpl = {};\n",
"\n",
"mpl.get_websocket_type = function() {\n",
" if (typeof(WebSocket) !== 'undefined') {\n",
" return WebSocket;\n",
" } else if (typeof(MozWebSocket) !== 'undefined') {\n",
" return MozWebSocket;\n",
" } else {\n",
" alert('Your browser does not have WebSocket support.' +\n",
" 'Please try Chrome, Safari or Firefox ≥ 6. ' +\n",
" 'Firefox 4 and 5 are also supported but you ' +\n",
" 'have to enable WebSockets in about:config.');\n",
" };\n",
"}\n",
"\n",
"mpl.figure = function(figure_id, websocket, ondownload, parent_element) {\n",
" this.id = figure_id;\n",
"\n",
" this.ws = websocket;\n",
"\n",
" this.supports_binary = (this.ws.binaryType != undefined);\n",
"\n",
" if (!this.supports_binary) {\n",
" var warnings = document.getElementById(\"mpl-warnings\");\n",
" if (warnings) {\n",
" warnings.style.display = 'block';\n",
" warnings.textContent = (\n",
" \"This browser does not support binary websocket messages. \" +\n",
" \"Performance may be slow.\");\n",
" }\n",
" }\n",
"\n",
" this.imageObj = new Image();\n",
"\n",
" this.context = undefined;\n",
" this.message = undefined;\n",
" this.canvas = undefined;\n",
" this.rubberband_canvas = undefined;\n",
" this.rubberband_context = undefined;\n",
" this.format_dropdown = undefined;\n",
"\n",
" this.image_mode = 'full';\n",
"\n",
" this.root = $('');\n",
" this._root_extra_style(this.root)\n",
" this.root.attr('style', 'display: inline-block');\n",
"\n",
" $(parent_element).append(this.root);\n",
"\n",
" this._init_header(this);\n",
" this._init_canvas(this);\n",
" this._init_toolbar(this);\n",
"\n",
" var fig = this;\n",
"\n",
" this.waiting = false;\n",
"\n",
" this.ws.onopen = function () {\n",
" fig.send_message(\"supports_binary\", {value: fig.supports_binary});\n",
" fig.send_message(\"send_image_mode\", {});\n",
" fig.send_message(\"refresh\", {});\n",
" }\n",
"\n",
" this.imageObj.onload = function() {\n",
" if (fig.image_mode == 'full') {\n",
" // Full images could contain transparency (where diff images\n",
" // almost always do), so we need to clear the canvas so that\n",
" // there is no ghosting.\n",
" fig.context.clearRect(0, 0, fig.canvas.width, fig.canvas.height);\n",
" }\n",
" fig.context.drawImage(fig.imageObj, 0, 0);\n",
" };\n",
"\n",
" this.imageObj.onunload = function() {\n",
" this.ws.close();\n",
" }\n",
"\n",
" this.ws.onmessage = this._make_on_message_function(this);\n",
"\n",
" this.ondownload = ondownload;\n",
"}\n",
"\n",
"mpl.figure.prototype._init_header = function() {\n",
" var titlebar = $(\n",
" '');\n",
" var titletext = $(\n",
" '');\n",
" titlebar.append(titletext)\n",
" this.root.append(titlebar);\n",
" this.header = titletext[0];\n",
"}\n",
"\n",
"\n",
"\n",
"mpl.figure.prototype._canvas_extra_style = function(canvas_div) {\n",
"\n",
"}\n",
"\n",
"\n",
"mpl.figure.prototype._root_extra_style = function(canvas_div) {\n",
"\n",
"}\n",
"\n",
"mpl.figure.prototype._init_canvas = function() {\n",
" var fig = this;\n",
"\n",
" var canvas_div = $('');\n",
"\n",
" canvas_div.attr('style', 'position: relative; clear: both; outline: 0');\n",
"\n",
" function canvas_keyboard_event(event) {\n",
" return fig.key_event(event, event['data']);\n",
" }\n",
"\n",
" canvas_div.keydown('key_press', canvas_keyboard_event);\n",
" canvas_div.keyup('key_release', canvas_keyboard_event);\n",
" this.canvas_div = canvas_div\n",
" this._canvas_extra_style(canvas_div)\n",
" this.root.append(canvas_div);\n",
"\n",
" var canvas = $('');\n",
" canvas.addClass('mpl-canvas');\n",
" canvas.attr('style', \"left: 0; top: 0; z-index: 0; outline: 0\")\n",
"\n",
" this.canvas = canvas[0];\n",
" this.context = canvas[0].getContext(\"2d\");\n",
"\n",
" var rubberband = $('');\n",
" rubberband.attr('style', \"position: absolute; left: 0; top: 0; z-index: 1;\")\n",
"\n",
" var pass_mouse_events = true;\n",
"\n",
" canvas_div.resizable({\n",
" start: function(event, ui) {\n",
" pass_mouse_events = false;\n",
" },\n",
" resize: function(event, ui) {\n",
" fig.request_resize(ui.size.width, ui.size.height);\n",
" },\n",
" stop: function(event, ui) {\n",
" pass_mouse_events = true;\n",
" fig.request_resize(ui.size.width, ui.size.height);\n",
" },\n",
" });\n",
"\n",
" function mouse_event_fn(event) {\n",
" if (pass_mouse_events)\n",
" return fig.mouse_event(event, event['data']);\n",
" }\n",
"\n",
" rubberband.mousedown('button_press', mouse_event_fn);\n",
" rubberband.mouseup('button_release', mouse_event_fn);\n",
" // Throttle sequential mouse events to 1 every 20ms.\n",
" rubberband.mousemove('motion_notify', mouse_event_fn);\n",
"\n",
" rubberband.mouseenter('figure_enter', mouse_event_fn);\n",
" rubberband.mouseleave('figure_leave', mouse_event_fn);\n",
"\n",
" canvas_div.on(\"wheel\", function (event) {\n",
" event = event.originalEvent;\n",
" event['data'] = 'scroll'\n",
" if (event.deltaY < 0) {\n",
" event.step = 1;\n",
" } else {\n",
" event.step = -1;\n",
" }\n",
" mouse_event_fn(event);\n",
" });\n",
"\n",
" canvas_div.append(canvas);\n",
" canvas_div.append(rubberband);\n",
"\n",
" this.rubberband = rubberband;\n",
" this.rubberband_canvas = rubberband[0];\n",
" this.rubberband_context = rubberband[0].getContext(\"2d\");\n",
" this.rubberband_context.strokeStyle = \"#000000\";\n",
"\n",
" this._resize_canvas = function(width, height) {\n",
" // Keep the size of the canvas, canvas container, and rubber band\n",
" // canvas in synch.\n",
" canvas_div.css('width', width)\n",
" canvas_div.css('height', height)\n",
"\n",
" canvas.attr('width', width);\n",
" canvas.attr('height', height);\n",
"\n",
" rubberband.attr('width', width);\n",
" rubberband.attr('height', height);\n",
" }\n",
"\n",
" // Set the figure to an initial 600x600px, this will subsequently be updated\n",
" // upon first draw.\n",
" this._resize_canvas(600, 600);\n",
"\n",
" // Disable right mouse context menu.\n",
" $(this.rubberband_canvas).bind(\"contextmenu\",function(e){\n",
" return false;\n",
" });\n",
"\n",
" function set_focus () {\n",
" canvas.focus();\n",
" canvas_div.focus();\n",
" }\n",
"\n",
" window.setTimeout(set_focus, 100);\n",
"}\n",
"\n",
"mpl.figure.prototype._init_toolbar = function() {\n",
" var fig = this;\n",
"\n",
" var nav_element = $('')\n",
" nav_element.attr('style', 'width: 100%');\n",
" this.root.append(nav_element);\n",
"\n",
" // Define a callback function for later on.\n",
" function toolbar_event(event) {\n",
" return fig.toolbar_button_onclick(event['data']);\n",
" }\n",
" function toolbar_mouse_event(event) {\n",
" return fig.toolbar_button_onmouseover(event['data']);\n",
" }\n",
"\n",
" for(var toolbar_ind in mpl.toolbar_items) {\n",
" var name = mpl.toolbar_items[toolbar_ind][0];\n",
" var tooltip = mpl.toolbar_items[toolbar_ind][1];\n",
" var image = mpl.toolbar_items[toolbar_ind][2];\n",
" var method_name = mpl.toolbar_items[toolbar_ind][3];\n",
"\n",
" if (!name) {\n",
" // put a spacer in here.\n",
" continue;\n",
" }\n",
" var button = $('');\n",
" button.addClass('ui-button ui-widget ui-state-default ui-corner-all ' +\n",
" 'ui-button-icon-only');\n",
" button.attr('role', 'button');\n",
" button.attr('aria-disabled', 'false');\n",
" button.click(method_name, toolbar_event);\n",
" button.mouseover(tooltip, toolbar_mouse_event);\n",
"\n",
" var icon_img = $('');\n",
" icon_img.addClass('ui-button-icon-primary ui-icon');\n",
" icon_img.addClass(image);\n",
" icon_img.addClass('ui-corner-all');\n",
"\n",
" var tooltip_span = $('');\n",
" tooltip_span.addClass('ui-button-text');\n",
" tooltip_span.html(tooltip);\n",
"\n",
" button.append(icon_img);\n",
" button.append(tooltip_span);\n",
"\n",
" nav_element.append(button);\n",
" }\n",
"\n",
" var fmt_picker_span = $('');\n",
"\n",
" var fmt_picker = $('');\n",
" fmt_picker.addClass('mpl-toolbar-option ui-widget ui-widget-content');\n",
" fmt_picker_span.append(fmt_picker);\n",
" nav_element.append(fmt_picker_span);\n",
" this.format_dropdown = fmt_picker[0];\n",
"\n",
" for (var ind in mpl.extensions) {\n",
" var fmt = mpl.extensions[ind];\n",
" var option = $(\n",
" '', {selected: fmt === mpl.default_extension}).html(fmt);\n",
" fmt_picker.append(option)\n",
" }\n",
"\n",
" // Add hover states to the ui-buttons\n",
" $( \".ui-button\" ).hover(\n",
" function() { $(this).addClass(\"ui-state-hover\");},\n",
" function() { $(this).removeClass(\"ui-state-hover\");}\n",
" );\n",
"\n",
" var status_bar = $('');\n",
" nav_element.append(status_bar);\n",
" this.message = status_bar[0];\n",
"}\n",
"\n",
"mpl.figure.prototype.request_resize = function(x_pixels, y_pixels) {\n",
" // Request matplotlib to resize the figure. Matplotlib will then trigger a resize in the client,\n",
" // which will in turn request a refresh of the image.\n",
" this.send_message('resize', {'width': x_pixels, 'height': y_pixels});\n",
"}\n",
"\n",
"mpl.figure.prototype.send_message = function(type, properties) {\n",
" properties['type'] = type;\n",
" properties['figure_id'] = this.id;\n",
" this.ws.send(JSON.stringify(properties));\n",
"}\n",
"\n",
"mpl.figure.prototype.send_draw_message = function() {\n",
" if (!this.waiting) {\n",
" this.waiting = true;\n",
" this.ws.send(JSON.stringify({type: \"draw\", figure_id: this.id}));\n",
" }\n",
"}\n",
"\n",
"\n",
"mpl.figure.prototype.handle_save = function(fig, msg) {\n",
" var format_dropdown = fig.format_dropdown;\n",
" var format = format_dropdown.options[format_dropdown.selectedIndex].value;\n",
" fig.ondownload(fig, format);\n",
"}\n",
"\n",
"\n",
"mpl.figure.prototype.handle_resize = function(fig, msg) {\n",
" var size = msg['size'];\n",
" if (size[0] != fig.canvas.width || size[1] != fig.canvas.height) {\n",
" fig._resize_canvas(size[0], size[1]);\n",
" fig.send_message(\"refresh\", {});\n",
" };\n",
"}\n",
"\n",
"mpl.figure.prototype.handle_rubberband = function(fig, msg) {\n",
" var x0 = msg['x0'];\n",
" var y0 = fig.canvas.height - msg['y0'];\n",
" var x1 = msg['x1'];\n",
" var y1 = fig.canvas.height - msg['y1'];\n",
" x0 = Math.floor(x0) + 0.5;\n",
" y0 = Math.floor(y0) + 0.5;\n",
" x1 = Math.floor(x1) + 0.5;\n",
" y1 = Math.floor(y1) + 0.5;\n",
" var min_x = Math.min(x0, x1);\n",
" var min_y = Math.min(y0, y1);\n",
" var width = Math.abs(x1 - x0);\n",
" var height = Math.abs(y1 - y0);\n",
"\n",
" fig.rubberband_context.clearRect(\n",
" 0, 0, fig.canvas.width, fig.canvas.height);\n",
"\n",
" fig.rubberband_context.strokeRect(min_x, min_y, width, height);\n",
"}\n",
"\n",
"mpl.figure.prototype.handle_figure_label = function(fig, msg) {\n",
" // Updates the figure title.\n",
" fig.header.textContent = msg['label'];\n",
"}\n",
"\n",
"mpl.figure.prototype.handle_cursor = function(fig, msg) {\n",
" var cursor = msg['cursor'];\n",
" switch(cursor)\n",
" {\n",
" case 0:\n",
" cursor = 'pointer';\n",
" break;\n",
" case 1:\n",
" cursor = 'default';\n",
" break;\n",
" case 2:\n",
" cursor = 'crosshair';\n",
" break;\n",
" case 3:\n",
" cursor = 'move';\n",
" break;\n",
" }\n",
" fig.rubberband_canvas.style.cursor = cursor;\n",
"}\n",
"\n",
"mpl.figure.prototype.handle_message = function(fig, msg) {\n",
" fig.message.textContent = msg['message'];\n",
"}\n",
"\n",
"mpl.figure.prototype.handle_draw = function(fig, msg) {\n",
" // Request the server to send over a new figure.\n",
" fig.send_draw_message();\n",
"}\n",
"\n",
"mpl.figure.prototype.handle_image_mode = function(fig, msg) {\n",
" fig.image_mode = msg['mode'];\n",
"}\n",
"\n",
"mpl.figure.prototype.updated_canvas_event = function() {\n",
" // Called whenever the canvas gets updated.\n",
" this.send_message(\"ack\", {});\n",
"}\n",
"\n",
"// A function to construct a web socket function for onmessage handling.\n",
"// Called in the figure constructor.\n",
"mpl.figure.prototype._make_on_message_function = function(fig) {\n",
" return function socket_on_message(evt) {\n",
" if (evt.data instanceof Blob) {\n",
" /* FIXME: We get \"Resource interpreted as Image but\n",
" * transferred with MIME type text/plain:\" errors on\n",
" * Chrome. But how to set the MIME type? It doesn't seem\n",
" * to be part of the websocket stream */\n",
" evt.data.type = \"image/png\";\n",
"\n",
" /* Free the memory for the previous frames */\n",
" if (fig.imageObj.src) {\n",
" (window.URL || window.webkitURL).revokeObjectURL(\n",
" fig.imageObj.src);\n",
" }\n",
"\n",
" fig.imageObj.src = (window.URL || window.webkitURL).createObjectURL(\n",
" evt.data);\n",
" fig.updated_canvas_event();\n",
" fig.waiting = false;\n",
" return;\n",
" }\n",
" else if (typeof evt.data === 'string' && evt.data.slice(0, 21) == \"data:image/png;base64\") {\n",
" fig.imageObj.src = evt.data;\n",
" fig.updated_canvas_event();\n",
" fig.waiting = false;\n",
" return;\n",
" }\n",
"\n",
" var msg = JSON.parse(evt.data);\n",
" var msg_type = msg['type'];\n",
"\n",
" // Call the \"handle_{type}\" callback, which takes\n",
" // the figure and JSON message as its only arguments.\n",
" try {\n",
" var callback = fig[\"handle_\" + msg_type];\n",
" } catch (e) {\n",
" console.log(\"No handler for the '\" + msg_type + \"' message type: \", msg);\n",
" return;\n",
" }\n",
"\n",
" if (callback) {\n",
" try {\n",
" // console.log(\"Handling '\" + msg_type + \"' message: \", msg);\n",
" callback(fig, msg);\n",
" } catch (e) {\n",
" console.log(\"Exception inside the 'handler_\" + msg_type + \"' callback:\", e, e.stack, msg);\n",
" }\n",
" }\n",
" };\n",
"}\n",
"\n",
"// from http://stackoverflow.com/questions/1114465/getting-mouse-location-in-canvas\n",
"mpl.findpos = function(e) {\n",
" //this section is from http://www.quirksmode.org/js/events_properties.html\n",
" var targ;\n",
" if (!e)\n",
" e = window.event;\n",
" if (e.target)\n",
" targ = e.target;\n",
" else if (e.srcElement)\n",
" targ = e.srcElement;\n",
" if (targ.nodeType == 3) // defeat Safari bug\n",
" targ = targ.parentNode;\n",
"\n",
" // jQuery normalizes the pageX and pageY\n",
" // pageX,Y are the mouse positions relative to the document\n",
" // offset() returns the position of the element relative to the document\n",
" var x = e.pageX - $(targ).offset().left;\n",
" var y = e.pageY - $(targ).offset().top;\n",
"\n",
" return {\"x\": x, \"y\": y};\n",
"};\n",
"\n",
"/*\n",
" * return a copy of an object with only non-object keys\n",
" * we need this to avoid circular references\n",
" * http://stackoverflow.com/a/24161582/3208463\n",
" */\n",
"function simpleKeys (original) {\n",
" return Object.keys(original).reduce(function (obj, key) {\n",
" if (typeof original[key] !== 'object')\n",
" obj[key] = original[key]\n",
" return obj;\n",
" }, {});\n",
"}\n",
"\n",
"mpl.figure.prototype.mouse_event = function(event, name) {\n",
" var canvas_pos = mpl.findpos(event)\n",
"\n",
" if (name === 'button_press')\n",
" {\n",
" this.canvas.focus();\n",
" this.canvas_div.focus();\n",
" }\n",
"\n",
" var x = canvas_pos.x;\n",
" var y = canvas_pos.y;\n",
"\n",
" this.send_message(name, {x: x, y: y, button: event.button,\n",
" step: event.step,\n",
" guiEvent: simpleKeys(event)});\n",
"\n",
" /* This prevents the web browser from automatically changing to\n",
" * the text insertion cursor when the button is pressed. We want\n",
" * to control all of the cursor setting manually through the\n",
" * 'cursor' event from matplotlib */\n",
" event.preventDefault();\n",
" return false;\n",
"}\n",
"\n",
"mpl.figure.prototype._key_event_extra = function(event, name) {\n",
" // Handle any extra behaviour associated with a key event\n",
"}\n",
"\n",
"mpl.figure.prototype.key_event = function(event, name) {\n",
"\n",
" // Prevent repeat events\n",
" if (name == 'key_press')\n",
" {\n",
" if (event.which === this._key)\n",
" return;\n",
" else\n",
" this._key = event.which;\n",
" }\n",
" if (name == 'key_release')\n",
" this._key = null;\n",
"\n",
" var value = '';\n",
" if (event.ctrlKey && event.which != 17)\n",
" value += \"ctrl+\";\n",
" if (event.altKey && event.which != 18)\n",
" value += \"alt+\";\n",
" if (event.shiftKey && event.which != 16)\n",
" value += \"shift+\";\n",
"\n",
" value += 'k';\n",
" value += event.which.toString();\n",
"\n",
" this._key_event_extra(event, name);\n",
"\n",
" this.send_message(name, {key: value,\n",
" guiEvent: simpleKeys(event)});\n",
" return false;\n",
"}\n",
"\n",
"mpl.figure.prototype.toolbar_button_onclick = function(name) {\n",
" if (name == 'download') {\n",
" this.handle_save(this, null);\n",
" } else {\n",
" this.send_message(\"toolbar_button\", {name: name});\n",
" }\n",
"};\n",
"\n",
"mpl.figure.prototype.toolbar_button_onmouseover = function(tooltip) {\n",
" this.message.textContent = tooltip;\n",
"};\n",
"mpl.toolbar_items = [[\"Home\", \"Reset original view\", \"fa fa-home icon-home\", \"home\"], [\"Back\", \"Back to previous view\", \"fa fa-arrow-left icon-arrow-left\", \"back\"], [\"Forward\", \"Forward to next view\", \"fa fa-arrow-right icon-arrow-right\", \"forward\"], [\"\", \"\", \"\", \"\"], [\"Pan\", \"Pan axes with left mouse, zoom with right\", \"fa fa-arrows icon-move\", \"pan\"], [\"Zoom\", \"Zoom to rectangle\", \"fa fa-square-o icon-check-empty\", \"zoom\"], [\"\", \"\", \"\", \"\"], [\"Download\", \"Download plot\", \"fa fa-floppy-o icon-save\", \"download\"]];\n",
"\n",
"mpl.extensions = [\"eps\", \"jpeg\", \"pdf\", \"png\", \"ps\", \"raw\", \"svg\", \"tif\"];\n",
"\n",
"mpl.default_extension = \"png\";var comm_websocket_adapter = function(comm) {\n",
" // Create a \"websocket\"-like object which calls the given IPython comm\n",
" // object with the appropriate methods. Currently this is a non binary\n",
" // socket, so there is still some room for performance tuning.\n",
" var ws = {};\n",
"\n",
" ws.close = function() {\n",
" comm.close()\n",
" };\n",
" ws.send = function(m) {\n",
" //console.log('sending', m);\n",
" comm.send(m);\n",
" };\n",
" // Register the callback with on_msg.\n",
" comm.on_msg(function(msg) {\n",
" //console.log('receiving', msg['content']['data'], msg);\n",
" // Pass the mpl event to the overriden (by mpl) onmessage function.\n",
" ws.onmessage(msg['content']['data'])\n",
" });\n",
" return ws;\n",
"}\n",
"\n",
"mpl.mpl_figure_comm = function(comm, msg) {\n",
" // This is the function which gets called when the mpl process\n",
" // starts-up an IPython Comm through the \"matplotlib\" channel.\n",
"\n",
" var id = msg.content.data.id;\n",
" // Get hold of the div created by the display call when the Comm\n",
" // socket was opened in Python.\n",
" var element = $(\"#\" + id);\n",
" var ws_proxy = comm_websocket_adapter(comm)\n",
"\n",
" function ondownload(figure, format) {\n",
" window.open(figure.imageObj.src);\n",
" }\n",
"\n",
" var fig = new mpl.figure(id, ws_proxy,\n",
" ondownload,\n",
" element.get(0));\n",
"\n",
" // Call onopen now - mpl needs it, as it is assuming we've passed it a real\n",
" // web socket which is closed, not our websocket->open comm proxy.\n",
" ws_proxy.onopen();\n",
"\n",
" fig.parent_element = element.get(0);\n",
" fig.cell_info = mpl.find_output_cell(\"\");\n",
" if (!fig.cell_info) {\n",
" console.error(\"Failed to find cell for figure\", id, fig);\n",
" return;\n",
" }\n",
"\n",
" var output_index = fig.cell_info[2]\n",
" var cell = fig.cell_info[0];\n",
"\n",
"};\n",
"\n",
"mpl.figure.prototype.handle_close = function(fig, msg) {\n",
" fig.root.unbind('remove')\n",
"\n",
" // Update the output cell to use the data from the current canvas.\n",
" fig.push_to_output();\n",
" var dataURL = fig.canvas.toDataURL();\n",
" // Re-enable the keyboard manager in IPython - without this line, in FF,\n",
" // the notebook keyboard shortcuts fail.\n",
" IPython.keyboard_manager.enable()\n",
" $(fig.parent_element).html('');\n",
" fig.close_ws(fig, msg);\n",
"}\n",
"\n",
"mpl.figure.prototype.close_ws = function(fig, msg){\n",
" fig.send_message('closing', msg);\n",
" // fig.ws.close()\n",
"}\n",
"\n",
"mpl.figure.prototype.push_to_output = function(remove_interactive) {\n",
" // Turn the data on the canvas into data in the output cell.\n",
" var dataURL = this.canvas.toDataURL();\n",
" this.cell_info[1]['text/html'] = '';\n",
"}\n",
"\n",
"mpl.figure.prototype.updated_canvas_event = function() {\n",
" // Tell IPython that the notebook contents must change.\n",
" IPython.notebook.set_dirty(true);\n",
" this.send_message(\"ack\", {});\n",
" var fig = this;\n",
" // Wait a second, then push the new image to the DOM so\n",
" // that it is saved nicely (might be nice to debounce this).\n",
" setTimeout(function () { fig.push_to_output() }, 1000);\n",
"}\n",
"\n",
"mpl.figure.prototype._init_toolbar = function() {\n",
" var fig = this;\n",
"\n",
" var nav_element = $('')\n",
" nav_element.attr('style', 'width: 100%');\n",
" this.root.append(nav_element);\n",
"\n",
" // Define a callback function for later on.\n",
" function toolbar_event(event) {\n",
" return fig.toolbar_button_onclick(event['data']);\n",
" }\n",
" function toolbar_mouse_event(event) {\n",
" return fig.toolbar_button_onmouseover(event['data']);\n",
" }\n",
"\n",
" for(var toolbar_ind in mpl.toolbar_items){\n",
" var name = mpl.toolbar_items[toolbar_ind][0];\n",
" var tooltip = mpl.toolbar_items[toolbar_ind][1];\n",
" var image = mpl.toolbar_items[toolbar_ind][2];\n",
" var method_name = mpl.toolbar_items[toolbar_ind][3];\n",
"\n",
" if (!name) { continue; };\n",
"\n",
" var button = $('');\n",
" button.click(method_name, toolbar_event);\n",
" button.mouseover(tooltip, toolbar_mouse_event);\n",
" nav_element.append(button);\n",
" }\n",
"\n",
" // Add the status bar.\n",
" var status_bar = $('');\n",
" nav_element.append(status_bar);\n",
" this.message = status_bar[0];\n",
"\n",
" // Add the close button to the window.\n",
" var buttongrp = $('');\n",
" var button = $('');\n",
" button.click(function (evt) { fig.handle_close(fig, {}); } );\n",
" button.mouseover('Stop Interaction', toolbar_mouse_event);\n",
" buttongrp.append(button);\n",
" var titlebar = this.root.find($('.ui-dialog-titlebar'));\n",
" titlebar.prepend(buttongrp);\n",
"}\n",
"\n",
"mpl.figure.prototype._root_extra_style = function(el){\n",
" var fig = this\n",
" el.on(\"remove\", function(){\n",
"\tfig.close_ws(fig, {});\n",
" });\n",
"}\n",
"\n",
"mpl.figure.prototype._canvas_extra_style = function(el){\n",
" // this is important to make the div 'focusable\n",
" el.attr('tabindex', 0)\n",
" // reach out to IPython and tell the keyboard manager to turn it's self\n",
" // off when our div gets focus\n",
"\n",
" // location in version 3\n",
" if (IPython.notebook.keyboard_manager) {\n",
" IPython.notebook.keyboard_manager.register_events(el);\n",
" }\n",
" else {\n",
" // location in version 2\n",
" IPython.keyboard_manager.register_events(el);\n",
" }\n",
"\n",
"}\n",
"\n",
"mpl.figure.prototype._key_event_extra = function(event, name) {\n",
" var manager = IPython.notebook.keyboard_manager;\n",
" if (!manager)\n",
" manager = IPython.keyboard_manager;\n",
"\n",
" // Check for shift+enter\n",
" if (event.shiftKey && event.which == 13) {\n",
" this.canvas_div.blur();\n",
" // select the cell after this one\n",
" var index = IPython.notebook.find_cell_index(this.cell_info[0]);\n",
" IPython.notebook.select(index + 1);\n",
" }\n",
"}\n",
"\n",
"mpl.figure.prototype.handle_save = function(fig, msg) {\n",
" fig.ondownload(fig, null);\n",
"}\n",
"\n",
"\n",
"mpl.find_output_cell = function(html_output) {\n",
" // Return the cell and output element which can be found *uniquely* in the notebook.\n",
" // Note - this is a bit hacky, but it is done because the \"notebook_saving.Notebook\"\n",
" // IPython event is triggered only after the cells have been serialised, which for\n",
" // our purposes (turning an active figure into a static one), is too late.\n",
" var cells = IPython.notebook.get_cells();\n",
" var ncells = cells.length;\n",
" for (var i=0; i= 3 moved mimebundle to data attribute of output\n",
" data = data.data;\n",
" }\n",
" if (data['text/html'] == html_output) {\n",
" return [cell, data, j];\n",
" }\n",
" }\n",
" }\n",
" }\n",
"}\n",
"\n",
"// Register the function which deals with the matplotlib target/channel.\n",
"// The kernel may be null if the page has been refreshed.\n",
"if (IPython.notebook.kernel != null) {\n",
" IPython.notebook.kernel.comm_manager.register_target('matplotlib', mpl.mpl_figure_comm);\n",
"}\n"
],
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
""
],
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
""
]
},
"execution_count": 31,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"plt.xlabel('x [m]')\n",
"plt.ylabel('y [m]')\n",
"plt.title(\"Contours (%d in total) of the head in layer %d with inactive section\" % (nc, layer))\n",
"plt.contour(xm, ym, Phi[layer], nc)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"* The white area is the aquifer part that was defined as inactive (impervious).\n",
"* The red trough is the well location.\n",
"* The two flanks and the front sides are closed (no FQ and no fixed head).\n",
"* The head at the back side is prescribed and maintained at zero."
]
},
{
"cell_type": "markdown",
"metadata": {
"collapsed": true
},
"source": [
"## Conclusion"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"In this chapter we have developed, from scratch, a full 3D finite difference model for which fixed heads, fixed flows and inactive subareas are be prescribed. While developing we also computed a concrete example of which the head contours were finally shown.\n",
"\n",
"To turn this model into a general Python function that can compute arbitrary 3D steady-state models of this kind, we only have to gather the developed lines and put them under a function definition, add some help strings (a doc string) and add error checking for convenience, while the data for this specific case, like the grid dimensions (`x`,`y`,`z`) the conductivities (`kx`,`ky`,`kz`), the prescribed flows `FQ` and initial heads `HI` plus the boundary array `IBOUND` that tells which cells have fixed heads and which are inactive, are to be passed as in user-given arguments at the function call like so:\n",
"\n",
"Phi = fdm3(x,y,z,kx,ky,kz,FQ,IH,IBOUND) # function that solves arbitrary 3D steady steate finite difference model\n",
"\n",
"We do this in the next section."
]
}
],
"metadata": {
"anaconda-cloud": {},
"kernelspec": {
"display_name": "Python [conda root]",
"language": "python",
"name": "conda-root-py"
},
"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.5.2"
}
},
"nbformat": 4,
"nbformat_minor": 1
}