diff --git a/util/qbox_angle.py b/util/qbox_angle.py
index fce817c..9e01e61 100755
--- a/util/qbox_angle.py
+++ b/util/qbox_angle.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # qbox_angle.py
 # extract angle defined by three atoms from Qbox output
 # use: qbox_angle.py name1 name2 name3 file.r
@@ -7,7 +7,7 @@
 import math
 
 if len(sys.argv) != 5:
-  print "use: ",sys.argv[0]," name1 name2 name3 file.r"
+  print("use: ",sys.argv[0]," name1 name2 name3 file.r")
   sys.exit()
 
 name1 = ""
@@ -77,7 +77,7 @@ def endElement(self, name):
       c = sp
       c = max(-1.0,min(1.0,sp))
       a = (180.0/math.pi)*math.acos(c)
-      print '%.4f' % a
+      print('%.4f' % a)
     elif name == "position":
       self.readPos1 = 0
       self.readPos2 = 0
diff --git a/util/qbox_distance.py b/util/qbox_distance.py
index fa023f7..35fadfa 100755
--- a/util/qbox_distance.py
+++ b/util/qbox_distance.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # qbox_distance.py
 # extract distance between two atoms from Qbox output
 # use: qbox_distance.py name1 name2 file.r
@@ -7,7 +7,7 @@
 import math
 
 if len(sys.argv) != 4:
-  print "use: ",sys.argv[0]," name1 name2 file.r"
+  print("use: ",sys.argv[0]," name1 name2 file.r")
   sys.exit()
 
 name1 = ""
@@ -54,9 +54,9 @@ def endElement(self, name):
       r2z = float(pos2[2])
       #print "r1: ",self.r1x,self.r1y,self.r1z
       #print "r2: ",self.r2x,self.r2y,self.r2z
-      print '%.4f' % math.sqrt((r1x-r2x)**2+
+      print('%.4f' % math.sqrt((r1x-r2x)**2+
                                (r1y-r2y)**2+
-                               (r1z-r2z)**2)
+                               (r1z-r2z)**2))
     elif name == "position":
       self.readPos1 = 0
       self.readPos2 = 0
diff --git a/util/qbox_dos.py b/util/qbox_dos.py
index 893c0d0..5a1c7b5 100755
--- a/util/qbox_dos.py
+++ b/util/qbox_dos.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # qbox_dos.py: extract electronic DOS from Qbox output
 # generate DOS plot in gnuplot format
 # use: qbox_dos.py [-last] [mu] emin emax [mu] width file.r
@@ -16,7 +16,7 @@
 
 argc=len(sys.argv)
 if (argc < 5) or (argc > 7) :
-  print "use: ",sys.argv[0]," [-last] [mu] emin emax width file.r"
+  print("use: ",sys.argv[0]," [-last] [mu] emin emax width file.r")
   sys.exit()
 
 # default chemical potential mu=0.0
@@ -29,7 +29,7 @@
   lastonly = True
   iarg += 1
   if (argc < 6):
-    print "use: ",sys.argv[0]," [-last] [mu] emin emax width file.r"
+    print("use: ",sys.argv[0]," [-last] [mu] emin emax width file.r")
     sys.exit()
 
 # check for mu argument
@@ -37,7 +37,7 @@
   mu = float(sys.argv[iarg])
   iarg += 1
   if (argc < 6):
-    print "use: ",sys.argv[0]," [-last] [mu] emin emax width file.r"
+    print("use: ",sys.argv[0]," [-last] [mu] emin emax width file.r")
     sys.exit()
 
 emin = float(sys.argv[iarg])
@@ -108,17 +108,17 @@ def accumulate_dos(self):
           self.dos_dn[j] += gauss(float(self.e[i])-self.mu-ej, width ) * self.weight
 
   def print_dos(self):
-    print "# ",infile," mu=",self.mu," spin=0 width=",width
+    print("# ",infile," mu=",self.mu," spin=0 width=",width)
     for j in range(ndos):
       ej = emin + j * de
-      print ej, self.dos_up[j]
+      print(ej, self.dos_up[j])
     if self.nspin == 2:
-      print
-      print
-      print "# ",infile," mu=",mu," spin=1 width=",width
+      print()
+      print()
+      print("# ",infile," mu=",mu," spin=1 width=",width)
       for j in range(ndos):
         ej = emin + j * de
-        print ej, self.dos_dn[j]
+        print(ej, self.dos_dn[j])
 
 parser = xml.sax.make_parser()
 handler = QboxOutputHandler()
