Coral Symbiosis Carbon Flow: A Numerical Model Study Spanning Cellular to Ecosystem Levels

Corals rely on a symbiotic relationship with algae (zooxanthellae), which reside in the host tissue and play a critical role for host metabolism through photosynthesis, respiration, carbon translocation, and calcification. These processes affect coral reefs on different scales from cellular to organismal and ecosystem levels. A process-based dynamic model was developed and coupled with a one-dimensional (1-D) biogeochemical model to describe coral photosynthesis, respiration, and carbon translocation at the cellular level, calcification and ion transport in different coral polyp components (i.e., coelenteron, calcifying fluid) at the organismal level; and the exchange of material between corals and the ambient seawater at the ecosystem level. Major processes controlling the carbon budget in internal symbiosis were identified. For the symbiont, photosynthesis is the primary carbon source and translocation to the host is the major sink. For the host, most of the carbon translocated from the symbiont is lost through mucus leakage. In the host dissolved inorganic carbon (DIC) pool, most of the carbon is obtained from the surrounding seawater through uptake; photosynthesis and calcification are the major sinks of DIC. Based on a series of scenario studies, the model produced increase of photosynthesis rate with decline of calcification rate under higher air pCO2 and associated carbonate chemistry variabilities in different polyp components. The model results support the hypothesis that elevated pCO2 stimulates photosynthesis, resulting in a reduced supply of DIC to calcification. Such coupled models allow the exploration of process-based mechanisms, complementing laboratory and field studies.


INTRODUCTION
Coral reefs are unique coastal ecosystems that are largely distributed in shallow tropical and subtropical seas. Coral reefs rival their terrestrial counterparts in both biodiversity and productivity, and have been described as the version of the ocean of a tropical rain forest (Morrison et al., 2013). The coral symbiosis system includes calcifying coral animals and their algal symbionts, which together are able to thrive in nutrient-deficient oceanic waters Frankowiak et al., 2016). The main players in the symbiosis are the heterotrophic animal, scleractinian ("reef-building") corals, and unicellular dinoflagellate algae of the genus Symbiodinium ("zooxanthellae"). Scleractinian corals from colonies of polyps are linked by a common gastrovascular system, with living tissues of a few millimeters in thickness calcified by calcium carbonate (CaCO 3 ). This coralzooxanthellae symbiosis involves nutrient recycling and allows algal photosynthesis, which provides organic carbon that fuels animal metabolism and growth, and as a result, leads to evolutionary success (Frankowiak et al., 2016).
Coral reef ecosystems are experiencing increasing pressure from various anthropogenic perturbations worldwide, namely, ocean warming and acidification (Hoegh-Guldberg et al., 2007). The increasing concentrations in air and seawater CO 2 , and a corresponding decrease in pH have lowered coral calcification rates (Andersson et al., 2015). The observed and projected decrease in corals can have a global influence by altering the oceanic and atmospheric CO 2 processes (Eyre et al., 2018;Kline et al., 2019). Understanding the tightly coupled feedback between water column chemistry and benthic coral reef metabolic processes is one of the challenges for projecting coral reef responses to ocean acidification at the ecosystem level (Shaw et al., 2015).
Photosynthesis/respiration and calcification/calcium carbonate dissolution are the two main important biochemical processes that influence the carbon cycle in the coral reef ecosystem. Heterotrophy plays a significant role in coral skeletal growth by enhancing the calcification and supply of metabolic inorganic carbon (Furla et al., 2000). In reef waters, numerous symbiotic coral species can absorb dissolved organic carbon (DOC) (Hoegh-Guldberg and Williamson, 1999) released by benthic alga photosynthesis and decomposition of POC (Zlotnik and Dubinsky, 1989;Mueller et al., 2014). Additionally, healthy corals are typically net producers of DOC via the release of mucus and/or dissolved organic materials that account for losses of 5-45% of photosynthetically fixed C (Crossland et al., 1980). Among these metabolic processes experienced by coral, photosynthesis and respiration drive the internal inorganic and organic carbon fluxes between the symbiont and host coral animals at the cellular level. Meanwhile, ion transport among different coral polyp components, host animal feeding and coral calcification occurs at the organismal level. The ecosystem level processes include the exchange of material between the coral community and the ambient seawater and also air-sea CO 2 interactions at the sea surface (Gattuso et al., 1999). Considering the different types and scales of metabolic processes in a coral reef community, it is essential to develop a better approach to quantify the associated energy fluxes, and assess their feedbacks to environment change.
Given that many different processes co-exist in coral reef ecosystems, numerical modeling is an effective way to investigate processes at each level of the ecosystem and to link these processes, allowing the simulation of coral reef systems under various environmental conditions. Earlier studies have attempted to quantify the symbiotic system mass and energy budgets using isotopic tracers with in situ measurements of carbon budget in photosynthesis, respiration, growth, and mucus and skeletal production Muscatine, 1990;Falkowski et al., 1993;Pupier et al., 2019). Cunning et al. (2017) used a dynamic bioenergetic model for the coral-symbiont system and showed that the balance of autotrophic and heterotrophic nitrogen sources affects the steady state of the symbiotic system. Numerical models have also been used to understand the physiological processes leading to thermal coral bleaching by studying photoinhibition and oxidative stress in zooxanthellae at the cellular level (Gustafsson et al., 2013;Gustafsson et al., 2014). Baird et al. (2018) nested process-based coral-symbiont models within biogeochemical/ecosystem models that considered temperature-mediated light-driven oxidative stress with the mechanistic representation of environmental forcing. Nakamura et al. (2013) used a coral polyp model to study the response of coral photosynthesis, respiration, and calcification at the organismal level. This coral polyp model was then coupled with a three-dimensional (3-D) hydrodynamic model to evaluate reef-scale coral calcification responses to ocean acidification and sea-level rise (Nakamura et al., 2018). Calcification models based on the coral polyp model (organismal level) have been developed to describe ion transport through epithelial cells via the transcellular pathway and intercellular transport through the paracellular pathway (Hohn and Merico, 2012;Hohn and Merico, 2015). These models have been used to scrutinize the day-night CO 2 system changes in the coral tissue, coelenteron, and calcifying fluid and have confirmed the importance of the different ion pathways under the scenario of ocean acidification. While many models have been developed, none have integrated the complexity across all scales in the coral reef system, which is essential to diagnose the stability of coral symbiosis.
The overarching goal of this research is to develop a model capable of simulating energy pathway at different scales in coral reef ecosystems (i.e., cellular, organismal, and ecosystem levels). We use this coupled model for the simulations of coral photosynthesis, respiration, calcification, heterotrophy feeding, mucus release, and associated carbon fluxes in symbiont cells and host cells. We also can check the carbonate chemistry (dissolved inorganic carbon (DIC), total alkalinity (TA), pH) statues in polyp components. The most important objective of this study is to evaluate the coral symbiotic system feedbacks for different atmospheric pCO 2 scenarios. To meet these objectives, this study includes the following components: 1. Build a process-based cellular module that describes coralsymbiont metabolism activities at the cellular level, 2. Build a coral polyp module to calculate ion transfer in the polyp, calcification, and flux exchanges among different coral polyp components, 3. Couple the above modules into ocean biogeochemical model, where the environmental parameters (i.e., light, temperature, currents, air pCO 2 ) and coral autotrophy and heterotrophy energy sources (i.e., DIC, DIN, and DOC) are provided. By doing so, we could examine the coral symbiosis feedback for different air pCO 2 through sensitivity test runs. Figure 1 provides the basic framework of our model system which presents the cell-based coral symbiosis module (Section 2.1), coral polyp module (Section 2.2), and coupling with the ocean hydro-biogeochemical module (Section 2.3). Our modeling approach is based on coupling the three different scaled modules, each of which simulates a different component of the coral system. The detailed explanations of each module are provided below. The model state variable initial conditions are listed in Table 1, with set of parameter values listed in Table 2.

