Open access peer-reviewed chapter

A TLM Formulation Based on Fractional Derivatives for Dispersive Cole-Cole Media

Written By

Mohammed Kanjaa, Otman El Mrabet and Mohsine Khalladi

Reviewed: February 3rd, 2021 Published: May 5th, 2021

DOI: 10.5772/intechopen.96378

Chapter metrics overview

157 Chapter Downloads

View Full Metrics

Abstract

An auxiliary differential equation (ADE) transmission line method (TLM) is proposed for broadband modeling of electromagnetic (EM) wave propagation in biological tissues with the Cole-Cole dispersion Model. The fractional derivative problem is surmounted by assuming a linear behavior of the polarization current when the time discretization is short enough. The polarization current density is approached using Lagrange extrapolation polynomial and the fractional derivation is obtained according to Riemann definition of a fractional α-order derivative. Reflection coefficients at an air/muscle and air/fat tissues interfaces simulated in a 1-D domain are found to be in good agreement with those obtained from the analytic model over a broad frequency range, demonstrating the validity of the proposed approach.

Keywords

  • computational electromagnetics
  • numerical methods
  • transmission line matrix
  • biological system modeling
  • Cole-Cole medium
  • fractional derivative

1. Introduction

In the last two decades there has been a growing interest in the interaction between biological tissues and electromagnetic field at microwave frequencies. New promising applications of this technology in biomedical engineering like microwave imaging [1, 2], minimal invasive cancer therapies as thermal ablations [3], ultra-wide band temperature dependent dielectric spectroscopy [4] and EM dosimetry [5], rely heavily on an accurate mathematical model of the response of these tissues to an external electromagnetic field. Numerical resolution of the propagation problem within these tissues requires a robust mathematical model and the previous incorporation of the dielectric data that at all the working frequencies.

A first model of the time response in a time varying electric field of biological tissues was formulated by Debye [6] through a time decaying polarization current jt=Ptresulting from the time variation of the polarization vector Pt, the permittivity for this model is given in Eq. (1) and the corresponding argand diagram is depicted in Figure 1.

Figure 1.

Argand diagram for a unique pole Debye model obtained by representingEq. (1)in the complexe plane.

εrωεr=εrωεrjεrω=Δεp1+jωτpE1

This model even if it fits the experimental results in liquids it loses its accuracy when applied over a large band of frequency or in the presence of more than one type of polar molecule. This non-Debye relaxation is attributed to the existence of different relaxation processes [7] each with its own relaxation time τand its amplitude Δε. A more accurate approximation to the behavior of this category of dielectrics is given by the sum of the individual relaxation processes leading to a multipole model, the permittivity is given as the sum of the pterms corresponding to the ppoles each with its own amplitude. In the case of the biological tissue there is a multiple contribution of each relaxation process resulting in a broadening of the relaxation zone. Cole-Cole [8] proposed an Argand diagram in which in which εrthe imaginary part of the complex permittivity is plotted as a function of its real part εr, following the empirical relation:

εrωεr=εrωεrjεrω=p=1pΔεp1+jωτpαpE2

where εris the optical relative permittivity, Δεp=εrsεris the amplitude of the pthpole, εsis the static relative permittivity, ωis the angular frequency ,τpis the relaxation time, αpthe parameter that indicates the broadening of the dispersion for this pole, which in the Cole-Cole model must satisfy 0<αp<1, in the case of αp=1the model is simplified to a Debye dispersion problem and j=1(Figure 2).

Figure 2.

Argand diagram for a simple pole Cole-Cole model obtained by representingEq. (2)in the complexe plane.

