Contents

Debugging

Aims:

  • Understand the importance of debugging

  • Learn some debugging techniques

  • Learn the basics of pdb

Software prerequisites

Contact details:

  • Dr. Rhodri Nelson

  • Room 4.85 RSM building

  • email: rnelson@ic.ac.uk

  • Teams: @Nelson, Rhodri B in #ACSE1 or #General, or DM me.

Notes

  • Much of the material presented in today’s lecture has been adapted from SciPy’s debugging lecture found here.

  • Scripts for todays lecture can be found in the files folder within Lecture10.

When working on technical software development projects it’s often the case that you’ll spend as much (if not more) time debugging, compared to writing, ‘actual’ code. This statement alone should make one thing abundantly clear, debugging is important and in order to increase your productivity you need to get good at it.

Over the course of time as you develop/encounter code with bugs, fix them and repeat this process countless times these are skills you will be ‘forced’ to acquire naturally. But as with everything, when starting out it doesn’t hurt to get a few tips from sources who have spent too much time searching for that elusive typo in the middle of a million lines of code. Thus, this lecture aims to provide you with a few tips and tricks to help you on your journey to becoming a debugging guru.

Learn what works in what language

This lecture focuses on debugging in Python. While the general skill set is transferable to any coding language, some of the tools are not. When working in an unfamiliar language it is therefore important to quickly find the right tools to suit your debugging style - for example we’ll first look at a couple of manual debugging strategies utilizing print statements and IPython, how would we transfer similar techniques to c++?

Avoiding bugs

  • We all, even the great Fabio Luporini, write buggy code.

  • Write your code with testing and debugging in mind.

  • Keep It Simple, Stupid (KISS):

    • What is the simplest thing that could possibly work?

  • Don’t Repeat Yourself (DRY):

    • Every piece of knowledge must have a single, unambiguous, authoritative representation within a system.

    • Constants, algorithms, etc.

  • Try to limit interdependencies of your code. (Loose Coupling).

  • Give your variables, functions and modules meaningful names (not cryptic/mathematics names).

There have been some infamous ‘bugs’ throughout history, here are a few:

The world’s first computer bug?

On September 9th, 1947 at 3:45 p.m., Grace Murray Hopper logged in her log book the first computer bug! Her log read: “First actual case of bug being found”. The term “bug” in computer science is now, of course, not taken literally. It is used to describe a flaw or failure in a computer program that causes it to produce an unexpected result or crash.

The Y2K bug

The year 2000 problem spawned fears of worldwide economic collapse, nuclear war, and led to an industry of consultants providing last-minute fixes.

Scary stuff, but clearly a bit over exaggerated. Nevertheless, the bug was real and billions of dollars were spend world wide fixing the bug. Indeed, not all systems were fixed in time and in Spain some parking meters failed, the French meteorological institute published the weather for 19100 and in Australia some ticket validation machines failed.

The Dhahran Missile

But on a more serious note, a bug in the software of an MIM-104 Patriot caused its system clock to drift by one third of a second over a period of one hundred hours - resulting in failure to locate and intercept an incoming Iraqi Al Hussein missile, which then struck Dharan barracks, Saudi Arabia (February 25, 1991), killing 28 Americans.

And here are some of the more ‘interesting’ bugs you may encounter while developing software:

The Heisenbug

From Wikipedia:

In computer programming jargon, a heisenbug is a software bug that seems to disappear or alter its behavior when one attempts to study it. The term is a pun on the name of Werner Heisenberg, the physicist who first asserted the observer effect of quantum mechanics, which states that the act of observing a system inevitably alters its state. In electronics the traditional term is probe effect, where attaching a test probe to a device changes its behavior.

Mandelbug

Mandelbugs do not change their properties or vanish like Heisenbugs. Instead, they are so unusual and complex that there is typically no practical solution to fix them. In fact, Mandelbugs may be so nondeterministic that many scholarly developers believe that fractal mathematics inventor Benoit Mandelbrot developed the Mandelbug to drive developers insane.

Bohrbug

A software bug which manifests reliably under a well-defined, but possibly unknown, set of conditions.

Schroedinbug

A software bug which manifests only when somebody debugging it finds out that it shouldn’t work at all.

pyflakes: fast static analysis

They are several static analysis tools in Python (most of which you’re already familiar with); to name a few:

The simplest of these tools is probably pyflakes:

  • Fast, simple

  • Detects syntax errors, missing imports, typos on names.

Another good recommendation is the flake8 tool which is a combination of pyflakes and pep8. Thus, in addition to the types of errors that pyflakes catches, flake8 detects violations of the recommendation in PEP8 style guide.