Coral-Symbiont Module
In theoretical models for coral symbiosis, a coral host provides a protective environment for zooxanthellae within its tissues, and there are dynamic fluxes of carbon and nitrogen between zooxanthella, coral host, and the environment through autotrophic and heterotrophic processes ( Figure 1). Of the total zooxanthellae photosynthetic carbon and bind nitrogen, some is used for algal respiration, some is used for maintenance, and some is used for algal cell growth; the rest are passed to the host for its carbon and energy needs, namely, respiration, growth, and production of mucus and skeletal material (Muscatine, 1990). In turn, the symbiont receives an influx of nutrient-rich compounds derived from the metabolic processes of the host (Falkowski et al., 1993). In addition to the internal cycle described above, there are other more important processes associated with energy sources for the coral host: heterotrophic feeding and obtaining DIC and DIN from the surrounding water via the advection from the coral polyp mouth or diffusion through the cell walls (Goiran et al., 1996;Tambuttéet al., 1996).
To describe the internal carbon cycle in the symbiosis, a basic cell-based framework of carbon and nitrogen fluxes in the symbiont system was adapted from Gustafsson et al. (2013). The coral-symbiont module constructed herein employs the cell as the basic unit and describes the carbon flows within and between the symbiont and host. It defines the specifications of carbon and nitrogen functions in energy reserve and functional pools for the symbiont and host respectively. The functional pool is defined as the minimum energy required for an organism to perform vital functions and structural maintenance. The energy reserve pool is determined as the energy required to perform metabolic processes and reserve organic carbon. In the model, C F S and C F H denote the carbon content of the functional pool, whereas C R S and C R H represented the reserve carbon pool. Corals regulate algal growth and activity by limiting their access to nitrogen . Therefore, in the model, nitrogen is assumed to be limited. Cellular N is contained FIGURE 1 | Schematic showing the fluxes of the coupled hydrodynamic, biogeochemical, and coral module. For the cell-base coral module, the symbols corresponding to these processes shown in Tables 3, 4 are listed for clarity. only in the symbiont and host functional pools, the carbon content is assumed to bind to nitrogen at the stoichiometric ratios (IRF FS ,IRF FH ). We also imposed a fixed stoichiometry for translocated carbon and nitrogen (IR T ). In this study, we focus on the discussions of carbon budget, the issues of nitrogen fluxes and budget will be addressed in a future paper. Equations for symbiosis model considering symbiont and coral host are listed in Tables 3 and 4 respectively.

