Added BiCGSTAB for the iterative solution of unsymmetric linear system Ax = RHS using the Preconditioned BiConjugate Gradient Stabilized method. Currently only supports diagonal preconditioners, and is optimized for the case of the identity matrix as a preconditioner.