Modeling and Observations of the Effects of the Alfvén Velocity Profile on the Ionospheric Alfvén Resonator

We have modeled the individual harmonic frequencies of the Ionospheric Alfvén Resonator (IAR) at Eskdalemuir by solving a one‐dimensional wave equation and using non‐uniform modeled Alfvén velocity profiles. By comparing the results of the modeling alongside harmonics obtained from the Eskdalemuir, UK, data set from 2013 to 2021, the effects of the non‐uniformity of the Alfvén velocity profile on the IAR are considered. We calculated the offset between the fundamental frequency and the harmonic frequency separation and found that this is not constant. From this parameter, we infer that the lower boundary condition of the electric field of the IAR is closest to a node, which agrees with previous studies. We compare the results of the non‐uniform model with previous uniform models and evaluate their interpretations and the implications for the lower boundary condition.


Introduction
Alfvén waves propagating along Earth's magnetic field lines are partially reflected when the Alfvén velocity gradient reaches a maximum due to changes in the plasma mass density in the ionosphere.The upper boundary is likely to be above the F region peak, and the lower boundary is toward the bottom of the ionosphere, where the Alfvén velocity maximizes.The resonance which occurs is termed the Ionospheric Alfvén Resonator (IAR) (Polyakov, 1976).
The first confirmed observation of IAR was at mid-latitudes (Belyaev et al., 1989), and since then IAR have been observed in magnetometer data across a range of latitudes.They have also been observed in sounding rocket (Cohen et al., 2013) and satellite (Hirano et al., 2005;Miles et al., 2018;Pakhotin et al., 2018) data.IAR typically occur over the range of 1-10 Hz but can extend up to 30 Hz (Beggan & Musur, 2018).The harmonic frequency separation (Δf) is the average frequency spacing between the harmonics and is often used in literature to characterize the frequency of the IAR (e.g., Belyaev et al., 1989;Nosé et al., 2017;Yahnin et al., 2003).The fundamental frequency is often obscured by low frequency magnetic activity in the data, and only its higher harmonics are observed, hence Δf is useful in describing IAR frequencies.
IAR have been detected in induction coil magnetometer data at Eskdalemuir Geophysical Observatory, UK (55.3°N, −3.2°E) (Beggan, 2014).We have previously completed a climatological study of IAR Δf for 9 years of data from Eskdalemuir covering 2013 to 2021 (Hodnett et al., 2023) looking at behavior patterns and occurrences of IAR as compared to the predictions from a global ionospheric model.In that study, we compared the measured Δf with modeled Δf values derived from empirical methods.To model Δf, the Alfvén velocity was required.We used the International Reference Ionosphere (IRI) (Bilitza et al., 2017) to estimate the plasma mass density and the International Geomagnetic Reference Field (IGRF) (Alken et al., 2021) to trace the magnetic field lines in order to obtain Alfvén velocity profiles.We also explored other electron density models.We found that the model 10.1029/2023JA032308 2 of 13 followed the general trends of the data well, with greater Δf toward midnight and during the winter time, when the plasma mass density is lower.
As well as studies regarding Δf, the frequencies of the individual harmonics of the IAR have previously been examined.For example, Prikner et al. (2014) used EISCAT measurements to model the IAR in a non-uniform ionosphere, and by applying a full-wave model (Prikner et al., 2007), modeled the first three harmonics and the corresponding altitude profile of the magnetic and electric fields.These harmonics were modeled for an event on 7 March 2001 and compared with magnetometer data from auroral latitudes.The modeled and observed harmonic frequencies were comparable and hence the altitude profiles of the modeled wave harmonics could be studied.
Previous studies which consider a uniform Alfvén velocity profile (Hebden et al., 2005;Potapov et al., 2022) have modeled the IAR by using the ratio of the harmonic frequencies for two cases of standing waves.Figure 1 shows this uniform model, where the lower boundary has (a) a node and (b) an antinode in the electric field.For each, the first three harmonics are plotted.The magnetic field B is directed toward Earth.In the node case (a), the frequencies scale as f, 3f, 5f etc.In the antinode case (b), the frequencies scale as f, 2f, 3f etc.By comparing measured IAR harmonic frequency values to this uniform model, Hebden et al. (2005) and Potapov et al. (2022) inferred information about the lower boundary condition of the IAR reflector.
The interactions between Alfvén waves and plasma has been modeled in a dipole geometry (Streltsov & Lotko, 2008).Lysak et al. (2013) also modeled the IAR in a dipolar geometry with non-uniform conductivity profiles.Lysak (1993) considered the non-uniform geometry of the IAR and found that the IAR frequency separations scale the same way in a more general geometry.
The frequency spacing between IAR harmonics is not, however, uniform and the spacing of the harmonics and their individual frequencies provide insight into the non-uniformity of the ionosphere.In this study, we have modeled the harmonics of the IAR at Eskdalemuir using a non-uniform Alfvén velocity profile by solving a one-dimensional wave equation.The modeled harmonics are then compared to the harmonics obtained from the Eskdalemuir data set, which consists of frequencies extracted from spectrograms between 2013 and 2021.From modeling and observations, we have investigated how the non-uniformity of the ionosphere affects the frequencies of the individual IAR harmonics.
In Section 2, we describe how uniform models (shown in Figure 1) have previously been used to explore the effects that a node or antinode lower boundary condition has on the frequencies of the IAR harmonics.In Section 3, we use non-uniform Alfvén velocity profiles to model the individual harmonics and in Section 4 we demonstrate how the individual harmonics are extracted from the spectrograms.To analyze the results of the model alongside the data, we start by calculating the offset of the fundamental frequency from Δf in Section 5 and then compare the results with the uniform Alfvén velocity profile analysis of Hebden et al. (2005) and Potapov et al. (2022) in Section 6.We conclude in Section 7.

