Numerical simulation methods play an important role in understanding the mechanical behavior of structures in nano and micro scales. Compared with experimental measurements in extremely small scales, computational modeling approaches are considered to have superior abilities to capture delicate structures’ reactions under complex loading and restrained environments. In conventional atomistic simulation approaches, computation of polyatomic structures with large amounts of atoms has obvious size limitations. Recently, the mechanical behavior of protein microtubules has attracted much attention from researchers. As a typical kind of polyatomic bio-structures, a single long microtubule contains up to billions of different types of atoms. In this study, mechanics of microtubules is explored to reveal novel insights. It is a real challenge to develop a more practical theory to consider this kind of polyatomic structure with both result accuracy and computing efficiency. The present research involves a bridging-scales technique based on intrinsic interatomic potential and a continuum description method, based on which the overall mechanical performance of microtubules is investigated. The systematic work begins with evaluation of interatomic potential by using a homogenization technique; large numbers of different types of atoms are replaced by a product of volume densities and the occupied space volumes. The potential energy stored between the basic subunit of microtubules, tubulin dimmers, is obtained from a mutual definite integral process between pair bodies. Without tracing every single atom, deformation of macromolecules components is determined by the proposed fictitious bond connecting central points of neighboring bodies. Unlike traditional continuum mechanics, material properties in this bridging-scales approach are described in an atomistic-continuum way. The material constants are not fundamental prerequisites but deliverable parameters under a minimum system energy state. A constitutive relationship is derived by choosing a representative unit cell in the established lattice structure formed by interconnected fictitious bonds. The orthotropic elastic properties of microtubules are predicted by equating the deformation energy of a lattice structure with that of an equivalent volume of the continuum, and then minimizing system energy in a selected representative unit cell. Further, a mesh-free theoretical and numerical framework based on a higher-order Cauchy-Born rule under higher-order gradients continuity has been specifically set up. This simulation scheme is generally applicable and can be employed to study the overall mechanical behavior of microtubules. Evaluation of strain energy in higher-order gradients continuum approach depends on both the first- and second-order deformation gradients and is determined by deformation of fictitious bonds. A higher-order Cauchy-Born rule is used to build connections between deformation of fictitious bonds in a lattice structure and higher-order continuum deformation gradients in the reciprocal continuum counterpart. Specifically, the second-order deformation gradients describe inhomogeneous deformation. Therefore, while the second-order deformation gradients are considered, approximation of fictitious bond vectors is largely enhanced; the derived outcomes are more reasonable. The development of a mesh-free theoretical and computational framework based on the established atomistic-continuum constitutive relationship is an important contribution of this research. This methodology has its own merits in conducting interpolation of displacements in continuity. While nodal displacements are the only unknown variables, it possesses intrinsic nonlocal properties in mesh-free approximation and automatically satisfies the requirement of the whole higher-order theoretical system. The specific procedure involved in this research is advanced, and has been actualized in the study of mechanical properties of microtubules in the following works. With the proposed methodology, elastic properties, transverse and longitudinal buckling and post-buckling behavior, vibration modes, natural frequencies and dynamic responses are numerically simulated by computer programming. Case studies are conducted. Response of microtubules under hydrostatic pressure is first determined. In this case, an undeformed microtubule can be viewed as having been formed by rolling up a lattice sheet into a cylindrical shape. This sheet is built of fictitious bonds and the rolling process is appropriately written as a set of equations but it is not a rigid transformation of the configuration. An energy optimization proceeds simultaneously to attain a minimum energy state. The microtubule deforms uniformly along axial direction under hydrostatic pressure and the critical value of buckling load and post-buckling structure transition are well simulated. After that, a strip model is proposed for modeling long and curved microtubules, in which the circumferential interactions can be ignored. Critical global buckling loads, vibration modes and frequencies and dynamic responses are further predicted. The influence of different mesh-free parameters, such as total number of mesh-free nodes, weight functions and the dimensionless size of the support domain, are tested. The overall performance of the mesh-free numerical scheme, such as convergence, accuracy and efficiency, are evaluated. Obtained results in all these cases are examined and compared with previous experiments and theoretical results reported by others. Simulation outcomes are proved to be in good agreement with available publications, which make it possible to explore and unveil subtle mechanical responses that are difficult to be measured by experiments in microscales. Finally, this work proposes a cellular automata mechanical approach for microtubules which is suitable for parallel computing required for assessing microtubule interwoven materials and polyatomic large systems. By means of parallelization, the solution algorithms are allocated to parallel processors or multiple personal computer (PC) clusters for parallel processing. By considering the intrinsic atomic interaction, the atomistic-based CA algorithm provides the same degree of accuracy as the atomic-based method, but at a much faster computing speed. This work generates a comprehensive understanding of microtubules and provides useful and valuable references for micro-mechanical research and engineering practices.