Hydraulic Tomography for Estimating the Diffusivity of Heterogeneous Aquifers Based on Groundwater Response to Tidal Fluctuation in an Artificial Island in Taiwan

1Department of Safety, Health and Environmental Engineering, National Yunlin University of Science and Technology, Douliou, Yunlin, Taiwan 2Research Center for Soil andWater Resources andNatural Disaster Prevention, National Yunlin University of Science and Technology, Douliou, Yunlin, Taiwan 3Doctoral Program, Graduate School of Engineering Science and Technology, National Yunlin University of Science and Technology, Douliou, Yunlin, Taiwan 4Department of Hydrology and Atmospheric Sciences, University of Arizona, Tucson, AZ, USA 5Department of Environment, AECOM Taiwan Corporation, Kaohsiung, Taiwan


Introduction
The efficient planning of groundwater resources is necessary, and aquifer hydrogeological parameters provide valuable information when addressing groundwater resource management issues.The transmissivity (), hydraulic conductivity (), storage coefficient (), and hydraulic diffusivity () are essential parameters for controlling the groundwater flow in aquifers [1].Additionally, a correct hydraulic parameterization of the aquifers has a direct impact for an accurate description of conservative and reactive transport in the subsurface [2,3].
For estimating the hydraulic parameters from heterogeneous aquifers, Meier et al. [4] assumed that  is heterogeneous and  is homogeneous in a synthetic aquifer.They collected the late time drawdown data results of observation wells in pumping tests and used Jacob's method to estimate the  and  results.They found a strong spatial variability of the estimated  result.Nevertheless, it often becomes problematic when determining the hydrogeological parameters of heterogeneous aquifers using analytical solutions of traditional homogeneous hypotheses [5].
For that reason and in order to obtain the heterogeneous hydrogeological parameters of an aquifer, Huang et al. [6] obtained the drawdown data of an aquifer in a real site using pumping tests.Those data were then used in numerical methods to estimate the heterogeneous hydrogeological parameters.Due to the geological heterogeneity, the fluid flow has the effect of preferential flow [7].Russian et al. [8] presented a relatively simple multicontinuum approach that could be used to link the scaling of the discharge of the frequency transfer function (FTF) to a stochastic description of the catchment heterogeneity in the fractured aquifers.Pedretti et al. [9] focused on the scale dependence of hydraulic parameters in heterogeneous fractured aquifers based on the concept of transfer functions (TF).Their results showed that the scale dependence of  was independent from the adopted formulation (single or dual-continuum), while  was more sensitive to the presence of multiple continua.What is more, other relative researches used the hydraulic tomography (HT) to prove the heterogeneity of aquifers in real fields [6,[10][11][12].
However, the aforementioned literature focuses on inland studies to estimate aquifers with a heterogeneous hydrogeological distribution field.Hodgkinson et al. [13] analyzed the geological heterogeneity of a back-barrier sand island using the geophysical method.Other research investigations characterized the influence of lithological heterogeneity in groundwater systems on island atolls [14][15][16].Such studies were based on the conditions of naturally occurring island aquifers with heterogeneous hydrogeology.Regarding the characteristics of aquifers in artificial islands, several studies investigated the groundwater behavior using analytical solutions [17,18].Furthermore, Li et al. [19] estimated the hydrogeological parameters using semi-numerical simulations in homogeneous aquifers.
As the groundwater level responses to tidal fluctuation, the tidal methods have been widely used as a cost-effective way to assess major hydrogeological parameters in coastal aquifers [20][21][22][23].Based on this feature, Gelhar [24] derived a formula for hydrogeological parameters and frequency (reciprocal of time) in unconfined aquifers.According to this, the natural tidal fluctuation in homogeneous aquifers can be used in the spectral analysis for determining  in coastal aquifers [25].
While the previously mentioned studies focused on the homogeneous hydrogeological parameters under uniform aquifers, the purpose of this study was to prove the existence of hydrogeological heterogeneity in artificial islands, as well as to develop a method for estimating the heterogeneous  field using the HT, when the groundwater level response is affected by the tidal fluctuation.The aquifer of an artificial island (Yunlin Offshore Industrial Park) in Taiwan was used as our case study.

Site Description.
For this research, we selected the artificial island that is located at the mid-west coast of Taiwan (Figure 1(a)).It is approximately 8 km long, 3 km wide, and it expands in a territory of 22.55 km 2 .Its creation occurred due to the artificial land reclamation (7 m under the surface level).Its major material is silty sand [27].For this study, 55 observation wells were used (Figure 1(a) circles).The average depth of all the 55 observation wells is 10 m.Most screen intervals are opened at the depth ranging from 1 to 10 m.Therefore, the aquifer is the phreatic aquifer.The long-term data (11/07/2013∼05/27/2014) of the rainfall, tidal, and groundwater levels are showed in Figure 1(b).It can be observed that January 2014 has no precipitation, whereas May 2014 presents the highest precipitation occurrence.Furthermore, the average tidal level was 0.4 m and the maximum tidal range was 4.23 m for the same period.The tide occurred every day and it was divided into two periods (semidiurnal tide), as it is shown in the shortterm data (01/01/2014∼01/03/2014) of Figure 1(b).It should be noted that this study ignored the influence of seawater intrusion.The detailed long-term groundwater levels of the 55 observation wells can be found in the Supplementary Materials (available here).This study used 55 piezometers that were equipped with pressure transducers for automatic water level measurements at 1-hr intervals (Formosa Petrochemical Corporation provided 45 items of original groundwater level data-observation wells numbers ob01 to ob45.The groundwater level data from the remaining 10 observation wells were collected using the Solinst Model 3001 Levelogger).We selected 5 observation wells (ob37, ob44, ob33, ob05, and ob17 located east, south, west, north, and middle, resp.)out of the total 55, to illustrate the long-term groundwater level of the artificial island, as it is shown in Figure 1(b).The highest rise of the groundwater level was observed during the rainfall period (May 2014).In order to avoid the precipitation influence, this research focused on the nonrainfall period (00:00 a.m., 01/01/2014∼11:00 a.m., 01/21/2014).
The slug test can be used to obtain the hydrogeological parameters of an aquifer.For the artificial island, the slug test results of 10 observation wells were obtained from a previous study [26], as it is shown in Table 1.It was found that the  range was within 0.211∼1.849m/day, and the average  was 0.866 m/day.However, the specific storage (  ) was not investigated in this study [26].According to Freeze and Cherry [28],   can be calculated by   = (  + ), where  represents weight density,   represents compressibility of a porous medium,  represents compressibility of water, and    represents porosity.Therefore,  of water is 9.8 + 3 (N/m 3 ); the average   of sand and  of water are 1.0−8 (m 2 /N) and 4.4 − 10 (m 2 /N), respectively [28, p. 55].For the artificial island, Chien and Lin [27] found that  is 0.41.Based on those analyses,   is 9.98 − 5 m −1 and the spatially averaged hydraulic diffusivity () is 3.61 + 2 m 2 /hr.

Hydrogeological Analysis of Homogeneous Aquifers.
Previous researches [24,25] estimated  using analytical solutions which assumed 1D homogeneous aquifers.According to Gelhar [24], the spectrum analysis is a technique for estimating the hydrogeological parameters by means of spatially distributed model based on the Dupuit theorem.This relates groundwater level response to stream fluctuation through a spectral representation of time series.Our study adopted and extended the approach by Gelhar [24] with estimating  in the artificial island using two different scenarios.One was the same solution for the 1D homogeneous aquifer (1D aquifer), as presented by Gelhar [24], and shown in Figure 2(a).The other one was the 2D axisymmetric homogeneous aquifer (2D aquifer), which is shown in Figure 2(b).The Gelhar solution [24] uses the power spectrum ratio of groundwater fluctuation and tidal fluctuation with the distance ratio of the observation well and aquifer to estimate .However, does this always apply in aquifers of artificial islands?
We considered the linearized form of the classical Dupuit approximation [29].It is given by Gelhar [24] as ).This study ignored the accretion ( = 0), and the solution of the Fourier transform of (1) that is given by Gelhar [24] as where   and  ℎℎ are the spectral density functions of the tidal fluctuation and groundwater fluctuation, respectively,  [L 2 /T] is the hydraulic diffusivity, Ω is the dimensionless frequency, and  [rad/T] is the angular frequency.The ratio of the power spectrum refers to the groundwater fluctuation response of an observation well that is affected by tidal fluctuation;  can be estimated from (2).However, (2) determines the length of the aquifer (L), while Gelhar [24] determines  from the impermeable boundary.The aquifer impermeable boundaries of an island are difficult to define.In the real aquifer, we estimated a constant equivalent length ( equiv ) representing the length of the aquifer.
This study used the software VSAFT2 [30] to simulate the aforementioned scenarios in 1D and 2D aquifers.The first scenario of this research was the same as the 1D aquifer presented by Gelhar [24].The grid size length and width were 10 m each. was 3.61+2 m 2 /hr as previously mentioned.We assumed 7 different L (100 m, 250 m, 500 m, 750 m, 1,000 m, 2,500 m, and 5,000 m) to estimate  equiv .For each L, the observation wells were placed in every 10 m.Then, we created a sin wave as the boundary condition (to represent the tidal fluctuation) with a 12-hr period and a tidal range of 2 m.The boundary behavior was a periodically varying head boundary condition on one end and an impermeable boundary on the other (Figure 2(a)).The total simulation time was 256 hr and the time step was 1 hr.
The simulation was executed based on the 7 different  of the 1D aquifer.From each L, we obtained 256 simulation groundwater level datasets in every observation well.Next, we used the Fast Fourier Transform (FFT) to analyze the spectrums of the groundwater level of each observation well and tidal fluctuation.Then, we applied (2) in the MATLAB built-in function (least squares fitting) to estimate  equiv for the 7 different  based on the aforementioned spectrums and the known conditions ( and /).When  is a fixed condition in the 1D aquifer, a set of  equiv for different aquifer lengths is expected to be obtained.
The second scenario was based on the 2D aquifer (Figure 2(b)).The setting conditions were the same as in the 1D aquifer (e.g., grid size, boundary condition, locations of observation wells, and L).However, because the 2D aquifer is an axisymmetric aquifer, the assumed  was the distance from the boundary to the center of a circle.Furthermore, the N direction was selected to estimate  equiv from the 7 different .We used the estimated  equiv method to the 2D aquifer in order to investigate the feasibility of the Gelhar solution [24] in the 2D aquifer.

