Overview of Computational Modeling in Nano/Micro Scaled Thin Films Mechanical Properties and Its Applications

Chang-Chun Lee and Pei-Chen Huang

1 Department of Power Mechanical Engineering,National Tsing Hua University,No.101,Section 2,Kuang-Fu Road,Hsinchu 30013,Taiwan.

Abstract:This research reviews the application of computational mechanics on the properties of nano/micro scaled thin films,in which the application of different computational methods is included.The concept and fundamental theories of concerned applications,material behavior estimations,interfacial delamination behavior,strain engineering,and multilevel modeling are thoroughly discussed.Moreover,an example of an interfacial adhesion estimation is presented to systematically estimate the related mechanical reliability issue in the microelectronic industry.The presented results show that the peeled mode fracture is the dominant delamination behavior of layered material system,with high stiffness along the bonding interface.However,the shear mode fracture being dominated as the polymer cover plate with low moduli is considered.The occurrence of crack advance is also significantly influenced by the interfacial crack length and applied loading.Therefore,this paper could serve as a guideline of several engineering cases with the assistance of computational mechanics.

Keywords:Computational mechanics,mechanical responses,nano/micro thin films.

1 Introduction

Over the past half-century,the pivotal concepts and theories of solid mechanics have been maturely developed to describe the mechanical responses of mesoscale structures.From the viewpoint of traditional solid mechanics,the major concern when developing analytic solutions is investigating the mechanical responses of bulk-scaled materials.However,the abundance of novel nano/microprocess and corresponding material structures with unique mechanical behavior implies the limitation of traditional solid mechanics.The aforementioned solution is insufficient for comprehensively considering the particular and complex behaviors of nano/microscale materials with narrow structural dimensions.On this basis,several various computation-based methods have been developed to describe and investigate the material behavior of nano/microscale thin films.Among all the computational methods,the atomistic continuum method (ACM),molecular dynamics (MD) method,and finite element analysis (FEA) method are the most widely used approaches adopted in the applications of computational mechanics,mainly in the investigation of material characteristics,mechanical responses,equivalent material approach,interfacial fracture behavior,and multiscale simulation approaches.

This review aims to explore the literature on the mechanical behavior and applications of various one-,two-,and three-dimensional nano/micro scale thin films,including graphene,carbon nanotube,microtubule,nanowire,electronic packaging,and flexible electronics.

2 Application of computational mechanics on properties estimation of nano-scaled thin films

2.1 Fundamental of molecular potential energy system

In the viewpoint of molecular systems,the mechanical behavior is determined by the connection between each atom.The forces between individual atoms that influences material properties are primarily described.The total energy can be described via forcefield method,which can be expressed as follows:

where Esys,Er,Eθ,EΦ,Eω,EvdW,Eelecare the energy of the total system,the stretching energy,the bending angle energy,the torsional energy in the bonding surface,the torsional energy out of the bonding surface,the energy of the Van der Waals force,and the energy of electrostatic force,respectively.These energies should be correctly considered in accordance with the lattice and crystal structure of the concerned material.In the following section,this paper will introduce the applications on material behavior estimation using ACM and MD methods.

2.2 Application of ACM-based approach

