The 3D Princeton Ocean Model is described in detail by Blumberg and Mellor (1987) and the Celtic/Irish Sea setup used for this study is detailed by Robins et al. (2013). In summary, the model is based upon a parallelised and sigma-coordinate (terrain-following vertical layers) representation of the primitive equations that describe the conservation of momentum, thermal energy and continuity. Variables within the model (e.g. velocity, elevation, turbulence and temperature) are solved using finite-difference discretisation on an orthogonal Arakawa-C grid in the horizontal plane, which provides an accurate representation of the velocity around complex coastlines such as within the Irish Sea, due to the staggered computation of elevation and velocity nodes.

The model domain is 280 km by 708 km, with a horizontal cell size of 1/30° (longitude) by 1/60° (latitude) resulting in a resolution of approximately 2 × 2 km. In the vertical plane, there are 20 equally segmented sigma layers, with an average and maximum distance between the surface layer and the next layer of 4.3 and 9.6 m respectively, these distances change spatially following the underwater terrain (Robins et al., 2013). Accurate simulations of the long-term circulation are a major challenge to ocean modellers and confidence in the predictions of the model can only be achieved through extensive validation of a range of variables (Robins et al., 2013). The year 1990 was chosen for validation based upon it being a mean representation of the long term (1989-1998) climate (Neill et al., 2010) and availability of temperature records to evaluate the models ability to simulate stratification that can be important for larval transport.

The simulation consisted of a six-month, three-dimensional, barotropic (tide-driven and wind-driven circulation) and baroclinic (density-driven circulation caused by surface solar heating) simulation for the period 1^st April – 30^th September 1990. The model was forced at the open ocean boundaries of the domain with six tidal constituents including the dominant lunar and solar cycles (Robins et al., 2013). These two constituents were validated by comparing the simulated amplitudes and phases of the lunar and solar elevations with records from 18 coastal stations with records of known values (Robins et al., 2013). Both baroclinic and barotropic components of circulation have the potential to play an important role in larval transport and as a result, the main oceanographic factors that transport larvae are the residual currents and tidal velocities. Oceanographic conditions such as, stratification, gyres, tidal range and depth-averaged residual currents for the period 1^st April – 30^th September 1990 are presented in Robins et al. (2013). Wind effects also contribute to the residuals, although due to the stochasticity of the wind over several weeks these are small compared to baroclinic components and are confined to the surface layer.

The ocean model simulated temperature change within the domain. The temperature data was outputted every two days; therefore, the temperature was interpolated between the two nearest recordings to the particle tracking model time step. The oceanographic model was run, and the outputs saved every 30 minutes to be used ‘offline’; this meant that the model only had to be run once and allowed the computing time of the subsequent PTM to be drastically reduced.

A Langragian 3D PTM was developed in MATLAB to simulate the individual movement of particles spatially and temporally. These are based on advection (movement by simulated 3D currents from the ocean model), sub-grid-scale turbulent mixing, and the individual particle vertical migration parameterised on scallop behaviour. Turbulent mixing was simulated in the hydrodynamic model grid scale (i.e., >2 km). However, sub-2 km scale turbulence exists in reality and, so, was parameterised within the PTM (Robins et al., 2013). This was based upon random displacement models (random walks). For the random displacement, the longitudinal change in position (Δχ, Equation 1a), the latitudinal change in position (Δy, Equation 1b) and the depth change (Δz, Equation 1c) over each time step (every save) of the hydrodynamic model of 1800 seconds (Δt, 30 minutes) can be expressed as:

Where R is a random number in the range of 0 – 1 and r is the standard deviation of R cos(2πR). The horizontal diffusivities Kx and Ky (m² s^-1) were set to (u/100)² and (v/100)², respectively.

Simulated larvae were released from locations that are known to have substantial populations of scallops identified from fishing data and stock assessment surveys (ICES, 2021); these were the waters surrounding the Isle of Man, a small population to the northeast of the Isle of Man (hereafter, referred to as Point of Ayre), North Wales, Llyn Peninsula, Tremadog Bay, Cardigan Bay, Tuskar, Celtic Sea and North Cornwall ( Figure 1 ).

Within the model the initial egg size for each larva was randomly assigned within the range 64.2 to 72.2 µm with a mean of 68.8 µm, source: Iroise Sea, France, (Cochard and Devauchelle, 1993; Paulet et al., 1995) to allow for variability among individuals.