diff --git a/util/qbox_eig.py b/util/qbox_eig.py
index ee27846..5d0bad2 100755
--- a/util/qbox_eig.py
+++ b/util/qbox_eig.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # qbox_eig.py: extract eigenvalues from Qbox output
 # use: qbox_eig.py kpoint n ispin file.r
 # extracts eigenvalue n at (ispin,kpoint)
@@ -10,8 +10,8 @@
 
 argc=len(sys.argv)
 if ( not ( argc in [3,4,6,7] ) ):
-  print "use: ",sys.argv[0]," [ispin] [kx ky kz] n file.r"
-  print " ispin = 0..1, n = 1..neig"
+  print("use: ",sys.argv[0]," [ispin] [kx ky kz] n file.r")
+  print(" ispin = 0..1, n = 1..neig")
   sys.exit()
 
 if argc == 7:
@@ -46,7 +46,7 @@
   n = int(sys.argv[1])
   infile = sys.argv[2]
 
-print "# ",infile," ispin=",ispin, " n=", n, " k=", kx, ky, kz
+print("# ",infile," ispin=",ispin, " n=", n, " k=", kx, ky, kz)
 
 # Qbox output handler to extract and process data
 class QboxOutputHandler(xml.sax.handler.ContentHandler):
@@ -81,9 +81,9 @@ def endElement(self, name):
   def print_eig(self):
     self.e = self.buffer.split()
     if n > int(self.n):
-      print "n>neig: neig=", self.n
+      print("n>neig: neig=", self.n)
     else:
-      print self.e[n-1]
+      print(self.e[n-1])
 
 parser = xml.sax.make_parser()
 handler = QboxOutputHandler()
diff --git a/util/qbox_force.py b/util/qbox_force.py
index 4ec5edd..de3d7f2 100755
--- a/util/qbox_force.py
+++ b/util/qbox_force.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # qbox_force.py
 # extract force of an atom from Qbox output
 # use: qbox_force.py atom_name file.r
@@ -7,7 +7,7 @@
 import math
 
 if len(sys.argv) != 3:
-  print "use: ",sys.argv[0]," atom_name file.r"
+  print("use: ",sys.argv[0]," atom_name file.r")
   sys.exit()
 
 # Qbox output handler to extract and process <atomset>
@@ -36,7 +36,7 @@ def endElement(self, name):
         fx = float(force[0])
         fy = float(force[1])
         fz = float(force[2])
-        print '%.8f'%fx,'%.8f'%fy,'%.8f'%fz
+        print('%.8f'%fx,'%.8f'%fy,'%.8f'%fz)
         self.readForce = 0
 
 atom_name = sys.argv[1]
diff --git a/util/qbox_maxforce.py b/util/qbox_maxforce.py
index c1278f3..bb2d01d 100755
--- a/util/qbox_maxforce.py
+++ b/util/qbox_maxforce.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # Copyright 2016 The Regents of the University of California
 # This file is part of Qbox
 #
@@ -9,10 +9,11 @@
 import os.path
 import xml.sax
 import sys
-import urllib2
+from urllib.error import HTTPError, URLError
+from urllib.request import urlopen
 
 def usage():
-  print "use: ",sys.argv[0]," {file|URL}"
+  print("use: ",sys.argv[0]," {file|URL}")
   sys.exit()
 
 argc=len(sys.argv)
@@ -76,9 +77,9 @@ def endElement(self, name):
         if ( fz*fz > fzmax*fzmax ):
           fzmax = fz
           z_name = self.atomname[i]
-      print '%10.3e'%fxmax,'%-8s'%x_name,\
+      print('%10.3e'%fxmax,'%-8s'%x_name,\
             '%10.3e'%fymax,'%-8s'%y_name,\
-            '%10.3e'%fzmax,'%-8s'%z_name
+            '%10.3e'%fzmax,'%-8s'%z_name)
       self.inAtomset = 0
 
 parser = xml.sax.make_parser()
@@ -96,14 +97,20 @@ def endElement(self, name):
 else:
   # attempt to open as a URL
   try:
-    f = urllib2.urlopen(input_source)
-    s = f.read(8192)
-    while ( s !="" ):
-      parser.feed(s)
-      s = f.read(8192)
-    f.close()
-  except (ValueError,urllib2.HTTPError) as e:
-    print e
-    sys.exit()
+     with urlopen(input_source, timeout=10) as f:
+        s = f.read(8192)
+        while ( s !="" ):
+           parser.feed(s)
+           s = f.read(8192)
+        f.close()
+  except HTTPError as error:
+      print(error.status, error.reason)
+      sys.exit()
+  except URLError as error:
+      print(error.reason)
+      sys.exit()
+  except TimeoutError:
+      print("Request timed out")
+      sys.exit()
 
 parser.reset()
