A method for predicting the propagation characteristics of low-dispersion very low frequency radio waves

Abstract

The present application belongs to the field of radio wave propagation theory and electromagnetic field numerical calculation technology, and discloses a low-dispersion very low frequency radio wave propagation characteristic prediction method.The present application method aims at the defects of large numerical dispersion, low calculation efficiency and high resource occupation of the existing FDTD method in long-distance, complex layered medium scene, and realizes high-precision and high-efficiency prediction of the very low frequency radio wave propagation characteristic in seawater-air-ionosphere medium by combining HIE-FDTD method and M(2,4)FDTD method through three-layer processing of the whole calculation area.The present application method solves the contradiction between reducing numerical dispersion error and reducing computer consumption of the traditional FDTD method, can greatly reduce computer memory resource occupation and improve calculation speed while ensuring prediction accuracy, is suitable for electromagnetic propagation simulation of large-size, complex medium and terrain, and has obvious advantages in phase delay prediction of very low frequency ground wave propagation.

Description

Technical Field

[0001] This invention belongs to the field of radio wave propagation theory and electromagnetic field numerical calculation technology, specifically relating to a method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves. Background Technology

[0002] Very low frequency (VLF) radio waves, due to their long propagation distance, strong penetration, and excellent anti-interference performance, have irreplaceable engineering application value in underwater communication, global navigation, and ionospheric detection. Their propagation scenarios involve complex layered media of seawater, air, and the ionosphere, each with significantly different characteristics: seawater is a highly conductive medium, and its electromagnetic parameters drastically affect propagation loss; the air layer is a relatively homogeneous medium and the main channel for long-distance propagation of radio waves; the ionosphere is a plasma medium with strong dispersion and anisotropy, and its parameters fluctuate dynamically due to factors such as solar activity and geomagnetic variations. Therefore, predicting the propagation characteristics of VLF radio waves in complex layered media requires simultaneously addressing three core needs: accurate simulation of highly conductive media, low-dispersion propagation in homogeneous media, and handling of dispersion and anisotropy in plasma media. This presents a significant technical challenge in the field of electromagnetic field numerical computation.

[0003] The Finite-Difference Time-Domain (FDTD) method, starting directly from Maxwell's curl equations, can sample electromagnetic fields in the time domain and comprehensively capture the physical processes of electromagnetic propagation in complex media, making it the mainstream method for very low frequency (VLF) radio wave propagation prediction. However, traditional FDTD methods have two major drawbacks: First, significant numerical dispersion errors occur, and the second-order spatial difference scheme easily leads to deviations between the numerical wavenumber and the physical wavenumber. The error accumulates and amplifies with the propagation distance, severely affecting the prediction accuracy of key parameters such as phase delay. Second, due to the Courant-Friedrichs-Lewy (CFL) stability condition, the time step is determined by the minimum spatial step. To ensure the simulation accuracy of highly conductive seawater and complex ionospheres, fine mesh partitioning is required, resulting in a surge in iterations and low computational efficiency, making it difficult to balance high accuracy and high efficiency.

[0004] To address the aforementioned shortcomings, scholars both domestically and internationally have proposed two types of improved FDTD methods: one type is the low-dispersion high-order FDTD method, represented by the modified second-order-in-time and fourth-order-in-space FDTD (M(2,4)FDTD) method. This method significantly reduces the numerical dispersion error in all directions by discretizing the spatial partial derivatives using the fourth-order difference and introducing a correction factor, while maintaining high accuracy even with coarse mesh partitioning. However, this method is still limited by the stability condition of CFL, and its high-order difference scheme is difficult to directly adapt to the dispersion characteristics of non-uniform media and plasmas, and its accuracy is prone to decrease at the seawater-air-ionospheric interface. Another type is the unconditionally stable FDTD method, represented by the Hybrid Implicit-Explicit FDTD (HIE-FDTD) method. By splitting the magnetic field components along different directions and using a hybrid explicit-implicit iterative solution, it breaks through the limitation of the time step by the CFL stability condition and can significantly improve computational efficiency. However, the HIE-FDTD method still uses a second-order spatial difference scheme, which has insufficient numerical dispersion suppression capability and lacks specific adaptation to the anisotropy and dispersion characteristics of ionospheric plasma, making it difficult to meet the high-precision prediction requirements for very low frequency electromagnetic wavelength distance propagation.

[0005] Currently, no single improvement method can simultaneously adapt to the propagation characteristics of complex layered media such as seawater, air, and the ionosphere. Traditional FDTD methods struggle to balance accuracy and efficiency. The M(2,4)FDTD method exhibits low dispersion but is limited by stability conditions and is unsuitable for complex media. The HIE-FDTD method is unconditionally stable but suffers from large dispersion errors and struggles to handle plasma characteristics. Therefore, a highly efficient improved FDTD method is needed, combining low dispersion, unconditional stability, and adaptability to complex layered media, to achieve high-precision and rapid prediction of propagation delay, field strength, and other characteristics of very low frequency radio waves, thereby providing technical support for related engineering applications. Summary of the Invention