Spatial Distribution of Groundwater Level Affected by
Tidal Fluctuation.This analysis investigated the spatial distribution based on the data obtained from groundwater levels of the 55 observation wells and tidal fluctuation in the frequency and time domains.Specifically, 484 groundwater and tidal level datasets were obtained for each observation well, as well as the tidal station during the nonrainfall period (00:00 a.m., 01/01/2014∼11:00 a.m., 01/21/2014).By following this analysis, it was possible to detect the existence of heterogeneous characteristics in the real aquifer.This study used two approaches to analyze the spatial distribution of groundwater levels of the 55 observation wells affected by tidal fluctuation in the frequency domain: (1) the continuous wavelet transform (CWT) and (2) the wavelet coherence (WC).
The wavelet transform can be used to decompose a time series over a time-scale space.It provides a visualization of power distribution along time and frequency.A detailed mathematical analysis of wavelet transform can be found in Walker [31].The CWT is widely applied in hydrology sciences [32][33][34][35][36]. Therefore, we used the CWT analysis to obtain the power spectrums of the groundwater levels of the 55 observation wells and tidal fluctuation.According to Torrence and Compo [37], we selected the Morlet function of the "mother" wavelet of CWT because it is a Gaussian function and provides accurate localization in the frequency domain.The tidal level is a regular fluctuation; therefore, we obtained the maximum spectrum on a fixed frequency with different time using the CWT approach.The average spectrum of the tidal level was calculated by averaging the maximum spectrum at different times and fixed frequencies.Similarly, we calculated the average spectrum of the groundwater level of each observation well for the same tidal frequency at different times.Then, these 55 average spectrums were divided by the average spectrum of the tidal level in order to obtain the tidal efficiency (TE) results of the 55 observation wells, as well as the TE ranges that were between 0 and 1. Last, the spatial distribution of the TE results was obtained using the kriging method (linear variogram function) of the software Surfer.
The WC analysis investigates the correlation of two signals in the frequency domain.It also analyzes the correlation between river stage and groundwater levels [36].We used the WC to analyze the correlation between the groundwater levels of the 55 observation wells and tidal fluctuation, as well as to calculate the average WC results at different times by fixed frequency (same frequency as with the CWT approach).Last, the spatial distribution of the WC results was obtained using the kriging method.Both CWT and WC analyses were based on the cwt and wcoherence of the MATLAB built-in function.
In the time domain, the correlation and time lag between the groundwater levels of the 55 observation wells and tidal fluctuation were investigated using the cross-correlation analysis.Prior to that, this study removed the seasonal trend of the groundwater and tidal levels.Thus, the tidal perturbation (ℎ  ) was obtained by the linear regression of the tidal level.The groundwater perturbation (ℎ  ) of each observation well was obtained by the polynomial regression of the groundwater level.ℎ  and ℎ  were evaluated using the following equation: where Cor() is the cross-correlation and L [T] is time lag or sampling point lag; ℎ  [L] and ℎ  [L] represent the time series of the tidal and groundwater perturbations, respectively; ℎ  and ℎ  are the mean perturbations; and  is the total number of data in time series.Because the tides have periodic fluctuation, the crosscorrelation between the tidal and groundwater levels at lag times has a periodic rise and fall.For this study, the periodic rise of the first maximum cross-correlation for all the 55 observation wells was used to obtain the spatial distributions of the first maximum cross-correlation and its time lag using the kriging method.The time lag refers to the time elapse of the groundwater level affected by tidal fluctuation at a given location lag.