The difficulty in the numerical implementation of such a model arises from the αpparameter which is a noninteger. The consequence is a fractional order differentiation in the time domain which is more challenging compared to the Debye time domain solution. To address this problem, authors in [9] used the Letnikov fractional derivative by introducing a Stirling asymptotic formula and a recursive relation for the polarization vector. The computational demands of the FDTD scheme were considerably reduced, but nevertheless it required the storage of a large number of previous values of the polarization vector. In [10], authors used the Z transform to formulate the frequency dependence between the electric flux density and the electric field, the fractional derivative was approximated by using a polynomial method. Rekanos et al. [11] used the Pade approximation, where the fractional power term is approximated using a Pade technique. All the latter approximations were made in the context of the FDTD method.

The TLM method, even if it is flexible and wide band and being a time domain method, cannot deal with the dispersive aspect of the Cole-Cole medium directly. In [12] an approach based on a convolution product between the susceptibility and the electric field and the temporal behavior is deduced by a DFT and a nonrecursive summation leading to a considerable computational cost and a nonnegligible error at high frequencies. The causality principle is used to justify the minor dependence of the recent susceptibility values on previous ones, but the problem of the fractional Differentegration wasn’t addressed.

In this work an ADE-TLM algorithm to model the Cole-Cole dispersion. The polarization current density is approached using an extrapolation with Lagrange polynomial method and the fractional derivation is obtained using the Riemann definition of the α-order derivative. The auxiliary differential equation is used to establish the update equation of the polarization currents which are included later in the general structure of the SCN-TLM node.

Advertisement

2. Formulation of the method

In the Cole-Cole model the relationship between the polarization current related to the pthpole and the electric field is given by:

Jpω=ε0Δεp1+jωτpαpEωE3

where jωτpαpis the power law function of frequency which in time domain results in a fractional derivative of order α:

Jpt+ταpDαpJpt=ε0ΔεpEttE4

One could consider an analytic solution for this equation, but due to the fractional derivative this can only be done for particular values of αp, hence a numerical solution is necessary with a previous discretization of the temporal values of the vectors.

In order to obtain the discretized expression of Jpt, the polarization current vector can be approached using a Lagrange interpolation in the time interval ti<t<tjwith cardinal functions:

Lit=j=1;ijnttjtitjE5

the interpolated expression of Jptis then obtained:

Jpt=i=1nJiLitE6

and by applying (5) in the time interval nΔt<t<n+1Δtwe obtain the cardinal functions:

Lnt=tn+1ΔtnΔtn+1ΔtLn+1t=tnΔtn+1ΔtnΔtE7

By substituting (7) in (6) the polynomial extrapolation of the polarization current density between time steps nΔtand n+1Δtcan be found in a straightforward manner as follow:

Jpt=tn+1ΔtnΔtn+1ΔtJpn+tnΔtn+1ΔtnΔtJpn+1E8

which after simplification can be rewritten as:

Jpt=Jpn+1JpnΔtt+n+1JpnnJpn+1E9

Furthermore, by substituting (8) into (3) we derive the auxiliary differential equation with a fractional order derivative given by:

Jpn+1+Jpn2+ταpJpn+1JpnΔtDαpt+Dαpn+1JpnnJpn+1=ε0ΔεpEtE10

The generalization of the derivation operator to arbitrary non-integer orders, has been subject to intensive research from mathematicians [13, 14] and in electromagnetism [15]. One of the definitions in [16] for the fractional derivative, symbolized by the operator Dxαis given for the power series with fractional exponents by considering a function hxdefined as a power series hx=i=1nAixaν+i/nfor ν>1and na positive integer, so each term xapcould have a noninteger exponent p=ν+i/n, then according to Reimann [15] a fractional derivative of order αfor this term is given by:

Dxαxap=dαxapdxaα=Γp+1Γpα+1xapαforx>aE11

where in our case αis a noninteger number [16].

Hence by applying the fractional derivation operator given in (11) to (10) and for p=1and a=0we obtain the update equation for the polarization current:

Jpn+1=ϕpψpJpn+ζpψpEn+1EnE12

where the update parameters ψp, ϕpand ζpare:

ψp=12+τpαpΔtΓ2Γ2αpt1αpnτpαpΓ1αptαpE13
ϕp=12τpαpΔtΓ2Γ2αpt1αp+n+1τpαpΓ1αptαpE14
ζp=ε0ΔεpΔtE15

To get the update equation of the electric field components we start from the Maxwell-Ampere equation, and by including the conductivity term:

×H=ε0εEt+σ0E+p=1pJpE16

where σ0is the static ionic conductivity, Eq. (16) formulated at time step n+12and by approaching Jn+12by the average of its values at time steps nand n+1

Jpn+12=Jpn+1+Jpn2E17

Finally, the updated electric field is given by:

Ex,y,zn+1=Ex,y,zn12σ0ΔtDp1ϕpψpΔtDJpx,y,zn+2ΔtD×Hx,y,zn+12E18
D=2ε0ε+pζpΔtψp+σ0ΔtE19
Advertisement

3. The TLM Formalisme

The TLM method is a numerical technique based on the discretization of the computational domain according to Huygens principle as an alternative to the Maxwell equations used in the FDTD method [17]. In this method the simulation domain is discretized in cells where a series of uniform transmission lines, parallel and series stubs and additional sources are used to take account of the real characteristics of the propagation through the medium as given by Maxwell equations. Therefore, instead of electric and magnetic field components the electromagnetic field is represented by voltage and current waves, propagating through the unit cell circuit referred to as symmetrical condensed node (SCN). The relationship between the electromagnetic field and the voltage and current waves at the time step nΔtand at the center of the node are formulated as follows [18].

As can be seen in Figure 3(b), the SCN consists of interconnected transmission lines. This structure models a unit cell of the propagation medium as in Figure 3(a). Each face of the cell corresponds to two ports orthogonal to each other and labeled from 1 to 12, these stubs model the propagation through free space and have therefore its characteristic impedance Z0=ε0μ0. Three additional open circuit (13,14,15) stubs must be added to the node to account for the dispersive behavior of the medium. The SCN nodes are connected to each other to form the simulation domain given in Figure 4 for a 1-D propagation.

Figure 3.

Structure of the TLM symmetrical condensed node, (a) unit volume of thedielectric material with dimensions Δx, Δy,and Δz, (b) equivalent SCN.

Figure 4.

Structure of the 1-D SCN grid used in this work to model the simulation domain.

The relationship between the electromagnetic field and the voltage and current waves at the time step nΔtand at the center of the node are formulated as follows [18]:

Ex,y,zn=Vx,y,znΔlHx,y,zn=Ix,y,znΔlE20

In this algorithm as in [11, 19] the dispersive properties of the medium are accounted for by adding voltage sources sVto each node, these sources constitute the counterparts of the terms in the FDTD formulation that have no equivalents in the TLM formalism. Therefore the voltage at the center of the node, as deduced from the conservation laws [20], becomes:

Vx,y,zn+1=Vx,y,zn+14+Yocx,y,zsVx,y,zn+1+sVx,y,zn+44+Yocx,y,zΔlΔtε0×Hx,y,zn+12E21

where Yocx,y,zindicates the normalized admittance of the open circuit stub added to the node.

When expressed in terms of incident voltages on the corresponding stubs and accounting for the voltage sources injected in the stubs (16,17,18) to complete the model of the Cole-Cole medium:

Vx,y,zn+1=24+Yox,y,zV1,3,5i+V2,4,6i+V9,8,7i+V12,11,10i+Yoc,x,y,zV13,14,15i+12sV16,17,18n+1E22

where the subscript iindicates the incident pulses on the indicated stubs.

In the TLM formalism the update equation issued from the ADE in the Eq. (12) becomes:

Jpx,y,zn+1=ϕpψpJpx,y,zn+ζpψpΔlVx,y,zn+1Vx,y,znE23

and the electric field update equation is also formulated in the TLM formalism as:

Vx,y,zn+1=Vx,y,zn12σ0ΔtDp1ϕpψpΔtΔlDJpx,y,zn+2ΔtΔlD×Hx,y,zn+12E24

The analogy between Eqs. (22) and (24) the expression of the normalized admittance of the stub added to the SCN node is obtained straightforward:

Yoc=4D2ε01E25

and the update equation for voltage sources at the center of the node:

sVx,y,zn+1=sVx,y,zn4σ0Δtε0Vx,y,znp21ϕpψpΔtΔlε0Jpx,y,zE26

that can be simplified to

sVx,y,zn+1=sVx,y,zn+C1Vx,y,zn+pC2pJpx,y,zE27

with update constants:

C1=4σ0Δtε0E28
C2p=21ϕpψpΔtΔlε0E29

A final step to establish the TLM formulation is to express the scattering process on the capacitive stubs:

V13,14,15rn=Vx,y,znV13,14,15inE30
Advertisement

4. Simulation and results

In order to verify the new proposed algorithm, and for the sake of simplicity a one dimensional problem is simulated, where a Cole-Cole medium (fat or muscle) occupies the region z0, while the rest of space is air, 2 dimensional and 3 dimensional propagation and reflection models can be naturally deducted by ensuring the connexion between TLM cells in the yaxis using ports 1,5,7 and 12, to account for the propagation along the xaxis ports 3,6,10 and 11 must be connected to adjacent cells as in Figure 5.

Figure 5.

2-D ADE-TLM simulation of an incident and reflected pulse on the air/muscle, the tissue is modeled by a 4 pole Cole-Cole model, the interface is located at the TLM cellz=100, at time steps 160 (a) and 218 (b).