Symbiont Cell
The symbiont produces fixed carbon through photosynthesis based on light and inorganic nutrient. The amount of light coming from the photosynthesis available radiation (PAR) at the sea surface is modified substantially relative to different kinds of scattering. One of the important factors is by the phytoplankton in the water column. The total photon flux density available for coral at the model bottom grid after phytoplankton attenuation is expressed as I, with the unit of µmol photon m −2 day −1 . It is also attenuated by coral skeleton and self-shading by surrounding symbionts following an amplification factor A (Eq. S-1), which indicates the influence of symbiont density on the downwelling irradiance. We assumed that with the increase of symbiont population, there are limited spaces in the host tissue, and S/S max is used as symbiont density function to specify the ability of packing symbionts in the host tissue. The total light available to symbiont is according to the photosynthesis-irradiance relationship (Eq. S-3), where P max S is assumed to be the amount of C needed for symbiont maximum growth rate m max S along with cell respiration (d s ) and cost of biosynthesis (g s ). HCO − 3 is considered in our model as the source of inorganic carbon for the symbiont photosynthesis (Goiran et al., 1996).
The photosynthesis organic C enters the symbiont cell reserve pool first, where part of them will be used to bind to the newly attained N for cell functional need. The V N S,max (Eq. S-5) denotes the maximum newly N uptake rate in the condition of the maximum growth rate of symbiont. There is considerable diel variability of maximal N assimilation based on the information from the literatures (Domotor and D'Elia, 1984;McAuley and Smith, 1995), thus V N S,max is limited in dark condition by the factor a. The actual N uptake follows Michaelis-Menten style equation, and considers that there is enough amount of C bind to these N at the Redfield ratio (V C S , Eq. S-6). Part of the C in the reserve pool will be lost due to mortality where a parameter m s is used to specify a natural mortality rate (M R S,C , Eq. S-7). Symbiont cells lost through natural mortality were either re-absorbed by the host (M S out Eq. S-9), or expelled from the symbiosis system (M C out , Eq. S-10). The excess of fixed carbon is used to produce new symbiont biomass (m S C R S ), and surplus, is translocated to the host reserve pool until C R H reaches C R H,max (T c , Eq. S-11). For the translocation rate, the photosynthesis carbon translocated to the host is limited by the size of the host reserve pool (1 −  biosynthesis (R S,C , Eq. S-13), mortality (M F S,C , Eq. S-8), translocates to the host (T CF , Eq. S-12), and forms new biosynthesis (m S C F S ). We use a constant C thres S to set up the half-saturation constant of respiration and cell maintenance, and describe that biosynthesis of proteins and nucleic acids is limited by the energy leftover in the functional pool.
In the model, the cell division is also considered. To achieve a change in the symbiont population, a new entry is assumed to occur only when there are enough N and C to support respiration and number of maintenance. When the structural biomass in a cell function pool reaches maximum size (C F S,max ), a cell starts to divide according to the growth rate, the population of symbiont cells increases (Eq. S-16). If there is not enough light or nutrients for symbiont photosynthesis, the cell population size decreases at the rate of natural mortality (Eq. S-17). The symbiont growth rate changes with carbon concentration in the function pool (Eq. S-18).

Host Cell
In this module we represent host heterotrophic feeding of DOC and DON by a Holling-type s-shaped curve with a maximum feeding rate and half-saturation constant k DOC (Eq. H-9). DOC concentration is calculated from the biogeochemical model in the ambient water. DOC is bounded with DON and depends on Redfield ratio mass ratio ( 16 * 14 106 * 12 = 0:176). Under the Nlimitation assumption, the N availability for synthesis of proteins and nucleic acids for functional host cell is bounded with certain amount of C. This amount of C comes from the host reserve pool. We used V C H as the N-limited C translocation rate from the host reserve pool to the functional pool. The C sources in the host reserve pool include heterotrophic DOC, translocated carbon from the symbiont, and re-ingestion of dead symbiont cells (Eq. S-9). The C loss includes the C bind to N in host functional pool (V C H ) and excrete outside of the polyp in the form of mucus to the ambient seawater (Eq. H-7). Carbon in the host functional pool is consumed for host biosynthesis and respiration (Eq. H-3), and part is excreted as mucus, which acts as DOC source to the environment water (Eq. H-6). C R H,max (Eq. H-8) denotes the carbon size limitation in the host reserve pool. As C R H approaches C R H,max , the C re-ingestion and translocation from the symbiont to the host will subsided (Eqs. S-9 and S-11). C R H,max is given as the maximum leftover C in the host functional pool after meeting the C bound for the host H,max will change accordingly as C F H changes.

Coral Polyp Module
A coral polyp consists of a bag-like invagination that opens to the surrounding seawater at the mouth ( Figure 2). The polyp module simplifies the coral polyp system into four components: ambient seawater, coral tissue, coelenteron, and calcifying fluid (Hohn and Merico, 2012;Nakamura et al., 2013). Fluxes across these boundary layers connect the four parts with different chemical components, driving the dynamics of the polyp model ( Figure 2). The associated mass balance terms are listed in Table 5.

Host Tissue
The host is assumed to be in direct contact with the external environment. We considered the exchange of inorganic form of C between the coral tissues and the ambient water. DIC supply is suggested as saturated at ambient water (Goiran et al., 1996), so DIC uptake from the seawater is parameterized with Michaelis-Menten kinetics, with a maximum uptake rate and a halfsaturation constant (Eq. H-1). The DIC in the host tissue then can be used by symbiont photosynthesis (Eq. S-3). The maintenance respiration products from symbiont and host in the form of DIC go back to the host tissue DIC pool (Eq. H-3). So the balance of these source and sink terms determine the change of DIC in the host tissue and also the equilibration with DIC in ambient seawater. We assumed the leftover carbon in the host tissue DIC pool can be available to transport to the coelenteron and calcifying fluid through transcellular transport (Furla et al., 2000). It is assumed to be light limited, where E I was defined as the switch of light conditions in the simulation following Hohn

Symbol Description Equation Unit
A Coral self-shading factor 1:26 + 1:39e (−6:48 S Smax ) (S-1) - The rate of C bind to the newly attained N V N S,max and Merico (2012). When I >0, E I = 1, while when I = 0, E I = 0 (Eqs. P-1 and P-2). For other carbon pathways, CO 2 was assumed to diffuse freely between the seawater and coral tissue (Eq.P-3, F CO 2 sea ), between the coral tissue and coelenteron through the oral endoderm (Eq.P-4, F CO 2 coe ), and between the host tissue and calcifying fluid through the aboral endoderm (Eq.P-5, F CO 2 cal ). So the mass balance of DIC in the host tissue (DIC tiss ) is determined by DIC uptake, the internal cycle of symbiosis metabolic processes (i.e., photosynthesis, respiration, and calcification), and the equilibria diffusion-associated removal or addition of DIC (Eq. P-23). DIN fluxes have the similar pathway as DIC, except they are not involved into calcification processes. About the DIC variance associated pH and TA, the coral tissue pH is known to be higher under light conditions than in the dark (Venn et al., 2009), implying that photosynthesis and calcification processes can regulate the coral host pH. However, considering the complexity of intracellular cell adjustments under environmental change, pH and TA in the host tissue were assumed to be constant with no flux equations specified. We also calculated the fluxes influencing DIN in the host tissue. We assumed that the coral animal acquires DIN from the surrounding seawater represented as ammonium Hoegh-Guldberg and Williamson, 1999). Ammonium diffuses into the host tissue pool and is taken up by symbiont are described by Michaelis-Menten style equation (Eq. H-2). The DIN in the host tissue (DIC tiss ) is under the balance of NH 4 uptake from ambient seawater by the host, host respiration, symbiont respiration, and loss by the symbiont uptake (Eq. P-24).

