# «OPTIMIZATION AND MONITORING OF GEOLOGICAL CARBON STORAGE OPERATIONS A DISSERTATION SUBMITTED TO THE PROGRAM IN EARTH, ENERGY AND ENVIRONMENTAL ...»

For reasonable values of Ns (in this work we use Ns = 500), the ‘approximate’ ˆ history match deﬁned by Equations 4.10 and 4.11 does not result in values of S that are close to those achieved by performing a detailed history match as described in Section 4.3.1. Thus, these m∗ are more likely to represent points on the tail of the j true posterior PDF. Nonetheless, this approach does provide some guidance regarding sensor placement and data weighting, and we observe it to lead to history matched

## 4.3. DATA ASSIMILATION AND SENSOR PLACEMENT

models with smaller errors than are observed using default placements and data weight parameters. In future work, it will be useful to assess the impact of using other (less approximate but still eﬃcient) approaches for determining m∗. jwhere wP, wBTT, wBHP ∈ [0, 1] are the data weight parameters to be optimized and θP, θBTT and θBHP are scaling parameters.

These scaling parameters are determined such that, given a ‘typical’ initial guess for m(ξ), the overall contributions from each of the three error components (sensor pressures, CO2 breakthrough times at sensors, injector BHPs) at the start of the minimization are of equal magnitude. In other words, the θ parameters scale the three basic error contributions so they are of comparable magnitude, and the w parameters (with w ∈ [0, 1]) then optimally weight the three error components. The terms I1, I2 and I3 designate identity matrices of appropriate dimension.

The speciﬁc values used in this work for the scalings are θP = 1, θBTT = 10−2 and θBHP = 1.6. The magnitude of θBTT is small relative to θP and θBHP because there is only a single measurement for breakthrough time at each sensor location, and we require this contribution to the misﬁt to be of the same magnitude as the contributions from the sensor and injection pressure data (note that, in our runs, pressure data are measured every 50 days during the injection period and every 500 days during the equilibration period). Breakthrough time for a particular sensor block is deﬁned as the time when the CO2 (gas) saturation in that block ﬁrst exceeds 0.01. The accuracy

## 70 CHAPTER 4. CLOSED-LOOP AQUIFER MANAGEMENT

of breakthrough time is limited by the size of the time step.To ﬁnd the optimal sensor placements and data weight parameters that minimize Equation 4.8, we apply the PSO method described in Section 2.2.2. The time consuming portion of this optimization is the initial simulation of the Ns prior models.

Once this is accomplished, the optimization proceeds quickly as no additional ﬂow simulations are required.

This completes the description of the key computational components required for our closed-loop aquifer management procedure.

**4.4 A-priori optimization results**

In this section we present optimization results for the placement of four horizontal CO2 injection wells and their respective time-varying injection rates (these rates are determined at initial time and then every ﬁve years over the 30-year injection period). We ﬁrst assess deterministic cases (known geology) with no brine cycling in Section 4.4.1, followed by cases with brine cycling in Section 4.4.2. Optimization with uncertain geology is considered in Section 4.4.3. In all cases, PSO is applied for the optimizations.

In this work, we will refer to ﬁve ‘true’ geological models, denoted TR1,..., TR5.

These synthetic models, generated using SGeMS, are consistent with the aquifer model speciﬁcations provided in Section 2.3. Figure 4.3 displays the three-dimensional porosity ﬁelds for two of these models. Detailed optimization results will be presented for model TR1. Results for the other four models will be provided in terms of key summary statistics.

The results for the deterministic cases (Sections 4.4.1 and 4.4.2) correspond to results presented earlier in Chapter 3, and our general ﬁndings are analogous. Here, however, we use diﬀerent geological models, consider a diﬀerent objective function, apply a diﬀerent optimization algorithm, and make extensive use of parallel computation. These results are also used for comparison to cases with uncertain geology.

4.4. A-PRIORI OPTIMIZATION RESULTS

** Figure 4.3: Porosity ﬁelds for two of the ﬁve ‘true’ geological models considered in this work 4.**

