Acoustoelasticity　has　significant　implications　in　seismic　exploration,　enabling　the　inference　of　geological　structures　and　the　properties　of　oil　and　gas　reservoirs　through　rock　physical　acoustoelastic　parameters,　as　well　as　elucidating　the　propagation　of　elastic　waves　in　subsurface　rocks.In　traditional　acoustoporoelastic　numerical　simulations,　third-order　elastic　constants　are　derived　from　the　Taylor　expansion　of　the　strain　energy　function.However,　under　conditions　of　high　effective　stress,　the　Taylor　approximation　can　lead　to　divergent　elastic　wave　velocities.To　address　this　limitation,　we　propose　employing　the　Padé　approximation　for　the　expansion　of　the　strain　energy　function,　as　it　provides　superior　accuracy　compared　to　the　Taylor　approximation.By　calibrating　the　Padécoefficients　using　experimental　data　for　various　rock　types,　the　resultant　acoustoelastic　constants　exhibit　a　reasonable　theoretical　limit　on　elastic　wave　velocities　as　effective　stress　increases.This　method　allows　for　a　more　precise　characterization　of　velocity　variations　with　stress,　especially　under　high-pressure　conditions.The　acoustoporoelasticity　numerical　simulation　is　mainly　limited　to　the　2D　cases,　so　we　consider　confining　prestressed　mode　in　the　3D　cases.　©　2024　Institute　of　Physics　Publishing.　All　rights　reserved.