Uniform Models of IAR Alfvén Velocity
For a uniform Alfvén velocity, V A , the frequency f of a harmonic of number N is given by Equation 1, where L is the IAR cavity size in which the wave resonates, and Φ is the phase factor, which relates the offset of the harmonic frequency f and Δf, where Δf = V A /2L (Polyakov & Rapoport, 1981).
It is accepted that the top boundary of the IAR cavity has an antinode in electric field of the Alfvén wave (e.g., Lysak, 1991); however, it is less certain what the conditions at the lower boundary are.By considering the two uniform cases presented in Figure 1, Hebden et al. (2005) used the ratios of the frequencies to determine Φ for case (a), where the lower boundary is a node, and case (b), where the lower boundary is an antinode.
As Equation 1 can be written in the form y = mx + c that is,: ) of the line.When the lower boundary is a node, the frequencies of the harmonics scale as f, 3f, 5f etc. and so Φ is 0.5.When the lower boundary is an antinode, the frequencies scale as f, 2f, 3f, and Φ is 1.Depending on which harmonic numbers are used for this calculation, then for the uniform node case, Φ will be equal to half integer values, 0.5, 1.5, 2.5 etc., and for the antinode case, Φ will be equal to integer values, 1, 2, 3 etc.Hence, the phase factor Φ is equal to the number of half wavelengths in the fundamental frequency, and is representative of the fraction of quarter wavelengths needed in addition to the integer number of half wavelengths for resonance to occur (Hebden et al., 2005).
Hebden et al. ( 2005) calculated Φ from 10 days of observations of IAR using data from the Sodankylä magnetometer (at high latitude in Finland), and found Φ to range between 0.5 and 1.0.The average value of Φ was 0.61.They suggested the lower boundary lay between being a node and an antinode-though generally closer to a node-and that further work was needed to determine the relationship between Φ and the conditions at the lower boundary.Potapov et al. (2022) have also examined the properties of the lower boundary, by assuming a uniform ionosphere and calculating the ratio between the harmonic frequencies of the IAR, for the first five harmonics.Like Hebden et al. (2005), Potapov et al. (2022) considered the two uniform cases in Figure 1.They presented the ratios of the frequencies in a triangular matrix, where the first row represents the ratio of the frequencies of the higher harmonics with the fundamental, the second row is the ratio between the higher harmonics and the first harmonic, and so on.
For the case where the lower boundary has a node in the electric field, the triangular matrix is given by matrix   1  , and for the case where the lower boundary is an antinode, the matrix is given by   2  .The ij notation used by Potapov et al. (2022) represents the rows i and the columns j.
To draw comparisons between the models and observations, the determinants (det(T ij )) of matrices  Hebden et al. (2005), when calculating det(T ij ) from the data, this method requires the fundamental to have been correctly identified in the data for the ratios in the triangular matrix to be correctly calculated.If the fundamental is missed and the second harmonic is used instead, the determinants of matrices   1  and   2  become 3.7 and 3, respectively.If the third harmonic is used in place of the fundamental, the determinants become 2.6 and 2.3.2022) calculated the matrix T ij and the determinant, ensuring that the fundamental was used.They found that the average determinant calculated from all observations was 8.35, which is closest to 9, suggesting that the lower boundary of the IAR has a node in the electric field (and hence an antinode in the magnetic field).When comparing these observations and corresponding determinants with solar cycle, they found that the determinants were largest during the solar maximum.Potapov et al. (2022) suggested that the deviations in the values of the determinants from the uniform cases may arise from a non-uniform resonator cavity.

