In this paper, we present SILCA-Newton-Krylov, a new method for accurate, efficient and robust time-domain VLSI circuit simulation. Similar to SPICE, SILCA-Newton-Krylov uses time-difference and Newton-Raphson for solving nonlinear differential equations from circuit simulation. But different from SPICE, SILCA-Newton-Krylov explores a preconditioned flexible generalized minimal residual (FGMRES) method, instead of traditional LU factorization, to solve the system of linear equations in the inner iteration loop. Our key contribution is to introduce an effective and robust quasi-Newton preconditioning scheme to ensure both the robustness and efficiency of iterative methods. Quasi-Newton methods have been explored previously to reduce the number of LU factorization for circuit simulation, however, their converge rate and stability often deteriorate. In this paper, quasi-Newton methods are applied to construct the preconditioners instead of directly to construct and solve circuit equations. Specifically, a systematic method for adaptive time-step size control and a systematic method to generate piecewise weakly nonlinear (PWNL) definition of nonlinear devices are proposed to construct quasi-Newton preconditioners so that the total number of LU factorization for preconditioning is minimized during the entire time-domain simulation. With the PWNL definition, the preconditioner is kept constant if all nonlinear devices reside in their present operating PWNL regions. When nonlinear devices switch their operating PWNL regions, the low-rank update technique is applied to update the preconditioner efficiently rather than to perform new LU factorization. To reduce further the preconditioning cost, we test incomplete LU preconditioners derived from the factorized full L and U matrices, as well as incomplete LU preconditioners followed by an FGMRES preconditioner. All these techniques combined lead to a new time-domain circuit simulation method, named SILCA-Newton-Krylov, which has been implemented into SPICE3. Experimental results on a collection of analog and digital circuits have shown that SILCA-Newton-Krylov is as robust and accurate as SPICE. SILCA-Newton-Krylov is especially attractive for simulating circuits with a massive amount of parasitic RLC elements for post-layout verification. For example, for a nonlinear circuit with power/ground networks with tens of thousand elements, SILCA-Newton-Krylov is shown to yield the SPICE-like accuracy and over 20X overall CPU time speedup over SPICE3, and furthermore the speedup increases with the size of a circuit.