Cracking in Nuclear Graphite

Cracking in Nuclear Graphite

Advanced gas cooled reactors (AGRs) have provided a major contribution to the security of the UK’s electricity supply for four decades.  It is important that they continue to operate safely until a new generation of nuclear power stations is built or until alternative baseload power generation technologies become viable. The cores of these reactors are composed of graphite bricks and, as the reactors age, the bricks are subject to oxidation and cracking which may ultimately determine the time when they will have to be decommissioned.  In this article we describe some of the work that has been undertaken over the last ten years to use mathematical modelling to help support the continued safe operation of these reactors.

1 Background

There are seven UK AGR power stations (14 reactors) operated by EDF Energy.  Hunterston B and Hinkley Point B were the first to be commissioned in 1976.

Figure 1: The lattice of graphite bricks (left) and individual bricks (right).

Graphite cores in AGRs consist of around 3,000 interlocking bricks, as shown in Figure 1.  There is a regular lattice with a radial keying system that allows for expansion and contraction.  There are vertical channels (composed of typically 12 layers of graphite bricks) into which fuel stringers and control rods are inserted.  Figure 1 also shows two graphite bricks, the larger of which is a fuel brick that is 900mm tall and 460mm in diameter; the fuel is placed in the brick bore.

Information on the state of the core is obtained at routine inspections, when the reactor is offline (an outage), that take place at least every 18 months. Measurements are made of the shape of selected channels and cameras enable any cracks that can be seen at the bore to be identified. There are fuel handling constraints that limit the number of channels that can be inspected at one time and about 10% of the core is inspected in this way at each outage.

Various crack morphologies have been seen, but those that have the potential to distort the core are full-height axial cracks, and these are the focus of attention here. Distortion of the core could eventually mean that free movement of the fuel and control rods into and out of the reactor cannot be guaranteed; this may eventually limit the reactor safe operational lifetime.

2 Early-life brick cracking

Bore-initiated cracking

When graphite is irradiated in a reactor it initially shrinks, but then later expands.  Graphite close to the fuel brick bore experiences a higher radiation dose than that at the outside of the brick so that initially stresses are compressive on the outside of the brick and in tension on the inside – see Figure 2. Later in life the directions of these stresses reverse.

Figure 2: Early-life stresses in a graphite brick.

The early-life stresses at the bores have led to some cracking of bricks at some of the AGRs. These early-life cracks are termed bore-initiated cracks.

Statistical modelling

There is a lot of variability in the system at different scales and deterministic models cannot predict which bricks will crack.  However, a non-homogeneous Poisson process can be used to model the brick cracking process and the general approach employed has been described [1, 2]. There are similarities with problems considered in reliability modelling, but key features of this problem are:

  • The bricks are not reparable; in many reliability modelling studies multiple failures may occur in sequence with the components being reparable.
  • When a brick is observed to have cracked, the exact time when this occurred is not known; in most reliability modelling studies component failure times are known.
  • The possibility exists that cracking rates could reach a peak and then decline, due to the stress reversal, rather than increasing or decreasing monotonically.
  • New data are obtained periodically, enabling models to be reviewed and updated.

The graphite bricks are assumed to age depending on an appropriate age measure.  An ensemble of models is considered with different underlying assumptions, but experience has shown that models that use an age measure related to radiation dose perform best: a modified version of the reactor burn-up (the total amount of energy produced by the reactor) is employed. The age measure is denoted by t and referred to simply as time. A given brick is assumed to be in one of three states:

  1. State 0, when it is uncracked.
  2. State 1, when it has one crack.
  3. State 2, when it has two (or more) cracks.

The probability of cracking at a given time is expressed in terms of a hazard function: h_0(t) is the cracking rate for transition from stated 0 to state 1 and h_1(t) is the cracking rate for transition from stated 1 to state 2.  If all the population of bricks under consideration have the same hazard functions with a constant ratio between the single and double cracking hazard functions after the cracking hazard commences, then:  

    \begin{equation*}h_0(t) = ug(t), \quad h_1(t) = \nu g(t).\end{equation*}

This simplification makes the calculations more straightforward; experience with modelling bore-initiated cracks has shown that little is gained from employing a more complex formulation. In any case, the number of doubly cracked bricks is small; for example, only four such cracks have been seen to date at the four reactors at Hunterston B and Hinkley Point B.

A subset of the population of bricks under consideration may be subject to a different risk from the rest of the population.  For example, bricks near the edge of the core experience enhanced stresses and are relatively more likely to crack. This is represented by introducing a factor r for this population (taking r = 1 for the main population). There can be multiplicative risk factors for different sub-populations.

The models that have been employed take the following general form:

    \begin{align*}g(t) &= r{\tau}^pq(\tau)I[\tau > 0],\\\tau &= t-t_0.\end{align*}

Here I is an indicator function which is unity when the argument is true and zero otherwise, and t_0 is the time when the risk is assumed to commence. The decay term q(\tau) is employed to enable the hazard function to reduce after the cracking risk peaks.  A simple linear form has been found to be satisfactory:

    \begin{equation*}q(\tau) = (1-\gamma \tau)I[\tau <1/\gamma].\end{equation*}