diff --git a/util/qbox_move_subsample.py b/util/qbox_move_subsample.py
index ce14182..567d395 100755
--- a/util/qbox_move_subsample.py
+++ b/util/qbox_move_subsample.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # Copyright 2018 The Regents of the University of California
 # This file is part of Qbox
 #
@@ -10,10 +10,11 @@
 import os.path
 import xml.sax
 import sys
-import urllib2
+from urllib.error import HTTPError, URLError
+from urllib.request import urlopen
 
 def usage():
-  print "use: ",sys.argv[0]," interval cmd {file|URL}"
+  print("use: ",sys.argv[0]," interval cmd {file|URL}")
   sys.exit()
 
 argc=len(sys.argv)
@@ -62,17 +63,17 @@ def endElement(self, name):
     elif name == "atomset":
       self.step += 1
       if ( self.step % interval == 0 ):
-        print "#",input_source,"iteration",self.step
+        print("#",input_source,"iteration",self.step)
         avec = self.cell_a.split()
         bvec = self.cell_b.split()
         cvec = self.cell_c.split()
-        print "set cell ",avec[0],avec[1],avec[2],\
+        print("set cell ",avec[0],avec[1],avec[2],\
           bvec[0],bvec[1],bvec[2],\
-          cvec[0],cvec[1],cvec[2]
+          cvec[0],cvec[1],cvec[2])
         for i in range(len(self.tau)):
-          print "move ",self.atomname[i]," to ",\
-            self.tau[i][0],self.tau[i][1],self.tau[i][2]
-        print cmd
+          print("move ",self.atomname[i]," to ",\
+            self.tau[i][0],self.tau[i][1],self.tau[i][2])
+        print(cmd)
       self.inAtomset = 0
 
 parser = xml.sax.make_parser()
@@ -90,14 +91,20 @@ def endElement(self, name):
 else:
   # attempt to open as a URL
   try:
-    f = urllib2.urlopen(input_source)
-    s = f.read(8192)
-    while ( s !="" ):
-      parser.feed(s)
-      s = f.read(8192)
-    f.close()
-  except (ValueError,urllib2.HTTPError) as e:
-    print e
-    sys.exit()
+     with urlopen(input_source, timeout=10) as f:
+        s = f.read(8192)
+        while ( s !="" ): 
+           parser.feed(s)
+           s = f.read(8192)
+        f.close()
+  except HTTPError as error:
+      print(error.status, error.reason)
+      sys.exit()
+  except URLError as error:
+      print(error.reason)
+      sys.exit()
+  except TimeoutError:
+      print("Request timed out")
+      sys.exit()
 
 parser.reset()
diff --git a/util/qbox_move_to.py b/util/qbox_move_to.py
index b657232..d5e6b5b 100755
--- a/util/qbox_move_to.py
+++ b/util/qbox_move_to.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # Copyright 2018 The Regents of the University of California
 # This file is part of Qbox
 #
@@ -9,10 +9,11 @@
 import os.path
 import xml.sax
 import sys
-import urllib2
+from urllib.error import HTTPError, URLError
+from urllib.request import urlopen
 
 def usage():
-  print "use: ",sys.argv[0]," [-iter i] {file|URL}"
+  print("use: ",sys.argv[0]," [-iter i] {file|URL}")
   sys.exit()
 
 argc=len(sys.argv)
@@ -72,18 +73,18 @@ def endElement(self, name):
       self.step += 1
       if ( self.step == iter ):
         if ( iter_option ):
-          print "#",input_source,"iteration",iter
+          print("#",input_source,"iteration",iter)
         else:
-          print "#",input_source
+          print("#",input_source)
         avec = self.cell_a.split()
         bvec = self.cell_b.split()
         cvec = self.cell_c.split()
-        print "set cell ",avec[0],avec[1],avec[2],\
+        print("set cell ",avec[0],avec[1],avec[2],\
           bvec[0],bvec[1],bvec[2],\
-          cvec[0],cvec[1],cvec[2]
+          cvec[0],cvec[1],cvec[2])
         for i in range(len(self.tau)):
-          print "move ",self.atomname[i]," to ",\
-            self.tau[i][0],self.tau[i][1],self.tau[i][2]
+          print("move ",self.atomname[i]," to ",\
+            self.tau[i][0],self.tau[i][1],self.tau[i][2])
       self.inAtomset = 0
       self.done = ( self.step >= iter )
 