Coelenteron
DIC and TA in the coelenteric fluid composition is controlled by the advective exchange of seawater through the mouth, transcellular diffusion with tissue, paracellular diffusion with seawater through the oral boundary layer, and calcifying fluid through the aboral endoderm layer. The paracellular diffusion fluxes between seawater and coelenteron were defined in Eqs. P-8 to P-9, where k p is the conductivity coefficient for the paracellular leakage process. The mouth of the polyp opens the coelenteron to the ambient seawater, the flux of advective exchange that will influence state variables between seawater and coelenteron was assumed to linearly depend on the concentration gradient between the two compartments with a constant exchange rate w. The exchange rate w is linearly proportional to the shear stress (t) of the 0.4 root (Hearn et al., 2001) (Eq. P-12). In calculation of t, we follow Eq. P-13 where r is water density, and C b is bottom drag coefficient. The calculated values of t is about 0.05 N m −2 . The H + transport by Ca-ATPase (F H ) also influences TA in the coelenteron (Eqs. P-25 and P-26).

Calcifying Fluid
As for the calcifying fluid, several processes influence the concentration of DIC and TA, namely, the diffusion of CO 2 , DIC transport from the host tissue to the calcifying fluid, paracellular diffusion of DIC and TA, and calcification.
Calcification is an important process in the calcifying fluid and is regulated by the cellular supply of ions, such as Ca 2+ and CO 2− 3 (Allemand et al., 2011). Corals also need energy to transport calcium, bicarbonate, and other ions to calcify their skeletons. Calcium ions are believed to be transported via a combination of calcium channels and an active calcium pump (Ca-ATPase) (Cohen and McConnaughey, 2003). The free energy from ATP hydrolysis is from host respiration. This is one of the ways the coral polyp model coupling with the coralsymbiont module, through the energy flux used for driven Ca-ATPase (Eq. P-18). In our model, the calcification formula was adapted from Nakamura et al. (2013). The Ca-ATPase transports  coe ) and is simplified to be a constant (DCa 2+ cons ), based on the assumption that the Ca 2+ difference is much smaller than the H + difference. The transport rate of H + through Ca-ATPase is estimated as the ratio of energy consumed by Ca-ATPase for 1 mol of H + , as per Eq. P-20. H + exchange is thought to be a crucial step in the coral biomineralization process and influences pH and TA in the coelenteron and calcifying fluid. H + in the calcifying fluid follows transcellular Ca 2+ /2H + transport (Eqs. P-26 and P-28).
The precipitation of aragonite is described by the equation k p (W a − 1) n p (Burton and Walter, 1990), where W a represents the concentration of Ca 2+ and CO 2− 3 dissolved in the calcifying fluid until the equilibrium concentration is reached (Eq. P-21). For our calculations, a constant value of [Ca 2+ ] = 10.4 mM is used (Nakamura et al., 2013). The CO 2− 3 ions involved in the calcification process are sourced from the calcifying fluid DIC pool. For our calculations, a constant value of [Ca 2+ ] = 10.4 mM was used (Nakamura et al., 2013). The CO 2− 3 ions involved in the calcification process are sourced from the host DIC pool. The calcification rate (G) was calculated following Eq. P-22. When 1 mol of CaCO 3 is produced (calcification), 1 mol of DIC cal and 2 mol of TA cal are consumed. The mass balance equations are shown in Eqs. P-27 and P-28.
Based on these flux terms, we calculated DIC and TA as the state variables in each of the polyp components, and their change with time associated with fluxes is derived using differential equations. A forward difference scheme is used for these equations. The units of all fluxes on the right-hand side of each equation are µmol C cm −2 s −1 , whereas divided by the thickness of each component (i.e., tissue thickness, coelenteron thickness, and bottom layer depth in the water column), and their units are consistent with the state variable on the left-hand side of the equation.

