new file mode 100644
--- /dev/null
+++ b/scripts/control/is_stabilizable.m
@@ -0,0 +1,96 @@
+# Copyright (C) 1993, 1994, 1995 John W. Eaton
+# 
+# This file is part of Octave.
+# 
+# Octave is free software; you can redistribute it and/or modify it
+# under the terms of the GNU General Public License as published by the
+# Free Software Foundation; either version 2, or (at your option) any
+# later version.
+# 
+# Octave is distributed in the hope that it will be useful, but WITHOUT
+# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+# FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+# for more details.
+# 
+# You should have received a copy of the GNU General Public License
+# along with Octave; see the file COPYING.  If not, write to the Free
+# Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
+
+function [retval,U] = is_stabilizable (a, b, tol)
+
+# Usage: [retval,U] = is_stabilizable (a {, b, tol})
+#
+# Returns retval = 1 if the system, a, is stabilizable, if the pair (a, b) is 
+# stabilizable, or 0 if not.
+#         U = orthogonal basis of controllable subspace.
+#
+# Controllable subspace is determined by applying Arnoldi iteration with
+# complete re-orthogonalization to obtain an orthogonal basis of the
+# Krylov subspace.
+#
+#   span ([b,a*b,...,a^   b]).
+#
+# tol is a roundoff paramter, set to 200*eps if omitted.
+#
+# See also: size, rows, columns, length, is_matrix, is_scalar, is_vector
+#     is_observable, is_stabilizable, is_detectable, krylov, krylovb
+
+# Written by A. S. Hodel (scotte@eng.auburn.edu) August, 1993.
+# Updated by A. S. Hodel (scotte@eng.auburn.edu) Aubust, 1995 to use krylovb 
+# Updated by John Ingram (ingraje@eng.auburn.edu) July, 1996 to accept systems
+# SYS_INTERNAL accesses members of system structure
+# $Revision: 1.1.1.1 $
+
+  if(nargin < 1)
+    usage("[retval,U] = is_stabilizable(a {, b ,tol})");
+  elseif(is_struct(a))
+    # sustem passed.
+    if(nargin == 2)
+      tol = b;          % get tolerance
+    elseif(nargin > 2)
+      usage("[retval,U] = is_stabilizable(sys{,tol})");
+    endif
+    sys = sysupdate(a,"ss");
+    a = sys.a;
+    b = sys.b;
+  else
+    # a,b arguments sent directly.
+    if(nargin > 3)
+      usage("[retval,U] = is_stabilizable(a {, b ,tol})");
+    endif
+  endif
+
+  if(exist("tol"))
+    [retval,U] = is_controllable(a,b,tol);
+  else
+    [retval,U] = is_controllable(a,b);
+    tol = 1e2*rows(b)*eps;
+  endif
+  
+  #disp("is_stabilzable: is_controllable returns")
+  #retval
+  #U
+  #disp("/is_stabilzable: is_controllable returns")
+
+  if( !retval & columns(U) > 0)
+    # now use an ordered Schur decomposition to get an orthogonal
+    # basis of the unstable subspace...
+    n = rows(a);
+    [ua,s] = schur(-(a+eye(n)*tol),'A');
+    k = sum( real(eig(a)) >= 0 );	# count unstable poles 
+
+    #disp("is_stabilizable: unstable poles found:")
+    #k
+    #s
+    #disp("/is_stabilizable: unstable poles found:")
+    
+    if( k > 0 )
+      ua = ua(:,1:k);
+      # now see if span(ua) is contained in span(U)
+      retval = (norm(ua - U*U'*ua) < tol);
+    else
+      retval = 1;			# all poles stable
+    endif
+  endif
+
+endfunction