Application of Experimental Design to Hydrogen Storage: Optimisation of Lignin-Derived Carbons. C

: Lignin is a signiﬁcant by-product of the paper pulping and biofuel industries. Upgrading lignin to a high-value product is essential for the economic viability of bioreﬁneries for bioethanol production and environmentally benign pulping processes. In this work, the feasibility of lignin-derived activated carbons for hydrogen storage was studied using a Design of Experiments methodology, for a time and cost-e ﬃ cient exploration of the synthesis process. Four factors (carbonisation temperature, activation temperature, carbonisation time, and activation time) were investigated simultaneously. Development of a mathematical model allowed the factors with the greatest impact to be identiﬁed using regression analysis for three responses: surface area, average pore size, and hydrogen uptake at 77 K and 1 bar. Maximising the surface area required activation conditions using the highest settings, however, a low carbonisation temperature was also revealed to be integral to prevent detrimental and excessive pore widening. A small pore size, vital for e ﬃ cient hydrogen uptake, could be achieved by using low carbonisation temperature but also low activation temperatures. An optimum was achieved using the lowest carbonisation conditions (350 ◦ C for 30 min) to retain a smaller pore size, followed by activation under the severest conditions (1000 ◦ C for 60 min) to maximise surface area and hydrogen uptake. These conditions yielded a material with a high surface area of 1400 m 2 g − 1 and hydrogen uptake of 1.9 wt.% at 77 K and 1 bar.