@@ -102,14 +103,20 @@ def endElement(self, name):
 else:
   # attempt to open as a URL
   try:
-    f = urllib2.urlopen(input_source)
-    s = f.read(8192)
-    while ( s !="" and not handler.done ):
-      parser.feed(s)
-      s = f.read(8192)
-    f.close()
-  except (ValueError,urllib2.HTTPError) as e:
-    print e
-    sys.exit()
+     with urlopen(input_source, timeout=10) as f:
+        s = f.read(8192)
+        while ( s !="" and not handler.done ):
+           parser.feed(s)
+           s = f.read(8192)
+        f.close()
+  except HTTPError as error:
+      print(error.status, error.reason)
+      sys.exit() 
+  except URLError as error:
+      print(error.reason) 
+      sys.exit() 
+  except TimeoutError:
+      print("Request timed out")
+      sys.exit() 
 
 parser.reset()
diff --git a/util/qbox_msd.py b/util/qbox_msd.py
index 8b83938..fb7d938 100755
--- a/util/qbox_msd.py
+++ b/util/qbox_msd.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # qbox_msd.py: compute mean-square displacement in an MD simulation
 # generate plot of <r^2>(t) in gnuplot format
 # use: qbox_msd.py species file.r [file.r ...]
@@ -8,7 +8,7 @@
 import math
 
 if len(sys.argv) < 3:
-  print "use: ",sys.argv[0]," species file [file ...]"
+  print("use: ",sys.argv[0]," species file [file ...]")
   sys.exit()
 
 species = sys.argv[1]
@@ -53,10 +53,10 @@ def endElement(self, name):
         dy = self.tau[i][1]-self.tau0[i][1]
         dz = self.tau[i][2]-self.tau0[i][2]
         disp2sum += dx*dx + dy*dy + dz*dz
-      print '%12.6f'%(disp2sum/len(self.tau))
+      print('%12.6f'%(disp2sum/len(self.tau)))
       self.step += 1
 
-print "# ",species
+print("# ",species)
 parser = xml.sax.make_parser()
 handler = QboxOutputHandler()
 parser.setContentHandler(handler)
@@ -64,5 +64,5 @@ def endElement(self, name):
   infile = sys.argv[i+2]
   parser.parse(infile)
 
-print
-print
+print()
+print()
diff --git a/util/qbox_position.py b/util/qbox_position.py
index f2f7e21..a02134f 100755
--- a/util/qbox_position.py
+++ b/util/qbox_position.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # qbox_position.py
 # extract position of an atom from Qbox output
 # use: qbox_position.py atom_name file.r
@@ -7,7 +7,7 @@
 import math
 
 if len(sys.argv) != 3:
-  print "use: ",sys.argv[0]," atom_name file.r"
+  print("use: ",sys.argv[0]," atom_name file.r")
   sys.exit()
 
 # Qbox output handler to extract and process <atomset>
@@ -36,7 +36,7 @@ def endElement(self, name):
         rx = float(pos[0])
         ry = float(pos[1])
         rz = float(pos[2])
-        print '%.8f'%rx,'%.8f'%ry,'%.8f'%rz
+        print('%.8f'%rx,'%.8f'%ry,'%.8f'%rz)
         self.readPos = 0
 
 atom_name = sys.argv[1]
diff --git a/util/qbox_repair_h2o.py b/util/qbox_repair_h2o.py
index 14c3144..e4588e5 100755
--- a/util/qbox_repair_h2o.py
+++ b/util/qbox_repair_h2o.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # qbox_repair_h2o.py: repair broken h2o molecules in a Qbox sys file
 # move hydrogen atoms across periodic boundaries to repair molecules
 # use: qbox_repair_h2o.py file.sys
@@ -39,9 +39,9 @@ def fold_in_ws(atom):
     a_cell = float(l[2])
     b_cell = float(l[6])
     c_cell = float(l[10])
-    print line
+    print(line)
   elif ( l[0] == "species" ):
-    print line
+    print(line)
   elif ( l[0] == "atom" ) & ( l[2] == "oxygen" ):
     olist.append([l[0],l[1],l[2],float(l[3]),float(l[4]),float(l[5])])
   elif ( l[0] == "atom" ) & ( l[2] == "hydrogen" ):
@@ -70,14 +70,14 @@ def fold_in_ws(atom):
             nearest_o = o;
   # print "min shift is: ",sx_min,sy_min,sz_min
   if ( sx_min != 0 ) | ( sy_min != 0 ) | ( sz_min != 0 ):