Modeling the IAR in a Non-Uniform Ionosphere
We have modeled the harmonics of the IAR by considering a non-uniform ionosphere, where the Alfvén velocity in the IAR cavity varies with altitude.The Alfvén velocity is given by Equation 4, where B is the magnetic field strength, ρ is the plasma mass density and μ 0 is the permeability of free space.
Alfvén velocity profiles can be modeled for the IAR at Eskdalemuir using the IRI (Bilitza et al., 2017) to obtain plasma mass density, and the IGRF (Alken et al., 2021) to obtain magnetic field strength.This is the same method used by Hodnett et al. (2023), where such Alfvén velocity profiles were employed to model Δf using a time-of-flight approach and compared with observations.Whilst capturing the general behavior of the frequency of the IAR, such an approach cannot provide information on any variability of the individual frequencies of the harmonics which may result from the non-uniformity in V A .Modeling the individual frequencies will also help to understand deviations from the uniform expected values of Φ and det(T ij ) which may be seen in the data, as introduced in Section 2.
As such, we have modeled the IAR harmonics (fundamental, second, third, fourth, and fifth) of the electric field at Eskdalemuir by solving a one-dimensional wave equation using the same non-uniform Alfvén velocity profiles as Hodnett et al. (2023).Previous modeling work has shown that the vertical structure of the ionosphere is the most important consideration for the IAR (Lysak et al., 2013).Therefore, to model the harmonics of the electric field, we used Equation 5, where ψ is the electric field, x is the position along the magnetic field and t is time.
The Alfvén velocity profile V A (x) is non-uniform and so dependant on x.Using separation of variables, with the temporal function being e −iωt , the wave equation can be written as a function of displacement only (Equation 6), where ω is the angular frequency of a given harmonic.
We solved Equation 6 computationally using the Runge-Kutta fourth order method (RK4).To do this, we refactor Equation 6 into two first order differential equations, where u is a variable: If the lower boundary of the electric field is a node, the initial condition is ψ(0) = 0.If the lower boundary of the electric field is an antinode, the initial condition is to maximize ψ, so ψ(0) = 1.Intermediate boundary conditions can also be modeled by using initial ψ(0) conditions between 0 and 1 (see Section 6).In both cases u(0) = 1.It is assumed that there is an antinode at the upper boundary.To solve the equation for all x, where the altitude range of the model is 100-2,000 km, we iterated through values of ω 2 and solved the equation with RK4 until the maximum at the upper boundary is found.Hence, the harmonics of the IAR can be modeled.Higher values of ω 2 are used to find higher harmonics.
Figure 2 illustrates the non-uniform electric field wave structures for the node (a) and antinode (b) case.When comparing with Figure 1, the effects of the non-uniformity are apparent, with the positions of the nodes and antinodes no longer being spaced evenly apart in each harmonic.
Next, we compare the modeled IAR frequencies with measured values from Eskdalemuir observatory.