Integrating pyflakes (or flake8) in your editor or IDE is highly recommended, it does yield productivity gains.

Debugging workflow

If you do have a non trivial bug, this is when debugging strategies kick in. There is no silver bullet. Yet, strategies help:

For debugging a given problem, the favorable situation is when the problem is isolated in a small number of lines of code, outside framework or application code, with short modify-run-fail cycles

  • Make it fail reliably. Find a test case that makes the code fail every time.

  • Divide and Conquer. Once you have a failing test case, isolate the failing code.

    • Which module.

    • Which function.

    • Which line of code.

=> isolate a small reproducible failure: a test case/minimum failing example

  • Change one thing at a time and re-run the failing test case.

  • Use the debugger to understand what is going wrong.

  • Take notes and be patient. It may take a while.

Note: Once you have gone through this process: isolated a tight piece of code reproducing the bug and fix the bug using this piece of code, add the corresponding code to your test suite.

Minimum failing examples

Given the occurrence of a particular bug, a minimum failing example (MFE) is the simplest piece of code that can reproduce an instance of the particular bug in question (which was probably discovered through building a larger more complex script).

When approaching someone for help debugging a problem, it is considered good ‘manners’ to provide them with an MFE. This will often drastically reduce the time taken for them to debug the problem. Further, when contributing to an open source project and reporting a problem (e.g. opening a bug report issue on Github), the developers will frequently ask ‘Please could you provide an MFE’.

Example:

This is a ‘real’ example. (Don’t worry about the details, this is just to illustrate and idea). I was once sent this long bit of code:

import numpy                   as np
import matplotlib.pyplot       as plot
import math                    as mt
import matplotlib.ticker       as mticker    
from   mpl_toolkits.axes_grid1 import make_axes_locatable
from   matplotlib              import cm
from   scipy.interpolate       import CubicSpline
from   scipy.interpolate       import interp1d
from   devito import *

from   examples.seismic        import TimeAxis
from   examples.seismic        import RickerSource
from   examples.seismic        import Receiver

nptx   = 861
nptz   = 311
x0     = 13650.
x1     = 56650. 
compx  = x1-x0
z0     = 0.
z1     = 9920.
compz  = z1-z0;
hxv    = (x1-x0)/(nptx-1)
hzv    = (z1-z0)/(nptz-1)

npmlx  = 30
npmlz  = 31

lx = npmlx*hxv
lz = npmlz*hzv

origin  = (x0,z0)
extent  = (compx,compz)
shape   = (nptx,nptz)
spacing = (hxv,hzv)

class d0domain(SubDomain):
    name = 'd0'
    def define(self, dimensions):
        x, z = dimensions
        return {x: ('middle', npmlx, npmlx), z: ('middle', 0, npmlz)}
d0_domain = d0domain()

class d1domain(SubDomain):
    name = 'd1'
    def define(self, dimensions):
        x, z = dimensions
        return {x: ('left',npmlx), z: z}
d1_domain = d1domain()

class d2domain(SubDomain):
    name = 'd2'
    def define(self, dimensions):
        x, z = dimensions
        return {x: ('right',npmlx), z: z}
d2_domain = d2domain()

class d3domain(SubDomain):
    name = 'd3'
    def define(self, dimensions):
        x, z = dimensions
        return {x: ('middle', npmlx, npmlx), z: ('right',nptz)}
d3_domain = d3domain()

grid = Grid(origin=origin,extent=extent,shape=shape,subdomains=(d0_domain,d1_domain,d2_domain,d3_domain),dtype=np.float64) 

nptxvel = 1407
nptzvel = 311
x0vel   = 0.        
x1vel   = 70300.     
z0vel   = 0.        
z1vel   = 9920.
hxvel   = (x1vel-x0vel)/(nptxvel-1)
hzvel   = (z1vel-z0vel)/(nptzvel-1)
Xvel    = np.linspace(x0vel,x1vel,nptxvel)
Zvel    = np.linspace(z0vel,z1vel,nptzvel)
velarq = np.fromfile('velocity_data/profilevel.bin',dtype='float32')
velarq = velarq.reshape((nptxvel,nptzvel))
fscale = 10**(-3) 
velarq = fscale*velarq

X0 = np.linspace(x0,x1,nptx)   
Z0 = np.linspace(z0,z1,nptz)
X1 = np.linspace((x0+0.5*hxv),(x1-0.5*hxv),nptx-1)
Z1 = np.linspace((z0+0.5*hzv),(z1-0.5*hzv),nptz-1)
v0 = np.zeros((nptx,nptz))
v1 = np.zeros((nptx-1,nptz-1))

