The increasing demand for wind energy to decarbonise electricity generation is pushing the design and manufacturing of wind turbine blades in a way that requires a deeper understanding of their structural response, such as the way damage in them initiates and propagates. One damage type of great relevance in composite structures is delaminations, specifically under fatigue loading to account for the entire life-in-service. However, the most advanced models developed by researchers have yet to generally become standard tools of the finite element programs used for structural analysis in the industry. The current work presents a user-defined cohesive element developed for ANSYS Mechanical APDL for modelling the growth under cyclic loads of multiple delaminations in the same structure. The fatigue propagation model earlier developed is expanded to allow for the introduction of delaminations due to fatigue initiation, such as in those cases where there are no initial pre-cracks from which the delamination may propagate. This ANSYS subroutine is evaluated against experimental results from propagation and initiation-dominated tests from literature. The fatigue propagation model is based on the Paris law and calculates the energy release rate through a J-integral contour in the cohesive zone. The implementation of this model is validated through the comparison against standardised fracture toughness experiments, showcasing its great accuracy in predicting the experimental Paris law compared to other models. The potential to model multiple cracks automatically in the structure is also exhibited by modelling modified standard Double Cantilever Beam specimens for quasi-static and fatigue loading. The fatigue initiation approach uses initiation S-N curve data to calculate the number of cycles to the introduction of delaminations, which can then grow according to the propagation model mentioned before. This approach is evaluated by comparing it to previous experimental results for initiation tests, such as Double Notched Shear specimens. These results show the capabilities of this ANSYS APDL subroutine as a tool in the modelling of initiation and propagation of delaminations in composite structures.