publicImage/gprMax
1
0
派生 0
镜像自地址 https://gitee.com/sunhf/gprMax.git 已同步 2025-08-07 04:56:51 +08:00
代码 工单 软件包 项目 版本发布 百科 动态

Corrected reading of Lorentz and Drude parameters.

浏览代码
这个提交包含在：
Craig Warren
2015-11-20 17:54:19 +00:00
父节点 0e2f50b987
当前提交 b241b2d24c
共有 6 个文件被更改，包括 30 次插入 和 26 次删除
下载 Patch 文件 下载 Diff 文件 展开所有文件 折叠所有文件

6
docs/source/features.rst
查看文件

@@ -16,7 +16,7 @@ The code has been completely re-written in Python/Cython. In the process a lot o
* Scripting in the input file.
* Built-in library of antenna models
* Anisotropic material modelling.
* Dispersive material modelling using multiple pole Debye, Lorenz or Drude formulations.
* Dispersive material modelling using multiple pole Debye, Lorentz or Drude formulations.
* Building heterogeneous objects using fractal distributions.
* Building objects with rough surfaces.
* Modelling soils with realistic dielectric and geometric properties.
@@ -28,7 +28,7 @@ New commands
* ``#python`` and ``#end_python`` are used to define blocks of the input file where Python code will be executed. This allows the user to use scripting directly in the input file.
* ``#material`` replaces ``#medium`` with a new syntax.
* ``#add_dispersion_debye`` is used to add Debye dispersive properties to a ``#material``.
* ``#add_dispersion_lorenz`` is used to add Lorenz dispersive properties to a ``#material``.
* ``#add_dispersion_lorentz`` is used to add Lorentz dispersive properties to a ``#material``.
* ``#add_dispersion_drude`` is used to add Drude dispersive properties to a ``#material``.
* ``#soil_peplinski`` is a soil mixing model that can be used with ``#fractal_box`` to generate soil(s) with more realistic dielectric and geometric properties.
* ``#cylindrical_sector`` (like a slice of pie shape) is a new object building command.
@@ -95,7 +95,7 @@ The input file has now been made scriptable by permitting blocks of Python code
Dispersive media
----------------
gprMax has always included the ability to represent dispersive materials using a single-pole Debye model. Many materials can be adequately represented using this approach for the typical frequency ranges associated with GPR. However, multi-pole Debye, Drude and Lorenz functions are often used to simulate the electric susceptibility of materials such as: water [PIE2009]_, human tissue [IRE2013]_, cold plasma [LI2013]_, gold [VIA2005]_, and soils [BER1998]_, [GIAK2012]_, [TEI1998]_. Electric susceptibility relates the polarization density to the electric field, and includes both the real and imaginary parts of the complex electric permittivity variation. In the new version of gprMax a recursive convolution based method is used to express dispersive properties as apparent current density sources [GIA2014]_. A major advantage of this implementation is that it creates an inclusive susceptibility function that holds, as special cases, Debye, Drude and Lorenz materials. For further details see the :ref:`material commands section <materials>`.
gprMax has always included the ability to represent dispersive materials using a single-pole Debye model. Many materials can be adequately represented using this approach for the typical frequency ranges associated with GPR. However, multi-pole Debye, Drude and Lorentz functions are often used to simulate the electric susceptibility of materials such as: water [PIE2009]_, human tissue [IRE2013]_, cold plasma [LI2013]_, gold [VIA2005]_, and soils [BER1998]_, [GIAK2012]_, [TEI1998]_. Electric susceptibility relates the polarization density to the electric field, and includes both the real and imaginary parts of the complex electric permittivity variation. In the new version of gprMax a recursive convolution based method is used to express dispersive properties as apparent current density sources [GIA2014]_. A major advantage of this implementation is that it creates an inclusive susceptibility function that holds, as special cases, Debye, Drude and Lorentz materials. For further details see the :ref:`material commands section <materials>`.
Realistic soils, heterogeneous objects and rough surfaces
---------------------------------------------------------