The PLD, identified from controlled laboratory experiments for P. maximus larvae from fertilisation to metamorphosis, can last 24 – 78 days depending on sea temperatures (Comely, 1972; Gruffydd and Beaumont, 1972; Beaumont and Barnes, 1992; Robert and Gérard, 1999). PLD is determined by temperature-dependent growth rates, with faster growth in warmer waters. However, within the domain where land meets the sea, and particularly in small embayment’s, the temperature was overpredicted by the model, reaching in excess of 25°C. As a result, where the temperature exceeded 17°C, the maximum attainable temperature in the domain was limited to 17°C within the PTM as identified by the National Oceanography Centre (2022). The average sea surface temperature and bottom temperature is 8°C on the 1^st April 1990 and is 15°C and 13.75°C on the 30^th September respectively (Robins et al., 2013). To simulate growth rates throughout the PLD, a novel approach was coded using growth rates ( Table 1 ), taken from Beaumont and Barnes (1992), at a given temperature and combining it with known sizes of each development stages at 16°C to estimate the size of the larvae at each stage within the larval development ( Table 2 ). This was important as the different larval stages (with different sizes) have a different behaviour in the water. Natural variability is expected in growth rates of individuals or from local food availability (Cochard and Gerard, 1987). Therefore, a random element was introduced to the growth rates (Equation 2a) as they were determined from laboratory studies and fed a stable diet. The larval length was calculated by adding the growth rate per iteration to the starting egg size (Equation 2b); this was calculated using the following equations:

Where L is the lowest growth rate for a given temperature ( Table 1 ), R is a random number in the range of 0 – 1 and Ra is the difference of the growth rate ( Table 1 ) of the two temperature points either side of the given temperature and ts is the number of seconds in a 24 hour period (86400). The growth was assumed to be linear (Beaumont and Barnes, 1992).

In the model the larvae were grouped into five behavioural classes ( Table 2 ) used to capture the development of behaviours in scallop larvae through their ontogeny. Vertical swimming behaviour was introduced into the model when the larvae became trochophores (77 µm) at a rate between 0.5 – 1.0 mm s^-1 towards the surface (Cragg, 1980). Once the trochophores developed into veligers (82 µm), the swimming speed was assigned within the range of 1.0 – 1.5 mm s^-1. Diel vertical migrations (DVMs) were introduced at the same time, sinking during the day and rising to the surface during the night. The depth below the surface of the DVMs that the scallop migrate to was limited to 10 m, as identified by Kaartveldt et al. (1987); Tremblay and Sinclair (1990) and Raby et al. (1994). Once the larvae became eyed veligers (184 µm) they begin to sink and stay just above the substrate. Many authors have described that sinking rates are quicker than swimming rates, but the former have not been quantified (Cragg, 1980; Gallager et al., 1996; Garland et al., 2002). Therefore, it was assumed that sinking rates were the same as the swimming rates for the PTM. All swimming rates were randomly assigned a value within the ranges described ( Table 2 ).

Estimates of the length at which larvae settle is described as between 225 – 258 μm (Gruffydd and Beaumont, 1972; Cochard and Devauchelle, 1993; Andersen et al., 2011). For the PTM the mean value of 241 μm was used as a settlement length threshold. Within the model the particles only settled when they entered the population area polygons delineated in Figure 1 . It has been found that the scallop larvae can delay their settlement typically by two weeks (Robert and Gérard, 1999; Robert and Nicholas, 2000); therefore, particles >241 μm were allowed 14 days to disperse to one of the defined populations. If they encountered a population boundary within this time, then they settled and if they had not entered a boundary after the 14 extra days then the larvae were deemed unsuccessful.

The particles have the potential to leave the model domain at open sea boundaries in the southwest in the Celtic Sea and in the northwest where the Irish Sea meets the North Atlantic ( Figure 1 ). As a result, larvae which encountered the domain boundaries were excluded from the analysis, this was limited to those released from the Celtic Sea and totalled 6.88% of particles released from this site. Larvae that reached a vertical boundary (sea surface or bed), were reflected back into the water column at their position on the previous iteration step. Similarly, larvae that encountered the coastline were reflected back to their previous position within the water column. This allowed for the maximum larval transport potential of each larval cohort to be investigated, following similar methods employed by North et al. (2006) and Robins et al. (2013).