Extracting IAR Frequencies From Magnetic Field Measurements
We extracted harmonics of the IAR observed at Eskdalemuir from nine years of almost continuous daily spectrograms (2013-2021).An example spectrogram for 15 September 2013 from 18:00 to 01:00 UT is shown in Figure 3.The y-axis shows the frequency, and the scale bar shows the intensity, where a lighter color represents more energetic signals.The IAR harmonics are the bright horizontal fringes, which increase in frequency until around 01:00 UT.Below 0.5 Hz, low frequency, higher intensity magnetic field activity is visible.Bright, vertical lines are indicative of lightning activity or occasional data loss.
The harmonics are identified by a neural network called a U-net (Ronneberger et al., 2015).We used the U-net developed by Marangio et al. (2020) and the method to extract frequencies is detailed in Hodnett et al. (2023).The input consists of grayscale spectrograms scaled to 256 × 256 pixels, where the time range is limited to "nighttime" periods (16:00-08:00 UT) and the frequency range is fixed between 0.5 and 6.5 Hz.The machine learning algorithm output consists of another 256 × 256 pixel image, where the value of each pixel ranges from 0 to 1 and represents the probability that it is part of an IAR harmonic.The temporal resolution is 3.75 min, given there are 16 pixels in each hour in the time range.
By finding the peaks above a threshold from the output image of the U-net, we can identify the individual frequencies of the harmonics.The U-net identified harmonics for 15 September 2013 are plotted in red in Figure 3.For the years 2013-2021, we have extracted the frequencies of the individual harmonics of the IAR at Eskdalemuir, for times where there are at least five harmonics visible.This ensures there are enough harmonics to calculate Δf and Φ.Overall, we find 12,488 time bins in the data where five harmonics have been identified by the U-net, across the 9-year data set.This is fewer than the number of points used in the study from Hodnett et al. (2023), as previously we used a less constrained approach and found times where at least four harmonics were identified to obtain Δf, rather than five.These instances where at least five harmonics were identified occurred over 644 separate days.From the Eskdalemuir data set, we can explore their individual frequencies alongside the inter-harmonic variation of the IAR harmonics, which are not constant as the ionosphere is not uniform.These results can then be compared with the uniform and non uniform models in order to interpret whether the lower boundary of the IAR has nodal or antinodal behavior.

