This study presents a novel three-dimensional (3D) Peridynamic (PD) model based on Peridynamic Differential Operator (PDDO) which is applicable to isotropic and orthotropic materials. It describes bond strain using the first-order derivative of displacement via the PDDO expression and establishes a direct mapping between bond force and bond strain using a micro-modulus matrix. Consequently, the bond force is also a function of the first-order displacement derivative, effectively reducing the constraints required in the governing equation. The concept of a bond-associated family is introduced to characterize the strain tensor of material points, generating pairwise equal and opposite force vectors within each bond which guarantee the symmetric stiffness matrix and numerical stability. The current method enables the direct imposition of natural and essential boundary conditions, making it easy for the failure analysis with complex geometric structures. Benefiting from the inherent advantages of PDDO, the proposed model is also suitable for irregular non-uniform discretization. The validity and accuracy of this approach are verified by predicting the elastic response of isotropic and orthotropic materials under uniaxial tension and pure shear deformation with uniform and non-uniform discretization. Moreover, stress and strain components can be computed via non-local PDDO expressions, enabling the direct application of stress-based failure criteria for damage initiation and propagation. The capability of this approach is demonstrated by considering a 3D structure with a crack subjected to different of boundary conditions.