C0x = np.zeros((nptx,nptzvel))
for j in range(nptzvel):
    x = Xvel
    z = velarq[0:nptxvel,j]
    #cs = CubicSpline(x,z)    
    cs = interp1d(x,z,kind='linear')
    xs = X0
    C0x[0:nptx,j] = cs(xs)
    
for i in range(nptx):
    x = Zvel
    z = C0x[i,0:nptzvel]
    #cs = CubicSpline(x,z)
    cs = interp1d(x,z,kind='linear')
    xs = Z0
    v0[i,0:nptz] = cs(xs)
          
C11x = np.zeros((nptx-1,nptzvel))
for j in range(nptzvel):
    x = Xvel
    z = velarq[0:nptxvel,j]
    cs = CubicSpline(x,z)
    xs = X1
    C11x[0:nptx-1,j] = cs(xs)
    
for i in range(nptx-1):
    x  = Zvel
    z  = C11x[i,0:nptzvel]
    #cs = CubicSpline(x,z)
    cs = interp1d(x,z,kind='linear')
    xs = Z1
    v1[i,0:nptz-1] = cs(xs)

t0 = 0.
tn = 20000.   
CFL = 0.4
vmax  = np.amax(v0) 
dtmax = np.float64((min(hxv,hzv)*CFL)/(vmax))
ntmax = int((tn-t0)/dtmax)+1
dt0   = np.float64((tn-t0)/ntmax)

time_range = TimeAxis(start=t0,stop=tn,num=ntmax+1)
nt         = time_range.num - 1

(hx,hz) = grid.spacing_map  
(x, z)  = grid.dimensions     
t       = grid.stepping_dim
dt      = grid.stepping_dim.spacing

f0     = 0.005
nfonte = 1
xposf  = 35150.
zposf  = 32.

src = RickerSource(name='src',grid=grid,f0=f0,npoint=nfonte,time_range=time_range,staggered=NODE,dtype=np.float64)
src.coordinates.data[:, 0] = xposf
src.coordinates.data[:, 1] = zposf

nrec   = nptx
nxpos  = np.linspace(x0,x1,nrec)
nzpos  = 32.

rec = Receiver(name='rec',grid=grid,npoint=nrec,time_range=time_range,staggered=NODE,dtype=np.float64)
rec.coordinates.data[:, 0] = nxpos
rec.coordinates.data[:, 1] = nzpos

u = TimeFunction(name="u",grid=grid,time_order=2,space_order=2,staggered=NODE,dtype=np.float64)

phi1 = TimeFunction(name="phi1",grid=grid,time_order=2,space_order=2,staggered=(x,z),dtype=np.float64)
phi2 = TimeFunction(name="phi2",grid=grid,time_order=2,space_order=2,staggered=(x,z),dtype=np.float64)

vel0 = Function(name="vel0",grid=grid,space_order=2,staggered=NODE,dtype=np.float64)
vel0.data[:,:] = v0[:,:]

vel1 = Function(name="vel1", grid=grid,space_order=2,staggered=(x,z),dtype=np.float64)
vel1.data[0:nptx-1,0:nptz-1] = v1

vel1.data[nptx-1,0:nptz-1] = vel1.data[nptx-2,0:nptz-1]
vel1.data[0:nptx,nptz-1]   = vel1.data[0:nptx,nptz-2]

src_term = src.inject(field=u.forward,expr=src*dt**2*vel0**2)

rec_term = rec.interpolate(expr=u)

x0pml  = x0 + npmlx*hxv 
x1pml  = x1 - npmlx*hxv 
z0pml  = z0            
z1pml  = z1 - npmlz*hzv 

def fdamp(x,z,i):
    
    quibar  = 0.05
          
    if(i==1):
        a = np.where(x<=x0pml,(np.abs(x-x0pml)/lx),np.where(x>=x1pml,(np.abs(x-x1pml)/lx),0.))
        fdamp = quibar*(a-(1./(2.*np.pi))*np.sin(2.*np.pi*a))
    if(i==2):
        a = np.where(z<=z0pml,(np.abs(z-z0pml)/lz),np.where(z>=z1pml,(np.abs(z-z1pml)/lz),0.))
        fdamp = quibar*(a-(1./(2.*np.pi))*np.sin(2.*np.pi*a))
      
    return fdamp