-    print "# current ",h[1]," at ", h[3],h[4],h[5]
-    print "# nearest O is at ", nearest_o[3],nearest_o[4],nearest_o[5]
-    print "# move ",h[1]," by ", sx_min, sy_min, sz_min
+    print("# current ",h[1]," at ", h[3],h[4],h[5])
+    print("# nearest O is at ", nearest_o[3],nearest_o[4],nearest_o[5])
+    print("# move ",h[1]," by ", sx_min, sy_min, sz_min)
     h[3] += sx_min
     h[4] += sy_min
     h[5] += sz_min
 
 for o in olist:
-  print o[0],o[1],o[2],'%10.5f'%o[3],'%10.5f'%o[4],'%10.5f'%o[5]
+  print(o[0],o[1],o[2],'%10.5f'%o[3],'%10.5f'%o[4],'%10.5f'%o[5])
 for h in hlist:
-  print h[0],h[1],h[2],'%10.5f'%h[3],'%10.5f'%h[4],'%10.5f'%h[5]
+  print(h[0],h[1],h[2],'%10.5f'%h[3],'%10.5f'%h[4],'%10.5f'%h[5])
diff --git a/util/qbox_torsion.py b/util/qbox_torsion.py
index e227445..b81f773 100755
--- a/util/qbox_torsion.py
+++ b/util/qbox_torsion.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # qbox_torsion.py
 # extract torsion angle defined by four atoms from Qbox output
 # use: qbox_torsion.py name1 name2 name3 name4 file.r
@@ -7,7 +7,7 @@
 import math
 
 if len(sys.argv) != 6:
-  print "use: ",sys.argv[0]," name1 name2 name3 name4 file.r"
+  print("use: ",sys.argv[0]," name1 name2 name3 name4 file.r")
   sys.exit()
 
 name1 = ""
@@ -117,7 +117,7 @@ def endElement(self, name):
         ss = max(min(we32,1.0),-1.0)
         a = (180.0/math.pi) * math.atan2(ss,cc)
 
-      print '%.4f' % a
+      print('%.4f' % a)
     elif name == "position":
       self.readPos1 = 0
       self.readPos2 = 0
diff --git a/util/qbox_velocity.py b/util/qbox_velocity.py
index e709988..0f011c6 100755
--- a/util/qbox_velocity.py
+++ b/util/qbox_velocity.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # qbox_velocity.py
 # extract velocity of an atom from Qbox output
 # use: qbox_velocity.py atom_name file.r
@@ -7,7 +7,7 @@
 import math
 
 if len(sys.argv) != 3:
-  print "use: ",sys.argv[0]," atom_name file.r"
+  print("use: ",sys.argv[0]," atom_name file.r")
   sys.exit()
 
 # Qbox output handler to extract and process <atomset>
@@ -36,7 +36,7 @@ def endElement(self, name):
         vx = float(velocity[0])
         vy = float(velocity[1])
         vz = float(velocity[2])
-        print '%.8f'%vx,'%.8f'%vy,'%.8f'%vz
+        print('%.8f'%vx,'%.8f'%vy,'%.8f'%vz)
         self.readVelocity = 0
 
 atom_name = sys.argv[1]
diff --git a/util/qbox_xyz.py b/util/qbox_xyz.py
index ae1930c..b09cfea 100755
--- a/util/qbox_xyz.py
+++ b/util/qbox_xyz.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # Copyright 2016 The Regents of the University of California
 # This file is part of Qbox
 #
@@ -10,10 +10,11 @@
 import os.path
 import xml.sax
 import sys
-import urllib2
+from urllib.error import HTTPError, URLError
+from urllib.request import urlopen
 
 def usage():
-  print "use: ",sys.argv[0]," [-first] {file|URL}"
+  print("use: ",sys.argv[0]," [-first] {file|URL}")
   sys.exit()
 
 argc=len(sys.argv)
@@ -75,11 +76,11 @@ def endElement(self, name):
       self.inPosition = 0
     elif name == "atomset":
       self.step += 1
-      print len(self.tau)
+      print(len(self.tau))
       avec = self.cell_a.split()
       bvec = self.cell_b.split()
       cvec = self.cell_c.split()
-      print self.step,\
+      print(self.step,\
       '%.6f'%(a0*float(avec[0])),\
       '%.6f'%(a0*float(avec[1])),\
       '%.6f'%(a0*float(avec[2])),\
@@ -88,11 +89,11 @@ def endElement(self, name):
       '%.6f'%(a0*float(bvec[2])),\
       '%.6f'%(a0*float(cvec[0])),\
       '%.6f'%(a0*float(cvec[1])),\
-      '%.6f'%(a0*float(cvec[2]))
+      '%.6f'%(a0*float(cvec[2])))
       for i in range(len(self.tau)):