The HT in Synthetic
Aquifers.This research investigated the estimated heterogeneous  field using the HT in the artificial island based on the groundwater level response affected by tidal fluctuation.Then, we evaluated the detailed spatial variations of the aquifer subsurface hydraulic properties using the HT, based on the successive linear estimator (SLE) [38][39][40].The HT has been evolved from the CAT (computerized axial tomography) scan concept of medical sciences and geophysics.The concept of HT involves stressing an aquifer by tidal fluctuation and collecting the head response of each observation well.The head response results of each observation well can be used to inverse the hydraulic parameters.This is a successful technique, as it can be used for artificial stress, such as pumping test [6,[10][11][12]41], or for natural stress, such as river fluctuation [36,42].Specifically for this study, we first tested a synthetic heterogeneous  field that is known (reference  field).Then, we simulated its groundwater levels based on the 55 observation wells responses induced by the known tidal fluctuation.According to the groundwater levels of the 55 observation wells and tidal fluctuation, we verified the estimated heterogeneous  field result using the HT.Although the reference  field was a known value in the synthetic aquifer, we investigated the feasibility of the HT to analyze the estimated heterogeneous  field in the synthetic aquifer.
Based on the study area range, a 2D horizontal domain of 70 × 70 square elements was built.Each element was 100 m × 100 m; the total element numbers were 2,267.As previously mentioned,  was 3.61 + 2 (m 2 /hr).The variance of ln  (i.e., the natural logarithm of ) was 0.5.The correlation scales were 1,500 m for both  and  directions.They were selected to be approximately half of the width of the study area.Additionally, they provided a description of the average size of the heterogeneity.Figure 3 shows the reference  field that has been generated by a spectral method random field generator [43].For this research, the boundary condition refers to a periodically varying head based on time series of the real tidal fluctuation.The initial condition was the distribution of the groundwater levels of the 55 observation wells using kriging method at 00:00 a.m., 01/01/2014.The simulation time lasted for 24 hr and the time step was 1 hr.We used the HT in the VSAFT2 to simulate the forward model and obtain the observed groundwater level of each observation well (24 datasets).The total observed groundwater level datasets of all the 55 observation wells were 1,320.
In the inverse model, we estimated the total number of the observation wells (the current observation well and all the previous ones; e.g., the 3rd observation well refers to the 3rd, 2nd, and 1st observation wells) that are required to obtain the optimal estimated heterogeneous  field.Specifically, we selected the observation wells using the grid method starting from the 1st observation well (observed groundwater levels of 24 datasets) and increasing it up to the 55th observation well (observed groundwater levels of 1,320 datasets).Then, we created a 2D horizontal domain with grids of 1,000 m × 1,000 m.The total number of grids was 30, and they were defined as G1, G2, . . ., G30.We started from the centered observation well and used the clockwise direction for the remaining observation wells (Figure 4).The observation wells were sequentially selected from 1 to 30 until they were 55 in total.If there were more than one observation wells in a grid, the one closer to the center of the grid was prioritized.If there were no observation wells in a grid, the grid was skipped.This process was repeated until all the observation wells were estimated.The sequential clockwise order is shown in Table 2. Based on our knowledge, this method for obtaining the optimal estimated heterogeneous  field is innovative for aquifer grid sampling in artificial islands; therefore, it is highlighted in this study.Sequential clockwise order of the observation wells ob14, ob16, ob55, ob13, ob51, ob52, ob03, ob02, ob07, ob47, ob34, ob32, ob24, ob25, ob04, ob41, ob21, ob10, ob01, ob53, ob40, ob22, ob05, ob38, ob46, ob43, ob20, ob23, ob15, ob17, ob42, ob31, ob12, ob37, ob45, ob08, ob33, ob48, ob50, ob36, ob29, ob39, ob26, ob28, ob49, ob09, ob54, ob44, ob18, ob35, ob06, ob27, ob11, ob19, ob30 The HT was applied from the 1st observation well to the 55th in order to investigate each estimated heterogeneous  field.Next, the mean absolute error (MAE), root mean square error (RMSE), and correlation coefficient (COR) between the estimated heterogeneous and reference  fields in different quantitative observation wells were calculated, so as to estimate the required number of observation wells for delivering the optimal estimated heterogeneous  field.The lower the MAE and RMSE and the higher the COR, the more consistent the estimates.We used the optimal estimated heterogeneous  field to simulate the forward model and obtain the predicted groundwater levels of the observation wells that were not used in the inverse model.Then, we calculated the error and correlation between the predicted and observed groundwater levels using the same observation wells, so as to validate the groundwater levels using the optimal estimated heterogeneous  field.