[0006] The purpose of this invention is to propose a method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves. This method overcomes the shortcomings of existing FDTD methods, such as large numerical dispersion, low computational efficiency, and high resource consumption in long-distance and complex layered media scenarios. Furthermore, it can achieve high-precision and high-efficiency prediction of the propagation characteristics of very low-frequency radio waves in seawater-air-ionospheric media.

[0007] To achieve the above objectives, the present invention adopts the following technical solution:

[0008] A method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves includes the following steps:

[0009] Step 1. Input the model file. The input model file includes the grid parameters of the computational domain, the electrical parameters of the radio wave propagation path, the parameters of the absorbing boundary PML, the parameters of the excitation source, the observation points, and the time settings.

[0010] Step 2. Set and initialize the parameters of the input model file;

[0011] Step 3. Add incentive sources;

[0012] Step 4. Update the upper and lower layer regions using the HIE-FDTD method. The electric field components in the directional direction are updated in the intermediate layer region using the M(2,4)FDTD method. The electric field component in the direction;

[0013] Step 5. Update the upper and lower layer regions using the HIE-FDTD method. The electric field components in the directional direction are updated in the intermediate layer region using the M(2,4)FDTD method. The electric field component in the direction;

[0014] Step 6. Update the magnetic field components of the upper and lower regions using the HIE-FDTD method, and update the magnetic field components of the middle region using the M(2,4)FDTD method;

[0015] Step 7. Update the incentive source;

[0016] Step 8. Determine if the current running time step is equal to the preset running time step; if yes, proceed to step 9; otherwise, update the current running time step. The value is Then proceed to step 4;

[0017] Step 9. Calculate the electromagnetic field components at the observation point and output the field strength at the sampling point.

[0018] Furthermore, based on the aforementioned method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves, this invention also proposes a computer device, which includes a memory and one or more processors.

[0019] The memory stores executable code, and when the processor executes the executable code, it implements the steps of the low-dispersion very low frequency radio wave propagation characteristic prediction method described above.

[0020] Furthermore, based on the aforementioned method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves, this invention also proposes a computer-readable storage medium storing a program thereon; when executed by a processor, this program is used to implement the steps of the aforementioned method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves.

[0021] The present invention has the following advantages:

[0022] As described above, this invention discloses a method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves. This method integrates the fourth-order higher-order difference scheme and correction factor of the M(2,4)FDTD method into the framework of the HIE-FDTD method. While significantly suppressing omnidirectional numerical dispersion errors, it inherits the hybrid explicit-implicit solution characteristics of the HIE-FDTD method, breaks through the CFL stability condition limitation, and can flexibly set the time step, effectively improving computational efficiency. This invention employs a three-layer partitioning calculation strategy, matching differentiated methods to the different medium characteristics of seawater, air, and the ionosphere. This leverages the advantages of low dispersion and high efficiency of the method while ensuring the simulation accuracy of complex layered regions, solving the problem that a single FDTD method is difficult to adapt to complex propagation scenarios in the seawater-air-ionosphere region. The overall solution significantly reduces computer memory usage and improves computing speed while ensuring prediction accuracy. It is suitable for very low frequency radio wave propagation simulation in electrically large, complex media and terrains. It can be widely used in fields such as very low frequency ground wave propagation, underwater communication, ionospheric detection and global navigation. It shows significant advantages in the prediction of very low frequency radio wave propagation phase delay and field strength, and has important engineering application value. Attached Figure Description

[0023] Figure 1 This is a flowchart of a method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves in an embodiment of the present invention.

[0024] Figure 2 This is a schematic diagram illustrating the principle of the low-dispersion very low frequency radio wave propagation characteristic prediction method in an embodiment of the present invention.

[0025] Figure 3 This is a schematic diagram of the grid arrangement and field distribution of the upper, middle, and lower layers of the computational region in an embodiment of the present invention.

[0026] Figure 4 This diagram illustrates a comparison of the field strength prediction results of the method of this invention and the traditional method under the same underwater source setup.

[0027] Figure 5 This is a schematic diagram comparing the attenuation results of very low frequency radio wave field strength at different frequencies and observation point depths with the method of this invention and the traditional method. Detailed Implementation

[0028] The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments:

[0029] Example 1

[0030] This invention specifically proposes an M(2,4)-HIE-FDTD method for predicting the propagation characteristics of low-dispersion very low frequency radio waves. It addresses the numerical simulation problem of very low frequency radio wave propagation in a complex layered medium of seawater-air-ionosphere by performing three-layer processing on the entire computational domain, combining the HIE-FDTD method and the M(2,4)FDTD method.

