Customise Debye material model

[Lauren Barr]

Hello,

First of all thanks for such a great package, gprMax has been incredibly useful for many of my projects!

The latest project I am working on involves studying some materials by extracting their real and imaginary permittivity as a function of frequency. I want to simulate these in gprMax and have only had limited success in fitting to my data using the built-in Debye model. Ideally I would use a double-Debye description of the material, but include relaxation parameters alpha and beta, following the equation here: https://en.wikipedia.org/wiki/Havriliak%E2%80%93Negami_relaxation

I think I need to adjust the `calculate_er` function in the `materials` subpackage, specifically this part:
if 'debye' in self.type:

    for pole in range(self.poles):

        er += self.deltaer[pole] / (1 + 1j * w * self.tau[pole])

And thought of modifying the last line to something like this, adding new attributes to the relevant classes, etc:

er += self.deltaer[pole] / (1 + (1j * w * self.tau[pole]) ** self.alpha[pole]) ** self.beta[pole]

Does this seem like a sensible approach, or can you see any obvious problems? Maybe you've thought of something like this before. Any advice or direction would be very appreciated before I get lost down a rabbit-hole of no return.

Many thanks,
Lauren

[Antonis Giannopoulos]

Hi Lauren,

If your material does not follow one of the classic or analytically available relaxation mechanisms that you can use directly (e.g. Debye, Lorentz, Drude)  then the only option is to try and fit a number of simple mechanisms like the Debye one two approximate the behaviour of your material in a given frequency range of interest. Outside that range the fit can be really bad but hopefully it will not affect you. 

The more Debye poles you use the better the fit will be but the more memory it will be needed for your model. If you look into our development branch there is work on this topic that slowly makes its way into gprMax under the Materials folder in the toolboxes. There are also papers on this. So, you can convert your material fitting process into an optimisation process for fitting the parameters of a number of Debye poles to match your experimental frequency data.

This has been done many times and it works. For example the models of echosorb absorbers available in gprMax have been done that way. Once you have the Debye poles that fit your data you just use them in gprMax as a normal material.

Sadly, you cannot add the alpha and beta parameters to generalise the relaxation function. This will be very nice to be able to do but not all alpha and beta combinations have closed form time domain representations that can be implemented in FDTD. So, the best option is to try and use more Debye poles to fit the experimental data. Also, our underlying formulation is very generic and maybe can fit other mechanisms but see our paper on how is done (number 4 in the publications list)

Changing the code as you suggest will not work as the formulation is in frequency domain and getting it back to the time domain there is no closed solution unless alpha and beta are very specific.

Hope this helps.

 Antonis

[Lauren Barr]

Hi Antonis,

Thanks so much for your response. I was worried there would be issues with implementation in the time domain and I'm glad I checked with you first! I will definitely look into expanding the number of poles to improve the fit, and take a look at those useful references.

Best wishes,

Lauren