Comparing Modeled and Observed Frequencies
Using the method outlined in Section 3, we modeled the first five harmonic frequencies of the IAR to correspond with the data set from 2013 to 2021.The model is run at hourly intervals, and so when choosing the time range of each day which had data, we modeled from the hour before the data started to the hour after it ended.We ran the model with the lower boundary condition both as a node and also as an antinode in the electric field; hence there are two models of harmonic frequencies that can be compared with the data.
We present example harmonics modeled for 20 June 2013, setting a node in the electric field at the lower boundary.The wave structure of the first three harmonics is presented in Figure 4 alongside the corresponding Alfvén velocity profile, for 23:00 UT.The values of ω 2 for the first five harmonics are presented in Table 1, alongside the resulting periods and frequencies.
The electric field ranges between −1 and 1 as it has been normalized to its  10.1029/2023JA032308 6 of 13 maximum value at the upper boundary.For this time, Δf between the first (fundamental) and second harmonic is 0.302, and 0.336 Hz between the second and third harmonic.These RK4 solutions can be compared to the time-of-flight calculation of Δf.The time-of-flight calculation leads to an average Δf of 0.320 Hz, which is similar to the frequency separations calculated by solving the wave equation.The mean Δf calculated from the RK4 model using all 5 harmonics is 0.316 Hz, and so gives a similar mean Δf result to the time-of-flight calculation.
The difference in Δf between different harmonics is due to the non-uniformity of the Alfvén velocity profile, as demonstrated in Figure 4a.
This demonstrates that the time-of-flight calculation of Δf accurately describes the separation of higher harmonics, but not the frequency of the fundamental, resulting in an offset of the harmonics and Δf, which was noted by Hodnett et al. (2023).This offset is not constant, which is a result of the non-uniformity of the ionosphere, which we next analyze.
The Runge-Kutta model was run for the hours 22:00-02:00 UT, with a node for the lower boundary condition.The first five modeled harmonics are shown in Figure 5 (plotted in blue).They are plotted over the spectrogram computed for Eskdalemuir channel 2 on 20 June 2013 from 20:00 to 04:00 UT.The harmonics identified by the U-net are plotted in red, and lie over the bright fringes (which are the IAR harmonics) in the spectrogram.
The harmonic number N is stated above each harmonic.As well as noting that the simple time-of-flight calculation results in an average Δf comparable to those returned by the RK4 solutions in Figure 4, in Figure 5 it is clear that the modeled (blue) harmonics lie over the harmonics in the spectrogram.
The U-net identified harmonics can be compared directly to the RK4 harmonics, in this example (20 June 2013 at 23:00 UT) for the 3rd, 4th, and 5th harmonic frequencies.For the data, to compare at the time 23:00 UT we have averaged the hour by calculating the mean of the data for that hour (as the data has a resolution of 16 points per hour).The harmonic frequencies for the 3rd, 4th, and 5th harmonic (N = 2, 3, 4) from the data are 0.817, 1.142, and 1.461 Hz (Δf = 0.322 Hz).In the model, they are 0.819, 1.143, and 1.446 Hz.We can see that for this time, the modeled harmonics are a good fit to the data.The means that the IRI is correctly capturing the characteristics of the ionosphere at Eskdalemuir during this time.
From the spectrogram in Figure 5, although not identified by the U-net analysis, by eye the fundamental appears to be visible from around 22:00-00:00,  slightly above the high intensity, low frequency activity, and below the first red line.The RK4 modeled fundamental fits well with this.Outside of this time range, the fundamental is not visible as it is obscured by the lower frequency activity.
In the uniform model, for a node at the lower boundary, from Equation 1 the fundamental frequency (N = 0), From the time-of-flight model, RK4 model and data, Δf ∼ 0.32 Hz at 23:00 UT.Using the uniform model, f 1 = 0.16 Hz.In the non-uniform RK4 model, f 1 = 0.18 Hz.This difference in uniform and non-uniform fundamental frequencies is a result of the non-uniformity of the Alfvén velocity profile, which changes the structure of the harmonics as seen in Figure 4.As such, the frequencies no longer scale as f, 3f, 5f for the node case in the model.
To explore this in more detail, we examined the offset of the fundamental frequency (f 1 ) with Δf in the non-uniform model and the data.We define the percentage difference between Δf and the fundamental frequency (f 1 ) as ζ: ζ is straightforward to calculate for both the node and antinode non-uniform models.However, in the Eskdalemuir data set we cannot be sure what the lowest harmonic number identified by the U-net is, given the fundamental is often obscured in the spectrogram by low frequency magnetic activity.To estimate f 1 from the data, Δf is subtracted from the lowest f identified by the U-net until f < Δf.This gives f 1 for the node case.If the lower boundary is an antinode, f 1 > Δf and so a Δf is added back, resulting in two fundamental frequencies estimated from the data for the node and antinode case.
The calculations of ζ for each UT from 16:00 to 01:00 are presented in Figure 6.The first column shows the data, the second shows the model and the third shows the difference between the two.The hours of 16:00, 17:00, and 01:00 UT have fewer data points than one quarter of the standard deviation of data points and so we have assessed   9) for Eskdalemuir from 2013 to 2021.The data for each UT from 16:00 to 01:00 UT has been averaged.There are no data points outside the time range that is plotted.The first row shows the node case, and the second the antinode case.Panels (a, d) show the data, with the black crosses representing the total number of data points before averaging in each UT (right hand y-axis).The total number of data points is 12,488.Green shaded bins correspond to UTs where there are a significant number of data points for analysis, with the number of data points being greater than the standard deviation of the total number of data points (σ = 860).Panels (b, e) show the RK4 models, for times where data are available.Panels (c, f) show the difference in ζ (model-data).The gray dashed line is plotted at 0. For each boxplot, the whiskers are 1.5× the interquartile range.The edges of the box are the first and third quartiles.The horizontal line in each box is the median, and any outliers are plotted as circles.
HODNETT ET AL.
10.1029/2023JA032308 8 of 13 that there is not enough data in those UT bins to conduct a thorough analysis.
The median difference between the modeled and observed ζ is −1.6% under the assumption of a node as the lower boundary condition and 10.3% under the assumption of an antinode lower boundary condition.Hence, from this analysis we conclude that the lower boundary of the electric field of the IAR is closer to being a node.
The modeled ζ values for the node case in Figure 6b increase with UT.The variation in ζ is due to the non-uniformity of the Alfvén velocity profile, and hence the increase in ζ with UT may correspond to an increase in non-uniformity with UT.To investigate this, for the UTs shown in Figure 6, we have taken the mean for all modeled Alfvén velocity profiles, and normalized them to their minimum value.These are plotted in Figure 7. From these normalized profiles, we can see that the shape of the profile changes over time.We compared the hours of 18:00 and 00:00 as these differ the most.We created two more profiles, which consisted of the topside of the 18:00 UT profile and the bottomside of the 00:00 UT profile (split at the minimum) and vice versa.From these four profiles, we solved the wave equation (Equation 5) and calculated ζ.We found that changing the topside and bottomside of the profile had similar effects, with approximately 40% of the change being attributed to the bottomside, and 60% being attributed to the topside.Therefore, it is likely that the change in ζ with time seen in the model is due to the changing shape of the Alfvén velocity profile (which becomes more non-uniform as the local time approaches midnight).In the data (Figure 6a), the same trend is less apparent, however it becomes clearer when disregarding the first few poorly sampled bins.There is also more spread in the data, suggesting a wider range of Alfvén velocity profile shapes.Therefore, the ζ parameter could be used to diagnose information about the shape of the Alfvén velocity profile of the ionosphere.