4.1 Optimization with known geology Our ﬁrst set of results is for geological model TR1 for CO2 storage with no brine cycling. There are a total of 30 optimization variables in this case. Recall that each PSO particle (we use 22 particles) is a potential solution; i.e., each particle deﬁnes a particular well conﬁguration and set of time-varying injection rates. Each PSO iteration thus requires 22 function evaluations (ﬂow simulations), which provide the objective function associated with each particle.

Results for three PSO optimization runs are shown in Figure 4.4. Since PSO is stochastic in nature, it is useful to run the optimization multiple times, assuming this is feasible computationally. Each run is limited to a maximum of 100 iterations, which corresponds to 2200 ﬂow simulations. As indicated above, using 22 computational cores, each optimization run requires about 5-8 hours of computation. All three runs appear to be essentially converged after 100 iterations. Recall that we use units for mobility of lbmol/[RB·cP], where 1 lbmol/[RB·cP] = 2.85×106 mol/[m3 ·Pa·s]. It is evident from Figure 4.4 that the results are similar in the three PSO runs; i.e., we consistently observe a reduction in the objective function (top layer CO2 mobility) from around 6000 to 3300 lbmol/[RB·cP]. The best run (shown in bold) will be considered in the subsequent discussion.

The default and optimal placements for the four wells, along with maps of the

## 72 CHAPTER 4. CLOSED-LOOP AQUIFER MANAGEMENT

corresponding time-averaged CO2 mobility in the top layer, are shown in Figure 4.5.The default conﬁguration entails simply placing the wells in a square pattern and injecting equal amounts of CO2 into each well. The optimization determines well placements and injection rates that take advantage of the interplay between geological heterogeneity and ﬂow physics such that the top layer CO2 mobility is minimized. In this case the total mobility of CO2 in the top layer is reduced from 6078 lbmol/[RB·cP] in the default case to 3272 lbmol/[RB·cP] in the optimal case.

In Figure 4.5(b), we see that, in the optimal case, all four injection wells are placed at the boundaries of the feasible region (inner square). As a result, the CO2 plume in the optimal case extends close to the boundary blocks (which lie outside the larger square). We note, however, that any CO2 that actually migrates into the boundary blocks is treated as if it were pure mobile CO2 in the top layer at the maximum mobility. This penalty is enforced to ensure that the optimization does not simply ‘drive’ CO2 into the boundary blocks. In the optimized results shown in Figure 4.5(b), less than 0.1% of the injected CO2 ﬂows into the boundary blocks. If we wished to maintain all CO2 within the aquifer domain, this requirement could be incorporated

** Figure 4.5: CO2 injection well placements with time-averaged (over 500 years) CO2 mobility in the top layer for geology TR1 into the optimization as a general constraint.**

The optimal injection rates for this example are shown in Figure 4.6. These results demonstrate that top layer CO2 mobility can be reduced by varying the injection rates of the four wells. Well 4 (denoted W4), and other wells to some extent, display a ‘bangbang’ nature, in which the well either injects at the maximum (0.4 in this case), or does not inject at all. As noted previously, it is diﬃcult to identify explicitly the complex relationships between local geological structures (or heterogeneities) and optimal well placements and rates due to the complex nature of ﬂow through heterogeneous media.

We expect that optimal well locations and injection rates will lead to early horizontal migration of CO2 and will additionally promote mixing to enhance dissolution and residual trapping. These eﬀects may be achieved by pulsing the injection rates, rather than using a constant rate.

** Figure 4.7 displays the distribution of CO2, among the three diﬀerent modes of trapping, within the overall model as a function of time.**

Results are shown for the default and optimized cases (for model TR1). At the end of the simulation time frame (500 years), we see a noticeable decrease in mobile CO2 fraction (from 0.294 in the default case to 0.214 in the optimized case), and a consequent increase in dissolved CO2 and residually trapped (immobile) CO2. Although the objective of the

## 74 CHAPTER 4. CLOSED-LOOP AQUIFER MANAGEMENT

0.3 Figure 4.6: Optimal injection rates for the four wells during the 30-year injection phase for geology TR1 optimization is to minimize mobile CO2 in the top layer of the model, these results demonstrate that optimum operation also acts to reduce mobile CO2 globally. As noted in Chapter 3, the increase in dissolution trapping likely occurs because the CO2 plume is more disconnected in the optimized case, which exposes it to a larger volume of unsaturated brine. The increase in residual trapping may be due in part to the fact that the plume in the optimized case displays a sharp decrease in saturation at the trailing edge, which acts to enhance hysteresis.