To address the physical limitation of traditional solid mechanics on the behavior estimation of nano-scaled thin films.The concept and approach of ACM are presented to link the relevance of traditional solid mechanics and distinct behavior of nano thin films.It should be noted that the mechanisms of replacing atomic bonds with truss,beam,and spring are totally different and needed to be considered carefully.Tsai et al.[Tsai,Tzeng and Tzou (2010)] showed the fracture characteristics of graphene sheet with hexagonal patterned atoms and described the bonded interaction between the atoms by the concept of the potential energy.A 0.08 J/m2energy release rate was obtained for the crack occurrence in graphene sheet architecture.Xiang et al.[Xiang and Liew (2011)] predicted the buckling behavior of microtubules and adopted a mesh-free framework based on the higher-order Cauchy-Born rule to describe the gradient continuum.The presented results showed that the critical compressive force was highly dependent on the length dimension of the concerned microtubules.A wide range of compressive forces ranging from 0.2 pN to 1.6 pN were predicted as the criteria for the buckling occurrence.Previous literature has revealed that two types of configuration,namely,zigzag and armchair [Shokrieh and Rafiee (2010)],is contained in the hexagon of graphene sheet.Moreover,the Young’s modulus of graphene sheet is estimated as 1.04 TPa,and the dependence of sheet thickness has also been observed.Phung-Van et al.[Phung-Van,Lieu,Nguyen-Xuan et al.(2017)] explored the distribution influence of carbon nanotube-reinforced composite nanoplates and investigated the corresponding mechanical responses under free vibration condition.They found that the size-dependent isogeometric behavior and the natural frequencies of the concerned material decreased with the increase of nonlocal parameters.Similar to the concept of composite mechanics,the representative volume element of graphene sheet is extracted to analyze the related strain energy under an applied cylindrical bending load [Odegard,Gates,Nicholson et al.(2002)].The strain energy significantly reduces with the increase of the applied bending radius.The axial buckling behavior of single-walled carbon nanotubes has been studied under different boundary conditions,such as simply supported,clamped-free,and clamped conditions [Ansari and Rouhi (2010)].The buckling curves show that the critical compressive buckling force tends to converge compared with the other two considered conditions when high aspect ratios are introduced.This phenomenon can be attributed to the chirality effect of clamped-free nanotubes dominating its buckling mechanism.A similar approach was also performed by Rouhi et al.[Rouhi and Ansari (2012)] in estimating the vibration and axial buckling responses of single-layer graphene sheets.The aforementioned studies focus on the concerned material system ignore the interaction effect between the nano/microscale material and related substrate or carrier structure.The collapsed adhesion of layered material system also plays an important role on the reliability issue.Accordingly,the collapsed adhesion between carbon nanotubes and the silicon substrate has been examined [Yuan and Wang (2018)].As the minimum potential energy principle and the energy-variational method are considered,the number effect of walls contained in multiwalled carbon nanotubes is systemically estimated.Results show that the cylindrical configuration becomes stable as the bending stiffness increases.

2.3 Application of MD-based approach

MD approach is another vital computational method that is widely adopted for depicting the potential energy system of nano thin films.This approach is beneficial in exploring the microscopic behavior of atoms and molecules.Some unique influence observed in nano-scaled structure,such as oriented defects,orientation dependence behavior,and grain boundary effect,have been examined using MD simulation.Lee et al.[Lee,Chang,Yang et al.(2013)] constructed a (5,0)/(8,0) single-walled carbon nanotube heterojunction to estimate the mechanical behavior of the concerned material.They performed MD-based tensile and compressive tests and estimated the Young’s modulus and yield stress by 0.2% strain linear regression.The authors indicated that the heterojunction of a single-walled carbon nanotube showed a slight difference on the extracted Young’s modulus;however,the yield stresses on tensile and compressive test conditions were significantly varied.These results could be attributed to the suffering of heterojunction structure from the additional shear force on the heterojunction region.In addition,a unique orientation-dependent behavior on the fracture properties of graphene was presented by Jhon et al.[Jhon,Jhon,Yeom et al.(2014)].The presented results indicated that nearly the same required stress magnitudes to failure occurrence for all different tensile orientations were observed.In other words,the quasi-isotropic elastic behavior for all tensile orientations is required for graphene sheet.Hajilar et al.[Hajilar and Shafei (2015)] studied the constituent properties of hydrated cement paste with distinct nanostructures.The simulated results provided by the MD simulation showed that the elastic properties of considered crystals were influenced by the pore effect contained in the C-S-H gel nanostructure.For nanostructure-reinforced composites,the temperature-dependent behavior of graphene/PMMA matrix with different graphene volume fractions was presented and compared with the rule of mixture results by Lin et al.[Lin,Xiang and Shen (2017)].The traditional method for the composite property estimation was unsuitable for the graphene-reinforced nanocomposites.Therefore,the graphene volume fraction significantly affected the sensitivity of temperature on the Young’s and shear moduli of the concerned graphene-reinforced composite.The buckling that occurs in the intramolecular junctions of single-walled carbon nanotubes also plays an important role.The strain rate effect was first estimated by Kang et al.[Kang,Li and Tang (2010)].The critical strain of buckling occurrence showed slight influence on the low strain rate condition but was sensitive to the high strain rate compressive loading.For silica materials,the temperature effect significantly influences the mechanical characteristics of ladder-like and helical-structured silica nanowire [Lee,Chen,Chen et al.(2014)].A comprehensive study has been performed to summarize the modeling and simulation technique of isolated carbon nanotubes [Raflee and Moghadam (2014)].The mechanical,thermal,buckling,and vibration characteristics of the material have been separately explored.The FEA-based accurate spring mass model enables the vibration characteristic estimation of circular and square graphene [Mahmoudinezhad and Ansari (2013)].Results show that the fundamental frequency tends to decrease with increase of graphene size.Many previous studies have considered the Young’s modulus of graphene at approximately 1 TPa [Shokrieh and Rafiee (2010);Zhang,Ma,Fan et al.(2014)].However,the grain boundary of graphene introduced in the fabrication process might influence the properties of the concerned graphene.The grain boundary-induced inflection observed in the graphene sheet significantly reduces the mechanical strength and Young’s moduli of graphene [Zhang,Duan,Zhang et al.(2013)].The mechanical characteristics of (6,0)/(8,0) coaxial and bias single-walled carbon nanotube heterojunctions can be used to obtain a high stiffness of bias heterojunction under tension loading compared with the coaxial heterojunction [Lee and Su (2013)].Lee et al.[Lee,Chang,Ju et al.(2011)] explored the local stress release mechanism of armchair singlewall zinc oxide (11,11) nanotube via phase transformation and demonstrated the structure rebuilding capability by the phenomenon of the boundary mismatch between different phases.They also investigated the properties of ultrathin zinc oxide with various diameters and cross-sectional shapes.Xiong et al.[Xiong and Cao (2015)] discussed the effectiveness of five potential functions on the description of single-layer MoS2behavior.The results indicated that the reactive empirical bond-order and constant valence force field potential were suitable for describing the mechanical behavior in elastic range.However,the Stillinger-Weber potential was reliable for characterizing the MoS2under an applied large deformation.Many researches have focused on the simulation-based studies because of the uncertainly and variation on nanoscale in situ experimental approach.Nevertheless,an in situ tensile testing approach is demonstrated to extract the fracture properties of graphene sheet,the measured value of fracture toughness is approximately 4.0 MPa*m1/2and shows acceptable variation compared with the MDsimulated results [Zhang,Ma,Fan et al.(2014)].