Terms in the Marine Biogeochemical Model
The biological part of our model represents the basic functions and biogeochemical processes in the lower trophic levels. The model is modified from the biogeochemical model consisting of the nitrogen cycle with carbonate chemistry (Fennel et al., 2006;Fennel and Wilkin, 2009;Xu et al., 2013). There are seven state variables: phytoplankton, zooplankton, nitrate, ammonium, small and large detritus, and chlorophyll. The time rate change of phytoplankton is influenced by the growth rate of phytoplankton, grazing by zooplankton, mortality, aggregation of phytoplankton to small and large detritus, and vertical sinking of the aggregates. Zooplankton was assumed to assimilate ingested phytoplankton while the remaining fraction is being transferred to the small detritus pool. Other zooplankton loss terms are excretion to ammonium and mortality. Phytoplankton and small detritus aggregate to large detritus. The nitrogen cycle includes phytoplankton uptake, zooplankton grazing, and microbial loop processes. Detritus remineralization feeds into the ammonium pool, and ammonium is subsequently nitrified to produce nitrate. Variables of the carbonate system (DIC, pCO 2 , carbonate ion concentration, bicarbonate ion concentration, pH and TA) determines the remaining four, but only DIC and TA behave conservatively with respect to mixing and temperature and pressure changes (Zeebe and Wolf-Gladrow, 2001). Thus, DIC and TA are included as active tracers (i.e., they were advected and diffused by the physical model). Local changes in DIC and TA occur due to primary production, respiratory processes, and gas exchange at the air-sea interface, whereas pCO 2 -the carbon variable relevant for the air-sea gas exchange of CO 2 -is calculated for the surface layer only. We also simulate DOC sources from phytoplankton excretion and decomposition of large detritus. The DOC sources include phytoplankton (E Phy ), zooplankton excretion (E Z ) and decomposition of large detritus (D LDet ) and lost by bacterial uptake (U DOC ). These processes were simplified from the work of Li et al., 2014), and are described in Table 6.

Model Run
First, we conducted a control run only for the coral symbiosis module and polyp module. In the control run, the environment variables were assumed consistent as listed in Table 1, only the light changes in day and night. There was no flux exchange between the coral system and the seawater. The model results in the control run were used to compare with observations collected in the laboratorial experiment, to make sure that our model is under reasonable assumption and parameterizations.
Based on the parameters set in the control run, we then designed a coupled run to integrate the coral modules with the hydro-biogeochemical model. Before coupling with the coral module, we did an idealized 1-D run for the hydrobiogeochemical model based on the Regional Ocean Modeling System (http://www.myroms.org, Haidvogel et al., 2008). This model was set up with 15 vertical layers and a maximum depth of 20 m, with timestep of 360 s. We assumed that there are no boundary inputs and horizontal advection in the 1-D simulation (Supplementary Figure 1). After reaching the steady state, the results were used as the initial input for the coupled model.
Coral reefs are typical benthic marine primary producers; therefore, to integrate the coral modules with the hydrobiogeochemical model, we assumed that two-way coupling took place at the bottom grid of the hydro-biogeochemical model (Figure 2). We set another idealized case with 15 vertical layers and grid size in both the x (dx) and y (dy) directions is 10,000 cm. Bottom grid height dz is set to 500 cm. The numerical simulation code was written in Fortran 90/95. The simulation was run with a time step of 0.001s. We used the temperature, salinity, light, velocity, nutrients, and diffusion coefficients from the output of the former idealized hydro-biogeochemical model run as the input for the coupled model. The hydrodynamic model was then forced only by the solar radiation with a diurnal cycle. The bulk formulae (Fairall et al., 2003) usually used to force the hydrofield in the ROMS system were closed in our experiment. So the temperature, salinity, and velocity were all constant over time. Vertical diffusion and advection exist between the upper water column and bottom benthic layer, where the coral reefs are located. We assumed the coral system as the point sources and sinks (Q src ) at the bottom layer of the biogeochemical model for the tracer variable DOC, DIC, DIN, and TA following Eq.
(1), where k = 1 indicates the bottom layer of the biogeochemical model, k = N indicates the surface layer, C C denotes the colony coverage, we assumed that there is 70% coverage in the mode grid. We use the variable DIC as an example here; the same formulas were applied to calculate DIN, TA, and DOC. The subscript "sea" denotes the DIC in the coral polyp model in the ambient seawater. The DOC, DIC, DIN, and TA are two-way coupled in the way that, when there are changes in DIC sea , TA sea , DIN sea , and DOC sea in the ambient water of the coral polyp model, these changes will influence the above biogeochemical model grids through vertical advection and diffusion.  In our model, pCO 2 in the seawater (pCO 2sw ) could be influenced by gas exchange and DIC biological uptake, including both the phytoplankton and coral systems. The gas exchange on the sea surface is modeled to be transported downward by vertical advection and diffusion and to influence DIC in the water column. These designs help us apply the model scenario experiments under different air pCO 2 of 340 and 800 ppm. The sea surface CO 2 gas exchange (CO 2_flux ) was calculated as follows: where k CO2 is the gas transfer velocity, CO 2_sol is the solubility coefficient of CO 2 calculated after Weiss, 1974), H is the depth of the water column. When air pCO 2 changed, it will influence the DIC concentration in the water column, and also the uptake of DIC in the ambient water by coral.

Control Run Results
In the control run, the environmental temperature is set as 26°C, and with salinity of 35 psu. The ambient water DIC sea , TA sea , DIN sea and DOC sea were set as the constants in Table 1, and there was no flux exchange between the coral system and the seawater, that means the right-hand side in Eqs. P-29-P-32 equal zero. The coral system variables ultimately reached a steady state after 500 days from the beginning of the model run. We only show a 24-h cycle for all the variables. The daily flux of photosynthetically fixed carbon in light-and shade-adapted measured by isotopic tracers in Muscatine et al. (1984) was used to compare with our modeled coral symbiosis carbon budget. We used the measurements in Al-Horani et al., 2003) to evaluate the modeled pH in different polyp components in the light and dark conditions. In addition, we also compared the measurements of DIC, TA, pH, in different coral polyp components, and W a value under natural light conditions (Cai et al., 2016).