Phase Factor and Determinants
In the uniform model presented by Hebden et al. (2005), as discussed in Section 2, Φ represents the fraction of quarter wavelengths in the wave structure in addition to the integer number of half wavelengths.However, in the non-uniform model, Φ no longer represents this as the wave structures are not uniform (as seen in Figure 2), and values of Φ not equal to 0.5 or 1 are expected even for purely node or antinode lower boundary conditions, respectively.Here we calculate Φ for our data set of harmonic frequencies measured at Eskdalemuir, following the technique of Hebden et al. (2005).We then compare this to the phase factor predicted from our modeling of the IAR under node and antinode lower boundary conditions.
To calculate Φ, we chose time bins where five harmonics were identified by the U-net, so that N in Equation 2 ranges from 0 to 4. Φ is calculated by plotting the harmonic number against the frequency of the harmonics, and finding the intercept/gradient (Equations 2 and 3).The Φ values are then corrected, because if the lowest frequency harmonic identified is not the fundamental, in the calculation of Φ for the following four harmonics, the value of Φ will be greater than 1.For example, if the second to sixth harmonic were identified, the phase factor will be its true value plus 1.Therefore, in our calculations, if the minimum value of Φ > 1 for a given day, then 1 is iteratively subtracted until the minimum is below 1.
The first column of Figure 8 is a demonstration of how Φ is calculated from the data.Panel (a) shows the harmonics identified by the U-net, for the day of data shown in the spectrogram in Figure 5. Panel (b) shows the initial Φ values, before they are corrected.As they are greater than 1, this means the fundamental has not been identified.This can be confirmed by looking at Figure 5 spectrogram with the overplotted U-Net identified and modeled harmonics.The Φ values are then corrected by subtracting 2 from each value, as shown in panel (c).The error in Φ is given by the green shaded area and is calculated by using an error of ±0.5 pixels in the U-net output images, 10.1029/2023JA032308 9 of 13 which corresponds to an error of ±0.025 Hz in the harmonic frequencies.The error can then be calculated for the intercept and gradient, and hence for Φ.
The same calculation can be made for both the node and antinode models of the IAR harmonics.No correction needs to be made for the node case as the fundamental is modeled.For the antinode case, Φ will be greater than 1 and so to compare with the data, we subtracted 1 from these values.Figure 8d shows the modeled harmonics for 20 June 2013, with a node for the lower boundary, and (e) shows the corresponding modeled Φ values.Panel (f) shows the difference between the modeled Φ and Φ calculated from the data.The difference between the model and data is small, with a mean difference of 0.040, hence in this example the lower boundary is likely to be a node.If this is the case, then because 2 has been subtracted from the uncorrected Φ values in the data (panel (b)), then by referring back to the spectrogram in Figure 5, we see that the modeled harmonics must be plotted over their corresponding harmonic in the data, that is, the lowest harmonic that has been identified (lowest red line in Figure 5) is the second harmonic.As such, we can confirm that the fundamental is visible in this spectrogram at 0.2 Hz, when the high activity noise reduces at around 21:45-23:00 UT.
For the complete analysis of the 9-year data set, we have calculated the phase factors from the Eskdalemuir data and performed the relevant correction.The resulting corrections subtracted from each data point for the 9-year data set are shown in Table 2.The maximum value subtracted was 6, however for the majority of the cases, only 1 was subtracted.For six data points, no correction was made, meaning that the fundamental was identified by the U-net.These six points correspond to two separate days (two points on one day and four on another).Out of the 644 days, after the corrections to Φ were made, 6 days remained which had some values greater than 1.This is because they also had values of Φ below 1 (Φ crossed 0 in these cases), and phase factors cannot be negative, so nothing was subtracted from any of the values.
The resulting phase factors for IAR observed at Eskdalemuir for 2013-2021 are presented in Figure 9.We have data from 16:00 to 01:00 UT as the U-net input images include the time range 16:00 to 08:00 UT, and IAR with at least five harmonics visible were not identified after 01:00-02:00 UT.The mean error calculated is ±0.036.The green shaded bins from 18:00 to 00:00 UT are the UTs where the total number of data points is greater than the standard deviation of the total number of data points across all UTs, and so we determined that there were enough data points to conduct analysis on them.Note that each UT bin contains data from that hour up to but not including the next hour.
The mean value for the corrected data Φ values for 18:00-00:00 UT bins (Figure 9, green shaded bins) is 0.54, the median is 0.55, the minimum value is 0.049 and the maximum is 1.12.When calculating the mean Φ for each hour, the values change, increasing from 18:00 to 20:00, decreasing from 20:00 to 22:00 and then decreasing again, which is shown in Figure 9.However, the change is small, and it is less than the standard deviation of the data set.As such, we have not observed any significant trends with time.The average data value of Φ is slightly bigger than the uniform case of 0.5, and is similar to the result of Hebden et al. (2005), who analyzed 10 days of data, which was 0.6.Hebden et al. (2005) interpreted their values of Φ as being close to a node, and that the deviation from 0.5 could be due to and intermediate lower boundary condition that lies between a node and an antinode.
The modeled phase factors for Eskdalemuir from 2013 to 2021 are plotted in panels (a) and (b) in Figure 10.Panel (a) shows the case where the model has been run with a node as the lower boundary condition, and panel (b) shows the case where the model has been run with an antinode as the lower boundary condition.The models were run for each hour that corresponded to a time where there was data.Note that both the y-axes for the top row of this figure have the same scale.
When calculating the modeled phase factors for the antinode case, all of the values were above 1.However, in order to compare with the data, we subtracted a correction factor of 1 from these values.The mean Φ for the node model in the time range 18:00-00:00 UT (Figure 10a) is 0.62, and the Φ values range from 0.51 to 0.71, which is slightly above the uniform node case of 0.5.The mean Φ for the antinode model in the time range 18:00-00:00 UT (Figure 10b) is 0.49, and the values range from 0.19 to 0.64.Φ has a much greater range in the antinode case, and overlaps with the node values, but in general the values are below the uniform node case of 0.5, and do not average to 1, which the uniform model predicts.The data derived value sits in between the modeled node and antinode phase factors.From this, although the values of the determinant from the data are closest to the antinode value of 5 for the uniform model, when non-uniformity is considered they are closest to the node values.This agrees with our previous analysis.Potapov et al. (2022) found that their values determined from data were closest to the node case but there was deviation which they suggested was likely to be due to non-uniformity, which again our modeling agrees with.
The phase factors calculated for the data points where the fundamental was identified ranged between 0.70 and 0.80.The modeled node Φ values were between 0.67 and 0.71 for these times, and the modeled antinode Φ values were between 0.41 and 0.50.From this we can also infer that the lower boundary is a node during these times.