This form for the decay term means that the risk of cracking eventually becomes zero. Less complex models are obtained for specific choices of the parameters in the general model (for example, by  setting \gamma = 0 to remove the decay term and by choosing p = 0 to specify a constant cracking risk after t = t_0).

Once the hazard functions are defined it is possible to calculate the various transition probabilities  between brick states and other quantities of interest in terms of the cumulative hazard functions defined by:

    \begin{equation*}H_0(t) = \int^t_{-\infty}h_0(\tau)\,\mathrm{d}\tau, \quad  H_1(t) = \int^t_{-\infty}h_1(\tau)\,\mathrm{d}\tau.      \end{equation*}

For example, the probability density function for the time to the first crack is h_0(t)\exp (-H_0(t)) which, in the absence of the decay term, can be written in the standard form for the Weibull density function:

    \begin{equation*}f(x) = \frac{\kappa}{\lambda}{\left(\frac{x}{\lambda}\right)}^{\kappa-1}\exp \left\{-{\left(\frac{x}{\lambda}\right)}^{\kappa}\right\}, \quad x>0,\end{equation*}

where \tau = x, p = \kappa -1 and ur = \kappa/\lambda^{\kappa}.

The log-likelihood function for the model can be expressed analytically. Various measures that represent the trade-off between the goodness of fit to the data (as indicated by the log-likelihood) and model complexity (as indicated by the number of fitted parameters) can be considered, and the AIC measure^1 (see, for example, Davison [3]) has been found to be suitable. Maximum likelihood parameter estimates can be used to provide a single set of model parameters for forward projections of core behaviour and Bayesian updating can be used to provide probability density functions for the parameters. The probability of transitions between brick states can be used to calculate the probability density functions for the numbers of transitions expected to be seen at an inspection.

Models are used prior to inspections to predict the number of cracks that are expected to be seen.   By subsequently comparing predictions with observations the best performing model(s) can be identified.

Example calculations

An example prediction for the number of bore-initiated cracks expected to be seen at an inspection is given in Table 1; this was for an inspection undertaken in 2011. Here the number of new singly cracked bricks (S) is shown on the vertical axis and the number of new doubly cracked bricks (D) on the horizontal axis; if more existing singly cracked bricks become doubly cracked than new singly cracked bricks appear, then S can be negative. For each possible outcome of the inspection a complementary cumulative probability value is given (expressed as a percentage) that indicates the number of possible observations that are calculated to be no more likely than the given observation. For example, the value shown in the cell S3-D2 indicates that the observation of three new singly cracked bricks and two new doubly cracked bricks is at least as likely as 6.44% of the possible outcomes.  The area in the table shaded light green indicates where the observations are expected to be 50% of the time, the area shaded light green or dark green 95% of the time and the area shaded light green, dark green or amber 99% of the time.  The actual observed number of cracks seen at the inspection (0 new double and 0 new singles) is indicated by the bold number.

Table 1: Example predictions for cumulative probability of outcomes for an inspection in 2011.

The use of this traffic light approach helps the understanding of which outcomes would be expected and which surprising at the inspection.  It is important to emphasise that to see all the predictions in the green zone over an extended period would actually be bad news for the model. This would indicate that the variability in the system was being overestimated; the occasional observation in the amber or even the red zones is to be expected. If the models have good predictive performance one would expect to obtain an approximately uniform distribution of the outcome values (this will not be exact because the outcomes are discrete). This has been confirmed over a large number of inspections, providing support for the predictive powers of the models.

There is good evidence that for the reactors that have operated longest, bore-initiated cracking has ceased in the bricks receiving the highest radiation dose (the best fit model parameter \gamma is positive and significantly different from zero).  This means that for these reactors it is possible to forecast with confidence the total numbers of bore-initiated cracks that will be seen over their lifetimes.

3 Brick cracking in later life

Keyway root cracking

From the time when the reactors were designed it had been expected that when the stresses reversed later on in the reactor’s life, so that the compressive stresses on the outside of the brick become tensile, some bricks would undergo cracking from the keyway roots (the base of the cut-outs on the outside of the bricks in Figure 1), with cracks progressing from the outside of the brick towards the bore; such cracks are referred to as keyway root cracks.  To date very few such cracks have been seen.  It is clear, however, that as evidence of the eventual cessation of bore-initiated cracking increases, keyway root cracking will become increasingly important.

As we have seen, for bore-initiated cracks we are now data rich: this is not currently the case for keyway root cracking. It has to be demonstrated that it is safe to operate the reactors between inspections – a different approach has to be employed in this data-poor period.

Stress analysis and brick shape modelling

Stress analysis calculations are employed to describe the expected evolution of brick shapes and keyway root cracking. Such calculations need a number of inputs, including mathematical models for graphite material properties, in particular for dimensional change (the shrinking and later expansion of the graphite as it is irradiated). Figure 3 shows the dimensional change rate data from a materials test reactor experiment plotted against average radiation dose, with ten example samples from a cubic power law model (the units of dose employed, SDU, refers to standard dose unit). The modelling covers the data well and illustrates the variability between different graphite samples.