Diurnal Variations in the Symbiosis
For the symbiont, the model simulated gross and net photosynthesis rate exhibited a diurnal variation with light intensity (Figures 3A, F). In the symbiont reserve pool, part of the carbon accumulated due to photosynthesis ( Figure 3D) was transferred to the symbiont functional pool, thus the symbiont functional pool carbon displayed high concentration in light ( Figure 3E). Some carbon in the functional pool are used for respiration, some are used for cell growth ( Figure 4C), and some are removed by loss of mortality. The left-over carbon was transferred to the host reserve pool ( Figure 3C). The symbiont respiration rate increased due to carbon increase in the symbiont functional pool, so it was higher during daytime ( Figure 3E). As for the host, two major carbon sources contributed to the host reserve pool: one was from the translocated carbon from the symbiont, and the other was from heterotrophic feeding of DOC from the ambient water. The 24-hour carbon cycle in the host reserve pool presented a strong diurnal signal, which was higher during the daytime due to translocation ( Figure 4A). The carbon in the host reserve pool was transferred to the host functional pool for the host respiration energy needs. Therefore, the carbon in the host functional pool accumulated during light ( Figure 4B), and the host respiration rate also increased in light ( Figure 4D). The respired carbon from the symbiont and host were recycled back to the DIC pool in the host. Part of DIC in the host DIC pool was reused by the symbiont for photosynthesis, and another part for calcification. The calcification rate increased with increasing aragonite saturation, which was associated with [CO 2− 3 ] in the calcifying fluid (Figures 4E, F).

Carbon Budget in the Symbiosis
Quantifying how much carbon is fixed by symbiont, how much is made available to the animal by translocation and then recycled back to the symbiont is very important to evaluate the health of the symbiotic association. To illustrate the daily energy budget in the symbiont-host internal cycle, we compared each term that contributes to the carbon budget in the symbiont and host cellular carbon pool (dC F S =dt + dC R S =dt, for details, please refer to Table 3, Eqs. S-13 and S-14) and also the inorganic carbon pool in the host tissue. The budget calculations were based on an integration over a diurnal cycle (24 h).
For the symbiont ( Figure 5A), photosynthesis was the major source of reduced carbon ( Figure 5A, black line), and translocation to host predominantly caused the loss of fixed carbon ( Figure 5A, purple line). Respiration also contributed to carbon loss, although its rate was an order of magnitude less than that of photosynthesis and translocation ( Figure 5A, red line). The photosynthesis rate is a function of irradiance. The gross photosynthesis rate is usually tens of times higher than the respiration rate in light conditions. The light respiration rate calculated herein was higher than the dark respiration rate. As for the daily carbon budget (Table 7), 21.50 µg C cm −2 of carbon was fixed through photosynthesis; approximately 30% (6.55 µg C cm −2 ) of this carbon was respired by the symbiont, 25% (5.23 µg C cm −2 ) was lost due to the symbiont mortality, 38% (8.23 µg C cm −2 ) was translocated to the coral reef host, and 6% (1.33 µg C cm −2 ) was used for the zooxanthellae growth.
Concerning the host, the translocated carbon ( Figure 5B, black line) contributed up to 79% of the carbon source (Table 7), whereas re-ingested carbon from dead symbiont cells accounted for 17% of the carbon source. The heterotrophic feeding rate revealed a diurnal cycle under the assumption that the coral polyp feeds and digests only during the night (Agostini et al., 2012) (Figure 5B, red line) and accounted for 4% of the reduced carbon income for the coral reef host. Approximately 84% of carbon was shown to be lost with the death of the host in the form of mucus, and 15% was used for host respiration.
The host tissue DIC pool is the inorganic carbon source for symbiont photosynthesis and calcification ( Figure 5C). As illustrated in Table 7, the DIC uptake by the coral host tissue from the surrounding seawater contributed up to 93% of the total DIC income in the host tissue. Carbon respired by the symbiont and host accounted for 5 and 2% of the DIC input, respectively.
Calcification and photosynthesis were the major processes of consuming DIC and accounted for 24 and 18% of DIC, respectively. Other DIC losses in the host tissue DIC pool originated from the diffusion of CO 2 to the calcifying fluid, seawater, and coelenteron, which were driven by the concentration gradient of CO 2 between the components. The CO 2 diffusion to the calcifying fluid contributed 4% of the DIC loss, whereas diffusion to the seawater and to the coelenteron contributed 2% of the DIC loss each.