The HT in Real Aquifers.
The HT was used in the real aquifer that has known conditions (e.g., boundary condition, initial condition, and simulation time).These conditions were the same as those with the synthetic aquifer.The observed groundwater levels were obtained from the historical records of the artificial island.
For the inverse model, that is, the use of the HT in the real aquifer, we followed the same procedures as those in the synthetic aquifer.We selected the observation wells according to the sequential clockwise order (Table 2).Then, we calculated the estimated heterogeneous  field error and correlation results by comparing 54 observation well pairs.Each pair included a current observation well and its following (e.g., 1st-2nd, 2nd-3rd, . . ., 54th-55th).By doing so, we estimated the required number of observation wells for delivering the optimal estimated heterogeneous  field.Next, we used the optimal estimated heterogeneous  field to simulate the forward model and obtain the predicted groundwater levels of the observation wells that were not used in the inverse model.We calculated the error and correlation by comparing the predicted and observed groundwater levels at the same observation wells in order to validate the optimal estimated heterogeneous  field.

Evaluation Criteria.
The MAE, RMSE, and COR (0 ≤ COR ≤ 1) were the performance statistics for evaluating (a) the similarity between the estimated heterogeneous and reference  fields in the synthetic aquifer; (b) the similarity between the estimated heterogeneous  fields of the observation well pairs in the real aquifer; (c) the similarity between the predicted and observed groundwater levels in the synthetic aquifer; and (d) the similarity between the predicted and observed groundwater levels in the real aquifer.

