This thesis project aims to advance the modelling of pressure gradient turbulent boundary layers (PG TBLs) and offer new insights into TBLs modelling. Previous analytical studies have explored various mathematical models, but this research introduces an extended unstacked Deep Operator Networks (DeepONets) architecture with double outputs and five branch parameters. The objective is to capture the mean velocity and Reynolds stress of turbulent boundary layers under pressure gradients. Numerical and experimental datasets of PG TBLs were accessed and utilized to train the DeepONets models. These models successfully predicted the mean velocity and Reynolds stress profiles using outer-scaled parameters. The DeepONets effectively learned the operator that describes the desired profiles based on input parameters, which correspond to the development of boundary layer thickness and pressure gradients. To identify the model with the best prediction performance, error statistics and distribution were examined across different configurations and dimensions. Furthermore, the individual and global sensitivity analyses revealed the relationship between input parameters and their influence on modelling PG TBLs with DeepONets.