A minimum of 10,000 particles were released from each site, as Robins et al. (2013) identified that this was sufficient to satisfactorily simulate larval transport in this region (i.e. the transport pattern remains similar for larger numbers of particles). The larvae were randomly assigned to starting positions within the coordinates of each delineated population ( Figure 1 ). Simulations were run for each site with all particles released on 1^st April 1990 at once at 12:00, based on a spring spawning of P. maximus (Baird, 1966). At the end of each simulation the percentage of larvae settling within the natal population boundary (retention), settling within a different population boundary (connectivity) and those which did not settle with one of the delineated boundaries (unsuccessful) were calculated. A sensitivity analysis on the date of release and the biological parameters was also undertaken on Cardigan Bay, see Section 2.6.

Overall, the simulation of velocity and temperature fields from the oceanographic model were in good agreement with observational records and other models developed for this region (Xing and Davies, 2001; Horsburgh and Hill, 2003; Hartnett et al., 2007). A tidal analysis conducted by Robins et al. (2013) showed that the simulated circulation compared well with 18 coastal tide gauge stations, giving root-mean square errors of approximately 14 cm (M2) and 6 cm (S2) in elevation and 10° and 12° in phase. Comparisons of the simulated seasonal mean SST with satellite data (National Oceanography Centre 2022), showed that vertically mixed regions were 1-2°C cooler in the model than the satellite data. Despite this the Irish Sea front to the east of the Isle of Man and the Celtic Sea front across the St George’s Channel were well simulated. It is important that these are accurately simulated within the model as they can act as barriers to larvae (Vera et al., 2022). In addition, the simulated temperatures compared well with observed monthly mean surface and bottom temperatures at Port Erin, Isle of Man (see Robins et al., 2013).

The depth profile of an individual particle parameterised on P. maximus larva, released from Cardigan Bay on the 1^st of April 1990 ( Figure 2 ), illustrates the simulated vertical swimming behaviour at different stages of the pelagic life, with migrations from the seabed to the surface once the larva developed into a trochophore after a PLD of 7 days ( Figure 2 , stage 2); exhibited DVM to 10 m when they became eyed-veligers at day 10 ( Figure 2 , stage 3) and at day 62 the pediveliger sank to the seafloor to just above the seabed ( Figure 2 , stage 4). The larvae are supposed to be next to the seabed surface, however, small differences of up to 10 m were seen for the larvae, most likely caused by the vertical velocities exceeding the directional swimming ability of the larvae, although the larvae only remained above the seabed for a short time. Despite their position above the seabed, the larvae were assumed to have settled, a common approach in larval modelling (Hartnett et al., 2007; Nicolle et al., 2013).

The mean PLD for all sites for the model was 95 days with the North Cornwall and the Llyn peninsula populations having the longest mean PLD at 109 and 107 days, respectively ( Figure 3A ). The smallest mean PLD was 83 days for North Wales. All sites except the Celtic Sea and Cornwall had larvae with PLDs across the range rather than aggregated around the mean. This is due to the larvae within the same location encountering areas of varying thermal properties that controls the growth rate. The mean larval transport distance of all sites combined (total distance travelled) was 2192 km with the longest mean transport distance of 3590 km for the Llyn Peninsula and the shortest for the Celtic Sea population, which is subject to low tidal velocities, at 1071 km ( Figure 3B ). In contrast, the pelagic larval distance for most of the populations tended to aggregate around the mean, with only a small number occasionally having a distance that was markedly different. However, larvae from within Cardigan Bay, Isle of Man and Tremadog Bay were subjected to varying tidal velocities such as strong tides around Anglesey, transporting them to the middle of the Irish Sea towards the Isle of Man. Interestingly, no relationship was observed between the PLD and distance travelled ( Figure 3 ).

The regression analysis of connectivity and retention against the population-averaged larval transport distance, for all simulations conducted (n = 41), suggests that the shorter the distance travelled the greater number of larvae that are retained ( Figure 4A , r² = 0.30). Conversely, no relationship between distance and connectivity was found ( Figure 4B , r²< 0.01).

The mean larval settlement success (connectivity and retention) from all populations was 41.3%; retention and total connectivity between populations accounted for 20.3% and 21.0%, respectively, with connectivity occurring primarily between neighbouring populations. Hence, on average, 58.7% of larvae did not find a successful settlement site as delineated in this study. North Wales had the greatest larval settlement success of 89.1%, this was due to retention and primarily connecting to the adjacent Isle of Man population ( Figure 5 ). Conversely, a high number of unsuccessful larvae were found from Point Ayre (99.5%) and North Cornwall populations (98.0%) ( Figure 5 ). Greatest retention was found within the Celtic Sea and Isle of Man populations at 69.5% and 67.11%, respectively, whilst the least retention was found at in the Point Ayre and the Llyn Peninsula populations at < 0.1% and 0.9%, respectively ( Figure 5 ). The number of populations that were connected varied; the Tremadog Bay population was the largest source of larvae, exporting to five different sites, whilst the Isle of Man population was the greatest sink for larvae with a total of four different sources ( Figure 5 ). The Celtic Sea and North Cornwall populations did not connect to any other sites, but larvae had the potential to be retained ( Figure 6 ). North Cornwall had a low retention of larvae (2.1%) and no external sources of larvae. Despite this, observations from Figure 6 demonstrate that North Cornwall had the potential for self-recruitment as most larvae were in close proximity to their release locations.