def generatemdamp():
    
    X0     = np.linspace(x0,x1,nptx)    
    Z0     = np.linspace(z0,z1,nptz)
    X0grid,Z0grid = np.meshgrid(X0,Z0)
    X1   = np.linspace((x0+0.5*hxv),(x1-0.5*hxv),nptx-1)
    Z1   = np.linspace((z0+0.5*hzv),(z1-0.5*hzv),nptz-1)
    X1grid,Z1grid = np.meshgrid(X1,Z1)
   
    D01 = np.zeros((nptx,nptz))
    D02 = np.zeros((nptx,nptz))
    D11 = np.zeros((nptx,nptz))
    D12 = np.zeros((nptx,nptz))
    
    D01 = np.transpose(fdamp(X0grid,Z0grid,1))
    D02 = np.transpose(fdamp(X0grid,Z0grid,2))
  
    D11 = np.transpose(fdamp(X1grid,Z1grid,1))
    D12 = np.transpose(fdamp(X1grid,Z1grid,2))
    
    return D01, D02, D11, D12

D01, D02, D11, D12 = generatemdamp();

dampx0 = Function(name="dampx0", grid=grid,space_order=2,staggered=NODE ,dtype=np.float64)
dampz0 = Function(name="dampz0", grid=grid,space_order=2,staggered=NODE ,dtype=np.float64)
dampx0.data[:,:] = D01
dampz0.data[:,:] = D02

dampx1 = Function(name="dampx1", grid=grid,space_order=2,staggered=(x,z),dtype=np.float64)
dampz1 = Function(name="dampz1", grid=grid,space_order=2,staggered=(x,z),dtype=np.float64)
dampx1.data[0:nptx-1,0:nptz-1] = D11
dampz1.data[0:nptx-1,0:nptz-1] = D12

dampx1.data[nptx-1,0:nptz-1]   = dampx1.data[nptx-2,0:nptz-1]
dampx1.data[0:nptx,nptz-1]     = dampx1.data[0:nptx,nptz-2]
dampz1.data[nptx-1,0:nptz-1]   = dampz1.data[nptx-2,0:nptz-1]
dampz1.data[0:nptx,nptz-1]     = dampz1.data[0:nptx,nptz-2]

# White Region
pde01   = Eq(u.dt2-u.laplace*vel0*vel0) 

# Blue Region
pde02a  = u.dt2   + (dampx0+dampz0)*u.dtc + (dampx0*dampz0)*u - u.laplace*vel0*vel0 
pde02b  = - (0.5/hx)*(phi1[t,x,z-1]+phi1[t,x,z]-phi1[t,x-1,z-1]-phi1[t,x-1,z])
pde02c  = - (0.5/hz)*(phi2[t,x-1,z]+phi2[t,x,z]-phi2[t,x-1,z-1]-phi2[t,x,z-1])
pde02   = Eq(pde02a + pde02b + pde02c)

pde10 = phi1.dt + dampx1*0.5*(phi1.forward+phi1)
a1    = u[t+1,x+1,z] + u[t+1,x+1,z+1] - u[t+1,x,z] - u[t+1,x,z+1] 
a2    = u[t,x+1,z]   + u[t,x+1,z+1]   - u[t,x,z]   - u[t,x,z+1] 
pde11 = -(dampz1-dampx1)*0.5*(0.5/hx)*(a1+a2)*vel1*vel1
pde1  = Eq(pde10+pde11)
                                                    
pde20 = phi2.dt + dampz1*0.5*(phi2.forward+phi2) 
b1    = u[t+1,x,z+1] + u[t+1,x+1,z+1] - u[t+1,x,z] - u[t+1,x+1,z] 
b2    = u[t,x,z+1]   + u[t,x+1,z+1]   - u[t,x,z]   - u[t,x+1,z] 
pde21 = -(dampx1-dampz1)*0.5*(0.5/hz)*(b1+b2)*vel1*vel1
pde2  = Eq(pde20+pde21)

stencil01 =  Eq(u.forward,solve(pde01,u.forward) ,subdomain = grid.subdomains['d0'])

subds = ['d1','d2','d3']

stencil02 = [Eq(u.forward,solve(pde02, u.forward),subdomain = grid.subdomains[subds[i]]) for i in range(0,len(subds))]
stencil1 = [Eq(phi1.forward, solve(pde1,phi1.forward),subdomain = grid.subdomains[subds[i]]) for i in range(0,len(subds))]
stencil2 = [Eq(phi2.forward, solve(pde2,phi2.forward),subdomain = grid.subdomains[subds[i]]) for i in range(0,len(subds))]

bc  = [Eq(u[t+1,0,z],0.),Eq(u[t+1,nptx-1,z],0.),Eq(u[t+1,x,nptz-1],0.),Eq(u[t+1,x,0],u[t+1,x,1])]

op = Operator([stencil01,stencil02] + src_term + bc + [stencil1,stencil2] + rec_term,subs=grid.spacing_map)

