Estimating rare event probabilities is a commonly encountered important problem in several engineering and scientific applications, most often observed in the form of probability of failure (PF) estimation or, alternatively and better sounding for the public, reliability estimation. In many practical applications, such as for structures, airplanes, mechanical equipment, and many more, failure probabilities are fortunately very low, from 10-4 to even 10-9 and less. Such estimations are of utmost importance for design choices, emergency preparedness, safety regulations, maintenance suggestions and more. Calculating such small numbers with accuracy however presents many numerical and mathematical challenges. To make matters worse, these estimations in realistic applications are usually based on high dimensional random spaces with numerous random variables and processes involved. A single simulation of such a model, or else a single model call, may also require several minutes to hours of computing time. As such, reducing the number of model calls is of great importance in these problems and one of the critical parameters that limits or prohibits use of several available techniques in the literature. This research is motivated by efficiently and precisely quantifying these probabilities, often encountered in reliability analysis of complex engineering systems, based on a developed framework termed Approximate Sampling Target with Postprocessing Adjustment (ASTPA), which herein is integrated with and supported by gradient-based Hamiltonian Markov Chain Monte Carlo (HMCMC) methods. Hamiltonian Markov Chain Monte Carlo sampling is characterized by much better scalability, faster mixing rates, is capable of generating samples with much weaker auto-correlation, even in complex high-dimensional parameter spaces, and has enjoyed broad-spectrum successes in most general settings. HMCMC adopts physical system dynamics, rather than a proposal probability distribution, and can be used to produce distant proposal samples for the integrated Metropolis step, thereby avoiding the slow exploration of the state space that results from the diffusive behavior of simple random-walk proposals. In this work, we aim to advance knowledge on Hamiltonian Markov Chain Monte Carlo methods, in general, with particular emphasis on its efficient utilization for rare event probability estimation in both Gaussian and Non-Gaussian spaces. This research also seeks to offer significant advancements in probabilistic inference and reliability predictions. Thus, in this context, we develop various Quasi-Newton based HMCMC schemes, which can sample very adeptly, particularly in difficult cases of high curvature, high-dimensionality and very small failure probabilities. The methodology is formally introduced, and the key theoretical aspects, and the underlying assumptions are discussed. Performance of the proposed methodology is then compared against state-of-the-art Subset Simulation in a series of challenging static and dynamic (time-dependent reliability) low- and high-dimensional benchmark problems. In the last phase of this work, with an aim to avoid using analytical gradients, within the proposed HMCMC-based framework, we investigate application of the Automatic Differentiation (AD) technique. In addition, to avoid use of gradients altogether and to improve the performance of the original SuS algorithm, we study the application of Quasi-Newton based HMCMC within the Subset Simulation framework. Various numerical examples are then presented to showcase the performance of the aforementioned approaches.