The biological parameters included in the model explained 51% of variance in the PLD and 29% of the average larval transport distance from analysis on Cardigan Bay, with the remainder of the variation in these two metrics attributed to oceanographic conditions ( Table 4 ). The release date contributed to 15% of the variation with a decrease in the PLD from 88.8 days when spawning was delayed from the 1^st of April to 67.1 days for the 1^st of June, due to increasing temperatures through the spring and summer. The exception from this trend was the simulation for the 1^st of July which had a slightly longer PLD ( Figure 7 ). The population-averaged transport distances increased from 1072 km in April to 2158 km in July with the exception of 1^st of April and the 14^th of June 1990 ( Figure 8 ). The change in the spawning date only contributed to 4% of the total variance within the model ( Table 4 ).

The proportion of larvae being retained larvae stayed broadly the same until June ( Figure 9A ) for each release date. This is due to larvae being transported away from the Cardigan Bay coast towards the Irish coastline driven by a thermal front that develops across the middle of Cardigan Bay. The total proportion of larvae connecting to other populations, primarily Tuskar, the Celtic Sea and to a lesser extent Tremadog Bay, also tended to be broadly the same at between 20 and 30% except for the larvae released on the 14^th of June ( Figure 9A ) which had 80% of the larvae connecting to other populations and 2% of larvae being retained in Cardigan Bay. In most instances where larvae were transported away from Cardigan Bay, no suitable settlement site was found. The change in the spawning date did not appear to drastically change distribution of the larvae, the larvae tended to stay within the area between south Wales and Ireland. The only exception to this occurred in the July simulation when large numbers of larvae were transported north towards the centre of the Irish Sea, with relatively few larvae traveling west towards the Tuskar site and then moving southwest towards the Celtic Sea site, along the thermal front as modelled in previous months leading up to July (see Supplementary Material for detailed distribution maps).

Of the variance presented in the larval transport distances, 22% can be explained by the DVM, which only explained 8% of the variance in the PLD ( Table 4 ). When no migrations occur and migrations up to 5 m occurred, most larvae travelled north to their final destination. When the DVM was 10 m or greater, the larval transport distance was much lower ( Figure 9B ) and they generally travelled towards the west and southwest along the Irish coast. Connectivity also increased with the depth of the DVM from 20% with no DVM to over 50% for DVMs to 20 and 25 m ( Figure 9B ).

Changes to larval size at settlement and changes to the time spawning could be delayed by, explained very little variance in the model for PLD (3% and 2% respectively) and larval transport distance (both <1%). Despite limited influence upon the PLD and transport distance, the retention of larvae increased when the delay period was extended but decreased when the size of settlement increased ( Figure 9D ). Conversely, the connectivity increased slightly with increasing settlement size ( Figure 9E ). The variation in growth rate had little effect contributing to <1% of the variance with little change in percentage of the retention and connectivity of the larvae ( Figure 9C ).

A random number was used within the model to simulate the natural variation in the environment. Observations from changing the random numbers in the parameters found differences in the distributions. Random numbers were used in assigning the individual larvae a starting position within a site, a starting egg size, the growth rate per iteration, the random walk and the swimming speed. However, the random number component explained <1% of the variance within the GLM.

