1. An optimal power flow (OPF) method of operating a power system, comprising the steps of: a) determining power flows for a base case operation and at least one contingency case operation of the power system, defined as an original discrete OPF problem; b) modifying the original discrete OPF problem to eliminate binary variables and introduce continuous relaxation of the original discrete OPF problem; c) introducing slack variables to transform inequalities present in the original discrete OPF problem, forming an “equality constrained” problem; d) applying an interior point method to produce a solution to the continuous relaxation; e) creating a reduced Newton system of equations defining the solution to the continuous relaxation; f) applying a chain rule process to the reduced Newton system of equations to create an irreducible set of linear equations defining the original discrete OPF problem; g) solving the irreducible set of linear equations and generating a solution to the original discrete OPF problem; and h) operating the power system according to the solution to the original discrete OPF problem; wherein in performing step a), the original discrete OPF problem is defined as follows: 
 min 
 q ( w )=½ x T Q x x+c x T x+ ½ t T Q t t+c t T t+ 
 ½ u * T Q u* u * +c u* T u * +½ u ** T Q u** u ** +c u** T u ** + 
 ½ u 0 T Q u0 u 0 +c u0 T u 0 
 s.t. x−At= 0 
 t−B 0 u 0 −B * u * −B ** u ** −H 0 z 0 =0 
 N z* z * =1 
 N z** z ** −M z z * =0 
 N u u ** −M u u * =0 
 0 ≦x≦d x 
 0 ≦t≦d t 
 D u * z * ≦u * ≦ D u * z * 
 D u ** z ** ≦u ** ≦ D u ** z ** 
 0 ≦u 0 ≦d u0 
 z * ε{0,1} n x * 
 z ** ε{0,1} n x ** where x stores bus-dependent variables; a subscript x indicates a bus-dependent variable at position x; t stores branch-dependent variables; u stores bus and branch control variables; Q represents reactive power constraints, c represents other constraints; the constraint x−At =0 represents flow balance for the power system; the constraint t−B 0 u 0 −B * u * −B ** u ** −H 0 z 0 =0 represents a discretized voltage law defined for each bus; the constraints N z* z * =1, N z** z ** −M z z * =0 and N u u ** −M u u * =0 model special order set-type 1 constraints; D stores upper bounds for x, t and u; D stores lower bounds for x, t and u; z * and z ** are used to test local variables within D and D ; and d x and d t are bounds for x and t, respectively.