3 Application of FEA-based approach on the design and mechanical behavior of structures with mixed nano/micro thin films

3.1 Interfacial delamination behavior of layered material system

On the basis of development of nano/microstructures and the corresponding manufacturing process,the multilayer stacking film architecture has been used extensively for producing semiconductors and smart devices.However,interfacial adhesive properties have become critically important because delamination generally occurs on the inherent poor bimaterial during fabrication.To investigate the adhesive behavior of the concerned dissimilar interface,several approaches have been developed to examine the adhesion and feasibility of delamination occurrence.The robust methodologies on interfacial fracture energy predictions using fracture-based FEA simulation integrated with four-point bending test (4PBT) experiment framework has been developed [Lee (2013)].Several approaches in the literature are described in detail in the following section.

3.1.1 J-integral estimation on interfacial adhesion property of dissimilar material system

J-integral approach is an energy-based method for calculating strain energy release rate (G-value) under the elastoplastic condition.Under the situation of pure elastic behavior,the estimated G-value enabled by J-integral approach is identical to the same magnitude of G-value.J-integral estimation can be expressed as follows:

whereWis the strain energy density per unit capacity;Tanduare the parameters of surface traction and displacement vectors alongsideΓcurve,respectively;Γis an arbitrary contour path around the crack tip;and ds is an infinitesimal section of the contour length alongΓ.Notably,the program of integral path is a major factor for the calculation of the interfacial adhesive property on the concerned interface included in the FEA model.A previous work revealed a robust approach on J-integral path determination with a rectangular shape [Lee (2013)].

3.1.2 Implantation of modified virtual crack closure (MVCCT) technique

The interfacial adhesive behavior might be complex because mixed fracture modes occur under intricate loading conditions.To address this issue,the MVCCT approach provides a superior and convenient solution for predicating the dominance of each fracture energy component.Under the situation of four-node element addressed in FEA simulation,the cracking energy of mode-I (GI) and mode-II (GII) can be expressed as follows:

where ΔAis the FEA element length multiplied by a unit width adopted in the simulation model;ndenotes the number of nodes alongside the crack front;andare the nodal forces along the vertical (Z-axis) and longitudinal (X-axis) directions,respectively,at node j1;andδzandδxare the corresponding displacements induced by the aforementioned forces between nodesi1andi2,respectively.This approach aids in determining the various cracking energy components of a dissimilar layered material system using a simple one-step procedure.

3.1.3 Peeling stress and stress intensity factor prediction enabled by crack tip opening displacements (CTOD) approach

The two aforementioned approaches use the concept of strain energy release rate to describe the interfacial adhesive properties of dissimilar films.From the viewpoint of designers and researchers,the criteria of peeling stress and stress intensity factor can be calculated as

where:

The relationship between stresses (vertical stress [σy] and in-plane shear stress [τxy]),relative displacements (δxand δy),and stress intensity factors (KI[opening mode] and KII[shearing mode]) can be explored.μiand νi(i=1,2) represent the shear modulus and Poisson’s ratio of individual materials,respectively.L and ε are the characteristic length and strength of oscillatory singularity,respectively.

The adhesions of several bilayer material systems,including SiO2/SiLK,TaN/SiLK,Ta/SiLK,and Si3N4/SiLK,are estimated by combining experimental work and simulation approaches.Simulated values of interfacial fracture energy are identical to the measured results from 4PBT.Therefore,the presented methodologies can precisely predict the interfacial fracture behavior for multilayer film systems [Lee (2013)].The cohesive zone model (CZM) is another widely used approach for investigating cracking resistance.The strain-dependent and mesh design effect is proposed to modify the traction-separation law [Tvergaard and Hutchinson (1996)].The advantages and limitations of the CZM method are summarized,and the prediction accuracy with different notch depths in threepoint bending architecture is compared [Elices,Guinea,Gomez et al.(2002)].Krull et al.[Krull and Yuan (2011)] performed atomistic simulations to provide suggestions on cohesive traction-separation law and observed a nonlinear traction-separation curve under hydrostatic stress.The criterion of nanoscale Cu/Si interface cracking occurrence in nano-cantilever test and modified 4PBT were presented by Yan et al.[Yan,Sumigawa,Shang et al.(2011)].Moreover,the peeling behavior of a bio-inspired nanoscale film on a rigid substrate with finite length is studied by theoretical and FEM approaches [Peng,Chen and Soh (2010)].The effects of grain morphology,elastic/plastic properties,residual stresses,and the form of the traction-separation law for the presented CZM model are predicted [Rezaei,Wulfinghoff and Reese (2017)].The mixed-mode fracture criterion of CTOD approaches is modeled;and the predicted values between experiment and prediction with various depth/width ratios are compared [Sutton,Deng,Ma et al.(2000)].The fracture behavior of BCC iron crystal under numerical tensile test is studied,and the amount of hydrogen that influences the crack propagation in BCC iron is explored [Telitchev and Vinogradov (2006)].The crack tip stress evolution of magnesium single crystal is estimated using large-scale atomic molecular massively parallel simulation [Tang,Kim,Jordon et al.(2011)].A special healing mechanism of cracking is demonstrated by atomistic simulations;several dislocations are observed around the crack tip under shear stress [Li,Fang,Liu et al.(2015)].Convergence issues in FEA modeling of cohesive interface is resolved by including a small viscous term in constitutive equations [Guo and Bower (2004)].The incremental crack length effect of adhesion strength on several low-k film system is examined by combining FEA simulation and 4PBT [Lee,Huang,Chang et al.(2009)].Lee et al.[Lee,Wei,Chian et al.(2015)] investigated the interfacial adhesion properties of encapsulated films in flexible displays.The presented results showed the feasibility of cracking advancing along the PSA/SiN interface or within the PSA film with the nearly cracking resistance.In another study,the adhesion behavior of the glass interposer is investigated via fracture-based FEA simulation,and the process-influenced high elastic modulus on adhered SiNxcoating is indicated to avoid the occurrence of crack propagation [Lee and Huang (2018)].The modified bilayer beam model is demonstrated to estimate the fracture properties of surface modification layers [Liu,Zhao and Zhang (2018)].

3.2 Strain engineering influence integrated with the layout pattern of the nano-scaled transistor

To meet the performance requirements of nanoscale transistor according to Moore’s law,strain engineering is used.Eventually,it has become the major technique for adjusting the device performance for various technology nodes.Stressors in strain engineering,including embedded source/drain (S/D),strain-relaxed buffer (SRB),and contact etch stop layer (CESL),are widely adopted in planar and three-dimensional device architecture.The effectiveness of S/D and SRB stressors is highly dependent on the extent of lattice mismatch mechanism between the stressor and the material and layout pattern of the concerned device.

The concentrated materials have the same crystal structure.On this basis,the lattice mismatch phenomenon can be described by the virtual thermal strain integrated with the estimation theory of Vegard’s law or the Sentaurus process [Eneman,Simoen,Verheyen et al.(2008)].The formula of Vegard’s law can be expressed as follows:

where aAand aBare the lattice parameters of the concentrated material in pure form,is the lattice constant of the concentrated compound combined with materials A and B,x represents the concentrated mole fraction of material B.Therefore,the corresponding lattice mismatch value f is estimated as follows:

To validate the accuracy and reliability of aforementioned lattice stress theory,a simulation work is performed to estimate the lattice mismatch induced stress magnitude of 20-nm Ge pMOSFET.The nano-scaled transistor is constructed under the assumption of plane strain condition and compared the simulated stress values with the stress data reported by Vincent et al.[Vincent,Shimura,Takeuchi et al.(2011)].Both Si- and Gebased pMOSFET architectures with S/D SiGe and S/D GeSn alloys,are simulated to contrast the longitudinal stress of concerned device channel and plotted in Fig.1.The simulated lattice stress is almost identical to the work proposed by previous literature.Hence,the accuracy of lattice stress theory is validated by these two S/D embedded stressor systems.To examine the lattice mismatch influence on device performance,Thomson et al.[Thompson,Armstrong,Auth et al.(2004)] advanced the concept and benefit of S/D embedded stressor on the p-type MOSFET.The stress-induced performance enhancement could be further quantified by using the piezoresistance behavior of each semiconductor material.An intrinsic-stressed CESL is also regarded as a promising technique for using the preferred stress status on device performance enhancement.The layout design effect and scalability of CESL with different intrinsic stress magnitude were investigated by Eneman et al.[Eneman,Verheyen,Keersgieter et al.(2007)].The results showed that the stress transfer efficiency of the designed CESL was sensitive to the spacer shape and poly-to-poly lengths.Furthermore,the authors presented the stress responses for each gate-first and -last procedure integrated with the different gate materials with a wide range of Young’s moduli.Gate removal during gatelast procedure further improves the stress transfer efficiency of S/D embedded stressors,such as Si1-xGexalloy.Thus,channel width dependence stress comparison between gatefirst and -last procedure of Si1-xGexS/D pMOSFET is explored by process-oriented simulation approach [Lee and Huang (2018)].The foundation effect of stretched gate and fully embedded SiGe array diffusion zones are adopted to analyze the performance of Sibased pMOSFETs [Lee and Hsieh (2018)].As an alternative solution to Si-based devices,Ge devices have an advantage of high mobility and compatible integration in current semiconductor process.Given the excellent characteristics of Ge,numerous studies have focused on the strain engineering and layout design of the Ge-based devices.Lee et al.[Lee and Huang (2018)] demonstrated the interaction effects of S/D GeSn alloy and CESL film on the performance of Ge-based pMOSFET.The results indicated that the gate-last procedure boosted the longitudinal channel stress but sacrificed the benefit on vertical stress because the stress relaxation mechanism occurred on the top surface of the channel during gate removal.In other study,the simulation approach of stress gradient on CESL significantly influences the stress distribution and magnitudes of Ge channel nMOSFET [Lee,Kuo,and Liu (2017)].On the basis of the piezoresistance behavior difference between Si and Ge materials,the biaxial tensile strain induced by CESL layer does not completely favor Ge-based nMOSFETs.On this basis,Lee et al.[Lee and Hsieh (2018)] proposed the superficial layout design of dummy gate to adjust the stress status of Ge channel and magnified the effectiveness of intrinsic stress coatings.The potential of groups III-V semiconductors with high native electron mobility is inspected as another solution for extending the roadmap for semiconductor devices.Similar to the strain technique adopted in group IV devices,many studies have focused on the effectiveness of S/D,SRB,and stressed CESL on device channels based on groups III-V semiconductors.Yeo et al.first demonstrated the strained In0.53Ga0.47As device featuring embedded In0.4Ga0.6As S/D [Chin,Gong,Liu et al.(2009)].The FEA simulated results showed that the lateral tensile and vertical compressive strains were introduced in the In0.53Ga0.47As channel by the in situ doped S/D In0.4Ga0.6As.Chang et al.examined the layout effects of device channel lengths and widths integrated with the InxGa1-xAs S/D stressor for n- and p-type MOSFETs [Lee,Chang,Sun et al.(2011);Chang,Liu and Ou-Yang (2012)].The change of band structure and the corresponding mobility variation induced by three orthogonal directional stress-strain were separately quantified.The SRB is another widely used lattice mismatch stressor.It tends to generate the bi-axial strain into the device channel rather than the uniaxial strain effect induced by embedded S/D.Kim et al.[Kim,Yokoyama,Taoka et al.(2008)] systemically investigated the concentration effect of InxGa1-xAs SRB on the In0.53Ga0.47As channel.Nearly 35% enhancement was estimated when the In0.6Ga0.4As-induced tensile strain was considered.In addition,the three-dimensional device architecture fin field-effect transistor is widely adopted in device technology nodes below 22 nm in accordance with high performances compared with the planar MOSFET.Similar to the study on planar MOSFET,several studies have focused on the performance influences of substrate material selection and gate-stack process [Xu,Ho,Choi et al.(2012);Eneman,Brunco,Witters et al.(2012)].The design rules of multiple stressor arrangement and relative enhancement in device mobility have also been explored [Xu,Ho,Choi et al.(2012)].The effect of substrate recess process for several device nodes beyond 14 nm is further estimated [Eneman,Brunco,Witters et al.(2012)].Increased recessed depth of the S/D region can slightly boost the induced device mobility.However,the effectiveness of SRB virtual substrate is weakened at that moment.These results prove the interaction effect between layout arrangements,and engineered stressors should be considered carefully.