Carbonate Chemistry in the Coral Polyp
The polyp model simulated the time series of the DIC and TA in the coelenteron and calcifying fluid of the polyp compartments. Calcification will directly influence carbonate chemistry in polyp component. It causes a carbon concentration gradient between the host tissue, coelenteron, and calcifying fluid, and drives the transcellular diffusion exchange, thereby affecting DIC and TA in the coelenteron and calcifying fluid. Calcification needs energy to transport Ca 2+ to the skeleton in exchange for H + , where the energy flux originates from host respiration. The model results showed increased respiration under light conditions ( Figures 3B,  D), thus, more ATP was generated during respiration and led to an increase in calcification in light ( Figure 4F). This suggests the coral biomineralization which has been identified by many coral scientists (Marshall, 1996;Gattuso et al., 1999). Due to calcification, DIC cal significantly decreased in the daytime ( Figure 6D). The enzyme Ca-ATPase removed protons from the calcifying fluid, thus TA cal and pH cal increased in the daytime ( Figures 6E, F). Based on the mass balance term that influences changes in DIC cal , calcification contributed most of the DIC loss in the calcifying fluid, the other carbon loss was CO 2 diffusion through the eukaryotic cell membranes from the calcifying fluid to the host tissue (F CO2_cal ). However, it was buffered by the paracellular transport of DIC from the coelenteron (F DIC_p_cal ). During the daytime, photosynthesis consumed CO 2 from the host tissue and created a CO 2 gradient that draws CO 2 from the calcifying fluid into the host tissue. In the dark, the photosynthesis rate reached zero, the ion transport became inactive, the calcification rate decreased, and DIC cal increased then.
In the coelenteron, DIC also displayed diurnal variations. DIC coe decreased in the daytime due to the paracellular diffusive transport caused by concentration gradients ( Figure 6A). When DIC in the calcifying fluid decreases due to calcification, the DIC concentration gradient generated DIC efflux from the coelenteron to the calcifying fluid. Although this efflux was buffered by  Muscatine et al. (1984). The shaded bars indicate bark periods.
advective flux from the seawater through the mouth and paracellular diffusive efflux from the seawater, DIC in the coelenteron still decreased to mass balance all the fluxes. The pH and TA in the coelenteron decreased due to the removal of protons from the calcifying fluid and their addition to the coelenteron ( Figures 6B, C), which corresponded to the increase of TA and pH in the calcifying fluid.

Coupled Run Results
In the coupled run, all conditions were the same as in the control run; however, the DIC, TA, NH4, and DOC values for the coral polyp model were no longer constant and followed Eqs. P-29-P-32. One of the coupled runs was under an air pCO 2 of 370 ppm. To test the effect of change air pCO 2 , we also did a coupled run under an air pCO 2 of 800 ppm. The change of air pCO 2 will affect the DIC in the seawater (Eqs. (2)-(3)), thus influence the coral DIC uptake from the seawater ( Table 4, Eq. H-1), and advection and/or diffusion exchanges of DIC with seawater ( Table 4, Eq. H-2).
Our numerical experiments to test the influence of pCO 2air on coral photosynthesis, calcification, and carbonate chemistry indicated higher photosynthesis rates and lower calcification rates when pCO 2air equaled 800 ppm (Figure 7, dashed line).
In our model, the saturation state of the surface ocean was controlled by pCO 2air . With higher pCO 2air , more CO 2 was dissolved in seawater, and the excess formed carbonic acid, which disassociated to form H + and HCO − 3 , resulting in an increase in [HCO − 3 ] and a decrease in [CO 2− 3 ]. As [HCO − 3 ] increased, the photosynthesis rate influenced by the HCO − 3 concentration ( Table 3, Eq. S-3) increased. This is consistent with the report that there are enhanced photosynthesis rates under higher pCO 2air in some corals (Herfort et al., 2008;Marubini et al., 2008). As proposed by Furla et al. (2000), with rising pCO 2air , more CO 2 diffuses to the calcification site, [HCO − 3 ] increased and decreases the [CO 2− 3 ], ultimately altering the overall calcification rate ( Figure 7B).
In the coupled run under different pCO 2air scenarios, the changes in DIC and TA in the polyp components are significant (Figure 8). DIC cal and DIC coe reveal the same pattern, both are opposite with the change of calcification curve ( Figures 8A, B), suggesting that calcification significantly contributes to the DIC variability in the calcifying fluid and then through the buffer fluxes due to the concentration gradient to influence the DIC concentration in coelenteron. The increase in DIC in the calcifying fluid under 800ppm pCO 2air condition are associated with decrease in calcification rate ( Figure 8A, black line). We associate the DIC increase in the coelenteron with the increase in diffusion/advection of DIC from the ambient seawater and calcifying fluid due to the concentration gradient.
TA in the calcifying fluid was under the balanced fluxes of calcification, proton transfer, and diffusion (Eq. P-28). The model shows higher magnitude of TA under 800 ppm pCO 2air in the calcifying fluid ( Figure 8C), likely associate with the decline in calcification rate and calcification-induced proton transfer out of the calcifying fluid site. The mode results suggest that calcification contributed more significantly to the TA variability in the calcifying fluid, since when 1 mol of CaCO 3 is produced, 2 mol of TA cal will be consumed. In our simulation, TA in the coelenteron is also increased under 800 ppm pCO 2air ( Figure 8D, black line). Since TA in the coelenteron is under the balanced fluxes of proton transfer, advection, and diffusion, H + transport by Ca-ATPase (F H ) will directly influence TA in the coelenteron. The higher TA in the coelenteron could be associated with decrease in calcification due to ocean acidification and a resultantly decrease in TA consumption by calcification.