Conclusion
We have investigated how the non-uniformity of the Alfvén velocity profile in the ionosphere affects the harmonic frequencies of the IAR.We have extracted the individual harmonic frequencies over nine years of data from 2013 to 2021, using a U-net, where at least five harmonics were identified.We then modeled the individual frequencies of the IAR by solving a one-dimensional wave equation using non-uniform Alfvén velocity profiles modeled using the IRI and IGRF, for two lower boundary conditions of the electric fields, one being a node and one being an antinode.We calculated a quantity ζ, the percentage difference between the harmonic frequency separation and the fundamental frequency (Equation 9) for the data and both the node and antinode models.By analyzing this derived parameter, we find that the lower boundary condition is likely to be a node, which agrees with previous studies.We suggest that ζ may provide insight into the non-uniformity and shape of the Alfvén velocity profile of the ionosphere.Previous uniform models where the IAR frequencies scale as odd numbers (node case) or natural numbers (antinode case) used the so-called phase factor, Φ, to infer the lower boundary condition.For the node case, the uniform model leads to a phase factor of 0.5, and for the antinode case it is 1.Here the phase factor has been calculated from our extensive data set by using the IAR frequencies and corresponding harmonic number.For the Eskdalemuir data presented here, the mean value is 0.54.An analysis of the phase factor using the non-uniform model reveals that it is more strongly controlled by the non-uniformity of the ionosphere, rather than the lower boundary condition, and as such can take on a range of values in both the node and antinode case.