Figure 1:Lattice mismatch stress estimation of nano-scaled transistor with embedded Si and Ge-based S/D stressor

3.3 Multi-level modeling for the FEA architecture with significant scale mismatch

In the microelectronic industry,the layered system design is commonly used to achieve structural integrity under the requirement of structure minimization.In addition,standard FEA simulation is a widely used approach for evaluating the failure possibility based on the thermomechanical stress-strain induced by fabrication or operating conditions.However,the simulation of a structure with significant size difference remains a challenge because of the difficulty and complexity of mesh generation.Thus,the concept and multipoint constraint (MPC) and specified boundary condition (SBC) approach are proposed to handle the foregoing issue on FEA model construction.The MPC method aims to link the uncommon nodes with corresponding degrees of freedom in the FEA simulation.This method can effectively resolve the problem on multiscale model construction and overcome the excessive consumption of computing source.The detailed theories of the MPC method are described as follows.Under the situation of two substances A and B,node M1represents one of the nodes in substance A,and node M1* refers to the uncommon node projected by node M1in element e1in substance B.The corresponding MPC equation can be expressed as follows:

where uM1and u*M1are the displacement vectors of nodes M1and M1*,respectively;Ne1(M1*) is the shape function coefficient matrix at node M1*;and ue1is the displacement vector of element e1.The multiple linear connections of MPC equations are considered the specific boundary conditions and could be expressed as

where H is the shape function coefficient matrix of concerned elements,and u is the global displacement vector.On the basis of the displacement constraints described in the aforementioned characterization,the potential energy of the global system can be described as follows:

where λ denotes the Lagrange multiplier vector;and K and F refers to the global stiffness matrix and global nodal load vector,respectively.

The general approach and corresponding representative example of MPC method have been demonstrated [Ainsworth (2001)].Lee [Lee and Huang (2018)] developed an estimation methodology for interfacial adhesion on multilayer coatings that have a large mismatch extent to conquer the scale mismatch issue on FEA model construction among microelectronic products.The results showed that the second-order submodeling could significantly reduce the sensitivity of FEA mesh density on the estimated energy release rate.In addition,a decent element size ratios design combined with the MPC approach gained accurate prediction results.Both presented approaches could serve as a guideline to overcome the mismatch extent of FEA mesh among dissimilar coatings architectures.The multilevel submodeling technique,that is the SBC approach,has been widely used for transferring the global impact on the concerned local region.Chip package interaction (CPI) effect between flip-chip (FC) and three-dimensional integrated circuit (3D-IC) has attracted extensive attention on its influence on mechanical reliability.CPI-induced stress effect on interconnected system contained in 3D-IC has been previously studied [Lofrano,Vandevelde and Gonzalez (2014)].A twisting occurrence is observed at the edge of Cu pillar,which shows the possibility of failure.Mercado et al.[Mercado,Goldberg,Kuo et al.(2003)] integrated the interface-fracture-mechanics-based FEA simulation with the multilevel submodeling technique to estimate the interfacial stresses and cracking driving forces.The presented results examined the suitable combination of SiO2and low-k layer as dielectric layers in copper/low-k structure could effectively improve reliability of flip chip package.To restrain the near-chip-edge crack issue in FC architecture,the SBCbased simulation is performed to study the effect of crack stop structure design;moreover,the crack driving force for bimaterial interface located at various regions (including leadfree solder and Cu-pillar with and without an additional Cu layer) is explored [Auersperg,Vogel and Lehr (2010)].Another study shows that the die attach process is a critical process that affects the reliability of Cu/low-k architecture [Zhang,Wang,Im et al.(2012)].The delaminated length that depends on energy release rate and the relevant phase angle for interfaces located in each metal level are estimated [Liu,Shaw,Lane et al.(2007)].The crack advance or delamination between the Cu/low-k interfaces is successfully predicted via tied-released technique,and the presented methodology is validated by the experimental 4PBT apparatus [Lee,Lee and Yang (2010)].A coupling method is demonstrated to link the finite element and MD,and its validation on applicability and productivity is enabled by two given examples with three-dimensional cracks.The propagation and dislocation mechanism of cracks and dislocations could be simulated when only fewer degrees of freedom are needed as compared with the direct numerical simulation [Talebi,Silani and Rabczuk (2015)].On the other hand,a coarse grained model is also revealed to estimate the dynamic crack propagation behavior in view of atom level.The rectangular discretization is superimposed on the fine scale model and to approximate the crack path on the equivalent crack surface.Thus,an efficient multi-scaled modeling approach is proposed to discuss the fracture behavior of complex crack pattern [Budarapu,Gracie and Yang et al.(2014)].Quasi-static crack growth mechanism is further studied by the adaptive multiscale method.The presented approach adopted the phantom node method and molecular statics model to separately accomplish the modeling of crack and crack tip [Budarapu,Gracie and Bordas et al.(2014)].An open-resource software based on the concurrent and semi-concurrent multiscale methods are developed to construct the three-dimensional crack modeling and also to offer an efficient software interface of coupling various software solutions in various programming languages [Talebi,Silani and Bordas et al.(2014)].