u.data[:]     = 0.
phi1.data[:]  = 0.
phi2.data[:]  = 0.

op(time=nt,dt=dt0)

def graph2d(U):    
    plot.figure()
    fscale =  10**(-3)
    scale = np.amax(U[npmlx:-npmlx,0:-npmlz])/50.
    extent = [fscale*x0pml,fscale*x1pml,fscale*z1pml,fscale*z0pml]
    fig = plot.imshow(np.transpose(U[npmlx:-npmlx,0:-npmlz]),vmin=-scale, vmax=scale, cmap=cm.gray, extent=extent)
    plot.gca().xaxis.set_major_formatter(mticker.FormatStrFormatter('%d km'))
    plot.gca().yaxis.set_major_formatter(mticker.FormatStrFormatter('%d km'))
    plot.axis('equal')
    plot.title('Map - Acoustic Problem with PML - Devito')
    plot.grid()
    ax = plot.gca()
    divider = make_axes_locatable(ax)
    cax = divider.append_axes("right", size="5%", pad=0.05)
    cbar = plot.colorbar(fig, cax=cax, format='%.2e')
    cbar.set_label('Displacement [km]')
    plot.draw()
    plot.savefig('acoustic_map_pml_devito.png',dpi=100)
    plot.show()
    
graph2d(u.data[0,:,:])

def graph2drec(rec):    
        fig1 = plot.figure()
        fscale =  10**(-3)
        scale = np.amax(rec[:,npmlx:-npmlx])/50.
        extent = [fscale*x0pml,fscale*x1pml, fscale*tn, fscale*t0]
        fig = plot.imshow(rec[:,npmlx:-npmlx], vmin=-scale, vmax=scale, cmap=cm.gray, extent=extent)
        plot.gca().xaxis.set_major_formatter(mticker.FormatStrFormatter('%d km'))
        plot.gca().yaxis.set_major_formatter(mticker.FormatStrFormatter('%d s'))
        plot.axis('equal')
        plot.title('Receivers Signal Profile with PML - Devito')
        ax = plot.gca()
        divider = make_axes_locatable(ax)
        plot.savefig('rec_map_pml_devito.png')
        plot.show()

graph2drec(rec.data)

which was now producing an incorrect result, but a couple of software versions ago the result was correct.

The above piece of code is long and the offender has not been pin-pointed, thus making it difficult to deduce where the bug is coming from. After some work, I managed to produce an MFE:

import numpy as np
from devito import *
class d0domain(SubDomain):
    name = 'd0'
    def define(self, dimensions):
        x, y = dimensions
        return {x: ('middle', 2, 2), y: ('right', 10)}
sdom = d0domain()
grid = Grid(shape=(10, 10), extent=(10, 10), subdomains=(sdom))
u = Function(name='u', grid=grid)
eq = Eq(u, u+1, subdomain = grid.subdomains['d0'])
op = Operator(eq)
op.apply()

This piece of code is short and simple an allowed the bug to be tracked down and fixed in less than an hour!

Manual debugging

Here are a couple of ‘manual’ debugging techniques:

  • Print statements Adding print statements to your code can often be a quick and effective method of pin-pointing problematic lines of code, e.g.,

...
x = some_function(inputs)
print(x.proprty1, x.property2)
...

  • from IPython import embed; embed() When possible, utilizing this command can be an incredibly useful tool for debugging sections of code:

...
x = some_function(inputs)
from IPython import embed; embed()
...

I’ll now show you a live example!

Using the Python debugger

The python debugger, pdb: https://docs.python.org/library/pdb.html, allows you to inspect your code interactively.

Specifically it allows you to:

  • View the source code.

  • Walk up and down the call stack.

  • Inspect values of variables.

  • Modify values of variables.

  • Set breakpoints.

Invoking the debugger

Ways to launch the debugger:

  1. Postmortem, launch debugger after module errors.

  2. Launch the module with the debugger.

  3. Call the debugger inside the module

Postmortem

Example 1: You’re working in IPython and you get a traceback.

Here we debug the file index_error.py. When running it, an IndexError is raised. Type %debug and drop into the debugger.

In [1]: %run index_error.py

---------------------------------------------------------------------------

IndexError                                Traceback (most recent call last)

/home/varoquau/dev/scipy-lecture-notes/advanced/debugging/index_error.py in <module>()

      6

      7 if __name__ == '__main__':

----> 8     index_error()

      9



/home/varoquau/dev/scipy-lecture-notes/advanced/debugging/index_error.py in index_error()

      3 def index_error():

      4     lst = list('foobar')

----> 5     print lst[len(lst)]

      6

      7 if __name__ == '__main__':