Currently, scallop stocks are managed individually (ICES, 2021) with direct management in territorial waters (12 nautical miles) by regional bodies, and nationally for offshore waters within each country’s exclusive economic zone. The retention and connectivity from this study are magnitudes higher than expectations as no mortality was introduced into the model. However, the general dispersal patterns suggest that there is a single metapopulation within the Celtic and Irish Seas, except for North Cornwall which was reliant on self-recruiting. However, other scallop modelling had identified metapopulations within the English Channel and along the French coast, not within the scope of this study, which have the potential to act as a source for North Cornwall and the Celtic Sea populations (Le Goff et al., 2017; Nicolle et al., 2017). The four populations between the west coast of Wales and Ireland (Llyn Peninsula, Tremadog Bay, Cardigan Bay and Tuskar) appeared to be key to the metapopulation, with a large interchange of larvae which is consistent with other larval models for the region (Hartnett et al., 2007; Coscia et al., 2013). In particular, the exchange of larvae between these three central populations (Cardigan Bay, Tuskar and Tremadog Bay) suggests that they may be resilient to over-fishing as external sources of larvae may be able to replenish the populations (Grimm et al., 2003). Additionally, they also provide larvae for each other and for populations in the south (the Celtic Sea) and north (Isle of Man, North Wales and Point Ayre). The northern area of the Celtic Sea population is reliant on the import of larvae and to maintain the population abundance, the Tuskar and Cardigan Bay populations will need strong management if they are to provide a continual source of larvae. In addition, the central populations appear to be a vital source of larvae for the Llyn Peninsula population which has a low retention rate (0.9%), which acts as a steppingstone connecting the central populations to the northerly populations; but modelling by Hartnett et al. (2007) found retention may occur in this region. The North Wales population is also a crucial population, appearing to be reliant on self-recruitment with low external recruitment of larvae but a large percentage of its larvae connected to the Isle of Man population, an important fishing ground (ICES, 2021). It should also be noted that not all scallop populations were simulated or represented in the model, with populations identified from extensive fishing pressure with defined boundaries used (ICES, 2021). Other populations of scallops are known to occur in waters surrounding Northern Ireland (AFBI, 2012), Scotland (Dobby et al., 2012) and France (Nicolle et al., 2013) with international cooperation needed to provide accurate data. Additionally, populations exist on smaller scales such as within Skomer Marine Nature Reserve, but the area of these populations could not be defined or were too small in scale to be accurately modelled.

The study does not consider that in the natural environment larvae may not settle and recruit into the established populations or the consideration of mortality; it only reflects the probability for a larva to complete its larval lifespan and find a settlement substrate in the defined populations. To simulate the larval transport and define the patterns of dispersal a year with average conditions was chosen. To achieve realistic predictions for a given year for management measures, validated hydrodynamic models should be used. Particularly as scallops take approximately two years from fertilisation to become reproductively active adults (Devauchelle and Mingant, 1991). Therefore, the fishery can benefit if accurate oceanographic models for this period can be created annually to predict patterns of recruitment. This could provide information about when recruitment may be high or low, which could be proven by a coordinated stock assessment that would yield observed densities and distributions that would substantially improve the modelling. The implementation combined with microchemistry or genetic methods, could help validate patterns from mathematical models, particularly from external sources (Nicolle et al., 2013; Hold et al., 2021). For instance, in 1985 large numbers of larvae of P. maximus were trapped in an area with unsuitable settlement substrate by a whirlpool located in St. Brieuc Bay, France, which was caused by strong winds blowing in the opposite direction of the tidal currents. The poor recruitment of individuals from the larvae in 1985 was confirmed by a noticeable decline in the catch of adults three years later (Le Pennec et al., 2003). Therefore, understanding and predicting the larval transport can significantly improve fisheries management and conservation aims and objectives. Current management of scallops is conducted as though they are small unrelated stocks, but the results open up the possibility for the management of the fishery to be internationally integrated, so that the metapopulation can be maintained and to ensure resilience by recruitment from external sources.

The modelling conducted as part of this study provides the most detailed simulations of the connectivity of the well-known established P. maximus populations describing a large interconnected metapopulation in the Irish and Celtic Seas. The larval transport patterns, and hence connectivity networks, appeared to be driven by the hydrodynamics, but the proportion of larvae that were retained or dispersed to other populations appeared to be driven by the biological parameters. In particular, the sensitivity analysis demonstrated seasonal variability in spawning and the DVM may play an important role in the transport to geographically distinct areas. The PLD was overestimated in this study due to uncertainty in the starting larval size and final settlement size, and thus needs to be investigated further. However, a delay in settlement (longer PLD) had little effect on retention but suggest that modelled PLD that are longer than those that occur would overestimate the connectivity to other populations.

This study and others published in recent years (Tian et al., 2009; Nicolle et al., 2017; Hold et al., 2021) demonstrate that the current model of managing the fisheries is unsuitable, with populations reliant on retention and in some cases reliant on external recruitment. A strong effort should be made internationally to manage scallops based on their sensitivity to fishing pressure with key source populations given robust management measures.