Results and Discussion
3.1.Homogeneous Aquifers.Based on the estimated  equiv results of the 1D aquifer (Figures 5(a) and 5(b)), when  was from 500∼5,000 m, the error results showed a steady state (error = 1.72%∼0.94%).This proves the effectiveness of the estimated  equiv at the 1D aquifer.Regarding the estimated  equiv results of the 2D aquifer (Figures 6(a) and 6(b)), due to the grid limitations of the VSAFT2, we could only analyze 4 aquifer lengths (100 m, 250 m, 500 m, and 750 m) in which the error results did not show a steady state (the error dropped between the estimated  equiv and the theory L).By comparing the error results of 1D and 2D aquifers (Figures 5(b) and 6(b)), it can be observed that the latter presents higher error for all its aquifer lengths.For instance, at 500 m, the 1D and 2D aquifer errors were 1.72% and 9.51%, respectively.The reason for this condition is that the 2D aquifer has two degrees of freedom.Therefore, the Gelhar solution [24] cannot be applied to the 2D aquifer.
In line with our findings, the 2D scenario for a homogeneous island aquifer cannot be described by the Gelhar solution [24].In other words, heterogeneous characteristics may exist in artificial island aquifers.This is why we adopted a heterogeneous concept to describe the hydrogeological characteristics of the real aquifer.