IndexError: list index out of range



In [2]: %debug

> /home/varoquau/dev/scipy-lecture-notes/advanced/debugging/index_error.py(5)index_error()

      4     lst = list('foobar')

----> 5     print lst[len(lst)]

      6



ipdb> list

      1 """Small snippet to raise an IndexError."""

      2

      3 def index_error():

      4     lst = list('foobar')

----> 5     print lst[len(lst)]

      6

      7 if __name__ == '__main__':

      8     index_error()

      9



ipdb> len(lst)

6

ipdb> print(lst[len(lst)-1])

r

ipdb> quit



In [3]:

Example 2: Post-mortem debugging without IPython

In some situations you cannot use IPython, for instance to debug a script that wants to be called from the command line. In this case, you can call the script with python -m pdb script.py:

$ python -m pdb index_error.py
> /home/varoquau/dev/scipy-lecture-notes/advanced/optimizing/index_error.py(1)<module>()
-> """Small snippet to raise an IndexError."""
(Pdb) continue
Traceback (most recent call last):
File "/usr/lib/python2.6/pdb.py", line 1296, in main
    pdb._runscript(mainpyfile)
File "/usr/lib/python2.6/pdb.py", line 1215, in _runscript
    self.run(statement)
File "/usr/lib/python2.6/bdb.py", line 372, in run
    exec cmd in globals, locals
File "<string>", line 1, in <module>
File "index_error.py", line 8, in <module>
    index_error()
File "index_error.py", line 5, in index_error
    print lst[len(lst)]
IndexError: list index out of range
Uncaught exception. Entering post mortem debugging
Running 'cont' or 'step' will restart the program
> /home/varoquau/dev/scipy-lecture-notes/advanced/optimizing/index_error.py(5)index_error()
-> print(lst[len(lst)])
(Pdb)

Step-by-step execution

Example: You believe a bug exists in a module but are not sure where.

For instance we are trying to debug wiener_filtering.py. Indeed the code runs, but the filtering does not work well.

  • Run the script in IPython with the debugger using %run -d wiener_filtering.py :

In [1]: %run -d wiener_filtering.py

*** Blank or comment

*** Blank or comment

*** Blank or comment

Breakpoint 1 at /home/varoquau/dev/scipy-lecture-notes/advanced/optimizing/wiener_filtering.py:4

NOTE: Enter 'c' at the ipdb>  prompt to start your script.

> <string>(1)<module>()

  • Set a break point at line 34 using b 34:

ipdb> n

> /home/varoquau/dev/scipy-lecture-notes/advanced/optimizing/wiener_filtering.py(4)<module>()

      3

1---> 4 import numpy as np

      5 import scipy as sp



ipdb> b 34

Breakpoint 2 at /home/varoquau/dev/scipy-lecture-notes/advanced/optimizing/wiener_filtering.py:34

  • Continue execution to next breakpoint with c(ont(inue)):

ipdb> c