[0031] Specifically, this invention addresses the physical process of very low frequency (VLF) radio waves originating from an underwater source, propagating through seawater and the air layer, and ultimately interacting with the ionosphere. It rationally layers the computational domain. In the lowest water layer, the HIE-FDTD method is employed, whose implicit characteristics robustly handle the high conductivity of seawater. In the middle air layer, where electromagnetic waves travel long distances, it is the main region for the accumulation of numerical dispersion errors; therefore, the low-dispersion M(2,4)FDTD method is used to effectively suppress dispersion. In the uppermost ionosphere layer, the HIE-FDTD method is again employed to facilitate the introduction of plasma parameters to accurately describe ionospheric effects. Field quantity exchange occurs directly between layers through grid interfaces, as illustrated in the schematic diagram. Figure 2 As shown in the figure This indicates the location of the observation point.

[0032] This invention leverages the advantages of the HIE-FDTD method in handling complex dispersive media and the M(2,4)FDTD method in handling low dispersion, simplifying the algorithm structure and achieving high-precision, stable numerical prediction of the propagation characteristics of very low frequency (VLF) radio waves in complex three-layer media. By performing three-layer processing on the seawater-air-ionospheric composite computational domain involved in VLF radio wave propagation, and integrating the fine-grid HIE-FDTD method suitable for dispersive anisotropic media with the low-dispersion coarse-grid M(2,4)FDTD method, the method can predict the propagation delay, field strength, and other characteristics of VLF radio waves, improving computational speed and reducing computer memory usage while ensuring prediction accuracy.

[0033] The M(2,4)FDTD method is partly derived from the integral form of Maxwell's equations, including two types of integral loops: small loops and large loops, and introduces loop correction coefficients. and The HIE-FDTD method achieves omnidirectional low dispersion through fourth-order spatial difference and weighted optimization. Part of the method overcomes the Courant stability condition by splitting the magnetic field components and employing a hybrid explicit / implicit solution, enabling autonomous selection of the time step and efficient computation. The computational domain employs a three-layer mesh layout and field distribution design, as illustrated in the diagram below. Figure 3 As shown, the ratio of coarse to fine mesh between each layer is set to an odd number to facilitate the transfer of field quantities.

[0034] The method of the present invention will be described in detail below.

[0035] likeFigure 1 As shown, the method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves specifically includes the following steps:

[0036] Step 1. Input the model file. The input model file includes the grid parameters of the computational domain, the electrical parameters of the radio wave propagation path, the parameters of the absorbing boundary PML, the parameters of the excitation source, the observation points, and the time settings.

[0037] In this embodiment, the process of generating the model file includes layering the computational domain, setting the mesh subdivision accuracy for each domain, setting the propagation path electrical parameter model, and setting observation points. Step 1 specifically includes:

[0038] The computational region is divided into layers, with the size of the lower layer being [size missing]. The size of the intermediate layer area is The size of the upper region is ,in and Representing the lower layer regions respectively direction and Number of grids in the direction, and Representing the intermediate layer regions respectively direction and Number of grids in the direction, and Representing the upper region respectively direction and The number of grid cells in each direction.

[0039] Set the mesh subdivision accuracy for each region, and set the mesh spatial step size for the lower region to be [value missing]. , ,in and Representing the lower layer regions respectively direction and The grid size in the direction; let the grid space step size of the intermediate layer region be... , ,in and Representing the intermediate layer regions respectively direction and The grid size in the direction; let the grid space step size of the upper region be... , ,in and Representing the upper region respectively direction and Grid size in different directions; uniform time step used in all layers. .

[0040] The electrical parameters that determine the propagation path of radio waves include the relative permittivity in free space. Electrical conductivity magnetic permeability Relative permittivity of the ionosphere With conductivity Relative permittivity of seawater region With conductivity .

[0041] The parameters for setting the absorbing boundary PML include the number of absorbing boundary layers. Maximum coordinate stretching factor Maximum polarization loss factor Maximum conductivity ,in .

[0042] Time settings include the current running time step. Time step and preset running time steps ,in Right now .

[0043] Configure the parameters of the excitation sources, including the number of excitation sources. Location of the excitation source, type of excitation source, and width of the Gaussian pulse. Gaussian pulse delay Signal frequency and amplitude The types of excitation sources include soft sources and hard sources.

[0044] Set up observation points, including the number of observation points. ;set up The starting grid position of the direction is , The direction ends at the grid position. , The starting grid position of the direction is , The direction ends at the grid position. Output field types include Electric field components in the direction , Electric field components in the direction and magnetic field components .

[0045] In this embodiment, the field components specifically include electric field components and magnetic field components.

[0046] Among them, electric field components Includes lower layer area, middle layer area, and upper layer area. Electric field components in the direction , , Electric field components Includes lower layer area, middle layer area, and upper layer area. Electric field components in the direction , , Magnetic field components Including magnetic field components in the lower, middle, and upper regions. , , .

[0047] This step involves differentially partitioning the computational domain into three layers and uniformly setting the time step. This allows the method of the present invention to fully adapt to the differences in electromagnetic properties of the seawater-air-ionosphere layered media. While ensuring the simulation accuracy of each region, it lays the foundation for the efficient execution of subsequent hybrid methods and improves the overall computational stability and adaptability.

[0048] Step 2. Set and initialize the parameters of the input model file.