-        print self.atomname[i],'%.6f'%self.tau[i][0],\
+        print(self.atomname[i],'%.6f'%self.tau[i][0],\
                                '%.6f'%self.tau[i][1],\
-                               '%.6f'%self.tau[i][2]
+                               '%.6f'%self.tau[i][2])
       self.inAtomset = 0
       self.done_first = True
 
@@ -111,14 +112,20 @@ def endElement(self, name):
 else:
   # attempt to open as a URL
   try:
-    f = urllib2.urlopen(input_source)
-    s = f.read(8192)
-    while ( s !="" and not (first_only and handler.done_first) ):
-      parser.feed(s)
-      s = f.read(8192)
-    f.close()
-  except (ValueError,urllib2.HTTPError) as e:
-    print e
-    sys.exit()
+     with urlopen(input_source, timeout=10) as f:
+        s = f.read(8192)
+        while ( s !="" and not (first_only and handler.done_first) ):
+           parser.feed(s)
+           s = f.read(8192)
+        f.close()
+  except HTTPError as error:
+      print(error.status, error.reason)
+      sys.exit()
+  except URLError as error:
+      print(error.reason)
+      sys.exit()
+  except TimeoutError:
+      print("Request timed out")
+      sys.exit()
 
 parser.reset()
diff --git a/util/qso.py b/util/qso.py
index 588f596..291d015 100644
--- a/util/qso.py
+++ b/util/qso.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 import xml.sax
 # quantum-simulation.org (QSO) definitions
 class Atom:
diff --git a/util/qso2qbox.py b/util/qso2qbox.py
index 7e05362..d414daf 100755
--- a/util/qso2qbox.py
+++ b/util/qso2qbox.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # Convert <atomset> elements from quantum-simulation.org (QSO) format
 # to Qbox input file
 # use: qso2qbox.py [-last] {file|URL}
@@ -9,11 +9,12 @@
 import os.path
 import xml.sax
 import sys
-import urllib2
+from urllib.error import HTTPError, URLError
+from urllib.request import urlopen
 import datetime
 
 def usage():
-  print "use: ",sys.argv[0]," [-last] {file|URL}"
+  print("use: ",sys.argv[0]," [-last] {file|URL}")
   sys.exit()
 
 argc=len(sys.argv)
@@ -47,24 +48,31 @@ def usage():
 else:
   # attempt to open as a URL
   try:
-    f = urllib2.urlopen(input_source)
-    str = f.read(8192)
-    while ( str !="" and not (first_only and handler.done_first) ):
-      parser.feed(str)
-      str = f.read(8192)
-    f.close()
-  except (ValueError,urllib2.HTTPError) as e:
-    print e
-    sys.exit()
+     with urlopen(input_source, timeout=10) as f:
+        s = f.read(8192)
+        while ( s !="" and not (first_only and handler.done_first) ):
+           parser.feed(s)
+           s = f.read(8192)
+        f.close()
+  except HTTPError as error:
+      print(error.status, error.reason)
+      sys.exit()
+  except URLError as error:
+      print(error.reason)
+      sys.exit()
+  except TimeoutError:
+      print("Request timed out")
+      sys.exit()
+
 
 parser.reset()
 
 # write Qbox input file
 datestr=datetime.datetime.utcnow().isoformat()+'Z'
-print "# converted",datestr,"from",input_source
-print "set cell ",s.atoms.cell.a,s.atoms.cell.b,s.atoms.cell.c
+print("# converted",datestr,"from",input_source)
+print("set cell ",s.atoms.cell.a,s.atoms.cell.b,s.atoms.cell.c)
 for sp in s.atoms.species_list:
-  print "species",sp.name,sp.href
+  print("species",sp.name,sp.href)
 for a in s.atoms.atom_list:
-  print "atom",a.name,a.species,a.position[0],a.position[1],a.position[2],a.velocity[0],a.velocity[1],a.velocity[2]
+  print("atom",a.name,a.species,a.position[0],a.position[1],a.position[2],a.velocity[0],a.velocity[1],a.velocity[2])
 
diff --git a/util/qso2qe.py b/util/qso2qe.py
index 4536e53..4993bed 100755
--- a/util/qso2qe.py
+++ b/util/qso2qe.py
@@ -1,17 +1,18 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # Convert from quantum-simulation.org format to QuantumEspresso input format
 # use: env ECUT=ecut qso2qe.py [-last] {file|URL}
 
 from qso import *
 import sys
 import os.path