> /home/varoquau/dev/scipy-lecture-notes/advanced/optimizing/wiener_filtering.py(34)iterated_wiener()

     33     """

2--> 34     noisy_img = noisy_img

     35     denoised_img = local_mean(noisy_img, size=size)

  • Step into code with n(ext) and s(tep): next jumps to the next statement in the current execution context, while step will go across execution contexts, i.e. enable exploring inside function calls:

ipdb> s

> /home/varoquau/dev/scipy-lecture-notes/advanced/optimizing/wiener_filtering.py(35)iterated_wiener()

2    34     noisy_img = noisy_img

---> 35     denoised_img = local_mean(noisy_img, size=size)

     36     l_var = local_var(noisy_img, size=size)



ipdb> n

> /home/varoquau/dev/scipy-lecture-notes/advanced/optimizing/wiener_filtering.py(36)iterated_wiener()

     35     denoised_img = local_mean(noisy_img, size=size)

---> 36     l_var = local_var(noisy_img, size=size)

     37     for i in range(3):

  • Step a few lines and explore the local variables:

ipdb> n

> /home/varoquau/dev/scipy-lecture-notes/advanced/optimizing/wiener_filtering.py(37)iterated_wiener()

     36     l_var = local_var(noisy_img, size=size)

---> 37     for i in range(3):

     38         res = noisy_img - denoised_img

ipdb> print(l_var)

[[5868 5379 5316 ..., 5071 4799 5149]

 [5013  363  437 ...,  346  262 4355]

 [5379  410  344 ...,  392  604 3377]

 ...,

 [ 435  362  308 ...,  275  198 1632]

 [ 548  392  290 ...,  248  263 1653]

 [ 466  789  736 ..., 1835 1725 1940]]

ipdb> print(l_var.min())

0

Oh dear, nothing but integers, and 0 variation. Here is our bug, we are doing integer arithmetic.

Other ways of starting a debugger

Raising an exception as a poor man’s break point

If you find it tedious to note the line number to set a break point, you can simply raise an exception at the point that you want to inspect and use IPython’s %debug. Note that in this case you cannot step or continue the execution.

Debugging test failures using nosetests

You can run, e.g., nosetests --pdb my_test.py to drop in post-mortem debugging on exceptions, and nosetests --pdb-failure to inspect test failures using the debugger. (That that the default nose package does not work with parameterised tests).

In addition, you can use the IPython interface for the debugger in nose by installing the nose plugin ipdbplugin. You can than pass the --ipdb and --ipdb-failure options to nosetests.

Calling the debugger explicitly

Insert the following line where you want to drop in the debugger:

import pdb; pdb.set_trace()

Graphical debuggers and alternatives

  • pudb is a good semi-graphical debugger with a text user interface in the console.

  • The Visual Studio Code integrated development environment includes a debugging mode.

  • The Mu editor is a simple Python editor that includes a debugging mode.

Debugger commands and interaction

Command

Description

l(list)

Lists the code at the current position

u(p)

Walk up the call stack

d(own)

Walk down the call stack

n(ext)

Execute the next line (does not go down in new functions)

s(tep)

Execute the next statement (goes down in new functions)

bt

Print the call stack

a

Print the local variables

!command

Execute the given Python command (by opposition to pdb commands

Debugger commands are not Python code

You cannot name the variables the way you want. For instance, you cannot override the variables in the current frame with the same name: use different names than your local variable when typing code in the debugger.

Getting help when in the debugger

Type h or help to access the interactive help:

ipdb> help



Documented commands (type help <topic>):

========================================

EOF    bt         cont      enable  jump  pdef   r        tbreak   w

a      c          continue  exit    l     pdoc   restart  u        whatis

alias  cl         d         h       list  pinfo  return   unalias  where

args   clear      debug     help    n     pp     run      unt

b      commands   disable   ignore  next  q      s        until

break  condition  down      j       p     quit   step     up



Miscellaneous help topics:

==========================

exec  pdb



Undocumented commands:

======================

retval  rv

Debugging segmentation faults using gdb

Note:

The following example makes use of python3 compiled with debugging symbols. On Unix systems this can be installed via the package manager, e.g. sudo apt install python3-dbg python3-dev on Ubuntu systems. (I’m sorry, I did not get chance to check how this should be installed in Windows).

If you have a segmentation fault, you cannot debug it with pdb, as it crashes the Python interpreter before it can drop in the debugger. Similarly, if you have a bug in C code embedded in Python, pdb is useless. For this we turn to the gnu debugger, gdb, available on Linux, MacOS and Windows.

Before we start with gdb, let us add a few Python-specific tools to it. For this we add a few macros to our ~/.gdbinit. The optimal choice of macro depends on your Python version and your gdb version. A simplified version has been added in .gdbinit located in the files folder, but feel free to read DebuggingWithGdb.

NOTE: Some of the paths specified in .gdbinit are specific to my machine. You’ll need to modify these.

To debug with gdb the Python script segfault.py, we can run the script in gdb as follows

(gdb) run segfault.py
Starting program: /home/rhodri/debugging/bin/python segfault.py
[Thread debugging using libthread_db enabled]
Using host libthread_db library "/lib/x86_64-linux-gnu/libthread_db.so.1".
[New Thread 0x7ffff4732700 (LWP 43662)]
[New Thread 0x7ffff1f31700 (LWP 43663)]
[New Thread 0x7ffff1730700 (LWP 43664)]
[New Thread 0x7fffecf2f700 (LWP 43665)]
[New Thread 0x7fffec72e700 (LWP 43666)]
[New Thread 0x7fffe7f2d700 (LWP 43667)]
[New Thread 0x7fffe572c700 (LWP 43668)]
[New Thread 0x7fffe4f2b700 (LWP 43669)]
[New Thread 0x7fffe072a700 (LWP 43670)]
[New Thread 0x7fffddf29700 (LWP 43671)]
[New Thread 0x7fffdb728700 (LWP 43672)]
[New Thread 0x7fffd8f27700 (LWP 43673)]
[New Thread 0x7fffd6726700 (LWP 43674)]
[New Thread 0x7fffd3f25700 (LWP 43675)]
[New Thread 0x7fffd1724700 (LWP 43676)]

Thread 1 "python" received signal SIGSEGV, Segmentation fault.
0x00007ffff6e9a4c9 in _aligned_strided_to_contig_size8 () from /home/rhodri/debugging/lib/python3.8/site-packages/numpy/core/_multiarray_umath.cpython-38-x86_64-linux-gnu.so
(gdb)

Then, to get a back-trace we can utilize the py-bt command:

(gdb) py-bt
Traceback (most recent call first):
  <built-in method concatenate of module object at remote 0x7ffff744aef0>
  <built-in method implement_array_function of module object at remote 0x7ffff744aef0>
  File "<__array_function__ internals>", line 261, in concatenate
  File "/home/rhodri/debugging/lib/python3.8/site-packages/numpy/core/arrayprint.py", line 848, in _leading_trailing
    self.floatmode = floatmode
  File "/home/rhodri/debugging/lib/python3.8/site-packages/numpy/core/arrayprint.py", line 337, in _leading_trailing
    _leading_trailing(a, edgeitems, index + np.index_exp[ :edgeitems]),
  File "/home/rhodri/debugging/lib/python3.8/site-packages/numpy/core/arrayprint.py", line 489, in _array2string
    data = _leading_trailing(data, options['edgeitems'])
  File "/home/rhodri/debugging/lib/python3.8/site-packages/numpy/core/arrayprint.py", line 468, in wrapper
    return f(self, *args, **kwargs)
  File "/home/rhodri/debugging/lib/python3.8/site-packages/numpy/core/arrayprint.py", line 1204, in array2string
    
  File "/home/rhodri/debugging/lib/python3.8/site-packages/numpy/core/arrayprint.py", line 2002, in _array_str_implementation
  <built-in method print of module object at remote 0x7ffff761efb0>
  File "segfault.py", line 16, in print_big_array
    print(big_array[-10:])
  File "segfault.py", line 23, in <module>
    l.append(print_big_array(a))

And there we have it, our segmentation fault is being caused by line 23 in segfault.py, specifically the command l.append(print_big_array(a)).

Note

For a list of Python-specific commands defined in the .gdbinit, read the source of this file.

Debugging distributed (e.g. MPI) computations

Later in this course you will encounter the mpi4py package. Debugging MPI based programs can be a real pain and tools such as pdb are often less useful. Utilising manual debugging techniques in tandem with tools such as tmpi can however be an effective method.

Wrap up exercise

The following script is well documented and hopefully legible. It seeks to answer a problem of actual interest for numerical computing, but it does not work… Can you debug it?

Python source code of to_debug.py 

"""