3.4 Applications of equivalent material approaches on the FEA modeling

Generally,the mechanical properties of engineering materials can be extracted via various in situ material tests,such as tension,compressive,and indentation tests.However,such tests only show the material properties in bulk scale,whereas the evolutions of microstructure and heterogeneous behaviors are ignored.A formulation of fracture-based finite element modeling is developed for microstructural materials [Chang,Wang,Sluys et al.(2002)].The stress-strain relationship is modeled as a continuum with lattice microstructure.The applicability of measurement and simulation approaches on nanostructure materials is enabled by multiscale simulation [Gates,Odegard,Frankland et al.(2005)].A simulation strategy of micro/macro mechanical behavior is proposed,in which large deformations and elasto-viscoplastic mechanism are considered to describe the history-dependent material responses of the matrix material [Kouznetsova,Brekelmans and Baaijens (2001)].However,the aforementioned methodology is not fully suitable for resolving simulation works with significant scale mismatch among all components in an entire FEA model.The mesh mismatch issue is regarded as a challenge on these types of simulation study.Hence,equivalent material test approaches combined with FEA simulation have been developed for modeling multi-heterogeneous materials to an equivalent unit material model.This approach can be used to estimate the mechanical responses of continuous scaling microelectronics.The interfacial behavior of Cu/low-k interconnect system is examined by the equivalent material method integrated with fracture-based FEA simulation [Lee,Yang,Wu et al.(2013)].The sensitivity on layout design in through-silicon via contained interposer is studied [Lee,Chiu,Hsia et al.(2009)].The temperature-dependent behavior is further adopted to analyze the mechanical reliability of microbump array in 3D-IC architecture;results show an effective method in resolving the complexity on FEA modeling construction for the nonconcerned region [Lee,Tzeng and Huang (2015)].The plasticity effect on thermomechanical behavior of metals is examined by the small strain theory integrated with multiscale FEA modeling [Munjas,Canadija and Brni? (2018)].

Figure 2:Schematic diagram of OLED encapsulation architecture with various cover plat material selection:(a) Metal cover plate;(b) PI cover plate

4 Finite element model of OLED encapsulation architecture

An example about interfacial adhesion estimation described in the following section in view of the interfacial delamination is the important failure mode in structures within nano and micro scaled thin films.The OLED encapsulation architecture with metal and polymer (PI) cover plate and the corresponding structure dimensions are shown in Figs.2(a) and 2(b),respectively.Under the assumption of plane strain,a two-dimensional finite element model is constructed to explore the peeling behavior of concerned interface between OLED and indium tin oxide (ITO) film.In view of the boundary conditions of the presented FEA model,the longitudinal displacements of all nodes located at the symmetry axis of the FEA model are fixed to generate the structural symmetry.In addition,the vertical displacement of the bottom node on the symmetry axis is constrained to refrain the phenomenon of rigid body motion.Subsequently,a convextype bending radius of the curvature is applied to estimate the induced stress of OLED encapsulation under various operating conditions.Therefore,the interfacial peeling behavior of the concerned interface between OLED and ITO film is estimated by the CTOD approach and fracture-based FEA simulation.