As graphite bricks are irradiated their shape changes from a simple cylinder to a wheatsheaf, resulting from differences in dimensional change in different parts of the brick due to variations in the amount of irradiation. Various measures of the brick shape can be considered, but here we use the difference between the middle brick diameter and that at the top or bottom of the brick; this distance initially grows, reaches a peak and then starts to reduce again and is referred to as a lambda factor.

The calculated evolution of the lambda factors compares well with reactor data, and it has been confirmed that most of the variability seen in the calculated evolution of the brick shapes derives from variability in the dimensional change model parameters.

Figure 3: Example samples for the unstressed dimensional change model.

Although the evolution of brick shapes is well represented in the stress analysis modelling, the processes that subsequently lead to cracking are not well understood.  A failure criterion related to the maximum in plane principal (MIPP) stress has been employed, but there remains considerable uncertainty over the best criterion to use; future data from the reactors should help clarify this.

Figure 4: Calculated times when lambda factors return to zero and bricks crack at Hunterston B (HNB).

Figure 4 shows a scatter plot of the time when the relevant lambda factor returns to zero and the calculated brick failure times.  The correlation coefficient for each layer is about 0.80, indicative of a strong relationship between the chosen metric for bore shape evolution and the calculated risk of keyway root cracking. This correlation is central to current methods for modelling the future evolution of the system: we are data rich for brick shape evolution measurements but data poor for brick cracking observations.  Using this correlation helps statistical models to make predictions for the number of keyway root cracks that will be seen at inspections.

Modelling long-term evolution

Using the experience gained from modelling the evolution of bore-initiated cracks, the evolution of keyway root cracking can be modelled.  Key aspects of this modelling are:

  • Comparing stress analysis calculations for brick shapes with measurements demonstrates that the key features of system variability are adequately represented.
  • By using the relationship between brick shape evolution and cracking risk, the variation in brick cracking times can be estimated.
  • Simulations of the evolution of keyway root cracking can be undertaken using hazard functions suggested by the stress analysis calculations, conditioned on the observed numbers of such cracks seen to date.
  • The simulations provide an important contribution to assessing the period over which it is safe to operate the reactor before another inspection is made.

At a reactor inspection at Hunterston B in October 2015 three keyway root cracks were observed.  This outcome was not totally surprising; a probability of seeing three or more cracks of around 12% had been calculated. There is some suggestion, however, that the initiation of cracking may be slightly earlier than assumed in the hazard function prior distribution; this will be reassessed when more information becomes available.

Figure 5: Fitted brick shape curves plotted for each observed brick in Hunterston B Layer 4. The brick with a keyway root crack is highlighted in black.

Figure 5 shows an example of the fitted curves for the evolution of the relevant lambda factor in Layer 4 (one of the peak-rated layers) at this inspection. One keyway root crack was seen in this layer, and it can be seen that this brick appears to be ageing faster than most, consistent with the expected correlation between brick shape evolution and keyway root cracking risks.  The variability in the evolution of the lambda factor reflects the variability in brick properties and the environments seen by the bricks; this is in turn reflected in the extended period over which the risk of keyway root cracking will be experienced.

4 Discussion

The modelling of cracking of graphite bricks illustrates how, in a complex system like a nuclear reactor, it is impossible to develop mathematical models that can represent all features of system evolution, partly because of the large number of sources of variability at different scales.  However, all the information that is available from measurements on the system needs to be used in formulating models of how the system is evolving.  With the early-life bore-initiated cracks, initially only very simple models with associated large uncertainties in long-term forecasts could be justified. As more data became available more complex models could be employed and the total number of such cracks that will be present in the reactors by the end of their lives can now be estimated with confidence.

With limited direct information on late-life keyway root cracking, indirect information, such as that for brick shape evolution, needs to be employed in order to represent system variability correctly and to provide the best possible estimates of cracking rates.

At all stages of the process, making predictions of what will be seen ahead of reactor inspections provides confidence in model predictive power and helps to identify where model assumptions need to be modified.

Philip Maul CMath CSci FIMA
Quintessa Ltd

Peter Robinson
Quintessa Ltd

Jenny Burrow
Quintessa Ltd

Alex Bond
Quintessa Ltd


The authors would like to thank EDF Energy for permission to publish this article.


AIC is the Akaike information criterion that provides an indication of the relative goodness of fit of statistical models to a given set of data based on information theory.


  1. Maul, P.R. and Robinson, P.C. (2006) Cracking down, Nucl. Eng. Int., vol. 51, pp. 44–49.
  2. Maul, P.R., Robinson, P.C. and Northrop P.J. (2011)  Statistical modelling of graphite brick cracking in advanced gas cooled reactors, App. Stat., vol. 60, no. 3, pp. 339–353.
  3. Davison, A.C. (2003) Statistical Models,  Cambridge University Press, Cambridge.

Reproduced from Mathematics Today, June 2017

Download the article, Cracking in Nuclear Graphite (pdf)

Image credit: Hinkley Point Nuclear Power Station Somerset, UK by ©  Joseph Golby /


Leave a Reply

Your email address will not be published. Required fields are marked *