A script to compare different root-finding algorithms.



This version of the script is buggy and does not execute. It is your task

to find an fix these bugs.



The output of the script sould look like:



    Benching 1D root-finder optimizers from scipy.optimize:

                brenth:   604678 total function calls

                brentq:   594454 total function calls

                ridder:   778394 total function calls

                bisect:  2148380 total function calls

"""

from itertools import product



import numpy as np

from scipy import optimize



FUNCTIONS = (np.tan,  # Dilating map

             np.tanh, # Contracting map

             lambda x: x**3 + 1e-4*x, # Almost null gradient at the root

             lambda x: x+np.sin(2*x), # Non monotonous function

             lambda x: 1.1*x+np.sin(4*x), # Fonction with several local maxima

            )



OPTIMIZERS = (optimize.brenth, optimize.brentq, optimize.ridder,

              optimize.bisect)





def apply_optimizer(optimizer, func, a, b):

    """ Return the number of function calls given an root-finding optimizer, 

        a function and upper and lower bounds.

    """

    return optimizer(func, a, b, full_output=True)[1].function_calls,





def bench_optimizer(optimizer, param_grid):

    """ Find roots for all the functions, and upper and lower bounds

        given and return the total number of function calls.

    """

    return sum(apply_optimizer(optimizer, func, a, b)

               for func, a, b in param_grid)





def compare_optimizers(optimizers):

    """ Compare all the optimizers given on a grid of a few different

        functions all admitting a signle root in zero and a upper and

        lower bounds.

    """

    random_a = -1.3 + np.random.random(size=100)

    random_b =   .3 + np.random.random(size=100)

    param_grid = product(FUNCTIONS, random_a, random_b)

    print("Benching 1D root-finder optimizers from scipy.optimize:")

    for optimizer in OPTIMIZERS:

        print('% 20s: % 8i total function calls' % (

                    optimizer.__name__, 

                    bench_optimizer(optimizer, param_grid)

                ))





if __name__ == '__main__':

    compare_optimizers(OPTIMIZERS)

What Every Computer Scientist Should Know About Floating-Point Arithmetic Code profiling and optimisation

By ESE, Imperial College London
© Copyright 2020.