We present results from an extensive analytic and numerical study of a twodimensional model of a square array of ultrasmall Josephson junctions. We include the ultrasmall self and mutual capacitances of the junctions, for the same parameter ranges as those produced in the experiments. The model Hamiltonian studied includes the Josephson, E J , as well as the charging, E C , energies between superconducting islands. The corresponding quantum partition function is expressed in different calculationally convenient ways within its path-integral representation. The phase diagram is analytically studied using a WKB renormalization group (WKB-RG) plus a self-consistent harmonic approximation (SCHA) analysis, together with non-perturbative quantum Monte Carlo simulations. Most of the results presented here pertain to the superconductor to normal (S-N) region, although some results for the insulating to normal (I-N) region are also included. We find very good agreement between the WKB-RG and QMC results when compared to the experimental data. To fit the data, we only used the experimentally determined capacitances as fitting parameters. The WKB-RG analysis in the S-N region predicts a low temperature instability i.e. a Quantum Induced Transition (QUIT). We carefully analyze the possible existence of the QUIT via the QMC simulations and carry out a finite size analysis of T QU IT as a function of the magnitude of imaginary time axis L τ . We find that for some relatively large values of α = E C E J (1 ≤ α ≤ 2.25), the L τ → ∞ limit does appear to give a non-zero T QU IT , while for α ≥ 2.5, T QU IT = 0. We use the SCHA to analytically understand the L τ dependence of the QMC results with good agreement between them. Finally, we also carried out a WKB-RG analysis in the I-N region and found no evidence of a low temperature QUIT, up to lowest order in α −1 .