[0049] In this embodiment, step 2 specifically includes:

[0050] The electric field components of the entire region , Magnetic field components electric field component coefficient , , , Magnetic field component coefficient , intermediate variables , , , , , intermediate variable coefficients , , , , , , , , and the field strength at the sampling point Initialize to zero.

[0051] in, Indicates the direction of rotation about the axis in cylindrical coordinates. Indicates the geometric correction term in cylindrical coordinates; intermediate variable Intermediate variables including lower and upper regions and Intermediate variables Intermediate variables including lower and upper regions and Intermediate variables Intermediate variables including lower and upper regions and Intermediate variables Intermediate variables including lower and upper regions and .

[0052] Initialize the model parameters in all meshes, and set the free space relative permittivity. Initialized to 1, earth conductivity Initialized to 0, permeability Initialize to 1.

[0053] Based on the path information in the input model file, assign values ​​to the relative permittivity, conductivity, and permeability of each region's mesh. Set the parameters of the absorbing boundary PML. , , ,in , , Calculate using the following formula:

[0054] (1)

[0055] in, Locations of PML layers and non-PML sections; ; The thickness of the PML absorption boundary.

[0056] Take the integer value, and ; ; according to set up, , , For grid step size, , The wavelength of the excitation source.

[0057] This step completes the initialization of the total field quantity and the medium parameters, and realizes the precise configuration of the PML parameters and the three-layer medium parameters. It can ensure that the initial field is undisturbed, the electromagnetic parameters of the medium interface are continuous, avoid the accumulation of initial errors, and ensure the numerical stability and reliable calculation results in subsequent iterations.

[0058] Step 3. Add incentive sources.

[0059] In this embodiment, step 3 specifically includes:

[0060] Add an excitation source as a single-frequency sinusoidal source, with its current source excitation. Represented as:

[0061] (2)

[0062] in, .

[0063] This step supports various typical very low frequency excitations, such as single-frequency sinusoidal sources, and can flexibly adapt to the radiation source requirements of different engineering scenarios, ensuring accurate excitation loading and high waveform fidelity, and providing reliable input for accurate prediction of very low frequency radio wave propagation characteristics.

[0064] Step 4. Update the upper and lower fine mesh regions using the HIE-FDTD method. Electric field components in the direction , The M(2,4)FDTD method is used to update the coarse mesh region of the intermediate layer. Electric field components in the direction .

[0065] In this embodiment, step 4 specifically includes:

[0066] Step 4.1. Use the M(2,4)FDTD method to analyze the intermediate layer region. Electric field components in the direction To update, the specific update formula is as follows:

[0067] (3)

[0068] in, for In time ,space The value at; Indicates time difference; and Representing the calculation area respectively direction and Spatial position of direction; and They represent direction and Spatial difference in direction; Including intermediate variables , , , , , , electric field components , Magnetic field components ; ,in , , These are used to represent the lower layer, middle layer, and upper layer, respectively.

[0069] (4)

[0070] in, and They represent the electric field component coefficients, respectively. and The value at the m-th computational grid; and represents the relative permittivity and conductivity at the m-th computational grid, respectively.

[0071] when Right now hour, This is a singularity. Furthermore, as the radius decreases, the cell volume decreases, and according to stability requirements, the time step becomes very small. Therefore, in The surrounding area requires special handling. This problem can be solved using Ampere's law.

[0072] exist Ampere's law is used to process the electromagnetic field components in the intermediate layer region. Update:

[0073] (5)

[0074] Step 4.2. Use the HIE-FDTD method to analyze the lower layer region. Electric field components in the direction To update, the specific update formula is as follows:

[0075] (6)

[0076] in, and They represent the electric field component coefficients, respectively. and The value at the m-th computational grid.

[0077] (7)

[0078] exist Ampere's law is used to process the electromagnetic field components in the lower region. Update:

[0079] (8)

[0080] in:

[0081] (9)

[0082] in, express exist Direction first The values ​​at each grid point express exist Direction first The values ​​at each grid point Including intermediate variable coefficients , , , , , , , , .

[0083] (10)

[0084] in, express The PML coordinate stretching factor at the m-th grid in the direction. express The PML polarization loss factor at the m-th grid in the direction. express PML conductivity at the m-th grid point in the direction; It represents the dielectric constant in a vacuum.

[0085] (11)

[0086] (12)

[0087] in, express The PML coordinate stretching correction factor at the m-th grid in the direction. express The PML loss distribution coefficient at the m-th grid in the direction. and They are respectively:

[0088] (13)

[0089] in, This indicates the starting coordinates of the PML layer.

[0090] Step 4.3. Use the HIE-FDTD method to analyze the upper region. Electric field components in the direction To update, the specific update formula is as follows:

[0091] (14)

[0092] in:

[0093] (15)

[0094] in, For collision frequency, The frequency of the Drude model; That is, the relative permittivity at infinite frequency.

[0095] (16)

[0096] exist Ampere's law is used to update the electromagnetic field components in the upper region:

[0097] (17)

[0098] in:

[0099] (18)

[0100] (19)

[0101] This step employs a hybrid update strategy that uses the HIE-FDTD method for the upper and lower layers and the M(2,4)FDTD method for the middle layer. This strategy leverages the robustness of the HIE-FDTD method to handle highly conductive seawater and the dispersive ionosphere, while also suppressing long-distance numerical dispersion in the air region through the M(2,4)FDTD method. This approach improves computational efficiency while maintaining accuracy.

[0102] Step 5. Update the upper and lower fine mesh regions using the HIE-FDTD method. Electric field components in the direction , The M(2,4)FDTD method is used to update the coarse mesh region of the intermediate layer. Electric field components in the direction .

[0103] In this embodiment, step 5 specifically includes:

[0104] Step 5.1. Use the M(2,4)FDTD method to analyze the intermediate layer region. Electric field components in the direction Update the electric field components. It is necessary to first determine the magnetic field components of the upper region. and the magnetic field components of the lower region Magnetic field components transmitted to the intermediate layer region The specific update formula is as follows:

[0105] (20)

[0106] in, This is the loop coefficient.

[0107] Step 5.2. Use the HIE-FDTD method to analyze the lower layer region. Electric field components in the direction Update the electric field components. It is necessary to first separate the magnetic field components of the intermediate layer region. Magnetic field components transmitted to the lower region The specific update formula is as follows:

[0108] (twenty one)

[0109] in, Indicates the magnetic permeability in vacuum; This represents the relative permeability at the m-th grid.

[0110] (twenty two)

[0111] (twenty three)

[0112] in, express The PML coordinate stretching factor at the m-th grid in the direction. express The PML polarization loss factor at the m-th grid in the direction. express PML conductivity at the m-th grid in the direction.

[0113] (twenty four)

[0114] (25)

[0115] Step 5.3. Use the HIE-FDTD method to analyze the upper region. Electric field components in the direction To update, the specific update formula is as follows:

[0116] (26)

[0117] in:

[0118] (27)

[0119] (28)

[0120] (29)

[0121] (30)

[0122] This step continues to use the partitioned differential electric field update mechanism, and uses the HIE-FDTD method implicit iteration and the M(2,4)FDTD method to update the electric field. This can break through the CFL condition limitation and significantly improve computational efficiency while ensuring numerical stability.

[0123] Step 6. Update the magnetic field components of the upper and lower fine mesh regions using the HIE-FDTD method. , The magnetic field components of the intermediate coarse grid region are updated using the M(2,4)FDTD method. .

[0124] In this embodiment, step 6 specifically includes:

[0125] Step 6.1. Based on the model file defined in Step 1, use the HIE-FDTD method to analyze the magnetic field components in the lower region. To update, the specific update formula is as follows:

[0126] (31)

[0127] Step 6.2. Use the HIE-FDTD method to analyze the magnetic field components in the upper region. To update, the specific update formula is as follows:

[0128] (32)

[0129] Step 6.3. Use the M(2,4)FDTD method to analyze the magnetic field components in the intermediate layer region. To update, the specific update formula is as follows:

[0130] (33)

[0131] in, and They represent the magnetic field component coefficients, respectively. and The value at the m-th computational grid; This is the loop coefficient.

[0132] (34)

[0133] This step employs a hybrid approach to update the magnetic field, inheriting the unconditional stability advantage of the HIE-FDTD method and the low dispersion advantage of the M(2,4)FDTD method. This achieves self-consistency in electromagnetic field coupling iteration, ensuring continuous field quantities and energy conservation, and is suitable for electrically large-size very low-frequency propagation simulations.

[0134] Step 7. Update the incentive source.

[0135] In this embodiment, the process of updating the excitation sources requires updating all sources according to the types of excitation sources set in the model file in step 1. This step can update the excitation sources in real time and accurately, ensuring stable injection of source signals.

[0136] Step 8. Loop and check the termination condition: Check if the current running time step is equal to the preset running time step, i.e., check if the maximum number of iterations has been reached; if yes, go to step 9; otherwise, update the current running time step. The value is Then proceed to step 4.

[0137] This step enables automatic control of the iterative process, with simpler computational logic and lower resource consumption.

[0138] Step 9. Calculate and output observations at the observation point: Calculate the electromagnetic field components at the observation point and output the field strength at the sampling point.

[0139] In this embodiment, step 9 specifically includes:

[0140] When the excitation source is a single-frequency sinusoidal source, the peak amplitude of the vertical electric field at the observation point is extracted using the peak monitoring method. To extract the peak amplitude of the vertical electric field at the observation point. Normalized to standard radiation conditions, when the antenna radiated power is At that time, the corresponding standard field strength for:

[0141] (35)

[0142] in, Indicates the length of the radiating element. This refers to the peak amplitude of the vertical electric field at the observation point directly extracted from the simulation results; That is, the corresponding standard radiant power and radiation element length The theoretical field strength, This is the field strength at the sampling point.