The critical peeling stress of the concerned interface of 300 MPa is calculated by the CTOD approach.The vertical and in-plane shear stresses are extracted to evaluate the possibility of interfacial peeling occurrence.The induced vertical and in-plane shear stresses on OLED film around the crack tip under various applied bending radius of curvatures are shown in Fig.3.The initial crack length of 0.03 mm is located at the interface between OLED and ITO film.Critical vertical and in-plane shear stresses are introduced on the OLED encapsulation with metal cover plate.A critical radius of curvature is estimated as 3.8 mm on the basis of the stress magnitudes of the second and third nodes that exceed the critical peeling stress.The cracking advance might have occurred in this situation.

Figure 3:Vertical and in-plane shear stress variation induced on OLED film of OLED encapsulation with metal cover plate under several applied convex-type radius of curvatures

In comparison with OLED-induced stress,the induced stress on ITO film is extracted to estimate the possibility of failure occurrence (Fig.4).A1 GPa vertical normal stress is generated under a 1 mm applied radius of curvature.Notably,an initial brittle failure mechanism might have occurred on ITO film because the critical normal stress of ITO film is nearly 1.2 GPa.These results show the different failure modes of brittle conducting ITO film.The PI cover plate is used to improve the flexibility and mechanical reliability of flexible displays under operating conditions.This approach is a promising solution to replace the metal material as the cover plate.The induced stress of OLED film is plotted in Fig.5 considering the PI cover plate.The critical bending radius is further scaled at approximately 1.75-2 mm,which can be attributed to the superior flexibility of PI cover plate with low moduli.However,the induced stress on ITO film does not significantly decrease when the PI cover plate is considered (Fig.6).This behavior is ascribed to the relatively high modulus of the ITO film compared with the OLED film.Moreover,the shifting of the neutral axis plane is another critical parameter that significantly influences the induced stress of concerned films with different cover plate materials.

Figure 4:Stress magnitudes variation introduced on the ITO conducting film embedded in the OLED encapsulation architecture with metal cover plate

Figure 5:Stress curves of OLED layer contained in OLED encapsulation with PI cover plate design under fixed crack length and several applied convex-type radius of curvatures

Meanwhile,the dependence of initial crack length on induced stress is estimated.The corresponding stress curves of OLED and ITO films are shown in Figs.7 and 8,respectively.As the initial crack length extends to 0.07 mm,the in-plane shearing stress becomes dominated on OLED- and ITO-induced stresses.At this moment,the critical bending radius of curvature is estimated as 6.0 mm,which is considerably lower compared with the OLED encapsulation with 0.03 mm initial crack length.This mechanism can be attributed to the OLED encapsulation with long delaminated length,has low overall stiffness,and difficult to resist and the driving forces of cracking advance.

Figure 8:ITO film stress with PI cover plate when an initial crack length of 0.07 mm is acquired

5 Conclusions

This study provides an overview of several feasible applications based on the computational mechanics,such as material property estimations,interfacial adhesive behavior,strain engineering,and multilevel modeling,on simulation approaches.Several widely used computational methods,namely,ACM,MD,and FEA,for estimating the mechanical responses of nano/micro thin films are comprehensively included.For all research fields and corresponding applications,the conclusions are outlined as follows.

1.Property estimation of nano/microstructures:The evolution of potential energy system on concerned nanoscale materials is explored extensively.Moreover,the effects of grain boundary,heterojunction,and phase transformation on nanoscale material behavior are systemically investigated.

2.Interfacial adhesion property of dissimilar material system:The prediction accuracy of several approaches is demonstrated and validated using fracture-based simulation and experimental architecture.

3.Strain engineering of nanoscale transistor:The effectiveness of multiple stressors integrated with the layout design for different groups of semiconductor material (groups IV and III-V) in two-/three-dimensional device architectures are critically reviewed.

4.Multilevel modeling for FEA architecture:MPC- and SBC-based methods are introduced to estimate the mechanical responses of scaling structures with significant scale mismatch.

5.Applications of equivalent material approaches:The formulation and simulation approaches of equivalent material method are developed to investigate the difficulty of material characteristics and simulation modeling on nano/micro thin films.

Acknowledgement:The authors would like to thank the National Center for Highperformance Computing for supporting this research and the Ministry of Science and Technology of Taiwan,R.O.C.for providing financial support under contract numbers MOST 106-2221-E-007-126-MY3 and MOST 107-2622-E-007-010-CC3,respectively.