This paper extends the capacity of the current sand production models by eliminating the influence of artificial conditions and numerical mesh on localization and deformation response in the sanding model. Past studies indicate strong size effects when using classical elastoplastic models. To rectify this deficiency, a fracture energy regularization method is implemented in the numerical model. The model incorporates both the geomechanical aspects (e.g. rock elastoplastic deformation and rock disaggregation), as well as the transport aspects (e.g. the role of seepage oil rock deformation and solid release). The model employs a Mohr-Coulomb flow theory of elastoplasticity with friction hardening/cohesion softening. Emphasis is given on calibration procedure and validation of the enriched model through back analysis of triaxial and uniaxial compression tests. Next, the model is used to compare the numerical predictions with laboratory data on sand production. 'The comparison incorporates the stress and deformation, as well as the sand volume. The calibration study shows that friction hardening and cohesion softening call satisfactorily reproduce numerically the weak sandstone response to various loading conditions. Further, computation results of strain softening material illustrates that a fracture energy regularization strategy enables the model to exhibit mesh invariance of the energy dissipation.