DISCUSSION
In this paper, we presented a three-step approach to incorporate processes in the coral symbiont system on different scales (i.e., cellular, organismal, and ecosystem community) into a coupled model to explore different functions of coral symbiosis system in an integrative way. This process-based model specifically explained how the biological processes regulate the biomineralization mechanisms at the cellular level, the importance of ion transport for the calcification processes at the organismal level, and the influence of environmental changes (pCO 2air ) on the coral biochemical processes and carbonate chemistry. When compare with the lab-experiments Falkowski et al., 1984), in the control run under the steady state, the model results indicate good representation of diurnal variations of the critical processes, e.g., photosynthesis, respiration, and carbon translocation in the coral symbiosis; the modeled ion concentration and diurnal variations in different polyp components are broadly consistent with the microelectrodes measurements and boron isotope systematics (Al-Horani et al., 2003;McCulloch et al., 2012;Cai et al., 2016), for example, pH in the calcifying fluid keeps higher than that in the coelenteron, and there are enhanced pH cal when photosynthesis and calcification are highest. This suggests that the model is based on reasonable assumptions and parameterizations.
The model-simulated carbon budget of symbiont photosynthesis (32.25 µg C cm −2 ) is lower than the values published earlier of 57.6 ± 23.1 µg C cm −2 by Houlbrèque et al., 2003) and 164.3 ± 17.3 µg C cm −2 for fed coral by Houlbrèque et al. (2004). We believe that the lower values in our model simulation are associated with the low DIN uptake rate with a maximum value of approximately 1.2 µg N cm −2 d −1 and approximately 0.005 µg N cm −2 d −1 at the end of the model run. Our model showed about two times lower host respiration rate than the observed rate in Falkowski et al. (1984). One explanation could be the low heterotrophic feeding rate in our model (about 0.9 µg C cm −2 in a daily budget), which will result in less C in the host functional pool, thus lower metabolism demand.
Mucus loss in our model simulation was 12.21 µg C cm −2 and was the same order of magnitude as a measured DOC release of 9.6 ± 0.72 µg C cm −2 for Acropora corals by Naumann et al. (2010). As denoted in Figure 6E, the modelled W a was close to the measurement in Cai et al. (2016). The modeled pH and W a in the calcifying fluid were close to the range of that measured using boron isotope systematics in McCulloch et al. (2012) which found that the coral shows species-dependent pH range of 8.4-8.7, and The biomass values are in the units of µg C per unit surface area of the host tissue. with W a of 15-25. In the simulation results, the calcification rate gradually decreases with the decrease of W a , this is consistent with the previous studies that there are both linear and nonlinear relationships between calcification rate and W a (Gattuso et al., 1999;Marubini et al., 2008). In our model, we built a connection between photosynthesis and calcification under the assumption that photosynthesis and calcification shared a common DIC pool in the host cell. Therefore, our model could predict the change in photosynthesis, calcification rate, and polyp carbonate chemistry under different air pCO 2 conditions. Elevated air pCO 2 was shown to enhance the photosynthesis rate due to the increase in  in the ambient seawater for the zooxanthellae uptake. This further implies that higher air pCO 2 could stimulate higher rates of coral calcification due to increased energy produced by photosynthesis to drive Ca-ATPase. However, in our simulations, there was a decreased calcification rate under higher pCO 2 (Figure 7), suggesting that photosynthesis and calcification compete for the internal supply of the host DIC pool. Photosynthesis decreases the internal DIC pool in the host cell and results in less CO 2− 3 available for calcification. These results support the hypothesis that elevated pCO 2 stimulates photosynthesis, resulting in a reduced DIC supply to calcification (Hoegh-Guldberg and Smith, 1989;Marubini and Thake, 1999).
The model we developed is an idealized approach that only considered a 1-D idealized with constant wind stress, current speed, temperature, and salinity simulation. Through this study, we built a frame work by coupling the coral ecosystem with surrounding water column with phytoplankton biogeochemical processes, and can be used to test the feedback between two ecosystems using different air pCO 2 scenarios. However, in the model presented herein, one of the weaknesses is we did not consider the stability of the symbiosis. In a stable corals symbiosis, the density of zooxanthellae and its symbiotic partner remain relatively constant under a given set of environmental conditions, either the growth rates of the zooxanthellae and the host cells are comparable, if different, excess the extra algae may be expelled. Although we have used S/ S max to specify the ability of packing symbionts in the host tissue, we did not consider the surface area of the host tissue increases with size, as well as coral colonies grow by depositing new skeleton and budding additional polyps. Thus, the C R H and C F H may be accumulated. Another weakness is the constant temperature set up in the current version. The temperature forcing function can be updated to investigate whether the model can reproduce the seasonally fluctuating symbiont density in host tissue (Fitt et al., 2000). In the model, we only considered the increased concentration of dissolved CO 2 , which resulted in a decreased carbonate concentration and therefore, a decreased saturation state. Heat is absorbed (endothermic) during the formation of CaCO 3 precipitate, and solubility decreases with increasing temperature (Mucci, 1983). It is also necessary to incorporate the effects of ocean warming on rates of calcification into the model. The breakdown of symbiosis has been related with temperature (Lesser, 1997;Jones and Hoegh-Guldberg, 2001). There is no simple explanation for the cause of bleaching as it is a complicated process where several different mechanisms are likely involved in any given coral, depending on its physiological state and environmental conditions, however, experiments conducted under stable environmental conditions are limited in representing actual coral reef system responses to environmental change. In order to understand and predict the responses of coral reefs to changing environment, improved models consider combined changes in environmental conditions such as light, temperature, and nutrition are necessary.
In this study, based on sophisticated parameterization, we developed a coral symbiosis model system and validated the model accuracy by comparing the results with laboratorial observations. Our model can reconstruct the diurnal cycle of coral metabolism, and give the right response of coral carbonate chemistry under different air pCO 2 scenarios. Our validations are limited, since some of the parameters are based on experimental and observations of short-term changes in a single environmental condition, and it is a challenge to use such parameters to investigate the long-term response of corals to large-scale changes in environmental variables. Therefore, more processed based formulations need to be developed to overcome simplified assumptions of parameters in the future work.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author.

ACKNOWLEDGMENTS
We are grateful to the two anonymous reviewers for the constructive comments that improved the manuscript a lot. This work forms part of a research initiative of UNESCO/IOC-WESTPAC-Sub-Commission on "Coral Reefs under Climate and Anthropogenic Perturbations" (i.e. IOC/WESTPAC-CorReCAP Project). We thank Prof. Yeemin and Prof. Hong for their helpful comments.