In addition to increasing dissolution and residual trapping, we also expect the optimized conﬁguration and injection schedule to result in slower migration of CO2 into the top layer. Figure 4.8 shows that the optimization does indeed lead to a delay in mobile CO2 in the top layer. In particular, we see that at the end of the injection period (30 years), the mobility of CO2 in the top layer for the default case is almost 10 times higher than that in the optimized case.

4.4. A-PRIORI OPTIMIZATION RESULTS

We now expand the optimization to include brine cycling. The geology is still assumed to be known, and we again use model TR1. The additional brine cycling parameters that appear in the optimization deﬁne the timing, duration and pumping fractions for three brine cycling events. There are a total of 47 optimization variables in these cases.

We consider optimizations for scenarios involving 0.01 PV, 0.03 PV and 0.09 PV of

## 76 CHAPTER 4. CLOSED-LOOP AQUIFER MANAGEMENT

brine cycled, where PV refers to the pore volume of the core region of the model (the central 15 × 15 × 8 block region shown in Figure 2.6). As noted previously, these correspond to 0.0036, 0.0108 and 0.0324 pore volumes of the storage region (the entire 25 × 25 × 8 block region shown in Figure 2.6). For each value of PV, one PSO run is performed, in which we again set the maximum number of iterations to 100. Although improved results would presumably be obtained by performing additional PSO runs, the results in Figure 4.4 suggest that the variation between runs is relatively small.By considering diﬀerent brine cycling volumes, we can generate the Pareto front for this problem (analogous to Figure 3.13 in Chapter 3), which essentially represents an optimal trade-oﬀ curve for a bi-objective optimization. This curve, shown in Figure 4.9 (note that the dashed line simply connects the four optimized cases), indicates the minimum top layer CO2 mobility for a speciﬁed PV of brine cycling.

The curve can also be viewed as providing the minimum brine cycling PV required to achieve a speciﬁed top layer CO2 mobility. It is evident from Figure 4.9 that CO2 mobility can be reduced substantially through brine cycling, though large PVs of brine are required to achieve dramatic reductions. Speciﬁcally, the use of 0.03 PV of brine cycling reduces top layer CO2 mobility by a factor of 1.4 relative to the case with no brine cycling, while the use of 0.09 PV reduces it by a factor of 3.8. Because the cost of high-volume brine cycling will be considerable, the Pareto front in Figure 4.9 provides an indication of the trade-oﬀ between level of risk (top layer CO2 mobility) and cost (PV of brine cycling).

We now present detailed results for the 0.03 PV brine cycling case. Figure 4.10 shows the well placements and the top layer CO2 mobility for the default and optimized cases. The default case entails the same well locations and CO2 injection schedule as in the case with no brine cycling. The default model also includes 10-year brine cycling events at 75, 125 and 175 years, with equal brine injection for each well. The optimized case results in a top layer CO2 mobility of 2343 lbmol/[RB·cP], compared to 3712 lbmol/[RB·cP] in the default case. The reduction in top layer CO2 mobility from brine cycling (for both the default and optimized cases) is apparent by comparing Figure 4.10 to Figure 4.5. It is also evident that the optimized well locations are similar, but not identical, to those for the case without brine cycling.

4.4. A-PRIORI OPTIMIZATION RESULTS

** Figure 4.11 shows the CO2 injection rate fractions for the optimized case along with the optimized parameters associated with the three brine cycling events.**

We see in Figure 4.11(b) that the majority of the brine cycling ‘budget’ is used in the ﬁrst brine cycling event, which lasts from year 28 to year 30. This early brine cycling event acts to slow the migration of CO2 to the top layer and to increase residual trapping at early times. This is evident from Figure 4.12, which shows the diﬀerent modes of CO2 trapping, as a function of time, for the default and optimized cases.

The remaining brine cycling is delayed for as long as possible (the maximum start time for a brine cycling event is 200 years), to enable CO2 to accumulate in the top layer prior to brine injection. The optimized brine cycling solution (Figure 4.12(b)) displays a clear increase in late-time immobile CO2 relative to that in the default case (Figure 4.12(a)).

0.2 0.3 0.15 0.2 0.1