[0143] This step can quickly output the field strength results at the observation points, meeting the needs of direct engineering applications.

[0144] In addition, to verify the effectiveness of the method proposed in this invention, the following specific experiments are also provided:

[0145] Experiment 1 is a verification experiment comparing the accuracy of the M(2,4)-HIE-FDTD method with the standard FDTD method.

[0146] The M(2,4)-HIE-FDTD method is the method of this invention, while the standard FDTD method is the traditional method.

[0147] The size of the computation region is: Traditional methods direction and The mesh size in each direction is respectively and The time step is The method of this invention direction and The mesh size in each direction is respectively and The time step is Both the air layer and the ionosphere are ideal dielectrics with a relative permittivity of . electrical conductivity The seawater layer is a high-loss medium with a relative permittivity. electrical conductivity In the experiment, the signal sources for both methods were set at a depth of -3m in the seawater layer, and the extraction depth of the observation point was consistent with the source depth, i.e. .

[0148] Figure 4 This diagram illustrates a comparison of the field strength prediction results between the method of this invention and a traditional method under the same water source setting. Figure 4 It can be seen that the field strength attenuation curves of the two methods almost overlap, and the prediction results have no significant error. This proves that the method of the present invention, in ideal media environments such as air and the ionosphere, can maintain the same high accuracy as the traditional fine-grid FDTD method under coarse meshing by combining the low dispersion scheme of the M(2,4)FDTD method with the unconditionally stable method of HIE-FDTD, thus verifying the accuracy of the method of the present invention in the basic media environment.

[0149] Experiment 2 is a verification experiment of the propagation characteristics of very low frequency radio waves at different frequencies and observation point depths.

[0150] The computing area is expanded to a power-large size scene and set as follows: The method of this invention is used to predict the propagation characteristics of very low frequency radio waves. direction and The mesh size in each direction is maintained as follows: and The time step is The electromagnetic parameters remained consistent with those in Experiment 1. In the experiment, the signal source was set at a depth of -3m in the seawater layer, and two typical frequencies in the very low frequency band, 18.2kHz and 24.5kHz, were selected. Two different observation point depths were also set. =-3m and =-6m, and simulations were performed for comparison.

[0151] Figure 5 This diagram illustrates a comparison of the attenuation results of very low frequency (VLF) radio wave field strength calculated along the horizontal direction using the method of this invention at different frequencies and observation point depths. -3m and -6m represent the extraction depth of the field strength observation point within the seawater layer. Figure 5 The two curves located at the top correspond to observation points at a depth of -3m. Figure 5 The two curves located at the bottom correspond to observation points at a depth of -6m. (From...) Figure 5 It can be seen that, at the same observation point depth, the field strength attenuation of the 24.5kHz signal is less than that of the 18.2kHz signal, which conforms to the frequency loss law of very low frequency radio wave propagation. At the same operating frequency, the field strength value at an observation point depth of -6m is generally lower than the corresponding result at -3m, indicating that the deeper the observation point extraction depth, the more significant the cumulative absorption loss of radio waves by the highly conductive seawater medium, and the more significant the field strength attenuation. The field strength curve generated by the method of this invention exhibits a smooth periodic fluctuation characteristic with propagation distance, which is highly consistent with the waveguide interference effect formed by the seawater-air-ionospheric layer interface, proving that the method of this invention can accurately handle complex propagation physical processes including ideal dielectric layers.

[0152] The method of this invention overcomes the contradiction between accuracy and efficiency in the traditional FDTD method. Compared with the traditional FDTD method, the method of this invention achieves accurate prediction of very low frequency radio wave propagation characteristics even with coarse grids by unconditionally and stably fusing the low dispersion scheme of the M(2,4)FDTD method with the HIE-FDTD method. At the same time, relying on differentiated grid partitioning, it expands the electrically large computational region, improves computational efficiency without sacrificing accuracy, and has higher engineering applicability.

[0153] The method of this invention employs a fine-grid HIE-FDTD method at the bottom layer to accurately simulate the high conductivity of seawater and complex underwater environments. The middle layer uses a coarse-grid M(2,4) FDTD method to significantly reduce numerical dispersion errors while maintaining accuracy. The top layer uses a fine-grid HIE-FDTD method to effectively handle the dispersion and anisotropy characteristics of the ionospheric plasma medium. High-order difference schemes, correction factors, and inter-regional field transfer techniques are utilized to achieve high-precision and rapid prediction of electromagnetic propagation problems. This method resolves the contradiction between reducing numerical dispersion errors and lowering computer consumption in traditional FDTD methods. While maintaining prediction accuracy, it significantly reduces computer memory usage and increases computational speed, making it suitable for electromagnetic propagation simulations of electrically large scales, complex media, and terrains. It has a particularly significant advantage in predicting the phase delay of very low frequency ground wave propagation.

[0154] Example 2

[0155] This embodiment 2 describes a computer device that includes a memory and one or more processors.