Results of Spatial Distribution of Groundwater Level
Affected by Tidal Fluctuation.In the frequency domain, the TE results between the groundwater levels of the 55 observation wells and tidal fluctuation are shown in Figure 7(a).The TE range was from 2.43 − 6 to 5.88 − 2. It can be observed that the TE results of the ob03 and ob43 presented the highest values: 5.88 − 2 and 5.44 − 2, respectively.Moreover, in the southeast coast of the real aquifer, the groundwater level was significantly affected by tidal fluctuation.The TE results showed irregular distribution, a fact that proves the existence of heterogeneous characteristics in the real aquifer.Additionally, the WC results (Figure 7(b)) showed a range from 0.63 to 0.99.Specifically, the results of 8 observation wells (ob03, ob09, ob24, ob34, ob39, ob43, ob44, and ob52) presented values higher than 0.9.This occurred because the aforementioned observation wells were relatively close to the boundary of the artificial island, so they were susceptible to tidal fluctuation.
In the time domain, the first maximum cross-correlation results (Figure 7(c)) showed a range from 0.03 to 0.99.The first maximum cross-correlation results of 4 observation wells (ob03, ob24, ob39, and ob43) presented values higher than 0.9.Again, these observation wells were relatively close to the boundary of the artificial island, so they were susceptible to tidal fluctuation.The time lag results (Figure 7(d)) showed that the time lag of the observation wells near the boundary was relatively short.The results of spatial distributions of the TE, WC, first maximum cross-correlation, and time lag analyses presented many similarities, as it is shown in Figures 7(a), 7(b), 7(c), and 7(d).The area around the boundary of the aquifer was affected by the tide, whereas the middle area of the artificial island was not.Regarding the observation wells, the ob03 and ob43 were affected the most and were followed by ob09, ob24, ob34, ob39, ob44, and ob52.Due to the tidal fluctuation in the areas of the aforementioned observation wells, the pressure wave propagation was relatively fast and therefore, the results of the TE, WC, first maximum cross-correlation, and time lag analyses showed irregular spatial distributions.This also proves the existence of heterogeneous characteristics in the real aquifer.