-import urllib2
+from urllib.error import HTTPError, URLError
+from urllib.request import urlopen
 
 # Get Ecut from environment variable, default=25
 ecut = os.getenv("ECUT",25.0)
 
 def usage():
-  print "use: ",sys.argv[0]," [-last] {file|URL}"
+  print("use: ",sys.argv[0]," [-last] {file|URL}")
   sys.exit()
 
 argc=len(sys.argv)
@@ -45,39 +46,45 @@ def usage():
 else:
   # attempt to open as a URL
   try:
-    f = urllib2.urlopen(input_source)
-    str = f.read(8192)
-    while ( str !="" and not (first_only and handler.done_first) ):
-      parser.feed(str)
-      str = f.read(8192)
-    f.close()
-  except (ValueError,urllib2.HTTPError) as e:
-    print e
-    sys.exit()
+     with urlopen(input_source, timeout=10) as f:
+        s = f.read(8192)
+        while ( s !="" and not (first_only and handler.done_first) ):
+           parser.feed(s)
+           s = f.read(8192)
+        f.close()
+  except HTTPError as error:
+      print(error.status, error.reason)
+      sys.exit()
+  except URLError as error:
+      print(error.reason)
+      sys.exit()
+  except TimeoutError:
+      print("Request timed out")
+      sys.exit()
 
 parser.reset()
 
 # write QE input file
 
-print "&control"
-print "  calculation = 'scf'"
-print "  pseudo_dir = './'"
-print "/"
-print "&system"
-print "  ibrav=0"
-print "  nat=",len(s.atoms.atom_list),", ntyp=",len(s.atoms.species_list),","
-print "  ecutwfc=",ecut
-print "/"
-print "&electrons"
-print "/"
-print "ATOMIC_SPECIES"
+print("&control")
+print("  calculation = 'scf'")
+print("  pseudo_dir = './'")
+print("/")
+print("&system")
+print("  ibrav=0")
+print("  nat=",len(s.atoms.atom_list),", ntyp=",len(s.atoms.species_list),",")
+print("  ecutwfc=",ecut)
+print("/")
+print("&electrons")
+print("/")
+print("ATOMIC_SPECIES")
 for sp in s.atoms.species_list:
-  print sp.name,sp.mass,sp.href
-print "ATOMIC_POSITIONS {bohr}"
+  print(sp.name,sp.mass,sp.href)
+print("ATOMIC_POSITIONS {bohr}")
 for a in s.atoms.atom_list:
-  print a.species,a.position[0],a.position[1],a.position[2]
-print "CELL_PARAMETERS bohr"
-print s.atoms.cell.a
-print s.atoms.cell.b
-print s.atoms.cell.c
-print "K_POINTS {gamma}"
+  print(a.species,a.position[0],a.position[1],a.position[2])
+print("CELL_PARAMETERS bohr")
+print(s.atoms.cell.a)
+print(s.atoms.cell.b)
+print(s.atoms.cell.c)
+print("K_POINTS {gamma}")
diff --git a/util/qso_plot_species.py b/util/qso_plot_species.py
index 2169a0a..cfa2498 100755
--- a/util/qso_plot_species.py
+++ b/util/qso_plot_species.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # Copyright 2018 The Regents of the University of California
 # This file is part of Qbox
 #
@@ -9,10 +9,11 @@
 import os.path
 import xml.sax
 import sys
-import urllib2
+from urllib.error import HTTPError, URLError
+from urllib.request import urlopen
 
 def usage():
-  print "use: ",sys.argv[0]," {file|URL}"
+  print("use: ",sys.argv[0]," {file|URL}")
   sys.exit()
 
 argc=len(sys.argv)
@@ -65,24 +66,24 @@ def startElement(self, name, attributes):
       self.buffer = ""
     elif (name == "norm_conserving_semilocal_pseudopotential"):
       self.inNCSLP = True
-      print "# Norm-conserving semilocal pseudopotential"
+      print("# Norm-conserving semilocal pseudopotential")
     elif (name == "norm_conserving_pseudopotential"):
       self.inNCP = True
-      print "# Norm-conserving pseudopotential"
+      print("# Norm-conserving pseudopotential")
     elif (name == "projector"):
       if self.inNCP:
         self.l = int(attributes["l"])
         self.size = int(attributes["size"])
-        print "# projector l=",self.l," size=",self.size
+        print("# projector l=",self.l," size=",self.size)
       if self.inNCSLP:
         self.l = int(attributes["l"])
         self.i = int(attributes["i"])
         self.size = int(attributes["size"])