[0156] The memory stores executable code, which, when executed by the processor, is used to implement the steps of the low-dispersion very low frequency radio wave propagation characteristic prediction method in Embodiment 1 above.

[0157] In this embodiment, the computer device can be any device or apparatus with data processing capabilities, and will not be described in detail here.

[0158] Example 3

[0159] This embodiment 3 describes a computer-readable storage medium storing a program that, when executed by a processor, implements the steps of a method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves.

[0160] The computer-readable storage medium can be an internal storage unit of any device or apparatus with data processing capabilities, such as a hard disk or memory, or an external storage device of any device with data processing capabilities, such as a plug-in hard disk, smart media card (SMC), SD card, flash card, etc.

[0161] Of course, the above description is only a preferred embodiment of the present invention. The present invention is not limited to the above-described embodiments. It should be noted that any equivalent substitutions or obvious modifications made by those skilled in the art under the guidance of this specification fall within the scope of this specification and should be protected by the present invention.

Claims

1. A method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves, characterized in that, The steps include the following: Step 1. Input the model file. The input model file includes the grid parameters of the computational domain, the electrical parameters of the radio wave propagation path, the parameters of the absorbing boundary PML, the parameters of the excitation source, the observation points, and the time settings. Step 2. Set and initialize the parameters of the input model file; Step 3. Add incentive sources; Step 4. Update the upper and lower layer regions using the HIE-FDTD method. The electric field components in the directional direction are updated in the intermediate layer region using the M(2,4)FDTD method. The electric field component in the direction; Step 5. Update the upper and lower layer regions using the HIE-FDTD method. The electric field components in the directional direction are updated in the intermediate layer region using the M(2,4)FDTD method. The electric field component in the direction; Step 6. Update the magnetic field components of the upper and lower regions using the HIE-FDTD method, and update the magnetic field components of the middle region using the M(2,4)FDTD method; Step 7. Update the incentive source; Step 8. Determine if the current running time step is equal to the preset running time step; if yes, proceed to step 9; otherwise, update the current running time step. The value is Then proceed to step 4; Step 9. Calculate the electromagnetic field components at the observation point and output the field strength at the sampling point; Step 1 specifically involves: The computational region is divided into layers, with the size of the lower layer being [size missing]. The size of the intermediate layer area is The size of the upper region is ,in and Representing the lower layer regions respectively direction and Number of grids in the direction, and Representing the intermediate layer regions respectively direction and Number of grids in the direction, and Representing the upper region respectively direction and Number of grid cells in the direction; Set the mesh subdivision accuracy for each region, and set the mesh spatial step size for the lower region to be [value missing]. , ,in and Representing the lower layer regions respectively direction and The grid size in the direction; let the grid space step size of the intermediate layer region be... , ,in and Representing the intermediate layer regions respectively direction and The grid size in the direction; let the grid space step size of the upper region be... , ,in and Representing the upper region respectively direction and Grid size in different directions; uniform time step used in all layers. ; The electrical parameters that determine the propagation path of radio waves include the relative permittivity in free space. Electrical conductivity magnetic permeability Relative permittivity of the ionosphere With conductivity Relative permittivity of seawater region With conductivity ; The parameters for setting the absorbing boundary PML include the number of absorbing boundary layers. Maximum coordinate stretching factor Maximum polarization loss factor Maximum conductivity ,in ; Time settings include the current running time step. Time step and preset running time steps ,in Right now ; Configure the parameters of the excitation sources, including the number of excitation sources. Location of the excitation source, type of excitation source, and width of the Gaussian pulse. Gaussian pulse delay Signal frequency and amplitude The types of excitation sources include soft sources and hard sources; Set up observation points, including the number of observation points. ;set up The starting grid position of the direction is , The direction ends at the grid position. , The starting grid position of the direction is , The direction ends at the grid position. ; Output field quantity types include Electric field components in the direction , Electric field components in the direction and magnetic field components ; Among them, electric field components Includes lower layer area, middle layer area, and upper layer area. Electric field components in the direction , , Electric field components Includes lower layer area, middle layer area, and upper layer area. Electric field components in the direction , , Magnetic field components Including magnetic field components in the lower, middle, and upper regions. , , .

2. The method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves according to claim 1, characterized in that, Step 2 specifically involves: The electric field components of the entire region , Magnetic field components electric field component coefficient , , , Magnetic field component coefficient , intermediate variables , , , , , intermediate variable coefficients , , , , , , , , and the field strength at the sampling point Initialize to zero; in, Indicates the direction of rotation about the axis in cylindrical coordinates. Indicates the geometric correction term in cylindrical coordinates; intermediate variable Intermediate variables including lower and upper regions and Intermediate variables Intermediate variables including lower and upper regions and Intermediate variables Intermediate variables including lower and upper regions and Intermediate variables Intermediate variables including lower and upper regions and ; Initialize the model parameters in all meshes, and set the free space relative permittivity. Initialized to 1, earth conductivity Initialized to 0, permeability Initialize to 1; Based on the path information in the input model file, assign values ​​to the relative permittivity, conductivity, and permeability of each region's grid; set the parameters of the absorbing boundary PML. , , ,in , , Calculate using the following formula: (1) in, Locations of PML layers and non-PML sections; ; The thickness of the PML absorption boundary; Take the integer value, and ; ; according to set up, , , For grid step size, , The wavelength of the excitation source.