Introduction
Lignin is a key component of lignocellulosic biomass, forming part of the plant cell walls [1]. This natural, highly aromatic biopolymer is a major waste product with over 70 million metric tons produced every year, primarily from the paper pulping industry [2]. Currently, lignin has no high-value application and is predominantly burnt as fuel [2,3]. The future economic viability of biorefineries for bioethanol production and environmentally benign pulping processes will depend on utilising this waste lignin and upgrading to high-value products. The conversion of lignin to activated carbon has received increasing attention due to its high carbon content (approximately 60 wt.%) and wide availability [4][5][6]. The use of lignin as a carbon feedstock could reduce the impact of greenhouse gas emissions associated with rotting or burning, as demonstrated for agroforestry waste [7]. The highly aromatic structure and availability of lignin also make it a potentially sustainable and economical replacement for the crude oil-based phenolic resin precursors widely used to produce high-value activated carbon materials [8,9]. Activated carbons are promising materials for the physisorptive storage of hydrogen, due to their low cost, high surface area, and ease of synthesis in bulk [10]. Activated carbons typically exhibit gravimetric hydrogen capacities of between 0.5 and 3.3 wt.% at 77 K and 1 bar, depending on the feedstock and synthesis conditions [11][12][13]. The preparation of activated carbons requires pyrolysis of a carbonaceous feedstock to form a char with a rudimentary pore structure. A chemical or physical activation procedure is then performed to enhance this porosity. Alternatively, both carbonisation and activation processes can be performed in one-step using a chemical agent [4,[14][15][16]. Physical activation of the char using CO 2 was chosen as the focus of this work, due to the lower environmental impact of the process in comparison to chemical activation [17][18][19][20]. Although production of activated carbons using either a chemical or physical activation method is an energy intensive process, the use of a low-carbon electricity system can cause a significant reduction in the associated global warming potential [19]. In the case of chemical activation, the impregnation of the chemical agent has significant environmental impacts in addition to those of the energy demand [20]. Physical activation with CO 2 has previously been utilised to develop high surface area carbons (1000 to 3400 m 2 g −1 ) with gravimetric hydrogen capacities of up to 2.3 wt.% at 77 K and 1 bar [21][22][23].
It is well established that hydrogen capacity tends to increase with material surface area, due to an increase in available adsorption sites [24,25]. As a general rule, gravimetric uptake of H 2 at 77 K up to 1 bar increases by 1.1 wt.% for every 500 m 2 g −1 of surface area the material has [11]. However, it has been shown that greater hydrogen uptake capacities can be achieved in materials with a large surface area composed of micropores (diameter less than 20 Å). In particular, narrow micropores between 5 and 7 Å allow for more efficient packing of adsorbed hydrogen at 77 K [25][26][27][28] with solid-like densities reported in micropores at temperatures up to 100 K and pressures as low as 0.2 bar [29]. The pore size distribution of a carbon is highly dependent on the processing conditions [17,[30][31][32][33][34]. The impact of multiple factors (such as synthesis temperature, dwell time, and gas flow rate) means the investigation of how these parameters affect material properties, and their optimal levels are both challenging and time-intensive.
Design of Experiments (DoE) is a powerful tool which can be utilised for the optimisation of both carbonisation and activation processes simultaneously. Development of a mathematical model allows for the creation of a map of the material space, as well as the identification of interactions between the experimental factors, an outcome that cannot be achieved with the more traditional one-variable-at-a-time approach [35][36][37][38]. A fractional factorial design, where only half of the experimental runs of the full experimental design are required, can be employed for a time and cost-efficient exploration of the system whilst obtaining an increased level of information in comparison to traditional experimental approaches.
The focus of this work is to show how the combined use of a fractional factorial experimental design with high-throughput measurement methods is a valuable tool for rapid exploration of a new experimental system. An enhanced level of insight into the experimental system can be gained in a much shorter experimental time frame; identifying promising materials and evaluating the importance of chosen factors. A more detailed investigation into the most promising materials can then be performed. We examine the effect of four factors across both the carbonisation and activation processes at two levels; low and high. Lignin carbons were compared to an industrially produced phenolic-resin derived carbon, TE7, therefore the same physical activation agent (CO 2 ) was chosen for consistency. High-throughput screening of the carbon materials employs three responses; Brunauer-Emmett-Teller (BET) surface area, average pore width derived from small-angle x-ray scattering, and hydrogen uptake at 77 K up to 1 bar. A more in-depth gas sorption analysis was then performed on the most promising materials identified by the experimental design.
The use of DoE for activated carbon optimisation has been utilised by several authors [30,31,[39][40][41] primarily for water treatment. One of the few accounts of optimising lignin derived materials for C 2019, 5, 82 3 of 17 hydrogen storage was published by Cheng et al. [42] who produced materials from hydrolytic lignin using chemical and physical activation, exhibiting hydrogen uptake of 2.0 to 2.6 wt.% at 77 K and 1 bar; however, further published research remains limited. To the authors' knowledge, this is the first use of the DoE method for the optimisation of lignin materials for hydrogen storage.

Activated Carbon Preparation
Carbon materials were prepared from an organosolv lignin feedstock extracted from pine, supplied by MAST Carbon International (Basingstoke, UK). The lignin was used as provided, as a fine homogeneous powder. Carbonisation and activation of the lignin material was performed in a horizontal tube furnace (Lenton, Hope, UK, LTF 12/75/750) with a separate ceramic work tube and gas-tight end seals, fitted with a digital PID controller (Eurotherm, Worthing, UK, 3216CC). Mass flow controllers (Bronkhorst, Newmarket, UK, EL-FLOW) calibrated with a correlated flow meter (Cole Parmer, Vernon Hills, IL, USA, 32044-14) controlled the flow rate of gas into the furnace. PureShield argon (BOC, Guildford, UK, 99.998 % purity) with a flow rate of 500 mL min −1 was used as the inert furnace gas and a heating rate of 5 • C min −1 used throughout.
Lignin (3 g) was placed in an alumina combustion boat (20 mL capacity) and inserted into the centre of the furnace. The furnace was purged with argon for 30 min at room temperature, then heated to 105 • C and held for 30 min. The sample was then heated from 105 • C to the desired carbonisation temperature (between 350 and 900 • C), where it was held for the specified carbonisation dwell time (30 to 360 min) to produce a carbon char. The system was cooled to 350 • C at a rate of 1 • C min −1 , then without rate control to room temperature. The carbon char was ground with a pestle and mortar and passed through a 35-mesh (500 µm) sieve. The ground char (1 g) was loaded into the furnace and heated under argon, as described previously, to the desired activation temperature (800 to 1000 • C). The furnace gas was switched to carbon dioxide (BOC, Guilford, UK, 99.998% purity) at 100 mL min −1 for the specified activation time (20 to 60 min) and then cooled in an argon atmosphere, as specified above.

Experimental Design
The optimal carbonisation and activation conditions were studied with a fractional factorial design (Resolution IV) based on four factors: carbonisation temperature (CTemp), carbonisation dwell time (CTime), activation temperature (ATemp), and activation dwell time (ATime). These factors were selected based on previous literature study specifically identifying factors with the most significant influence on adsorptive properties. Each factor was studied at a high (+1), low (−1), and midpoint (0) level, outlined in Table 1. The ranges were determined from values in the literature, except for activation dwell time which was chosen to ensure a range of burn-off levels, an accepted measure of the extent of activation [14,43,44], were considered and was determined using separate thermogravimetric experiments (see Figures S1 and S2 in the Supplementary Information). The factorial design is illustrated graphically in Figure S3 in the Supplementary Information. The minimum number of experiments for a two-level, four-factor design is 2 4 = 16. However, from only half of these experiments, the main effects and two-order interactions between the factors can be obtained, allowing the study of four factors at two levels in just 8 experimental runs. Consequently, all factors were varied simultaneously over 11 experiments: 8 variations plus 3 repeats; the experimental matrix is shown in Table 2. Individual experiments, and their resulting carbons, are referred to throughout the following as Exp X, where X corresponds to the experiment number listed in Table 2. The effect of the combinations of the four factors on three responses was determined: BET surface area (S BET ), average pore width derived from small-angle scattering (w avg : SAXS), and the hydrogen uptake at 77 K and 1 bar. The responses were chosen to allow efficient high-throughput screening of materials. A narrow pore size is important for hydrogen uptake and gas sorption is frequently used to investigate the pore size distribution of nanoporous materials. However, the collection of a full isotherm is time-intensive and unsuitable for screening large numbers of samples. The BET area is a good predictor of trends in hydrogen uptake within the same group of materials and has been shown to be a useful tool for approximating the number of adsorption sites [24,25]. Several authors have also previously shown a good correlation between the average pore size derived from SAXS and gas sorption [45][46][47]. Consequently, the BET area is used in the experimental design as a screening tool and SAXS was chosen to determine the average pore width of carbon materials.
Data evaluation was performed using the software MODDE Pro Version 11.0.2.2309 (Umetrics, Umeå, Sweden). The influence of factors was determined by fitting of a second-order polynomial model using partial least squares regression analysis [31,36]: where y is the response, x terms are the factors, β 0 is a constant term, β terms are the model parameters and ε is the residual response variation not explained by the model [31,35,36]. Determining numerical values for model parameters (β) established the impact of each factor on the measured response, and how said factors combined in influencing the response overall to determine the optimal process conditions [35,36].

Activated Carbon Characterisation
The surface area of chars and the adsorptive properties of select activated carbons were characterised by gas sorption analysis using a 3Flex analyser from Micromeritics (Norcross, GA, USA). Nitrogen (BOC, Guildford, UK, 99.999%) at 77 K or carbon dioxide (BOC, Guildford, UK, 99.999%) at 273 K was used as the adsorptive gas. Chars and activated carbons (50-100 mg) were degassed at 250 • C or 300 • C, respectively, for 4 h at dynamic 10 −7 mbar vacuum. Isotherm data were collected up to 1 bar using filler rods from p/p o 10 −8 -1 for N 2 and p/p o 10 −5 -10 −2 for CO 2 . The p/p o range for CO 2 is limited in comparison to N 2 due to its higher saturation pressure at 273 K [48,49]. Helium (BOC, Guildford, UK, 99.999%) was used to calculate the free space following isotherm collection. The BET surface area was calculated according to the procedure from ISO 9277:2010 [50]. Due to the highly microporous nature of samples, the Rouquerol correction was applied with the range p/p o 0.01-0.03 and used to calculate the BET surface area. The performance of lignin-derived activated carbons and an industrially produced carbon TE7 (MAST Carbon International, Basingstoke, UK) was characterised by sorption of hydrogen (BOC, Guildford, UK, 99.999%) at 77 K. An identical procedure to the above was followed, except isotherms were collected using absolute pressure dosing from 5 to 900 mbar.
The pore size of activated carbons was screened using small-angle X-ray scattering (SAXS). SAXS measurements were performed at the University of Bristol using a Ganesha 300XL system (SAXSLAB, Copenhagen, Denmark) under vacuum to provide information on the average pore sizes. The instrument uses Cu K-α radiation (λ = 1.54 Å) and a two-dimensional Pilatus X-ray Detector. Activated carbon samples were loaded into 1.5 mm borosilicate-glass capillary tubes (Capillary Tube Supplies, Withiel, UK). Three 10 min measurements in different configurations were performed at room temperature on each sample to collect data over the range q = 0.004-2.6 Å −1 . Data processing and background subtractions were performed using the program SAXSGUI by SAXSLAB. The three SAXS patterns were merged post-data processing using PRIMUS [51].
The low q part of the SAXS curve is often approximated by a single size parameter, the radius of gyration R g [52]. For spherical pores with homogeneous density, the average pore width (w avg ) can be calculated by [53].
Following the low q range examined by the Guinier approximation, higher values of q contain information regarding the nanoparticle shape. The scattering intensity falls off, in a region termed the "power law regime", according to [52]: where n is the fractal dimension and is dependent on the dimensions of the scattering particle. The scattering pattern oscillates in a manner characteristic of the particle shape (or form), known as the form factor P(q) [53]. Theoretical SAXS curves can be fitted to the collected pattern in order to describe the size and shape of scattering objects, and the interactions between them [52]. Scattering in activated carbons was first described by Kalliat et al. [54] and the model was later modified by Gibaud et al. [55]. The scattering intensity I(q) is related to the total surface area of macro A and micropores B, plus a constant background C: To examine the surface morphology of lignin chars, scanning electron micrographs (SEMs) were collected at the University of Bristol (JEOL UK, London, UK, IT300) in secondary electron mode at an accelerating voltage of 10-15 keV and a working distance of 10 mm.

Influence of Factors on Lignin-Derived Carbon Properties
In-depth characterisation of the organosolv extracted lignin feedstock was reported in the authors' previous work [56]. The completed experimental design and measured responses are shown in Table 2. The raw data was evaluated using replicate plots ( Figure S4 in the Supplementary Information). The measured variation between the three replicates was much smaller than the overall variation in the data set, indicating a small replicate error, which will not affect the results of regression analysis. Normally distributed data is advantageous for regression analysis, improving the efficiency of data analysis and enhancing model validity [35]. To this end, a logarithmic transformation was applied to the surface area response. No transformation was required for the average pore size response, and transformation did not affect the hydrogen uptake response (data distributions before and after transformation are included in Figure S5 in the Supplementary Information).

Activated Carbon Surface Area
The surface area of each activated carbon was determined by application of the BET equation to partial nitrogen isotherms. Repeat measurements for selected samples indicated a small measurement uncertainty for BET surface area of up to ±40 m 2 g −1 (isotherms and BET fitting parameters are shown in Figure S6 and Table S1 in the Supplementary Information). Regression coefficients, showing the influence of factors on the activated carbon surface area, are illustrated in Figure 1a. Four factors were determined to influence the model validity: CTemp, ATemp, ATime, and the squared term ATemp*ATemp. Activation conditions were shown to have the greatest influence on the development of the BET surface area with the activation temperature and time having a positive effect i.e., increasing either led to an improved surface area. Increasing the activation temperature has been shown to cause an exponential increase in the rate of reaction between the carbon char and CO 2 [15,57,58] as higher temperatures lead to a larger release of volatile matter from the char. Similarly, increasing activation time allows for the further reaction of CO 2 with the char, also resulting in an increased emission of volatile matter. Both effects promote the creation of new porosity, leading to a higher burn off of the char and improving the surface area of activated carbons [17,18,30,31,34,[59][60][61]. The positive squared term ATemp*ATemp implies a nonlinear relationship between surface area and activation temperature. This relationship is discussed in detail in Section 3.1.2; in brief, it is likely related to the rate of gasification.
Carbonisation temperature is shown in Figure 1a to have a negative effect; increasing the temperature leads to a reduction in BET surface area. Previous literature [32,33,59] has attributed such a negative effect to the formation of an "intermediate melt" at high carbonisation temperatures; secondary devolatilisation and redeposition of high molecular weight volatiles which can obstruct the passage of activation gas into the rudimentary porous structure of the char, adversely affecting activated carbon burn-off and adsorptive properties. The presence of an intermediate melt was examined for a high-and low-temperature lignin char, carbonised at 900 • C or 350 • C for 30 min (see Table 2, Exp 5 & Exp 6), using scanning electron microscopy and gas sorption.
No evidence of an intermediate melt was observed in SEM micrographs ( Figure S7 in the Supplementary Information). N 2 BET surface area results of the 900 • C and 350 • C chars (see Figure S8 in the Supplementary Information) also showed no presence of an intermediate melt; the high-temperature char has a larger BET surface area (150 m 2 g −1 ) than the low-temperature char (<1 m 2 g −1 ). Diffusion limitations of N 2 into ultra-micropores (pores below 7 Å in width) at cryogenic temperatures are well-known phenomena [48,49]. Porosity inaccessible to N 2 was examined using CO 2 isotherms at 273 K ( Figure 2), which corroborate the BET surface area results: the high-temperature char C 2019, 5, 82 7 of 17 exhibits increased CO 2 uptake in comparison to the low-temperature char, with a maximum uptake of 70 cm 3 g −1 and 15 cm 3 g −1 , respectively, indicating the high-temperature char has an enhanced porous network, unhindered by an intermediate melt.
the development of the BET surface area with the activation temperature and time having a positive effect i.e., increasing either led to an improved surface area. Increasing the activation temperature has been shown to cause an exponential increase in the rate of reaction between the carbon char and CO2 [15,57,58] as higher temperatures lead to a larger release of volatile matter from the char. Similarly, increasing activation time allows for the further reaction of CO2 with the char, also resulting in an increased emission of volatile matter. Both effects promote the creation of new porosity, leading to a higher burn off of the char and improving the surface area of activated carbons [17,18,30,31,34,[59][60][61]. The positive squared term ATemp*ATemp implies a nonlinear relationship between surface area and activation temperature. This relationship is discussed in detail in Section 3.1.2; in brief, it is likely related to the rate of gasification. Carbonisation temperature is shown in Figure 1a to have a negative effect; increasing the temperature leads to a reduction in BET surface area. Previous literature [32,33,59] has attributed such a negative effect to the formation of an "intermediate melt" at high carbonisation temperatures; secondary devolatilisation and redeposition of high molecular weight volatiles which can obstruct the passage of activation gas into the rudimentary porous structure of the char, adversely affecting activated carbon burn-off and adsorptive properties. The presence of an intermediate melt was examined for a high-and low-temperature lignin char, carbonised at 900 °C or 350 °C for 30 min (see Table 2, Exp 5 & Exp 6), using scanning electron microscopy and gas sorption.
No evidence of an intermediate melt was observed in SEM micrographs ( Figure S7 in the Supplementary Information). N2 BET surface area results of the 900 °C and 350 °C chars (see Figure  S8 in the Supplementary Information) also showed no presence of an intermediate melt; the hightemperature char has a larger BET surface area (150 m 2 g -1 ) than the low-temperature char (<1 m 2 g −1 ). Diffusion limitations of N2 into ultra-micropores (pores below 7 Å in width) at cryogenic temperatures are well-known phenomena [48,49]. Porosity inaccessible to N2 was examined using CO2 isotherms at 273 K (Figure 2), which corroborate the BET surface area results: the hightemperature char exhibits increased CO2 uptake in comparison to the low-temperature char, with a maximum uptake of 70 cm 3 g −1 and 15 cm 3 g −1 , respectively, indicating the high-temperature char has an enhanced porous network, unhindered by an intermediate melt. The negative effect of carbonisation temperature, consequently, likely results from the widening of existing porosity in the high-temperature char. On activation, this leads to an increased formation of meso-and micropores, which have a reduced contribution to surface area compared to micropores. Other authors [62,63] have also noted materials carbonised at low temperatures yield improved surface area in comparison to chars prepared at high temperature. This was linked to a larger micropore volume in the low-temperature chars, and an increased meso-and macropore volume in the high-temperature chars.
Response-contour surface plots can be used to gain a better understanding of how the critical factors identified by regression analysis affect the modelled system by giving a visual representation. The surface area response-contour plot is shown in Figure 3. The model predicts a linear response The negative effect of carbonisation temperature, consequently, likely results from the widening of existing porosity in the high-temperature char. On activation, this leads to an increased formation of meso-and micropores, which have a reduced contribution to surface area compared to micropores. Other authors [62,63] have also noted materials carbonised at low temperatures yield improved surface area in comparison to chars prepared at high temperature. This was linked to a larger micropore volume in the low-temperature chars, and an increased meso-and macropore volume in the high-temperature chars.
Response-contour surface plots can be used to gain a better understanding of how the critical factors identified by regression analysis affect the modelled system by giving a visual representation. The surface area response-contour plot is shown in Figure 3. The model predicts a linear response between all factors and the surface area, which is reflected in the response plot. In summary, materials with the highest surface area are obtained using extended activation dwell times, low carbonisation temperatures, and high activation temperatures.

Activated Carbon Average Pore Size
The activated carbon average pore size was screened using SAXS; patterns and fitting parameters are included in Figure S9 and Table S2 in the Supplementary Information. Regression coefficients, illustrated in Figure 1b, indicate carbonisation and activation temperature have the greatest influence on average pore size. Four factors were determined to influence model validity; CTemp, ATemp, the squared term ATemp*ATemp, and the interaction term CTemp*ATemp. The positive influence of carbonisation temperature on activated carbon pore size is consistent with the negative effect this factor exhibits on the surface area (Section 3.1.1). Increased carbonisation temperature leads to an enhanced porous network in the resulting char, as demonstrated by CO2 isotherms (Figure 2) of the 350 °C and 900 °C lignin chars. This porous network is enhanced in the activation step, widening porosity and broadening the average pore size [32,33,63]. Conversely, because lowtemperature chars have a less extensive porous network, the pore widening effect is less dominant on activation, instead favouring the creation of new porosity and resulting in a smaller average pore size [30,33,[63][64][65]. A similar effect, well-reported in the literature, occurs at high activation temperatures: the processes of pore widening and pore wall collapse are more dominant than the development of new micropores at higher activation temperatures, leading to a broadening of the average pore width. At low activation temperatures, the reverse has been observed, with pore creation the dominant process leading to a narrower pore width [17,18,30,31,59,60,66].
For average activated carbon pore size, regression analysis has highlighted the two-factor interaction term CTemp*ATemp as significant. A two-factor interaction is when the impact of one factor is dependent on the level of a second factor [35]. The negative CTemp*ATemp interaction implies broadening of the average pore size is greater at both high carbonisation temperature and high activation temperature. The effect on the system is demonstrated by the response-contour surface plot in Figure 4; the largest pore widths are obtained when these parameters are at their highest levels. This interaction is likely related to the negative effect of carbonisation temperature on activated carbon surface area. High-temperature chars have a more extensive porous network and larger pores

Activated Carbon Average Pore Size
The activated carbon average pore size was screened using SAXS; patterns and fitting parameters are included in Figure S9 and Table S2 in the Supplementary Information. Regression coefficients, illustrated in Figure 1b, indicate carbonisation and activation temperature have the greatest influence on average pore size. Four factors were determined to influence model validity; CTemp, ATemp, the squared term ATemp*ATemp, and the interaction term CTemp*ATemp. The positive influence of carbonisation temperature on activated carbon pore size is consistent with the negative effect this factor exhibits on the surface area (Section 3.1.1). Increased carbonisation temperature leads to an enhanced porous network in the resulting char, as demonstrated by CO 2 isotherms (Figure 2) of the 350 • C and 900 • C lignin chars. This porous network is enhanced in the activation step, widening porosity and broadening the average pore size [32,33,63]. Conversely, because low-temperature chars have a less extensive porous network, the pore widening effect is less dominant on activation, instead favouring the creation of new porosity and resulting in a smaller average pore size [30,33,[63][64][65]. A similar effect, well-reported in the literature, occurs at high activation temperatures: the processes of pore widening and pore wall collapse are more dominant than the development of new micropores at higher activation temperatures, leading to a broadening of the average pore width. At low activation temperatures, the reverse has been observed, with pore creation the dominant process leading to a narrower pore width [17,18,30,31,59,60,66].
For average activated carbon pore size, regression analysis has highlighted the two-factor interaction term CTemp*ATemp as significant. A two-factor interaction is when the impact of one factor is dependent on the level of a second factor [35]. The negative CTemp*ATemp interaction implies broadening of the average pore size is greater at both high carbonisation temperature and high activation temperature. The effect on the system is demonstrated by the response-contour surface plot in Figure 4; the largest pore widths are obtained when these parameters are at their highest levels. This interaction is likely related to the negative effect of carbonisation temperature on activated carbon surface area. High-temperature chars have a more extensive porous network and larger pores than low-temperature chars. Activation increases this pore size further with higher temperatures contributing further to a broader pore size distribution [17].  The squared term ATemp*ATemp (first introduced in Section 3.1.1) in both the average pore size and BET surface area models exhibits large confidence intervals which cross the x-axis, indicating this term is statistically insignificant. However, removal of this term from either model is detrimental to the model validity, reducing the Q 2 value. This behaviour indicates the term is important; however, the current model does not entirely account for the observed behaviour. ATemp*ATemp implies a positive, nonlinear relationship for both responses with the activation temperature. This relationship is likely related to the rate of CO2 activation, which is temperature-dependent, and known to increase exponentially below 1500 °C [15, 57,58]. During activation, oxygen chemisorbs to active sites on the char surface. The reaction rate is proportional to the concentration of occupied sites, which is, in turn, dependent on the activation temperature and reactivity of the carbon surface [15,57]. Consequently, we theorise that burn-off, directly affecting pore size and surface area, could show the same exponential relationship as the reaction rate, although a more complex model is required to investigate this further.
The absence of carbonisation dwell time as a significant factor in Figure 1 indicates that CTime has little effect on either the activated carbon surface area or average pore width in this system. Previous research has highlighted the carbonisation dwell time as an important factor [30,33,63]. In this work, the irrelevance of CTime can be attributed to the small sample size (3 g), and relatively slow ramp rate (5 °C min −1 ) used throughout the experiments. The ramp rate is restricted to avoid cracking of the ceramic work tube in the furnace. Volatile emission occurs over a longer period at 5 °C min −1 in comparison to a rapid or flash pyrolysis (for example, at 20 °C min −1 ). Consequently, the ramp rate can be considered an extension of the dwell time. If complete emission of volatiles from the char (devolatilisation) occurs in this initial, slow ramping period, then the subsequent carbonisation dwell time becomes less significant.

Activated Carbon Hydrogen Uptake
The influence of factors on the hydrogen uptake at 77 K is illustrated in Figure 1c. Activation temperature and dwell time appear to have the greatest influence on activated carbon performance. Inclusion of the linear CTemp and squared ATemp*ATemp terms improves the model validity, indicating they are also significant. The influencing factors are akin to those affecting BET surface The squared term ATemp*ATemp (first introduced in Section 3.1.1) in both the average pore size and BET surface area models exhibits large confidence intervals which cross the x-axis, indicating this term is statistically insignificant. However, removal of this term from either model is detrimental to the model validity, reducing the Q 2 value. This behaviour indicates the term is important; however, the current model does not entirely account for the observed behaviour. ATemp*ATemp implies a positive, nonlinear relationship for both responses with the activation temperature. This relationship is likely related to the rate of CO 2 activation, which is temperature-dependent, and known to increase exponentially below 1500 • C [15, 57,58]. During activation, oxygen chemisorbs to active sites on the char surface. The reaction rate is proportional to the concentration of occupied sites, which is, in turn, dependent on the activation temperature and reactivity of the carbon surface [15,57]. Consequently, we theorise that burn-off, directly affecting pore size and surface area, could show the same exponential relationship as the reaction rate, although a more complex model is required to investigate this further.
The absence of carbonisation dwell time as a significant factor in Figure 1 indicates that CTime has little effect on either the activated carbon surface area or average pore width in this system. Previous research has highlighted the carbonisation dwell time as an important factor [30,33,63]. In this work, the irrelevance of CTime can be attributed to the small sample size (3 g), and relatively slow ramp rate (5 • C min −1 ) used throughout the experiments. The ramp rate is restricted to avoid cracking of the ceramic work tube in the furnace. Volatile emission occurs over a longer period at 5 • C min −1 in comparison to a rapid or flash pyrolysis (for example, at 20 • C min −1 ). Consequently, the ramp rate can be considered an extension of the dwell time. If complete emission of volatiles from the char (devolatilisation) occurs in this initial, slow ramping period, then the subsequent carbonisation dwell time becomes less significant.

Activated Carbon Hydrogen Uptake
The influence of factors on the hydrogen uptake at 77 K is illustrated in Figure 1c. Activation temperature and dwell time appear to have the greatest influence on activated carbon performance. Inclusion of the linear CTemp and squared ATemp*ATemp terms improves the model validity, indicating they are also significant. The influencing factors are akin to those affecting BET surface area, indicating for this system a large surface area is important for hydrogen uptake. The interplay between the three responses and the factors with the most significant influence (carbonisation and activation temperature) can be examined in Figure 5. Clearly, the hydrogen uptake and surface area responses exhibit analogous behaviour: maximum uptake is achieved in the material with the highest surface area. The adsorptive properties of the best-performing activated carbons (Exp 5 and Exp 8 in Table 2), with the largest BET surface area and H 2 uptake, identified from the experimental design are examined in more detail using N 2 sorption in Section 3.2.  Table 2), with the largest BET surface area and H2 uptake, identified from the experimental design are examined in more detail using N2 sorption in Section 3.2.

Analysis of Model Fit
Model validity was assessed using R 2 and Q 2 parameters. Measured between 0 and 1, the parameters should be high and preferably separated by no more than 0.3 [35]. The agreement between model and experimental values for all responses was examined and is shown in Figure 6, with the parameters of ANOVA analysis listed in Table 3. For all responses, the experimental points are clustered close to the line predicted by the model, resulting in high R 2 and Q 2 values. The regression coefficient plots in Figure 1 show several terms in the pore size and hydrogen uptake models have large confidence intervals that cross the x-axis, indicating they are statistically insignificant. However, removal of these terms is detrimental to the model, indicating that while important, the model does not entirely account for the observed behaviour, leading to lower Q 2 values. ANOVA analysis, however, indicates all models are statistically significant at the 95% confidence level (p < 0.05, for all models p ≤ 0.01), and the pore size model shows no lack of fit (p > 0.05, for this model p = 0.17).

Analysis of Model Fit
Model validity was assessed using R 2 and Q 2 parameters. Measured between 0 and 1, the parameters should be high and preferably separated by no more than 0.3 [35]. The agreement between model and experimental values for all responses was examined and is shown in Figure 6, with the parameters of ANOVA analysis listed in Table 3. For all responses, the experimental points are clustered close to the line predicted by the model, resulting in high R 2 and Q 2 values. The regression coefficient plots in Figure 1 show several terms in the pore size and hydrogen uptake models have large confidence intervals that cross the x-axis, indicating they are statistically insignificant. However, removal of these terms is detrimental to the model, indicating that while important, the model does not entirely account for the observed behaviour, leading to lower Q 2 values. ANOVA analysis, however, indicates all models are statistically significant at the 95% confidence level (p < 0.05, for all models p ≤ 0.01), and the pore size model shows no lack of fit (p > 0.05, for this model p = 0.17). of fit is artificial; with a true lack of fit, both values would be small. Instead, this is caused by the excellent replicates obtained for surface area (512-518 m 2 g −1 ) and hydrogen uptake (1.28-1.29 wt.%). Altering one of these values so there is more considerable variation between replicates (by 50 m 2 g −1 or 0.1 wt.% for example) yields models with no lack of fit (p = ~0.2).  Table 2. The summary of fit parameters R 2 and Q 2 for each plot are listed in Table 3.

Adsorptive Characteristics of Optimal Activated Carbon
The processing conditions for Exp. 5 and Exp. 8 were identified as producing the bestperforming materials; activated carbons exhibited the highest H2 uptake (1.9 wt.% and 1.6 wt.%, respectively, at 77 K and 1 bar) and BET surface area (1409 m 2 g −1 and 1055 m 2 g −1 , respectively). The adsorptive properties of these materials were examined in more detail using N2 sorption. Both carbons exhibited an isotherm with typical Type I behaviour; adsorption primarily occurred at p/po < 0.1 ( Figure S10 in the Supplementary Information).
The adsorptive properties of the carbon are displayed in Figure 7a and are at the higher end of those reported in the literature. Previous studies of lignin carbons from physical activation have reported BET surface areas of 46 to 1853 m 2 g −1 and pore volumes of 0.01 to 0.83 cm 3 g −1 [62,[67][68][69][70][71]. Both carbons are highly microporous, with micropores contributing 84 % (Exp. 5) and 93 % (Exp. 8) to the total pore volume. The average pore size predicted by gas sorption (wavg: GSA) is smaller than the value derived from SAXS analysis (wavg: SAXS); 6.6 Å and 9.3 Å respectively for Exp. 5 and 7.4 Å and 10.0 Å respectively for Exp. 8. This small variation is due to differences in the measurement method; SAXS measures all porosity, whereas gas sorption only considers pores open to nitrogen. For this reason, although SAXS has proved efficient for high throughput screening, the volumeweighted average pore size derived from gas sorption is more applicable to H2 uptake capacity, because it considers only accessible porosity.   Table 2. The summary of fit parameters R 2 and Q 2 for each plot are listed in Table 3. ANOVA results indicate the H 2 uptake and surface area models do exhibit a lack of fit (p < 0.05, for these models p ≤ 0.01). However, the high R 2 and Q 2 values for these models shows that this lack of fit is artificial; with a true lack of fit, both values would be small. Instead, this is caused by the excellent replicates obtained for surface area (512-518 m 2 g −1 ) and hydrogen uptake (1.28-1.29 wt.%). Altering one of these values so there is more considerable variation between replicates (by 50 m 2 g −1 or 0.1 wt.% for example) yields models with no lack of fit (p =~0.2).

Adsorptive Characteristics of Optimal Activated Carbon
The processing conditions for Exp. 5 and Exp. 8 were identified as producing the best-performing materials; activated carbons exhibited the highest H 2 uptake (1.9 wt.% and 1.6 wt.%, respectively, at 77 K and 1 bar) and BET surface area (1409 m 2 g −1 and 1055 m 2 g −1 , respectively). The adsorptive properties of these materials were examined in more detail using N 2 sorption. Both carbons exhibited an isotherm with typical Type I behaviour; adsorption primarily occurred at p/p o < 0.1 ( Figure S10 in the Supplementary Information).
The adsorptive properties of the carbon are displayed in Figure 7a and are at the higher end of those reported in the literature. Previous studies of lignin carbons from physical activation have reported BET surface areas of 46 to 1853 m 2 g −1 and pore volumes of 0.01 to 0.83 cm 3 g −1 [62,[67][68][69][70][71]. Both carbons are highly microporous, with micropores contributing 84 % (Exp. 5) and 93 % (Exp. 8) to the total pore volume. The average pore size predicted by gas sorption (w avg : GSA) is smaller than the value derived from SAXS analysis (w avg : SAXS); 6.6 Å and 9.3 Å respectively for Exp. 5 and 7.4 Å and 10.0 Å respectively for Exp. 8. This small variation is due to differences in the measurement method; SAXS measures all porosity, whereas gas sorption only considers pores open to nitrogen. For this reason, although SAXS has proved efficient for high throughput screening, the volume-weighted average pore size derived from gas sorption is more applicable to H 2 uptake capacity, because it considers only accessible porosity.
A pore size distribution, derived by the fitting of the 2D-NLDFT model to the N 2 isotherm, is shown in Figure 7b. The pore size distribution for both the Exp. 5 and Exp. 8 carbons is centred around pores 5-6 Å in diameter, with further, smaller contributions to pore volume from porosity approximately 8 Å and 12 Å in diameter. However, widening of porosity in the Exp. 8 carbon is evident. The pore size distribution is shifted to slightly wider pores, which appears to have adversely affected the adsorptive properties (illustrated in Figure 7a), resulting in reduced BET surface area and total pore volume in comparison to Exp. 5. A pore size distribution, derived by the fitting of the 2D-NLDFT model to the N2 isotherm, is shown in Figure 7b. The pore size distribution for both the Exp. 5 and Exp. 8 carbons is centred around pores 5-6 Å in diameter, with further, smaller contributions to pore volume from porosity approximately 8 Å and 12 Å in diameter. However, widening of porosity in the Exp. 8 carbon is evident. The pore size distribution is shifted to slightly wider pores, which appears to have adversely affected the adsorptive properties (illustrated in Figure 7a), resulting in reduced BET surface area and total pore volume in comparison to Exp. 5. , the total (Vtot) and micropore (Vmicro) pore volumes; the BET (SBET) and micropore (Smicro) surface areas; and the volume-weighted average pore size (wavg: GSA) derived from (b) the pore size distribution, which was calculated by fitting of the 2D-NLDFT model.

Hydrogen Uptake of Optimal Activated Carbons
Hydrogen isotherms for the Exp 5 and Exp 8 activated carbons are shown in Figure 8a. For comparison, a well-characterised and industrially produced activated carbon, TE7, was also investigated [29,72,73]. TE7 is derived from crude oil-based phenolic resin; there is great interest in using the highly aromatic lignin as a cheap and sustainable replacement feedstock. There is minimal hysteresis between the adsorption and desorption branches of all isotherms, indicating adsorption was completely reversible. The lignin carbons exhibit a gravimetric hydrogen uptake of 1.6-1.9 wt.% and TE7 exhibits an uptake of 2.0 wt.%, comparable to previously reported values for this material [73]. , the total (V tot ) and micropore (V micro ) pore volumes; the BET (S BET ) and micropore (S micro ) surface areas; and the volume-weighted average pore size (w avg : GSA) derived from (b) the pore size distribution, which was calculated by fitting of the 2D-NLDFT model.

Hydrogen Uptake of Optimal Activated Carbons
Hydrogen isotherms for the Exp 5 and Exp 8 activated carbons are shown in Figure 8a. For comparison, a well-characterised and industrially produced activated carbon, TE7, was also investigated [29,72,73]. TE7 is derived from crude oil-based phenolic resin; there is great interest in using the highly aromatic lignin as a cheap and sustainable replacement feedstock. There is minimal hysteresis between the adsorption and desorption branches of all isotherms, indicating adsorption was completely reversible. The lignin carbons exhibit a gravimetric hydrogen uptake of 1.6-1.9 wt.% and TE7 exhibits an uptake of 2.0 wt.%, comparable to previously reported values for this material [73].
The optimum pore size for efficient storage of hydrogen in micropores is between 5 and 7 Å [25][26][27][28]. For Exp. 5 and Exp. 8, micropores below 7 Å account for around 70 % and 60 % of the total pore volume, respectively. The increased H 2 uptake observed for Exp. 5 (1.9 wt.%) can be attributed to the increased surface area and quantity of micropores in the material. The efficiency of these narrow micropores can be seen in Figure 8b; the Exp. 5 carbon prepared in this work outperforms materials in the literature with a higher surface area. The hydrogen uptake is comparable to the uptake of hydrolytic lignin carbons prepared by Cheng et al. (2000 m 2 g −1 , 2.0 wt.%) [42] despite a lower surface area. Importantly, the hydrogen uptake of Exp. 5 is on par with the industrially produced TE7, indicating that lignin could be a promising replacement for phenolic resin precursors. This work has shown organosolv lignin can be upgraded to a high value activated carbon, producing a material with a high surface area and reasonable hydrogen uptake at 77 K up to 1 bar.  [21][22][23]) biomass-derived carbons from chemical activation (data from [12,13,22,[74][75][76][77][78]) other lignin-derived carbons from chemical or physical activation (data from [39,42]) carbide-derived carbons (data from [79]) and activated carbide-derived carbons (data from [79]). The x and y error bars represent the standard deviation between repeat BET surface area and H2 uptake measurements, respectively.
The optimum pore size for efficient storage of hydrogen in micropores is between 5 and 7 Å [25][26][27][28]. For Exp. 5 and Exp. 8, micropores below 7 Å account for around 70 % and 60 % of the total pore volume, respectively. The increased H2 uptake observed for Exp. 5 (1.9 wt.%) can be attributed to the increased surface area and quantity of micropores in the material. The efficiency of these narrow micropores can be seen in Figure 8b; the Exp. 5 carbon prepared in this work outperforms materials in the literature with a higher surface area. The hydrogen uptake is comparable to the uptake of hydrolytic lignin carbons prepared by Cheng et al. (2000 m 2 g −1 , 2.0 wt.%) [42] despite a lower surface area. Importantly, the hydrogen uptake of Exp. 5 is on par with the industrially produced TE7, indicating that lignin could be a promising replacement for phenolic resin precursors. This work has shown organosolv lignin can be upgraded to a high value activated carbon, producing a material with a high surface area and reasonable hydrogen uptake at 77 K up to 1 bar.
Carbons produced using a chemical activation procedure exhibit a higher surface area (1800-2500 m 2 g -1 ) and hydrogen uptake (1.6-3.3 wt.%) in comparison to the Exp. 5 and Exp. 8 lignin carbons [12,13,22,[74][75][76][77][78]. The best-performing biomass-or lignin-derived activated carbons are shown to the right-hand side of Figure 8b; materials exhibit BET surface areas over 3000 m 2 g -1 and large H2 uptake values of 2.3-3.3 wt.%. Despite their large hydrogen capacity, these materials exhibit a much broader pore size (from 8 to 40 Å) in comparison to the carbons in this work [13,21,22,42]. Consequently, the large hydrogen uptake can be attributed to their high BET surface area. Efficient hydrogen uptake, required to maximise both volumetric and gravimetric capacity, is reliant on increasing the BET surface area of a material whilst maintaining a narrow pore size of 5-7 Å [25]. The Design of Experiments approach has been shown to be a useful tool for exploring a new activated carbon process; able to handle multiple factors and resulting in valuable insights into the carbonisation and activation behaviour of a material. Future utilisation of this technique will be invaluable in enhancing the surface area of activated carbons whilst maintaining the narrow pore size required for hydrogen storage.  [21][22][23]) biomass-derived carbons from chemical activation (data from [12,13,22,[74][75][76][77][78]) other lignin-derived carbons from chemical or physical activation (data from [39,42]) carbide-derived carbons (data from [79]) and activated carbide-derived carbons (data from [79]). The x and y error bars represent the standard deviation between repeat BET surface area and H 2 uptake measurements, respectively.
Carbons produced using a chemical activation procedure exhibit a higher surface area (1800-2500 m 2 g −1 ) and hydrogen uptake (1.6-3.3 wt.%) in comparison to the Exp. 5 and Exp. 8 lignin carbons [12,13,22,[74][75][76][77][78]. The best-performing biomass-or lignin-derived activated carbons are shown to the right-hand side of Figure 8b; materials exhibit BET surface areas over 3000 m 2 g −1 and large H 2 uptake values of 2.3-3.3 wt.%. Despite their large hydrogen capacity, these materials exhibit a much broader pore size (from 8 to 40 Å) in comparison to the carbons in this work [13,21,22,42]. Consequently, the large hydrogen uptake can be attributed to their high BET surface area. Efficient hydrogen uptake, required to maximise both volumetric and gravimetric capacity, is reliant on increasing the BET surface area of a material whilst maintaining a narrow pore size of 5-7 Å [25]. The Design of Experiments approach has been shown to be a useful tool for exploring a new activated carbon process; able to handle multiple factors and resulting in valuable insights into the carbonisation and activation behaviour of a material. Future utilisation of this technique will be invaluable in enhancing the surface area of activated carbons whilst maintaining the narrow pore size required for hydrogen storage.

Conclusions
An investigation into the preparation of activated carbons from lignin utilising a Design of Experiments approach has provided valuable insight into the use of these materials for hydrogen storage. Regression analysis revealed that carbonisation and activation temperature had the most significant influence on the three measured responses; BET surface area, SAXS-derived average pore size, and hydrogen uptake at 77 K and 1 bar. The maximum surface area was obtained using the severest activation conditions (1000 • C for 60 min). However, a low carbonisation temperature (350 • C) was also important to minimise the effect of pore widening, which is detrimental to the material surface area. The lowest carbonisation and activation temperatures (350 • C and 800 • C, respectively) resulted in the smallest pore size. Maximum hydrogen uptake (1.9 wt.%) was obtained using the lowest carbonisation conditions (350 • C for 30 min), to retain a small pore size for efficient hydrogen packing, followed by activation under the severest conditions (1000 • C for 60 min) to maximise surface area (1400 m 2 g −1 ). This work has shown the potential application of organosolv lignin as a feedstock to produce high surface area activated carbons for hydrogen storage. The Design of Experiments method provided valuable insight into this synthetic process and will be invaluable for the further optimisation of these materials to enhance the material surface area while maintaining narrow micropores.