Results of the HT Used in Synthetic
Aquifers.The groundwater levels of all the 55 observation wells were used to investigate the estimated heterogeneous  field using the  HT in the synthetic aquifer.The error and correlation results between the estimated heterogeneous and reference  fields are shown in Figure 8(a).It can be observed that they presented a steady state from the 26th to the 55th observation well.Hence, the optimal estimated heterogeneous  field was obtained after reaching the 26th observation well (RMSE: 0.410 m 2 /hr and COR: 0.830).The arithmetic mean and variance of the total estimated heterogeneous  fields are shown in Figure 8(b).Furthermore, the results of the arithmetic mean and variance of the estimated heterogeneous  field in the 26th observation well were closer to the reference  field.
It should be noted that from the 26th observation well to the 55th, the arithmetic mean and variance presented a steady state.Therefore, in this study, the estimated heterogeneous  field of the 26th observation well was used as the optimal estimated heterogeneous  field.Figure 9(a) presents the optimal estimated heterogeneous  field.The high and low  patterns were generally similar to those of the reference  field (Figure 3).By comparing the  results between the optimal estimated heterogeneous and reference  fields at the same elements (2,267 datasets in total) in Figure 9(b), it can be observed that the  results are on or near the 45-degree line, showing a high correlation (COR: 0.830).Next, we compared the predicted groundwater levels of the optimal estimated heterogeneous  field and the observed groundwater levels of the reference  field.Both were obtained from the observation wells that were not used in the inverse model (the last 29 observation wells), as shown in Figure 10 (696 datasets in total).It can be observed that the groundwater level results are on or near the 45-degree line.Also, at the observation well ob35, the values of the predicted groundwater levels were higher than those of the observed groundwater level.This shows that the optimal estimated heterogeneous and reference  fields had few similarities in the region around the ob35.However, the error and correlation results of the groundwater levels of the last 29 observation wells were 0.086 m and 0.977, respectively.Hence, the optimal estimated heterogeneous  field presented in general the same hydrogeological characteristics as the reference  field.This proves that the HT was successfully used in our study for investigating the optimal estimated heterogeneous  field in the synthetic aquifer.

Results of the HT Used in Real
Aquifers.Next, we used the HT to the real tested scenario and analyzed the estimated heterogeneous  field error and correlation for all the observation well pairs (Figure 11(a)).According to Figure 11(a), the optimal estimated heterogeneous  field could be obtained after reaching the 50th-51st observation well pair.The RMSE and COR of this pair were 0.02 m 2 /hr and 1.000, respectively.The arithmetic mean and variance of the total estimated heterogeneous  fields are shown in Figure 11(b).The arithmetic mean and variance of the estimated heterogeneous  field from the 50th-51st observation well pair to the 54th-55th presented a steady state.Therefore, the estimated heterogeneous  field of the 50th observation well was used as the optimal estimated heterogeneous  field.
Figure 12(a) shows the results of the optimal estimated heterogeneous  field.It can be observed that the real aquifer has high  regions at middle and southeast, indicating that these high  regions are more sensitive to tidal fluctuation.Also, the high and low estimated  patterns of the regions that are close to the perimetric boundary of the artificial island in Figure 12(a) resemble the high and low first maximum crosscorrelation results in Figure 7(c).However, this phenomenon cannot be observed in the middle areas of the artificial island.Furthermore, the regions of the high estimated  results in Figure 12(a) have low uncertainty, which is the distribution of the residual variances of the estimated , shown in Figure 12(b).For that reason, the estimated results near the high  regions are more reliable.
Next, we compared the predicted groundwater levels of the optimal estimated heterogeneous  field and the observed groundwater levels of the historical records of the real aquifer.Both were obtained from the observation wells that were not used in the inverse model (the last 5 observation wells), as shown in Figure 13(a) (120 datasets in total).It can be observed that the groundwater level results are on or near the 45-degree line.The RMSE and COR results of the groundwater levels of the last 5 observation wells were 0.021 m and 0.997, respectively.Hence, the optimal estimated heterogeneous  field presented in general the same hydrogeological characteristics as the real aquifer.Also, this study used the optimal estimated heterogeneous  field to simulate the forward model and obtain the predicted groundwater levels of all the 55 observation wells from 3 different time periods, each of which lasted for 72 hr (01/08/2014∼ 01/10/2014, 01/11/2014∼01/13/2014, and 01/14/2014∼01/16/2014; 3,960 datasets/period).Then, we compared the predicted groundwater levels of the optimal estimated heterogeneous  field and the observed groundwater levels of the historical records of the real aquifer, as shown in Figures 13(b), 13(c), and 13(d).It can be observed that the correlation results were higher than 0.99.Again, this demonstrated that the optimal estimated heterogeneous  field presented in general the same hydrogeological characteristics as the real aquifer.Additionally, the HT sufficiently investigated the optimal estimated heterogeneous  field.