Figure 1 .
Figure 1.The first three harmonics of the electric field of the Alfvén wave, for a uniform Alfvén velocity.The magnetic field is shown by the arrow B. (a) Shows the case with a node in electric field at the lower boundary.(b) Shows the case with an antinode in electric field at the lower boundary.
Using observations from Mondy (mid-latitude in south-eastern Russia) for the years 2009-2019, Potapov et al. (2022) selected 12 days of data from each year and found the frequencies of the first five harmonics.Potapov et al. (

Figure 2 .
Figure 2. Wave structures for the electric field for the first three harmonics for a non-uniform Alfvén velocity.Panel (a) shows the case where the lower boundary is a node and (b) is the antinode case.

Figure 3 .
Figure 3.The spectrogram for the Eskdalemuir east-west (channel 2) coil for 15 September 2013 from 18:00 to 02:00 UT and 0-4 Hz.Harmonics identified by the U-net are plotted in red.

Figure 4 .
Figure 4. (a) Modeled Alfvén velocity profile and (b-d) first three harmonics for 20 June 2013 at 23:00 UT, for a node at the lower boundary.

Figure 5 .
Figure 5. Spectrogram for Eskdalemuir CH2 from 20:00 UT on 20 June 2013 to 04:00 UT on 21 June 2013.The U-net identified harmonics are plotted in red and the Runge-Kutta modeled harmonics are plotted in blue.The harmonic number for each harmonic is labeled.

Figure 6 .
Figure 6.The ζ value (see Equation9) for Eskdalemuir from 2013 to 2021.The data for each UT from 16:00 to 01:00 UT has been averaged.There are no data points outside the time range that is plotted.The first row shows the node case, and the second the antinode case.Panels(a, d)  show the data, with the black crosses representing the total number of data points before averaging in each UT (right hand y-axis).The total number of data points is 12,488.Green shaded bins correspond to UTs where there are a significant number of data points for analysis, with the number of data points being greater than the standard deviation of the total number of data points (σ = 860).Panels (b, e) show the RK4 models, for times where data are available.Panels (c, f) show the difference in ζ (model-data).The gray dashed line is plotted at 0. For each boxplot, the whiskers are 1.5× the interquartile range.The edges of the box are the first and third quartiles.The horizontal line in each box is the median, and any outliers are plotted as circles.

Figure 7 .
Figure 7. Mean Alfvén velocity profiles modeled by the International Reference Ionosphere and International Geomagnetic Reference Field, for UTs 18:00-01:00.They have been normalized by their minimum value.

Figure 8 .
Figure8.The left column shows the data over the time rage that the U-net identified Ionospheric Alfvén Resonator (approximately 22:00 on 20 June 2013-just after 01:00 UT).The right column shows the model, setting the lower boundary to be a node, from 22:00 to 02:00 UT to cover time where there is data.Panel (a) shows the harmonics which were extracted from the data using the U-net.(b) Shows the initial Φ values calculated from the data before any correction.(c) Shows the Φ values after they have been corrected.(d) Shows the modeled harmonics, (fundamental, second, third, fourth, and fifth), with a node at the lower boundary.From this, Φ is calculated, shown in panel (e).Panel (f) presents the difference between the modeled Φ and the corrected Φ from the data.The modeled Φ was interpolated to gain values at the same cadence as the data.The error in Φ is shown by the green shaded area.

Figure 9 .
Figure9.Phase factors calculated from the Eskdalemuir for 2013-2021.Time in UT is presented on the x-axis.Each UT contains the mean Φ for that hour for each day that we have data.The gray dashed line at 0.5 represents Φ for the node case for a uniform Alfvén velocity profile.The boxplot is plotted in the same way as Figure6.