-        print "# projector l=",self.l," i=",self.i," size=",self.size
+        print("# projector l=",self.l," i=",self.i," size=",self.size)
       self.buffer = ""
     elif (name == "local_potential"):
       self.size = int(attributes["size"])
-      print "# local potential, size=",self.size
+      print("# local potential, size=",self.size)
       self.buffer = ""
     elif (name == "radial_potential"):
       self.buffer = ""
@@ -92,21 +93,21 @@ def characters(self, data):
 
   def endElement(self, name):
     if (name == "symbol"):
-      print "# symbol:",self.buffer
+      print("# symbol:",self.buffer)
       self.buffer = ""
     elif name == "atomic_number":
-      print "# Z:",self.buffer
+      print("# Z:",self.buffer)
       self.buffer = ""
     elif name == "valence_charge":
       self.valence_charge=int(self.buffer)
-      print "# valence charge:",self.valence_charge
+      print("# valence charge:",self.valence_charge)
       self.buffer = ""
     elif (name == "mass"):
-      print "# mass:",self.buffer
+      print("# mass:",self.buffer)
       self.buffer = ""
     elif (name == "mesh_spacing"):
       self.mesh_spacing = float(self.buffer)
-      print "# mesh spacing:",self.mesh_spacing
+      print("# mesh spacing:",self.mesh_spacing)
       self.buffer = ""
     elif (name == "norm_conserving_semilocal_pseudopotential"):
       self.inNCSLP = False
@@ -117,9 +118,9 @@ def endElement(self, name):
       for i in range(len(self.p)):
         r = i * self.mesh_spacing
         val = float(self.p[i])
-        print '%.6f'%r, '%.10e'%val
-      print
-      print
+        print('%.6f'%r, '%.10e'%val)
+      print()
+      print()
       self.buffer = ""
     elif (name == "projector"):
       if self.inNCSLP:
@@ -127,19 +128,19 @@ def endElement(self, name):
         for i in range(len(self.p)):
           r = i * self.mesh_spacing
           val = float(self.p[i])
-          print '%.6f'%r, '%.10e'%val
-        print
-        print
+          print('%.6f'%r, '%.10e'%val)
+        print()
+        print()
         self.buffer = ""
     elif (name == "radial_potential"):
-      print "# radial potential"
+      print("# radial potential")
       self.p = self.buffer.split()
       for i in range(len(self.p)):
         r = i * self.mesh_spacing
         val = float(self.p[i])
-        print '%.6f'%r, '%.10e'%val
-      print
-      print
+        print('%.6f'%r, '%.10e'%val)
+      print()
+      print()
       self.buffer = ""
 
 parser = xml.sax.make_parser()
@@ -157,14 +158,20 @@ def endElement(self, name):
 else:
   # attempt to open as a URL
   try:
-    f = urllib2.urlopen(input_source)
-    s = f.read(8192)
-    while ( s !="" ):
-      parser.feed(s)
-      s = f.read(8192)
-    f.close()
-  except (ValueError,urllib2.HTTPError) as e:
-    print e
-    sys.exit()
+     with urlopen(input_source, timeout=10) as f:
+        s = f.read(8192)
+        while ( s !="" ):
+           parser.feed(s)
+           s = f.read(8192)
+        f.close()
+  except HTTPError as error:
+      print(error.status, error.reason)
+      sys.exit()
+  except URLError as error:
+      print(error.reason)
+      sys.exit()
+  except TimeoutError:
+      print("Request timed out")
+      sys.exit()
 
 parser.reset()
diff --git a/util/xyz2qbox.py b/util/xyz2qbox.py
index 3f5a032..98e7c87 100755
--- a/util/xyz2qbox.py
+++ b/util/xyz2qbox.py
@@ -1,4 +1,4 @@
-#!/usr/bin/env python2
+#!/usr/bin/env python3
 # generate a Qbox input file from an xyz file
 # use: xyz2qbox.py xyzfile
 import sys
@@ -7,7 +7,7 @@
 buf = line.split()
 nat = int(buf[0])
 line = f.readline()
-print "#",line,
+print("#",line)
 for i in range(nat):
   line = f.readline()
   buf = line.split()
@@ -15,4 +15,4 @@
   x = float(buf[1])/0.529177
   y = float(buf[2])/0.529177
   z = float(buf[3])/0.529177
-  print "atom ",name+str(i+1)," ",name+"_species",'%9.4f'%x,'%9.4f'%y,'%9.4f'%z
+  print("atom ",name+str(i+1)," ",name+"_species",'%9.4f'%x,'%9.4f'%y,'%9.4f'%z)