The incident wave is a derivative Gaussian pulse given by Einc=tt0τ2exp4π(tt02τ2êxwhere t0=200Δtand τ=220Δtpolarized in the xdirection, and propagates along the zaxis. In the TLM implementation, the considered 1-D simulation domain (see Figure 6) consists of a grid of 1000 cells interconnected on the z axis as in Figure 4, 500 of which were used to model the Cole-Cole medium and the remaining cells were used to model the air. The cell dimension Δlsatisfies the stability condition Δt=Δl/2c, the time step is set as Δt=0.125ps. The simulation domain is truncated by an unsplit PML layer [21] of 10 cells at the beginning and the end.

The fourth order Cole-Cole parameters for fat tissue as well as for muscle tissue are listed in Table 1 [9]. The values of the electric field are recorded at a point P located at the center of the SCN node 10 cells before the air-medium interface. Figure 7(a) shows the incident impulse at iteration 399 propagating in the air and Figure 7(b) at iteration 927 depicts the reflected and transmitted impulse on the air/muscle interface.

Tissueεσ0Δε1τ1(ps)α1Δε2τ2(ps)α2Δε3τ3(μs)α3Δε4τ4(ps)Δε4
Muscle2.50.03597.960.83515.920.93.3E4159.150.951.0E715.9150.99
Fat40.2507.230.97000353.680.91.2E6318.310.92.5E72.2741

Table 1.

Cole-Cole constants for tissues of muscle and fat.

Figure 6.

2-D ADE-TLM simulation of a reflected and transmitted pulse on the air/muscle , the tissue is modeled by a 4 pole Cole-Cole model, the interface is located at the TLM cell z = 100 , at (a) time step 305, and (b) at time step 446.

In Figures 8 and 9, the simulation results for both cases (fat and muscle) are compared to the theoretical reflection at the air/Cole-Cole medium interface which can be obtained by using the following equation [22]:

R=1εr1+εrE31

Figure 7.

Propagation of an incident derivative Gaussian pulse through Cole-Cole medium model of human muscle (a) at iteration 399 (incident pulse) (b) at iteration 677 (reflected and transmitted pulse).

Figure 8.

Reflection coefficient (dB) at normal incidence on an air/muscle interface.

A good agreement over the whole frequency band is observed. The slight discrepancy between the numerical and theoretical results can be ascribed to the approximation made to the polarization current value when performing the Lagrange extrapolation on each iteration. The propagation and reflection through a 2-D TLM grid modeling the incidence on an air/muscle interface is presented in Figures 5 and 10 at different time steps, the interface is located at z=100in a 10by 200cells grid, Figure 5(a) shows the Gaussian wave before the arrival on the interface while Figure 10(a) and (b) depicts it after reflection and transmission. A more precise simulation that targeted the Dielectric properties of pathological thyroid tissue types including adenoma, thyroiditis and the properties of the healthy thyroid based on the results of the experimental study presented in a previous published work in [2], and for the healthy thyroid tissue [23], in both works experimental results are used to extract the 2 poles Cole-Cole model parameters. Based upon this parameters we conducted a simulation for the 4 tissue types, the results of this simulation are presented in Figure 11(a) and (b) for the air/adenoma and the air/normal thyroidian interface respectively, Figure 12(a) and (b) depicts the reflection coefficient on an air/thyroyditis affected thyroidian tissue and an air/tumor respectively also in this case a perfect agreement over the working frequency is observed.

Figure 9.

Reflection coefficient (dB) at normal incidence on an air/fat interface.

Figure 10.

2-D ADE-TLM simulation of a reflected and transmitted pulse on the air/muscle, the tissue is modeled by a 4 pole Cole-Cole model, the interface is located at the TLM cellz=100, at time steps are 305 and 446.

Figure 11.

Reflection coefficient (dB) at normal incidence on an air/adenoma (a), air/Normal tyroidian tissue interface (b).

Figure 12.

Reflection coefficient (dB) at normal incidence on an air.

Advertisement

5. Conclusion

Numerical methods are an essential part of the modeling process, new propagation media with anomalous relaxation properties impose innovative modeling methods. An effective ADE-TLM formulation way to model electromagnetic wave propagation in biological tissues using Cole-Cole dispersion model. The fractional-order derivative in the Cole-Cole model is tackled by using a polynomial extrapolation of the polarization current. This approximation reduces considerably both the numerical cost of the simulation and its complexity since only two previous values of the Jfield are needed at each iteration. a Reimann derivation on the obtained polynomial extrapolation of the polarization current is then performed to solve the fractional order derivative. The presented results indicate the accuracy of our approach.

Advertisement

Conflict of interest

The authors declare no conflict of interest.

References

  1. 1. Martellosio A, Pasian M, Bozzi M, Perregrini A, Mazzanti F, Svelto PE, et al. L: 0.5-50 ghz dielectric characterisation of breast cancer tissues. Electronics Letters. 2015;51(13):974-975
  2. 2. Gavazzi S, Limone P, De Rosa G, Molinari F, and Vecchi G: Comparison of microwave dielectric properties of human normal, benign and malignant thyroid tissues obtained from surgeries: A preliminary study, Biomedical Physics Engineering Express, vol. 4, no. 4, p. 047 003, 2018
  3. 3. Ruvio G, Eaton-Evans J, Shahzad A, OHalloran M. Numerical evaluation of microwave thermal ablation to treat small adrenocortical masses. International Journal of RF and Microwave Computer-Aided Engineering. 2018;28(3):e21236
  4. 4. Ley S, Schilling S, Fiser O, Vrba J, Sachs J, Helbig M. Ultra-wideband temperature dependent dielectric spectroscopy of porcine tissue and blood in the microwave frequency range. Sensors. 2019;19(7):1707
  5. 5. Chakarothai J, Wake K, Watanabe S. Convergence of a single-frequency FDTD solution in numerical dosimetry. IEEE Transactions on Microwave Theory and Techniques. 2016;64(3):707-714
  6. 6. Debye P. Part i. dielectric constant. energy absorption in dielectrics with polar molecules. Transactions of the Faraday Society. 1934;30:679-684
  7. 7. Cole KS, Cole RH. Dispersion and absorption in dielectrics i. alternating current characteristics. The Journal of Chemical Physics. 1941;9(4):341-351
  8. 8. Rekanos IT, Yioultsis TV. Approximation of grunwald-letnikov fractional derivative for FDTD modeling of Cole-Cole media. IEEE Transactions on magnetics. 2014;50(2):181-184
  9. 9. Guo B, Li J, Zmuda H. A new FDTD formulation for wave propagation in biological media with Cole-Cole model. IEEE Microwave and Wireless Components Letters. 2006;16(12):633-635
  10. 10. Rekanos IT, Papadopoulos TG. An auxiliary differential equation method for FDTD modeling of wave propagation in Cole-Cole dispersive media. IEEE Transactions on Antennas and Propagation. 2010;58(11):3666-3674
  11. 11. Barba I. Cabeceira A C L, Panizo M, and Represa J: Modelling dispersive dielectrics in TLM method, International Journal of Numerical Modelling: Electronic Networks. Devices and Fields. 2001;14(1):15-30
  12. 12. Samko SG, Kilbas AA, Marichev OI, et al. Fractional integrals and derivatives. Gordon and Breach Science Publishers, Yverdon Yverdon-les-Bains, Switzerland, 1993. In: vol. 1993
  13. 13. Techniques in the fractional calculus,” in The Fractional Calculus, ser. Mathematics in Science and Engineering, K. B. Oldham and J. Spanier, Eds., vol. 111, Elsevier, 1974, pp. 133–160
  14. 14. N. Engheta: n the role of fractional calculus in electromagnetic theory,” IEEE Antennas and Propagation Magazine, vol. 39, no. 4, pp. 35–46, Aug. 1997
  15. 15. Engheta N. On fractional calculus and fractional multipoles in electromagnetism. IEEE Transactions on Antennas and Propagation. Apr. 1996;44(4):554-566
  16. 16. Wharmby AW, Bagley RL. The application of the fractional calculus model for dispersion and absorption in dielectrics in terahertz waves. International Journal of Engineering Science. 2015;93:1-12
  17. 17. Jin H, Vahldieck R. Direct derivations of TLM symmetrical condensed node and hybrid symmetrical condensed node from maxwell’s equations using centered differencing and averaging. IEEE Transactions on Microwave Theory and Techniques. Dec. 1994;42(12):2554-2561
  18. 18. Christopoulos C, The Transmission-Line Modeling (TLM) Method in Electromagnetics. Morgan Claypool, 2006
  19. 19. Cabeceira ACL, Barba I, Grande A, Represa J. A 2D-TLM model for electromagnetic wave propagation in chiral media. Microwave and Optical Technology Letters. 2005;46(2):180-182
  20. 20. Yaich MI, Khalladi M. The far-zone scattering calculation of frequency-dependent materials objects using the tlm method. IEEE Transactions on Antennas and Propagation. Nov. 2002;50(11):1605-1608
  21. 21. Yaich MI, Kanjaa M, Adraoui S, Mounirh K, Khalladi M. An unsplit formulation of the 3D-PML absorbing boundary conditions for TLM-method in time domain. In: 2018 6th International Conference on Multimedia Computing and Systems (ICMCS). May 2018. pp. 1-5
  22. 22. Taflove A, Hagness SC. Computational electrodynamics: the finite-difference time-domain method. Artech house. 2005
  23. 23. Gabriel S, Lau RWand Gabriel C1996b The dielectric properties of biological tissues: III. Parametric models for the dielectric spectrum of tissues Phys. Med. Biol. 41 227193

Written By

Mohammed Kanjaa, Otman El Mrabet and Mohsine Khalladi

Reviewed: February 3rd, 2021 Published: May 5th, 2021