22
docs/source/input.rst
查看文件

@@ -314,10 +314,10 @@ For example to create a model of water with a single Debye pole, :math:`\epsilon
* Temporal values associated with pole frequencies and relaxation times should always be greater than the time step :math:`\Delta t` used in the model.
#add_dispersion_lorenz:
#add_dispersion_lorentz:
-----------------------
Allows you to add dispersive properties to an already defined ``#material`` based on a multiple pole Lorenz formulation (see :ref:`capabilities` section). For example, the susceptability function for a single-pole Lorentz material is given by:
Allows you to add dispersive properties to an already defined ``#material`` based on a multiple pole Lorentz formulation (see :ref:`capabilities` section). For example, the susceptability function for a single-pole Lorentz material is given by:
.. math::
@@ -335,21 +335,21 @@ The syntax of the command is:
.. code-block:: none
#add_dispersion_lorenz: i1 f1 f2 f3 f4 f5 f6 ... str1
#add_dispersion_lorentz: i1 f1 f2 f3 f4 f5 f6 ... str1
* ``i1`` is the number of Lorenz poles.
* ``f1`` is the difference between the zero-frequency relative permittivity and the relative permittivity at infinite frequency, i.e. :math:`\Delta \epsilon_{rp1} = \epsilon_{rsp1} - \epsilon_{r \infty p1}` , for the first Lorenz pole.
* ``f2`` is the frequency (Hertz), :math:`\omega_{p1}`, for the first Lorenz pole.
* ``f3`` is the damping coefficient (Hertz), :math:`\delta_{p1}`, for the first Lorenz pole.
* ``f4`` is the difference between the zero-frequency relative permittivity and the relative permittivity at infinite frequency, i.e. :math:`\Delta \epsilon_{rp2} = \epsilon_{rsp2} - \epsilon_{r \infty p2}` , for the second Lorenz pole.
* ``f5`` is the frequency (Hertz), :math:`\omega_{p2}`, for the second Lorenz pole.
* ``f6`` is the damping coefficient (Hertz), :math:`\delta_{p2}`, for the second Lorenz pole.
* ``i1`` is the number of Lorentz poles.
* ``f1`` is the difference between the zero-frequency relative permittivity and the relative permittivity at infinite frequency, i.e. :math:`\Delta \epsilon_{rp1} = \epsilon_{rsp1} - \epsilon_{r \infty p1}` , for the first Lorentz pole.
* ``f2`` is the frequency (Hertz), :math:`\omega_{p1}`, for the first Lorentz pole.
* ``f3`` is the damping coefficient (Hertz), :math:`\delta_{p1}`, for the first Lorentz pole.
* ``f4`` is the difference between the zero-frequency relative permittivity and the relative permittivity at infinite frequency, i.e. :math:`\Delta \epsilon_{rp2} = \epsilon_{rsp2} - \epsilon_{r \infty p2}` , for the second Lorentz pole.
* ``f5`` is the frequency (Hertz), :math:`\omega_{p2}`, for the second Lorentz pole.
* ``f6`` is the damping coefficient (Hertz), :math:`\delta_{p2}`, for the second Lorentz pole.
* ...
* ``str1`` identifies the material to add the dispersive properties to.
.. note::
* You can continue to add triplets of values for :math:`\Delta \epsilon_{rp}`, :math:`\omega_p` and :math:`\delta_p` for as many Lorenz poles as you have specified with ``i1``.
* You can continue to add triplets of values for :math:`\Delta \epsilon_{rp}`, :math:`\omega_p` and :math:`\delta_p` for as many Lorentz poles as you have specified with ``i1``.
* The relative permittivity in the ``#material`` command should be given as the relative permittivity at infinite frequency, i.e. :math:`\epsilon_{r \infty}`.
* Temporal values associated with pole frequencies and relaxation times should always be greater than the time step :math:`\Delta t` used in the model.

2
gprMax/input_cmds_file.py
查看文件

@@ -130,7 +130,7 @@ def check_cmd_names(processedlines):
singlecmds = dict.fromkeys(['#domain', '#dx_dy_dz', '#time_window', '#title', '#messages', '#num_threads', '#time_step_stability_factor', '#time_step_limit_type', '#pml_cells', '#excitation_file', '#src_steps', '#rx_steps'], 'None')
# Commands that there can be multiple instances of in a model - these will be lists within the dictionary
multiplecmds = {key: [] for key in ['#geometry_view', '#material', '#soil_peplinski', '#add_dispersion_debye', '#add_dispersion_lorenz', '#add_dispersion_drude', '#waveform', '#voltage_source', '#hertzian_dipole', '#magnetic_dipole', '#rx', '#rx_box', '#snapshot', '#pml_cfs']}
multiplecmds = {key: [] for key in ['#geometry_view', '#material', '#soil_peplinski', '#add_dispersion_debye', '#add_dispersion_lorentz', '#add_dispersion_drude', '#waveform', '#voltage_source', '#hertzian_dipole', '#magnetic_dipole', '#rx', '#rx_box', '#snapshot', '#pml_cfs']}
# Geometry object building commands that there can be multiple instances of in a model - these will be lists within the dictionary
geometrycmds = ['#edge', '#plate', '#triangle', '#box', '#sphere', '#cylinder', '#cylindrical_sector', '#fractal_box', '#add_surface_roughness', '#add_surface_water', '#add_grass']

12
gprMax/input_cmds_multiuse.py
查看文件

@@ -432,7 +432,7 @@ def process_multicmds(multicmds, G):
if G.messages:
print('Debye-type disperion added to {} with delta_epsr={}, and tau={} secs created.'.format(material.ID, ','.join('%4.2f' % deltaer for deltaer in material.deltaer), ','.join('%4.3e' % tau for tau in material.tau)))
cmdname = '#add_dispersion_lorenz'
cmdname = '#add_dispersion_lorentz'
if multicmds[cmdname] != 'None':
for cmdinstance in multicmds[cmdname]:
tmp = cmdinstance.split()
@@ -452,7 +452,7 @@ def process_multicmds(multicmds, G):
raise CmdInputError("'" + cmdname + ': ' + ' '.join(tmp) + "'" + ' material(s) {} do not exist'.format(notfound))
for material in materials:
material.type = 'lorenz'
material.type = 'lorentz'
material.poles = poles
material.average = False
for pole in range(1, 3 * poles, 3):
@@ -461,12 +461,12 @@ def process_multicmds(multicmds, G):
material.tau.append(float(tmp[pole + 1]))
material.alpha.append(float(tmp[pole + 2]))
else:
raise CmdInputError("'" + cmdname + ': ' + ' '.join(tmp) + "'" + ' requires positive values for the permittivity difference and relaxation times, and relaxation times that are greater than the time step for the model.')
raise CmdInputError("'" + cmdname + ': ' + ' '.join(tmp) + "'" + ' requires positive values for the permittivity difference and frequencies, and associated times that are greater than the time step for the model.')
if material.poles > Material.maxpoles:
Material.maxpoles = material.poles
if G.messages:
print('Lorenz-type disperion added to {} with delta_epsr={}, tau={} secs, and alpha={} created.'.format(material.ID, ','.join('%4.2f' % deltaer for deltaer in material.deltaer), ','.join('%4.3e' % tau for tau in material.tau), ','.join('%4.3e' % alpha for alpha in material.alpha)))
print('Lorentz-type disperion added to {} with delta_epsr={}, omega={} secs, and gamma={} created.'.format(material.ID, ','.join('%4.2f' % deltaer for deltaer in material.deltaer), ','.join('%4.3e' % tau for tau in material.tau), ','.join('%4.3e' % alpha for alpha in material.alpha)))
cmdname = '#add_dispersion_drude'
@@ -497,12 +497,12 @@ def process_multicmds(multicmds, G):
material.tau.append(float(tmp[pole ]))
material.alpha.append(float(tmp[pole + 1]))
else:
raise CmdInputError("'" + cmdname + ': ' + ' '.join(tmp) + "'" + ' requires positive values for the relaxation times, and relaxation times that are greater than the time step for the model.')
raise CmdInputError("'" + cmdname + ': ' + ' '.join(tmp) + "'" + ' requires positive values for the frequencies, and associated times that are greater than the time step for the model.')
if material.poles > Material.maxpoles:
Material.maxpoles = material.poles
if G.messages:
print('Drude-type disperion added to {} with tau1={} secs, and tau2={} secs created.'.format(material.ID, ','.join('%4.3e' % tau for tau in material.tau), ','.join('%4.3e' % alpha for alpha in material.alpha)))
print('Drude-type disperion added to {} with omega={} secs, and gamma={} secs created.'.format(material.ID, ','.join('%4.3e' % tau for tau in material.tau), ','.join('%4.3e' % alpha for alpha in material.alpha)))
cmdname = '#soil_peplinski'

12
gprMax/materials.py
查看文件

@@ -28,7 +28,7 @@ class Material():
maxpoles = 0
# Types of material
types = ['standard', 'debye', 'lorenz', 'drude']
types = ['standard', 'debye', 'lorentz', 'drude']
# Properties of water from: http://dx.doi.org/10.1109/TGRS.2006.873208
waterer = 80.1
@@ -102,12 +102,16 @@ class Material():
if self.type == 'debye':
self.w[x] = self.deltaer[x] / self.tau[x]
self.q[x] = -1 / self.tau[x]
elif self.type == 'lorenz':
wp2 = (2 * np.pi * (1 / self.tau[x])) * (2 * np.pi * (1 / self.tau[x]))
elif self.type == 'lorentz':
# tau for Lorentz materials are pole frequencies
# alpha for Lorentz materials are the damping coefficients
wp2 = (2 * np.pi * self.tau[x]) * (2 * np.pi * (1 / self.tau[x]))
self.w[x] = -(wp2 * self.deltaer[x]) * j / np.sqrt(wp2 - (self.alpha[x] * self.alpha[x]))
self.q[x] = -self.alpha[x] + np.sqrt(wp2 - (self.alpha[x] * self.alpha[x])) * j
elif self.type == 'drude':
wp2 = (2 * np.pi * (1 / self.tau[x])) * (2 * np.pi * (1 / self.tau[x]))
# tau for Drude materials are pole frequencies
# alpha for Drude materials are the inverse of relaxation times
wp2 = (2 * np.pi * self.tau[x]) * (2 * np.pi * self.tau[x])
self.se += wp2 / self.alpha[x]
self.w[x] = - (wp2 / self.alpha[x])
self.q[x] = - self.alpha[x]

2
tools/inputfile_new2old.py
查看文件

@@ -43,7 +43,7 @@ newfile += '_oldstyle'
print("Attempting to convert inputfile '{}' to use old syntax...\n".format(inputfile))
newcommands = ['#add_dispersion_lorenz', '#add_dispersion_drude', '#fractal_box', '#soil_peplinski', '#python', '#end_python', '#add_surface_roughness', '#add_surface_water', '#add_grass', '#magnetic_dipole', '#pml_cfs', '#cylindrical_sector', '#time_step_limit_type']
newcommands = ['#add_dispersion_lorentz', '#add_dispersion_drude', '#fractal_box', '#soil_peplinski', '#python', '#end_python', '#add_surface_roughness', '#add_surface_water', '#add_grass', '#magnetic_dipole', '#pml_cfs', '#cylindrical_sector', '#time_step_limit_type']
materials = {}
debyes = []
waveforms = []