Conclusions
This study proved the effectiveness of the HT for estimating  of heterogeneous aquifers based on groundwater level responses that are affected by tidal fluctuation in the artificial island.According to our research results, the following should be highlighted: (1) The spectrum analysis was used to estimate  equive of 1D and 2D aquifers at different .We found that the error result of the estimated  equive of the 2D aquifer was higher than that of the 1D aquifer by 7% in the same distance (L).This occurred due to the fact that the 2D aquifer has two degrees of freedom.Therefore, the 2D scenario for a homogeneous island aquifer cannot be described by the Gelhar solution [24].
(2) The CWT, WC, cross-correlation, and time lag analyses were used to investigate the spatial distributions of the groundwater level affected by tidal fluctuation.We found that the pressure wave propagation was relatively fast in the southeast coast of the real aquifer, due to the tidal fluctuation.Additionally, the spatial distribution results of TE, WC, first maximum crosscorrelation, and time lag analyses were irregular.This proves the existence of heterogeneous hydrogeological characteristics in the artificial island.
(3) We investigated the estimated heterogeneous  field in the synthetic aquifer using the HT.The results showed that the optimal estimated heterogeneous  field could be obtained by the 26th observation well; the RMSE and COR were 0.410 m 2 /hr and 0.830, respectively.By comparing the optimal estimated  heterogeneous and reference  fields, we found similarities in both high and low  regions.According to the predicted groundwater level results of the 29 last observation wells, the RMSE and COR were 0.086 m and 0.977, respectively.Hence, it has been demonstrated that the HT can sufficiently investigate the optimal estimated heterogeneous  field in synthetic aquifers.
(4) We investigated the estimated heterogeneous  field in the real tested scenario using the HT.The results showed that the optimal estimated heterogeneous  field could be obtained by the 50th observation well; the RMSE and COR were 0.020 m 2 /hr and 1.000, respectively.The optimal estimated heterogeneous  field and the first maximum cross-correlation distribution showed that the high correlation regions
where ℎ(, ) [L] is the hydraulic head, x [L] is the distance of the observation well from the tidal boundary,   [L −1 ] is specific storage, K [L/T] is the hydraulic conductivity, and  [L/T] represents accretion.An aquifer of finite length L [L] is connected to the tidal boundary (Figures 2(a) and 2(b)

Figure 4 :
Figure 4: Grid numbers and locations of the observation wells used for the inverse model.

Figure 5 :Figure 6 :
Figure 5: Results of 1D aquifer of the (a) scatter plot of the theory  versus the estimated  equiv and (b) error graph of the estimated  equiv .

Figure 7 :
Figure 7: Spatial distributions of the groundwater level affected by tidal fluctuation.(a) Tidal efficiency (TE) results based on the continuous wavelet transform (CWT) approach between the groundwater and tidal levels; (b) wavelet coherence (WC) results based on the wavelet coherence approach between the groundwater and tidal levels; (c) first maximum cross-correlation results between the groundwater and tidal levels; and (d) time lag results between the groundwater and tidal levels.

Figure 8 :Figure 9 :
Figure 8: Total number of used observation wells for investigating the estimated heterogeneous  field results in the synthetic aquifer.(a) Error and correlation graph of the estimated  versus the reference ; (b) the ln (mean) and ln (variance) graph of the estimated  and reference .

Figure 10 :
Figure10: Scatter plot of the predicted groundwater levels of the optimal estimated heterogeneous  field and the observed groundwater levels of the reference  field in the synthetic aquifer.

4 /hＬ 2 )Figure 11 :
Figure 11: Total number of used observation wells for investigating the estimated heterogeneous  field results in the real aquifer.(a) Error and correlation graph of the estimated  obtained from the current observation well and its following (observation well pair); (b) the ln (mean) and ln (variance) graph of the estimated  and initial .

Figure 12 :
Figure 12: Optimal estimated heterogeneous  field in the real aquifer (50th observation well used for the inverse model).(a) Optimal estimated heterogeneous  field; (b) uncertainty of  estimation.

Table 1 :
Slug test results of the 10 observation wells.

Table 2 :
Sequential clockwise order of the observation wells used for the inverse model.