3. The method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves according to claim 2, characterized in that, Step 3 specifically involves: Add an excitation source as a single-frequency sinusoidal source, with its current source excitation. Represented as: (2) in, .

4. The method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves according to claim 3, characterized in that, Step 4 specifically involves: The M(2,4)FDTD method was used to analyze the intermediate layer region. Electric field components in the direction To update, the specific update formula is as follows: (3) in, for In time ,space The value at; Indicates time difference; and Representing the calculation area respectively direction and Spatial position of direction; and They represent direction and Spatial difference in direction; Including intermediate variables , , , , , , electric field components , Magnetic field components ; ,in , , These are used to represent the lower layer region, the middle layer region, and the upper layer region, respectively. (4) in, and They represent the electric field component coefficients, respectively. and The value at the m-th computational grid; and Let represent the relative permittivity and conductivity at the m-th computational grid, respectively; exist Ampere's law is used to process the electromagnetic field components in the intermediate layer region. Update: (5) The HIE-FDTD method was used to analyze the lower-level region. Electric field components in the direction To update, the specific update formula is as follows: (6) in, and They represent the electric field component coefficients, respectively. and The value at the m-th computational grid; (7) exist Ampere's law is used to process the electromagnetic field components in the lower region. Update: (8) in: (9) in, express exist Direction first The values ​​at each grid point express exist Direction first The values ​​at each grid point Including intermediate variable coefficients , , , , , , , , ; (10) in, express The PML coordinate stretching factor at the m-th grid in the direction. express The PML polarization loss factor at the m-th grid in the direction. express PML conductivity at the m-th grid point in the direction; It represents the dielectric constant in a vacuum; (11) (12) in, express The PML coordinate stretching correction factor at the m-th grid in the direction. express The PML loss distribution coefficient at the m-th grid in the direction. and They are respectively: (13) in, Indicates the starting coordinates of the PML layer; The HIE-FDTD method was used to analyze the upper region. Electric field components in the direction To update, the specific update formula is as follows: (14) in: (15) in, For collision frequency, The frequency of the Drude model; (16) exist Ampere's law is used to update the electromagnetic field components in the upper region: (17) in: (18) (19)。 5. The method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves according to claim 4, characterized in that, Step 5 specifically involves: The M(2,4)FDTD method was used to analyze the intermediate layer region. Electric field components in the direction Update the electric field components. It is necessary to first determine the magnetic field components of the upper region. and the magnetic field components of the lower region Magnetic field components transmitted to the intermediate layer region The specific update formula is as follows: (20) in, Loop coefficient; The HIE-FDTD method was used to analyze the lower-level region. Electric field components in the direction Update the electric field components. It is necessary to first separate the magnetic field components of the intermediate layer region. Magnetic field components transmitted to the lower region The specific update formula is as follows: (21) in, Indicates the magnetic permeability in vacuum; This represents the relative permeability at the m-th grid. (22) (23) in, express The PML coordinate stretching factor at the m-th grid in the direction. express The PML polarization loss factor at the m-th grid in the direction. express PML conductivity at the m-th grid point in the direction; (24) (25) The HIE-FDTD method was used to analyze the upper region. Electric field components in the direction To update, the specific update formula is as follows: (26) in: (27) (28) (29) (30)。 6. The method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves according to claim 5, characterized in that, Step 6 specifically involves: The HIE-FDTD method was used to analyze the magnetic field components in the lower region. To update, the specific update formula is as follows: (31) The HIE-FDTD method was used to analyze the magnetic field components in the upper region. To update, the specific update formula is as follows: (32) The magnetic field components in the intermediate layer region were analyzed using the M(2,4)FDTD method. To update, the specific update formula is as follows: (33) in, and They represent the magnetic field component coefficients, respectively. and The value at the m-th computational grid; Loop coefficient; (34)。 7. The method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves according to claim 6, characterized in that, Step 9 specifically involves: The peak amplitude of the vertical electric field at the observation point is extracted using the peak monitoring method. ; The extracted vertical electric field peak amplitude at the observation point Normalized to standard radiation conditions, when the antenna radiated power is At that time, the corresponding standard field strength for: (35) in, Indicates the length of the radiating element. This is the field strength at the sampling point.

8. A computer device comprising a memory and one or more processors, wherein the memory stores executable code, characterized in that, When the processor executes the executable code, it implements the steps of the low-dispersion very low frequency radio wave propagation characteristic prediction method as described in any one of claims 1 to 7.

9. A computer-readable storage medium having a program stored thereon, characterized in that, When the program is executed by the processor, it implements the steps of the method for predicting the propagation characteristics of low-dispersion very low-frequency radio waves as described in any one of claims 1 to 7.

Images