Processing math: 75%
Research article

Modeling electromagnetism in and near composite material using two-scale behavior of the time-harmonic Maxwell equations

  • Received: 10 April 2016 Accepted: 28 March 2017 Published: 04 May 2017
  • The main purpose of this article is to study the two-scale behavior of the electromagnetic field in 3D in and near composite material. For this, time-harmonic Maxwell equations, for a conducting two-phase composite and the air above, are considered. Technique of two-scale convergence is used to obtain the homogenized problem.

    Citation: Canot Hélène, Frénod Emmanuel. Modeling electromagnetism in and near composite material using two-scale behavior of the time-harmonic Maxwell equations[J]. AIMS Mathematics, 2017, 2(2): 269-304. doi: 10.3934/Math.2017.2.269

    Related Papers:

    [1] Tariq Mahmood . The zero-energy limit and quasi-neutral limit of scaled Euler-Maxwell system and its corresponding limiting models. AIMS Mathematics, 2019, 4(3): 910-927. doi: 10.3934/math.2019.3.910
    [2] Geoffrey Beck, Akram Beni Hamad . Electromagnetic waves propagation in thin heterogenous coaxial cables. Comparison between 3D and 1D models. AIMS Mathematics, 2024, 9(4): 8981-9019. doi: 10.3934/math.2024438
    [3] Bin He . Developing a leap-frog meshless methods with radial basis functions for modeling of electromagnetic concentrator. AIMS Mathematics, 2022, 7(9): 17133-17149. doi: 10.3934/math.2022943
    [4] Maha M. Helmi, Ali M. Mubaraki, Rahmatullah Ibrahim Nuruddeen . Analysis of horizontally polarized shear waves on a highly inhomogeneous loaded bi-material plate. AIMS Mathematics, 2023, 8(1): 2117-2136. doi: 10.3934/math.2023109
    [5] Zuojun Ma . Pattern formation of a volume-filling chemotaxis model with a bistable source. AIMS Mathematics, 2024, 9(11): 30816-30837. doi: 10.3934/math.20241488
    [6] Gao Chang, Chunsheng Feng, Jianmeng He, Shi Shu . Stability analysis of a class of nonlinear magnetic diffusion equations and its fully implicit scheme. AIMS Mathematics, 2024, 9(8): 20843-20864. doi: 10.3934/math.20241014
    [7] Rui Hou, Yongqing Xu, Jinhua Fan, Yuanguo Zhu . Short time asymptotics for American maximum options with a dividend-paying asset. AIMS Mathematics, 2022, 7(8): 13977-13993. doi: 10.3934/math.2022772
    [8] Ling Zhu . Asymptotic expansion of a finite sum involving harmonic numbers. AIMS Mathematics, 2021, 6(3): 2756-2763. doi: 10.3934/math.2021168
    [9] Ali Hassani . Singular expansion of the wave kernel and harmonic sums on Riemannian symmetric spaces of the non-compact type. AIMS Mathematics, 2025, 10(3): 4775-4791. doi: 10.3934/math.2025219
    [10] Ibrahim-Elkhalil Ahmed, Ahmed E. Abouelregal, Doaa Atta, Meshari Alesemi . A fractional dual-phase-lag thermoelastic model for a solid half-space with changing thermophysical properties involving two-temperature and non-singular kernels. AIMS Mathematics, 2024, 9(3): 6964-6992. doi: 10.3934/math.2024340
  • The main purpose of this article is to study the two-scale behavior of the electromagnetic field in 3D in and near composite material. For this, time-harmonic Maxwell equations, for a conducting two-phase composite and the air above, are considered. Technique of two-scale convergence is used to obtain the homogenized problem.


    1. Introduction

    We are interested in the time-harmonic Maxwell equations in and near a composite material with boundary conditions modeling electromagnetic field radiated by an electromagnetic pulse (EMP). An electromagnetic pulse is a short burst of electromagnetic energy. It may be generated by a natural occurrence such like a lightning strike, meteoric EMP, EMP caused by geomagnetic Storm or nuclear EMP. This focuses on what happens over a period of time of a millisecond during the peak of the first return stroke. We study the electromagnetic pulse caused by this lightning strike. This is the first step of a larger study which goal is to understand the behavior of the electromagnetic field and its interaction with a composite material.

    EMP interference is generally damaging to electronic equipment. A lightning strike can damage physical objects such as aircraft structures, either through heating effects or disruptive effects of the very large magnetic field generated by the current. Structures and systems require some form of protection against lightning. Every commercial aircraft is struck by lightning at least once a year on average. Aircraft lightning protection is a major concern for aircraft manufacturers. Increasing its use of composite materials, up to 53% for the latest Airbus A350, and 50% for the Boeing B787, aircrafts offer increased vulnerability facing lightning. Earlier generation aircrafts, whose fuselages were predominantly composed of aluminum, behave like a Faraday cage and offer maximum protection for the internal equipment. Currently, in aircrafts, composite materials consisting of a resin enclosing carbon fibers have significant advantages in terms of weight gain and therefore fuel saving. Yet, because aluminium conducts 100 to 1000 times more than composite, we lose the Faraday effect. Modern aircrafts have seen also the increasing reliance on electronic avionics systems instead of mechanical controls and electromechanical instrumentation. For these reasons, aircraft manufacturers are very sensitive to lightning protection and pay special attention to aircraft certification through testing and analysis.

    There are two types of lightning strikes to aircraft: the first one is the interception by the aircraft of a lightning leader. The second one, which makes about 90% of the cases, is when the aircraft initiates the lightning discharge by emitting two leaders when it is found in the intense electric field region produced by a thundercloud, our approach applies in this case. When the aircraft flies through a cloud region where the atmospheric electric field is large enough, an ionized channel, called a positive leader, merges from the aircraft in the direction of the ambient electric field. Laroche et al. [15], at an altitude of 6000 m, observed an atmospheric electric field close to 50 kV/m inside the storm clouds, 100 kV/m to the ground. When upward leader connects with the downward negative leader of the cloud, a return stroke is produced and a bright return stroked wave travels from aircraft to cloud. The lightning return strokes radiate powerful electromagnetic fields which may cause damage to aircraft electronic equipment. Our work is devoted to the study of the electromagnetic waves propagation in the air and in the composite material. In this artificial periodic material, the electromagnetic field satisfies the Maxwell equations.

    We evaluate the electromagnetic field within and near a periodic structure when the period of this microstructure is small compared to the wavelength of the electromagnetic wave. Our model is composed by air above the composite fuselage and we study the behavior of the electromagnetic wave in the domain filled by the composite material, representing the skin aircraft, and the air. We build the 3D model, under simplifying assumptions, using linear time-harmonic Maxwell equations and constitutive relations for electric and magnetic fields. Composite materials consist of conducting carbon fibers, distributed as periodic inclusions in a matrix (epoxy resin). We impose a magnetic permeability μ0 uniform and an electrical permittivity ϵ=ϵ0ϵ, where ϵ is the relative permittivity depending of the medium. In the future, we will enrich this model by adding complexity and we will consider non uniform magnetic permeability and electrical permittivity.

    Now, we account for some characteristic values. In the first place we focus on the boundary conditions as we consider them as the source. Then, we use on the upper frontier, the magnetic field induced by the peak of the current of the first return stroke

    ¯Hd=I2πr, (1)

    with current intensity I=200 kA and r the radius of the lightning strike, this is the worst aggression that can suffer an aircraft, and we deduce a characteristic electric field ¯E=20 kV/m. In our model we consider that we have very conductive-but not perfect conductors-carbon fibers and an epoxy resin whose conduction depends on its doping rate. The conductivity of the air is non-linear. Air is a strong insulator [29] with conductivity of the order of 1014 S.m1 but beyond some electric solicitation, the air loses its insulating nature and locally becomes suddenly conductive. The ionization phenomenon is the only cause that can make the air conductor of electricity. The ionized channel becomes very conductive.

    Our mathematical context is periodic homogenization. We consider a microscopic scale ε, which represents the ratio between the diameter of the fiber and thickness of the composite material. So, we are trying to understand how the microscopic structure affects the macroscopic electromagnetic field behavior. Homogenization of Maxwell equations with periodically oscillating coefficients was studied in many papers. N. Wellander homogenized linear and non-linear Maxwell equations with perfect conducting boundary conditions using two-scale convergence in[26] and [27]. N. Wellander and B. Kristensson homogenized the full time-harmonic Maxwell equation with penetrable boundary conditions and at fixed frequency in [28]. The homogenized time-harmonic Maxwell equation for the scattering problem was done in F. Guenneau, S. Zolla and A. Nicolet [12]. Y. Amirat and V. Shelukhin perform two-scale homogenization time-harmonic Maxwell equations for a periodical structure in [5]. They calculate the effective dielectric ε and effective electric conductivity σ. They proved that homogenized Maxwell equations are different in low and high frequencies. The result obtained by two-scale convergence approach takes into account the characteristic sizes of skin thickness and wavelength around the material.

    On of the parameter we account for in our model: δ=1¯ω ¯σμ0, where ¯σ is the characteristic conductivity and ¯ω the order of the magnitude of the pulsation shares much with the definition of theoretical thickness skin δ=2¯ω¯σμ0. The thickness skin is the depth at which the surface current moves to a factor of e1. Indeed, at high frequency, the skin effect phenomenon appears because the current tends to concentrate at the periphery of the conductor. On the other side, at low frequencies the penetration depth is much greater than the thickness of the plate which means that a part of the electric field penetrates the composite plate. We use the theory of two-scale convergence introduced by G. Nguetseng [12] and developed by G. Allaire [3].

    The paper is organized as follows : in Section 2 we specify the geometry of the model and the dimensionless equations converting the problem into an equivalent one with which we work in the following sections. In Section 3 we perform the mathematical analysis of the model. In particular, we introduce the weak formulation of the problem for the electric field and we regularize it using divergence term. We establish the existence and uniqueness result for the regularized Maxwell equations thanks to Lax-Milgram Theorem. We conclude this section by estimate of the electric field. The last section is devoted to the homogenization of the problems for electric field using the two-scale convergence concept.


    2. Modeling

    This section is dedicated to the complete mathematical model we will study in this paper. First, we consider a problem that seems relevant with the perspective of propagation of the electromagnetic field in the air and in the skin of aircraft fuselage made of composite material. Secondly, we make a scaling of this model and finally we operate simplifications. If desired, the reader can go directly to the mathematical analysis knowing that the problem to be studied is given by (63), (68) equipped with boundary conditions (66), (67).


    2.1. Notations and setting of the problem

    We consider set ˜Ω={(˜x,˜y,˜z)R3, ˜y(¯L,d)} for ¯L and d two positive constants, with two open subsets ˜Ωa and ˜P (see Figure 1). The air fills ˜Ωa and we consider that the composite material, with two materials periodically distributed, stands in domain ˜P.

    Figure 1. The global domain.

    We assume that the thickness ¯L of the composite material is much smaller than its horizontal size. We denote by e the lateral size of the basic cell ˜Ye of the periodic microstructure of the material. The cell is composed of a carbon fiber in the resin. We define now more precisely the material, introducing:

    ˜P={(˜x,˜y,˜z)R3/¯L<˜y<0}, (2)

    which is the domain containing the material. Now we describe precisely the basic cell. For this we first introduce the following cylinder with square base:

    ˜Ze=[e2,e2]×[e,0]×R, (3)

    We consider α such that 0<α<1, and ˜Re=αe2. We set

    ˜De={(˜x,˜y)R2/(˜x2+(˜y+e2)2)<(˜Re)2}. (4)

    We define the cylinder containing the fiber as (see fig 1):

    ˜Ce=˜De×R. (5)

    Then the part of the basic cell containing the matrix is

    ˜YeR=˜Ze˜Ce, (6)

    and by definition, the basic cell ˜Ye is the couple

    (˜YeR,˜Ce). (7)

    The composite material results from a periodic extension of the basic cell. More precisely the part of the material that contains the carbon fibers is

    ˜Ωc=˜P{(ie,je,0)+˜Ce,iZ,jZ}, (8)

    where the intersection with ˜P limits the periodic extension to the area where stands the material. Set {(ie,je,0)+˜Ce,iZ,jZ} is a short notation for

    {(˜x,˜y,˜z)R3,iZ,jZ,(xb,yb,zb)˜Ce;˜x=xb+ie,˜y=yb+je,˜z=zb}. (9)

    In the same way the part of the material that contains the resin is

    ˜Ωr=˜P{(ie,je,0)+˜YeR}, (10)

    or equivalently

    ˜Ωr=˜P{(ie,je,0)+˜Ze˜Ce}=(R×(¯L,0)×R)˜Ωc. (11)

    So the geometrical model of our composite material is couple (˜Ωc,˜Ωr). Now, it remains to set the domain that contains the air:

    ˜Ωa={(˜x,˜y,˜z)/0˜y<d}. (12)

    We consider that d is of the same order as ¯L and we introduce the upper frontier ˜Γd={(˜x,˜y,˜z)/˜y=d} of domain ˜Ω. On this frontier we will consider that the electric field and magnetic field are given. We also introduce the lower frontier ˜ΓL={(˜x,˜y,˜z)/˜y=¯L} with those definitions we have ˜Ωa˜P=, ˜Ωc˜Ωr=, ˜P=ΩrΩc, ˜Ω=Ωa˜P=ΩaΩrΩc, and for any (˜x,˜y,˜z)˜Ω=˜Γd˜ΓL and, we write ˜n, the unit vector, orthogonal to ˜Ω and pointing outside ˜Ω. We have:

    ˜n=e2on˜Γd˜n=e2on˜ΓL. (13)

    In the following we need to describe what happens at the interfaces between resin and carbon fibers, and resin and air. So we define Γra={(˜x,˜y,˜z) /  ˜y=0} and Γcr the boundary of the set defined by (9).


    2.2. Maxwell equations

    In ˜Ω, we now write a PDE model that has to do with electromagnetic waves radiated from return stroke. We are well aware that the model we write is a simplified one. Nonetheless, it seems to be well dimensioned for our problem which consists in making homogenization. It is well known see Maxwell [17] the propagation of the electromagnetic field is described by the Maxwell equations which write:

    ˜Dt+טH=˜J, (14)
    ˜Bt+טE=0, (15)
    ˜D=˜ρ, (16)
    ˜B = 0, (17)

    in RטΩ.

    In (14)-(17), ∇ and ∇ are the curl and divergence operators. ˜E(t,x,y,z) is the electric field, ˜H(t,x,y,z) the magnetic field, ˜D(t,x,y,z) the electric induction, ˜B(t,x,y,z) the magnetic induction and ˜ρ(t,x,y,z) is the charges density see T. Abboud and I. Terrasse [1].

    System of Maxwell equations ((14)-(17)) is completed by the constitutive laws which are given in R×˜Ω by:

    ˜D=ϵ0ϵ˜E, (18)
    ˜B=μ0˜H. (19)

    where μ0 and ϵ0 are the permeability and permittivity of free space. ϵ is the relative permittivity of the domains defined by

    ϵ|~Ωa=1,ϵ|~Ωr=ϵr,ϵ|~Ωc=ϵc, (20)

    where ϵr and ϵc are positives constants. In order to account for energy transfer between the electromagnetic compartment and the propagation of the electric charges, we take for granted the Ohmic law, in RטΩ

    ˜J=σ˜E, (21)

    where σ is the electric conductivity. Its value depends on the location:

    σ|~Ωa=σa,σ|~Ωr=σr,σ|~Ωc=σc, (22)

    where ~Ωa, ~Ωr and ~Ωc were defined in (12), (10) and (8).

    Figure 2. Left: The global microstructure in 2D. Right: Z-cell of the periodic structure..

    2.3. Boundary conditions

    For mathematical as well as physical reasons we have to set boundary conditions on ~Γd and ~ΓL. On ~Γd we will write conditions that translate that ˜E and ˜H are given by the source located in ˜y=d. The way we chose consists in setting:

    ˜Eטn=˜Edטn;˜Hטn=˜HdטnonR×~Γd, (23)

    where ˜Ed, ˜Hd are functions defined on ˜Γd for any tR. On ˜ΓL, we chose something simple, i.e:

    טEטn=0onR×~ΓL, (24)

    that translate that ˜E does not vary in the ˜y-direction near ˜ΓL.

    Problem (14)-(21) supplemented with (23) and (24), is considered as containing all physics we want to account for. In the following we will consider simplifications of it.


    2.4. Time-harmonic Maxwell equations

    The first simplification we make, consists in considering the harmonic version of the Maxwell equations (14)-(22). This simplification is used in many references studying electromagnetic phenomena and especially for lightning applications [14], in spite of the fact that it considers implicitly that every fields and currents are waves of the form, for all ˜ω  R:

    a(˜x,˜y,˜z)cos(˜ωt+ϕ(˜x,˜y,˜z))=e[a(˜x,˜y,˜z)expi˜ωtexpiϕ(˜x,˜y,˜z)], (25)

    where ˜ω is the pulsation, ϕ(˜x,˜y,˜z) is the phase shift of the wave and a(˜x,˜y,˜z) is its amplitude. In particular, it supposes ˜Ed, ˜Hd in (23) are of the form, for all ˜wR:

    ˜Ed(t,˜x,˜z)=e(˜Ed(˜x,˜z)expi˜ωt), (26)
    ˜Hd(t,˜x,˜z)=e(˜Hd(˜x,˜z)expi˜ωt), (27)

    where ˜Ed and ˜Hd take into account the amplitude and the phase shift of their corresponding fields. Taking (21) into account, the time-harmonic Maxwell equations, which describe the electromagnetic radiation, are written:

    טHi˜ωϵ0ϵ˜E=σ˜E,  Maxwell - Ampere equation (28)
    טE+i˜ωμ0˜H=0,  Maxwell - Faraday equation (29)
    (ϵ0ϵ˜E)=˜ρ, (30)
    (μ0˜H)=0, (31)

    where ˜E(t,˜x,˜y,˜z)=e(˜E(˜x,˜y,˜z)expi˜ωt) and ˜H(t,˜x,˜y,˜z)=e(˜H(˜x,˜y,˜z)expi˜ωt), (˜x,˜y,˜z)˜Ω. The magnetic field ˜H can be directly computed from the electric field ˜E

    ˜H=1iωμ0טE. (32)

    Now, for the electric approach, taking the curl of equation (32) yields an expression of טH in term of ×טE. Inserting טH in (28) we get the following equation for the electric field:

    ×טE+(˜ω2μ0ϵ0ϵ+i˜ωμ0σ)˜E=0in˜Ω. (33)

    Taking the divergence of the equation (28) yields the natural gauge condition:

    [(i˜ωϵ0ϵ+σ)˜E]=0in˜Ω. (34)

    Notice that i˜ωϵ0+σ is equal to i˜ωϵ0+σa in ~Ωa, to i˜ωϵ0ϵr+σr in ~Ωr and to i˜ωϵ0ϵc+σc in ~Ωc, those quantities being all nonzero. Then (34) is equivalent to:

    ˜E|Ωa=0  in  ~Ωa, ˜E|Ωr=0  in  ~Ωr,   ˜E|Ωc=0  in  ~Ωc. (35)

    with the transmission conditions

    (i˜ωϵ0+σa)˜E|~Ωa.˜n=(i˜ωϵ0ϵr+σr)˜E|~Ωr.˜n  on  ~Γra,(i˜ωϵ0ϵr+σr)˜E|~Ωr.˜n=(i˜ωϵ0ϵc+σc)˜E|~Ωc.˜n  on  ~Γcr. (36)

    Summarizing, we finally obtain the PDE model:

    ×טE+(˜ω2μ0ϵ0ϵ+i˜ωμ0σ)˜E=0in˜Ω. (37)

    According to the tangential trace of the Maxwell-Faraday equation (29) we obviously obtain that using boundary condition (23), is equivalent to using:

    טE×e2=i˜ωμ0˜Hd(˜x,˜z)×e2on~Γd (38)

    where ˜Hd is defined in (27) and where we used (13). And on ~ΓL we have the following boundary condition:

    טE×e2=0on~ΓL. (39)

    2.5. Scaling

    In this subsection we propose a rescaling of system ((37)-(39)), we will consider a set of characteristic sizes related to our problem. Physical factors are then rewritten using those values leading to a new set of dimensionless and unitless variables and fields in which the system is rewritten. The considered characteristic sizes are : ¯ω the characteristic pulsation, ¯σ the characteristic electric conductivity, ¯E the characteristic electric magnitude, ¯H the characteristic magnetic magnitude. We also use the already introduced thickness ¯L of the plate ˜P. We then introduce the dimensionless variables: x=(x,y,z) with x=˜x¯L, y=˜y¯L, z=˜z¯L and fields E, H and σ that are such that

    {E(ω,x)=1¯E˜E(¯ωω,¯Lx,¯Ly,¯Lz),H(ω,x)=1¯H˜H(¯ωω,¯Lx,¯Ly,¯Lz),σ(x)=1¯σ˜σ(¯Lx,¯Ly,¯Lz), (40)

    Taking (22) into account, σ also reads:

    σ(x)=σa¯σif0¯Lyd,σ(x)=σr¯σif(¯Lx,¯Ly,¯Lz)˜Ωr,σ(x)=σc¯σif(¯Lx,¯Ly,¯Lz)˜Ωc. (41)

    Doing this gives the status of units to the characteristic sizes. Since, for instance:

    Ex(ω,x)=¯L¯E˜E˜x(¯ωω,¯Lx,¯Ly,¯Lz), (42)

    using those dimensionless variables and fields and taking (41) into account, equation (37) gives:

    ¯E ××E(ω,x)(¯L2¯ω2c2ϵω2+i¯σ ¯ω ω¯L2μ0σ(x,ω))¯EE(ω,x,y,z)=0, (43)

    for any (ω,x) such that (¯ωω,¯Lx,¯Ly,¯Lz)˜Ω. Now we exhibit

    ¯λ=2πc¯ω, (44)

    which is the characteristic wave length and

    ¯δ=1¯ω ¯σμ0, (45)

    which is the characteristic skin thickness. Using those quantities equation (43) reads, for any (ω,x)Ω:

    ××E(ω,x)+(4π2¯L2¯λ2ω2+i¯L2¯δ2σa¯σω)E(ω,x)=0when0¯Lyd,××E(ω,x)+(4π2¯L2¯λ2ϵr ω2+i¯L2¯δ2σr¯σω)E(ω,x)=0when(¯Lx,¯Ly,¯Lz)˜Ωr,××E(ω,x)+(4π2¯L2¯λ2ϵc ω2+i¯L2¯δ2σc¯σω)E(ω,x)=0when(¯Lx,¯Ly,¯Lz)˜Ωc. (46)

    In the following expressions, ¯L¯λ and ¯L¯δ appearing in the equations above will be rewritten in terms of a small parameter ε.

    The boundary conditions are written

    ×E(ω,x)×e2=iω¯ωμ0¯L¯E~Hd(¯Lx,¯Lz)×e2  when  (¯Lx,¯Ly,¯Lz)~Γd,×E(ω,x)×e2=0  when  (¯Lx,¯Ly,¯Lz)~ΓL. (47)

    The characteristic thickness of the plate ¯L is about 103m and the size of the basic cell e is about 105m. Since e is much smaller than the thickness of the plate ¯L, it is pertinent to assume the ratio e¯L equals a small parameter ε:

    e¯L˜ 102=ε. (48)

    Then, in what concerns the characteristic pulsation ¯ω, in the tables below we consider several values. The lightning is seen as a low frequency phenomenon. Indeed, energy associated with radiation tracers and return stroke are mainly burn by low and very low frequencies (from 1 kHz to 300 kHz). Components of the frequency spectrum are however observed beyond 1GHz (see [16]). So, in the case when we want to catch low frequency ie we consider ¯ω=100 rad/s, (in our study we will consider ¯ω=106 rad/s), for medium frequency we set ¯ω=1010 rad/s and for high frequency phenomena ¯ω=1012 rad/s. Then, concerning the characteristic electric conductivity it seems to be reasonable to take for ¯σ the value of the effective electric conductivity of the composite material. Yet this choice implies to compute a coarse estimate of this effective conductivity at this level.

    For this we take into account that the composite material is composed of carbon fibers and epoxy resin. The resin can be doped, which increases strongly its conductivity, or not. The tables below summarize the cases when the resin is doped and also when the resin is not doped. We also account for the fact there is not only one effective electric conductivity but a first one in the fiber direction : the effective longitudinal electric conductivity (in cases 1, 2, 5 and 6 of the tables below), and a second effective electric conductivity, in the direction transverse to the fibers (considered in cases 3, 4, 7 and 8). In this context, we consider the basic model which is based on the electrical analogy and the law of mixtures. It corresponds to the Wiener limits: the harmonic average and the arithmetic average. The effective values are the extreme limits of the conductivity of the composite introduced by Wiener in 1912 see S. Berthier p 76 [7].

    The effective longitudinal electric conductivity corresponding of the upper Wiener limit is expressed by the equation:

    ¯σ=σlong=fc σc+(1fc) σr, (49)

    where fc=πα24 is the volume fraction of the carbon fiber.

    The effective transverse electric conductivity corresponding of the lower Wiener limit is expressed by

    ¯σ=σtrans=1fcσc+(1fc)σr. (50)

    For the computation, we take values close to reality. We consider composite materials with similar proportions of carbon and resin, this means that α is close to 12. When the resin is not doped σr1010 S.m1 is much smaller than σc40000 S.m1. Then, ¯σ=σlong is close to πα24 σcσc and ¯σ=σtrans is close to σr(1πα24)σr.

    Now, we express the electric conductivity of the air in terms of ¯σ, we consider two possibilities. The first one is relevant for a situation with a ionized channel. The second one of situation with a strong atmospheric electric field but without a ionized channel. In this situation air is not ionized and has a low conductivity. All possible situations are gathered in the tables below. Cases 5 to 8 are associated with the first situation with air conductivity σa being σlightning=4242 S.m1 for an ionized lightning channel see [13]. Cases 1 to 4, to the second one, with σa=1014 S.m1.

    All calculations of the different cases of the tables are detailed in Annex A. In our study we consider the case 6 for ω=106~rad.s1, which corresponds to the air ionized, a resin doped and the effective longitudinal electric conductivity of the carbon fibers.

    As it is well known the tables confirm that at high frequencies the thickness of the plate is much greater than the skin depth. This one depends on ¯σ and ¯ω and decreases strongly for high conductivity or high frequencies. For ¯ω=1010 rad.s1 and ¯σ=4104 S.m1, the effective conductivity in the direction of the carbon fibers, which the skin effect phenomenon appears. Indeed, for high frequencies ω=1012rad.s1 and when ¯σ is the effective conductivity is in direction of the carbon fibers i.e. in high conductivity, ¯δ=104 m. In low frequencies and low conductivity δ is large so the electromagnetic wave can penetrate the composite material. The high conductivity limits the penetration of the electromagnetic wave to a boundary layer whose depth is about ¯δ.

    Now, we will discuss on the values of ¯E and ¯ρ. It seems that the density of electrons in a ionized channel is about 1010 part.m3. Hence we take ¯ρ=1010. When the air is not ionized, the charge density is much smaller, and we choose: ¯ρ=1.

    For the boundary conditions, in the context of the case 6 and ¯ω=106 rad/s, we consider the peak of the current of the first return stroke. Then the magnetic field magnitude ¯H is ¯Hd given by (1).

    ×E(x,ω)×e2=iω¯ωμ0¯L¯E¯HdHd(x,z)×e2, (51)

    where ¯HdHd(x,z)=~Hd(¯Lx,¯Lz) and where ¯ωμ0¯L¯E¯Hd being of order 1, with the characteristic electric field ¯E=20 kV/m.

    From the physical spatial coordinates (˜x,˜y,˜z)˜Ω we define y=(ξ,ν,ζ) with ξ=˜xe,ν=˜ye,ζ=˜ze or equivalently ξ=xε,ν=yε,ζ=zε. And we now introduce Y, the basic cell. It is built from: Z=[12,12]×[1,0]×R and the set C=D×R with the disc D defined by:

    D={(ξ,ν)R2 /ξ2+(ν+12)2<R2}, (52)

    and R=α2. The set Ωc is then defined as:

    Ωc={(i,j,0)+C,iZ,jZ}. (53)

    We denote Yr as Yr=ZC and then the set

    Ωr={(i,j,0)+Yr,iZ,jZ}. (54)

    Then unit cell Y is defined as Y=(Yr,C). Finally, we define the domain Ωa:

    Ωa={y=(ξ,ν,ζ)  /  ν>0}. (55)

    Using this, we will give a new expression of the sets in which the variables range in equations (46). We see the following:

    (¯Lx,¯Ly,¯Lz)˜Ωr{(¯Lx,¯Ly,¯Lz)˜P,(¯Lex,¯Ley,¯Lez)Ωr, (56)

    i.e.

    (¯Lx,¯Ly,¯Lz)˜Ωr{(¯Lx,¯Ly,¯Lz)˜P,(xε,yε,zε)Ωr. (57)

    In the same way:

    (¯Lx,¯Ly,¯Lz)˜Ωc{(¯Lx,¯Ly,¯Lz)˜P,(xε,yε,zε)Ωc, (58)

    and:

    0¯Lyd{yR2¯Lyd, (59)

    or

    (¯Lx,¯Ly,¯Lz)˜Ωa{¯Lyd(xε,yε,zε)Ωa. (60)

    We define:

    Σε(y)=Σε(ξ,ν,ζ)={Σεain  Ωa,Σεrin  Ωr,Σεcin  Ωc, (61)

    where Σεa=σa¯σ¯L2¯δ2,Σεr=σr¯σ¯L2¯δ2 and Σεc=σc¯σ¯L2¯δ2 have their expressions in term of ε given from Tables above depending on the case we are interested in. The detail of this expressions are in appendix B. The model that we present is the case ω=106~rad.s1, η=5, Σεa=ε, Σεr=ε4 and Σεc=1.

    Defining also mapping

    ψε:R3R3(x,y,z)(xε,yε,zε), (62)

    we can set Ωεa as ψ1ε(Ωa)(R×[0,d¯L]×R), Ωεr as ψ1ε(Ωr)˜P and Ωεc as ψ1ε(Ωc)˜P. We also define the boundaries Γd={xR3,  y=d¯L} and ΓL={xR3,  y=¯L} and interfaces Γra={xR3,  y=0} and Γεcr=Ωc. Hence equation (46) reads:

    ××Eε+(ω2εηϵ+i ω σε(x,y,z))Eε=0   in  Ω, (63)

    where Ω=ΩεaΩεrΩεc={xR3,1<y<d¯L} does not depend on ε. Only its partition in Ωεa, Ωεr and Ωεc is ε-dependent where

    σε(x,y,z)=Σε(xε,yε,zε), (64)

    with Σε given by (61) and

    εη=4π2¯L2λ2, (65)

    with the value of η0 extracted from Tables, and where we replace E by Eε, to clearly state that it depends on ε.

    Equation (63) is provided with the following boundary conditions:

    ×Eε×e2=iωHd(x,z)×e2 on Γd, (66)

    coming from (51). And, coming from (47),

    ×Eε×e2=0  on  ΓL. (67)

    From (63) we can deduce the condition on the divergence of Eε which can be written in two ways. As previously in (34), (35) and (36) we obtain:

    [(ω2εηϵ+iωσε)Eε]=0  in  Ω, (68)

    which will be preferentially used with (63) and its second one is

    Eε|Ωεa=0  in  Ωεa,  Eε|Ωεr=0  in  Ωεr,  Eε|Ωεc=0  in  Ωεc, (69)

    with the transmission conditions on the interfaces Γra and Γεcr

    (ω2εη+iωΣεa) Eε|Ωεan|Ωεa=(ω2εηϵr+iωΣεr) Eε|Ωεrn|Ωεr  on  Γra,(ω2εηϵr+iωΣεr) Eε|Ωεrn|Ωεr=(ω2εηϵc+iωΣεc) Eε|Ωεcn|Ωεc  on  Γεcr. (70)

    Before treating mathematically the question we are interested in, we make a last simplification. Since it seems clear that physical relevant phenomena occur in the upper part of the plate. The boundary condition on the lower boundary of the plate has very little influence on the physics of what happens in the upper part, we consider that the lower boundary of Ω is located in y= in place of y=1, making the second boundary condition useless. Besides, as ¯L and d are of the same order it seems reasonable to set Γd={xR3,  y=1} and consequently

    {Ω={xR3,y<1}=ΩεaΩεrΩεc,  with,Ωεa=ψ1ε(Ωa),Ωεr=ψ1ε(Ωr),Ωεc=ψ1ε(Ωc), (71)

    with ψε defined in (62). We have that the border of Ω is Γd. In the following section we will establish existence and uniqueness results.


    3. Mathematical analysis of the models


    3.1. Preliminaries

    We are going to make precise the variational formulation. First of all, we need to introduce the following functional spaces. We have the standard function spaces L2(Ωε)=[L2(Ωε)]3

    H(curl,Ω)={uL2(Ω): ×uL2(Ω)},H(div,Ω)={uL2(Ω):  uL2(Ω)}, (72)

    with the usual norms:

    u2H(curl,Ω)=u2L2(Ω)+×u2L2(Ω),u2H(div,Ω)=u2L2(Ω)+u2L2(Ω). (73)

    They are well known Hilbert spaces.

    We use in this paper, the trace spaces H12(curl,Γd) and H12(div,Γd) defined by

    H12(curl,Γd)={uH12(Γd,R3), (nu)|Γd=0, curlΓduH12(Γd,R3)}, (74)
    H12(div,Γd)={uH12(Γd,R3), (nu)|Γd=0, divΓduH12(Γd,R3)} (75)

    where the surface divergence divΓdu and the surface rotation curlΓdu are defined by

    (divΓdu,V)L2(Γd)=(u,ΓdV)L2(Γd,R3),   VC1(Γd)curlΓdu=n(×u|Γd) (76)

    and the surface gradient ΓdV is defined by the orthogonal projection of on Γd, n denotes the outward unit vector normal to Γd.

    Finally we recall the trace theorems, see J.C Nédélec [19] for the demonstration, stating that the traces mappings

    γT:H(curl,Ω)H12(curl,Γd), that assigns any uH(curl,Ω) its tangential components n×(u×n), is continuous and surjective, that is:

    γT(u)H12(curl,Γd)CγTuH(curl,Ω), uH(curl,Ω)

    γt:H(curl,Ω)H12(div,Γd), that assigns any uH(curl,Ω) its tangential components u×n, is continuous and surjective:

    γt(u)H12(div,Γd)CγtuH(curl,Ω), uH(curl,Ω).

    Moreover, H12(div,Γd) is the dual of H12(curl,Γd) and one has the Green's formula:

    Ω(×uVu×V)dx=u×n,VTΓd (u,V)H(curl,Ω). (77)

    We define the next space:

    X(Ω)={uH(curl,Ω) | u|ΩεaL2(Ωεa),u|ΩεrL2(Ωεr), u|ΩεcL2(Ωεc)}. (78)

    Our variational space is:

    Xε(Ω)={uX(Ω) | (ω2εη+iωσε|Ωεa)u|Ωεae2=(ω2εηεr+iωσε|Ωεr)u|Ωεre2,           (ω2εηϵr+iωσε|Ωεr)u|Ωεrnε|Ωεr=(ω2εηεc+iωσε|Ωεc)u|Ωεcnε|Ωεc. (79)

    Finally

    Xε(Ω)={uX(Ω) | (ω2εη+iωΣεa)u|Ωεae2=(ω2εηεr+iωΣεr)u|Ωεre2,           (ω2εηϵr+iωΣεr)u|Ωεrnε|Ωεr=(ω2εηεc+iωΣεc)u|Ωεcnε|Ωεc}. (80)

    This space is equipped with the norm

    u2Xε(Ω)=u2L2(Ω)+u|Ωεa2L2(Ωεa)+u|Ωεr2L2(Ωεr)+u|Ωεc2L2(Ωεc)+×u2L2(Ω).

    3.2. Weak formulation

    Now, we introduce the variational formulation of our problem (63), (66) and (67) for the electric field. Integrating (63) over Ω and using the Green's formula and (66) we obtain

    Ω×EεׯV dx+Ωεa(ω2εη+iωΣεa)Eε¯V dx+Ωεc(ω2εηϵc+iωΣεc)Eε¯V dx+Ωεr(ω2εηϵr+iωΣεr)Eε¯V dx=Γd(×Eε×e2)¯VT dσ=ΓdiωHd×e2¯VT dσ (81)

    where ¯V is the complex conjugate of V and VT=(e2×V)×e2. We introduce the sesquilinear form depending on parameters η and ε:

    {For Eε,VXε(Ω),aε,η(Eε,V)=Ω×EεׯV dx+i=a,r,cΩεi(ω2εηϵi+iωΣεi) Eε¯V dx. (82)

    Hence, the weak formulation of (63), (66) and (67) that we will use is the following:

    {Find EεXε(Ω)  such as    VXε(Ω)  we have:aε,η(Eε,V)=iωΓdHd×e2¯VT dσ. (83)

    Integrating by parts in the variational formulation (81), we find the following transmission problem:

    {××Eε+(ω2εη+i ω Σεa)Eε=0      in  Ωεa,××Eε+(ω2εηϵr+i ω Σεr)Eε=0     in  Ωεr,××Eε+(ω2εηϵc+i ω Σεc)Eε=0     in  Ωεc.Eε|Ωεa×e2=Eε|Ωεr×n|Ωεr  on  Γra,Eε|Ωεr×n|Ωεr=Eε|Ωεc×n|Ωεc  on  Γεcr,×Eε|Ωεa×e2=×Eε|Ωεr×n|Ωεr  on  Γra,×Eε|Ωεr×n|Ωεr=×Eε|Ωεc×n|Ωεc  on  Γεcr, (84)

    where e2 is the unit outward normal to Ωεa, n|Ωεr is the unit outward normal to Ωεr and n|Ωεc is the unit outward normal to Ωεc. We refer to Annex C for the proof that transmission problem (84) is equivalent to ((63), (66), (67), (69)).


    3.3. Regularized Maxwell equations for the electric field

    The sesquilinear form aε,η is not coercive on Xε(Ω), so we regularize it adding terms involving the divergence of Eε in Ωεa, Ωεr and Ωεc. Thanks to the additional terms, existence and uniqueness of the regularized variational formulation solution will be established by the Lax-Milgram theory. Let s be an arbitrary positive number, we define the regularized formulation of problem (83):

    {Find EεXε(Ω)  such that for any  VXε(Ω)aε,ηR(Eε,V)=aε,η(Eε,V)+sΩεaEε¯V dx+sΩεrEε¯V dx+sΩεcEε¯V dx=iωΓdHd×e2¯VT dσ. (85)

    For any ε>0 and any η0, sesquilinear form aε,ηR(.,.) is continuous over Xε(Ω) thanks to the continuity conditions. We will show that it is also coercive. The following proposition was inspired by article [9] Lemma 1.1.

    Proposition 3.1. For any ε>0, for any η0 and for any s>0, there exists a positive constant ω0 which does not depend on ε and such that for all ω(0,ω0), there exists a positive constant C0 depending on εr,εc,s,ω but not on ε such that:

     EεXε(Ω),  [exp(iπ4) aε,ηR(Eε,Eε)]C0EεXε(Ω) (86)

    Proof. We have:

    [exp(iπ4) aε,ηR(Eε,Eε)]=aεR(Eε,Eε)Ωεaω2εη|Eε|2 dxΩεrω2εηϵr|Eε|2 dxΩεcω2εηϵc|Eε|2 dx. (87)

    with

    aεR(Eε,Eε)=Ω|×Eε|2 dx+sΩεa|Eε|2 dx                 +sΩεr|Eε|2 dx+sΩεc|Eε|2 dx                 +ΩεaωΣεa|Eε|2 dx+ΩεrωΣεr|Eε|2 dx                 +ΩεcωΣεc|Eε|2 dx. (88)

    We have the following estimate:

    |aεR(Eε,Eε)|min{1,ω,s}(×Eε2L2(Ω)+Eε2L2(Ωεa)+Eε2L2(Ωεr)                  +Eε2L2(Ωεc)+Eε2L2(Ω)). (89)

    Then we have:

     aεR(Eε,Eε)| min{1,ω,s}Eε2Xε(Ω). (90)

    Returning to formulation (86), for η0, since max(Σεa,Σεr,Σεc)>εη, inequality (86) is valid with C0=min{1,ω,s} as soon as ω2min{1,εr,εc}<min{1,ω,s} or ω<min{1,ω,s}min{1,εr,εc}. This ends the proof of Proposition 3.1.

    Thanks to Proposition 3.1 we can state the existence and uniqueness of the solution to regularized problem (85).

    Theorem 3.2. Under the assumptions of Proposition 3.1, there exists a unique solution Eε to regularized problem (85).

    proof The sesquilinear form aε,ηR is continuous, bounded, coercive thanks to Proposition 3.1 and the right hand side is continuous on Xε(Ω), then problem (85) has a unique solution in Xε(Ω) thanks to the Lax-Milgram Lemma.


    3.4. Existence, uniqueness and estimate

    Theorem 3.3. For any ε>0, for any η0, there exists a positive constant ω0 which does not depend on ε and such that for all ω(0,ω0), there exists a unique solution of (84) or ((63), (66), (67), (69)).

    proof We show that for an appropriate choice of s that Eε satisfies all equations (84) or ((63), (66), (67), (69)). It is obvious that any solution of (84) or of ((63), (66), (67), (69)) is also solution to (85). Indeed, since from (84) or from ((63), (66), (67), (69)) we have Eε|Ωεa=0, Eε|Ωεr=0, Eε|Ωεc=0, the additional terms sΩεaEε¯V dx + sΩεrEε¯V dx + sΩεcEε¯V dx cancel in (85).

    Uniqueness follows from that if Eε1 and Eε2 are two solutions to (63) with the boundary condition (67) their difference eε=Eε2Eε1 satisfies the problem (63) with (67). Then it comes

    Ω|×eε|2 dx+Ωεa(ω2εη+iωΣεa)|eε|2 dx                   +Ωεc(ω2εηϵc+iωΣεc)|eε|2 dx+Ωεr(ω2εηϵr+iωΣεr)|eε|2 dx                   =0. (91)

    Taking the imaginary part of the expression we get ΩεaωΣεa|eε|2 dx+ΩεcωΣεc|eε|2 dx+ΩεrωΣεr|eε|2 dx=0 and then eε=0.

    Let us consider the reciprocal assertion, according to the same proof of S. Hassani, S. Nicaise, A. Maghnouji in [23], we define H10(Ωεc,Δ) the subspace of ψH10(Ωεc) such that Δ(ψ)L2(Ωεc).

    Let Eεs be the solution of the regularized formulation (85). In (85) we take a test function V=ψ where ψH10(Ωεc,Δ), extended by zero outside Ωεc. We get:

    ΩεcsEεs(ψ) dx+Ωεc(ω2εηϵc+iωΣεc)Eεsψ dx=0. (92)

    By Green's formula, ψH10(Ωεc,Δ), we obtain:

    ΩεcEεs(Δψ+ω2εηϵciωΣεcs ψ)  dx=0. (93)

    Thus, if we choose s such that ω2εηϵciωΣεcs is not an eigenvalue of (Δdir,Ωεc): the Laplacian operator in Ωεc with Dirichlet condition on its boundary, then for all φL(Ωεc)2 there exists ψH10(Ωεc,Δ) solution of

    Δψ+ω2εηϵciωΣεcs ψ=φ. (94)

    Then, we conclude that

    Eεs|Ωεc=0. (95)

    A similar argument in Ωεa yields Eεs|Ωεa=0 for s such that ω2εηiωΣεas is not an eigenvalue of (Δdir,Ωεa). In the same way, we obtain in Ωεr, Eεs|Ωεr=0 with s such that ω2εηϵriωΣεrs is not an eigenvalue of (Δdir,Ωεr).

    Hence Eεs=0 in Ωεc, this cancels the additional term sΩεcEεs¯V dx in (85). In the same way, Eεs=0 in Ωεr and Eεs=0 in Ωεa cancel sΩεrEεs¯V dx and sΩεaEεs¯V dx in (85). So, (85) becomes (81). Applying Green's formula, we find (63).

    Theorem 3.4. Under the assumptions of Theorem 3.2, EεXε(Ω) solution of (85) satisfies

    EεXε(Ω)C (96)

    with C=CγtCγTC0HdH(curl,Ω)

    Proof. The sesquilinear form aε,ηR(Eε,V) is coercive, weak formulation (85) becomes:

    C0Eε2Xε(Ω)(exp(iπ4)aε,ηR(Eε,Eε))                   |exp(iπ4)aε,ηR(Eε,Eε)|=|aε,ηR(Eε,Eε)|                   |ΓdiωHd×e2EεT  dσ|                  Hd×e2H12(div,Γd)EεTH12(curl,Γd)                  CγtCγTHd×e2H(curl,Ω)EεH(curl,Ω) (97)

    where EεT=e2×(Eε×e2) and the continuous dependence of the trace norm with C=CγtCγTC0HdH(curl,Ω) gives:

    Eε2Xε(Ω)CEεH(curl,Ω)CEεXε(Ω). (98)

    4. Homogenization

    With the aim to obtain a convergence result for the problem (63), (66) and (67), we propose an approach based on two-scale convergence. This concept was introduced by G. Nguetseng [3, 4] and specified by G. Allaire [21, 22] which studied properties of the two-scale convergence. M. Neuss-Radu in [20] presented an extension of two-scale convergence method to the periodic surfaces. Many authors applied two-scale convergence approach D. Cionarescu and P. Donato [8]}, N. Crouseilles, E. Frénod, S. Hirstoaga and A. Mouton [10], Y. Amirat, K. Hamdache and A. Ziani [2] and also A. Back, E. Frénod [6]. This mathematical concept were applied to homogenize the time-harmonic Maxwell equations S. Ouchetto, O. Zouhdi and A. Bossavit [24], H.E. Pak [25].

    In our model, the parallel carbon cylinders are periodically distributed in direction x and z, as the material is homogenous in the y direction, we can consider that the material is periodic with a three directional cell of periodicity. In other words, introducing Z=[12,12]×[1,0]2, function Σε given by (61) is naturally periodic with respect to (ξ,ζ) with period [12,12]×[1,0] but it is also periodic with respect to y with period Z.

    Now, we review some basis definitions and results about two-scale convergence.


    4.1. Two-scale convergence

    We first define the function spaces

    {H#(curl,Z)={uH(curl,R3):  u  is Z-periodic}H#(div,Z)={uH(div,R3):  u  is Z-periodic} (99)

    and where H(curl,R3) and H(div,R3) are defined by (72) with Ωε replaced by R3. We introduce

    L2#(Z)={uL2(R3),u is Z-periodic}, (100)

    and

    H1#(Z)={uH1(R3),u is Z-periodic}, (101)

    where H1(R3) is the usual Sobolev space on R3. First, denoting by C0#(Z) the space of functions in C0(R3) and Z-periodic, C00(R3) the space of continuous functions over R3 with compact support, we have the following definitions:

    Definition 4.1. A sequence uε(x) in L2(Ω) two-scale converges to u0(x,y) L2(Ω,L2#(Z)) if for every V(x,y)C00(Ω,C0#(Z))

    limε0Ωuε(x)V(x,x/ε) dx=ΩZu0(x,y)V(x,y) dxdy. (102)

    Proposition 4.2. If uε(x) two-scale converges to u0(x,y)L2(Ω,L2#(Z)), we have for all v(x)C0(¯Ω) and all w(y)L2#(Z)

    limε0Ωuε(x)v(x)w(xε) dx=ΩZu0(x,y)v(x)w(y) dxdy. (103)

    Theorem 4.3. (Nguetseng) Let uε(x)L2(Ω). Suppose there exists a constant c>0 such that for all ε

    uεL2(Ω)c.

    Then there exists a subsequence of ε (still denoted ε) and u0(x,y)L2(Ω,L2#(Z)) such that:

    uε(x)u0(x,y). (104)

    Proposition 4.4. Let uε(x) be a sequence of functions in L2(Ω), which two-scale converges to a limit u0(x,y)L2(Ω,L2#(Z)). Then uε(x) converges also to u(x)=Zu0(x,y)dy in L2(Ω) weakly. Furthermore, we have

    limε0uεL2(Ω)u0L2(Ω×Y)uL2(Ω). (105)

    Remark 4.5.:-For any smooth function u(x,y), being Z-periodic in y, the associated sequence uε(x)=u(x,xε) two-scale converges to u(x,y).

    -Any sequence uε that converges strongly in L2(Ω) to a limit u(x), two-scale converges to the same limit u(x).

    -If uε admits an asymptotic expansion of the type uε(x)=u0(x,x/ε)+εu1(x,x/ε)+ε2u2(x,x/ε)+..., where the functions ui(x,y) are smooth and Z-periodic in y, two-scale convergence allows to identify the first term of the expansion u0(x,y) with the two-scale limit of uε and the two-scale limit of uε(x)u0(x,xε)ε with u1(x,y) see (Frénod, Raviart and Sonnendrucker [11]).

    Proposition 4.6. Let uε(x) in L2(Ω). Suppose there exists a constant c>0 such that for all ε

    uεL2(Ω)c.

    Up to a subsequence, uε(x) two-scale converges to u0(x,y)L2(Ω,L2#(Z)) such that:

    u0(x,y)=u(x)+~u0(x,y), (106)

    where ~u0(x,y)L2(Ω,L2#(Z)) satisfies

    Z~u0(x,y) dy=0, (107)

    and u(x)=Zu0(x,y) dy is a weak limit in L2(Ω).

    Proof. uε(x) is bounded in L2(Ω), then by application of Theorem 4.3, we get the first part of the proposition. Furthermore by defining ~u0 as

    ~u0(x,y)=u0(x,y)Zu0(x,y)dy, (108)

    we obtain the decomposition of u0.

    Defining x=(x;y;z), y=(ξ;ν;ζ), we have

    Proposition 4.7. Let uε(x) be bounded in H(curl,Ω). Then, up to a subsequence, there exists a function u1L2(Ω,H#(curl,Z)) such that

    ×uε(x)x×u0(x,y)+y×u1(x,y), (109)

    where u0 is given by Proposition 4.6.

    Proof. From Theorem 4.3, since uε and ×uε are bounded in L2(Ω) then, up to a subsequence, they two-scale converge to u0(x,y)L2(Ω,L2#(Z)) and η0(x,y)L2(Ω,L2#(Z)). So we have for all V(x,y)C00(Ω;C0#(Z)):

    limε0Ωuε(x)V(x,x/ε) dx=ΩZu0(x,y)V(x,y)dxdy, (110)
    limε0Ω×uε(x)V(x,x/ε) dx=ΩZη0(x,y)V(x,y)dxdy. (111)

    Next, by integration by parts, we have:

    Ω×uε(x)V(x,x/ε) dx=Ωuε(x)(x×V(x,x/ε)+1εy×V(x,x/ε)) dx. (112)

    If we choose a test function VC00(Ω,C0#(Z)) such that y×V=0, passing to the limit in the left-hand side (111) we get

    Ωx×uε(x)V(x,x/ε) dxΩZu0(x,y)x×V(x,y) dxdy                                          =ΩZx×u0(x,y)V(x,y) dxdy. (113)

    This means that with the difference between (111) and (113):

    ΩZ[η0(x,y)x×u0(x,y)]V(x,y) dxdy=0, (114)

    for all functions VC10(Ω) with y×V=0. It follows that function η0(x,y)x×u0(x,y) is orthogonal to functions with zero rotational in L2(Ω,H#(curl),Z). This implies that there exists a function u1L2(Ω,H#(curl,Z)) such that

    y×u1(x,y)=η0(x,y)x×u0(x,y). (115)

    Thus

    ×uε(x)x×u0(x,y)+y×u1(x,y). (116)

    Proposition 4.8. Let uε(x) be a bounded sequence in H(curl,Ω). Then a subsequence uε(x) can be extrated from ε such that, letting ε0

    uε(x)u(x)+yΦ(x,y). (117)

    where ΦL2(Ω,H1#(Z)) is a scalar-valued function and where uL2(Ω). And we have

    ×uε(x)x×u(x)  weaklyin  L2(Ω). (118)

    where u(x) is given by Proposition 4.6.

    proof. Proof of (117), for any V(x,y)C10(Ω,C1#(Z)), we have

    Ω×uε(x)V(x,xε) dx=Ωuε(x){x×V(x,xε)+1εy×V(x,xε)} dx. (119)

    Multiplying by ε we have

    εΩ×uε(x)V(x,xε) dx=Ωuε(x){εx×V(x,xε)+y×V(x,xε)} dx. (120)

    Taking the two-scale limit as ε0 we obtain

    0=ΩZu0(x,y)y×V(x,y) dxdy, (121)

    which implies that y×u0(x,y)=0. Thus u0(x,y) is a gradient with respect to the variable y for some scalar function Φ(x,y). And according to Proposition (4.6) u0(x,y) can be written as u0(x,y)=u(x)+yΦ(x,y), where u(x)=Zu0(x,y)dy for some scalar function Φ(x,y).

    Next, we choose a test function V(x)L2(Ω). Integration by parts yields:

    limε0Ω×uε(x)V(x) dx=limε0Ωuε(x)×V(x) dx=ΩZu0(x,y) dy×V(x) dx=Ω×u(x)V(x) dx. (122)

    These results are important properties of the two-scales convergence. We note that the usual concepts of convergence do not preserve information concerning the micro-scale of the function. However, the two-scale convergence preserves information on the micro-scale.


    4.2. Homogenized problem

    We will explore in this section the behavior of electromagnetic field Eε using the two-scale convergence to determine the homogenized problem. We place in the context of the case 6 with δ>L and ¯ω=106rad.s1, then we have η=5 and Σεa=ε, Σεr=ε4, Σεc=1 which gives the following equation:

    ××Eεω2ε5k(ϵ)Eε+iω[(1εC(xε)+ε41εR(xε))1{y<0}+ε1{y>0}]Eε=0, (123)

    where for a given set A, 1A stands for the characteristic function of A and where 1εA(x)=1A(xε), hence 1εC and 1εR are the characteristic functions of the sets filled by carbon fibers and by resin. And where k(ϵ)=(ϵc1εC(x)+ϵr1εR(x))1{y<0}+1{y>0}.

    Remark 4.9. We recall that ϵc and ϵr are respectively the relative permittivity of the carbon fibers and the resin. You should not confused with the microscopic scale ε.

    On this purpose, we have the following Theorem:

    Theorem 4.10. Under assumptions of Theorem 4.3, sequence Eε solution of (85) or (84) or ((63), (66), (67), (69)) converges to E(x)L2(Ω) which is the unique solution of the homogenized problem:

    {θ1x×x×E(x)+iωθ2E(x)=0  in Ω,θ1x×E(x)×e2=iωHd×e2  on Γd,x×E(x)×e2=0  on ΓL, (124)

    with θ1=ZIdyχ(y) dy and θ2=Z1C(y)(Idyχ(y)) dy.

    And where the scalar function χ is the unique solution, up to an additive constant in the Hilbert space of Z periodic functions H1#(Z), of the following boundary value problem

    {y(χ(y))=0  in ZΩC,[χn]=nj  on ΩC,[χ]=0  on ΩC, (125)

    where [f] is the jump across the surface of ΩC, nj, j={1,2,3} is the projection on the axis ej of the normal of ΩC.

    Proof. Step 1: Two-scale convergence. Due to the estimate (96), Eε is bounded in L2(Ω). Hence, up to a subsequence, Eε two-scale converges to E0(x,y) belonging to L2(Ω,L2#(Z)). That means for any V(x,y)C10(Ω,C1#(Z)), we have:

    limε0ΩEε(x)V(x,xε) dx=ΩZE0(x,y)V(x,y) dydx. (126)

    Step 2: Deduction of the constraint equation. We multiply the equation (123) by oscillating test function Vε(x)=V(x,xε) where V(x,y)C10(Ω,C1#(Z)):

    Ω×Eε(x)(x×Vε(x,xε)+1εy×Vε(x,xε))+[ω2ε5k(ϵ)    +iω((1εC(xε)+ε41εR(xε))1{y<0}+ε1{y>0})]EεVε(x,xε) dx    =iωΓdHd×e2(e2×V(x,1,z,ξ,1ε,ζ))×e2 dσ. (127)

    Integrating by parts, we get:

    ΩEε(x)(x×x×Vε(x,xε)+1εy×x×Vε(x,xε)+1εx×y×Vε(x,xε)+1ε2y×y×Vε(x,xε))+[ω2ε5k(ϵ)3+iω(1εC(xε)+ε41εR(xε))1{y<0}+ε1{y>0}]Eε(x)Vε(x,xε) dx=iωΓdHd×e2(e2×V(x,1,z,ξ,1ε,ζ))×e2 dσ. (128)

    Now we multiply (128) by ε2 and we pass to the two-scale limit, applying Theorem 4.3 we obtain:

    ΩZE0(x,y)(y×y×V(x,y)) dydx=0. (129)

    We deduce the constraint equation for the profile E0:

    y×y×E0(x,y)=0. (130)

    Step 3. Looking for the solutions to the constraint equation. Multiplying Equation (130) by E0 and integrating by parts over Z leads to

    Zy×y×E0(x,y)E0(x,y)  dy=Z|y×E0(x,y)|2  dy=0. (131)

    We deduce that equation (131) is equivalent to

    y×E0(x,y)=0. (132)

    Moreover a solution of (132) is also solution of (130). So (130) and (132) are equivalent.

    Hence, from Proposition (117) we conclude that E0(x,y) can be decomposed as

    E0(x,y)=E(x)+yΦ0(x,y). (133)

    Step 4. Equations for E(x) and Φ0(x,y). The divergence equation of (123) is multiplied with V(x,xε)=εv(x)ψ(xε), where vC10(Ω) and ψH1#(Z). Theorem 4.3 and integration by parts yields for all ψH1#(Z) and vC10(Ω)

    limε0Ω{ω2ε5k(ϵ)Eε(x)+iω[(1εC(xε)+ε41εR(xε))1{y<0}+ε1{y>0}]Eε(x)}εv(x)ψ(xε) dx=limε0Ω{ω2ε5k(ϵ)Eε(x)+iω[1εC(xε)+ε41εR(xε))1{y<0}+ε1{y>0}]Eε}(εv(x)ψ(xε)+v(x)yψ(xε)) dx=ΩZv(x)yψ(y)[iω1C(y)E0(x,y)] dydx=0, (134)

    from which it follows that

    y[iω1C(y)E0(x,y)]=0, (135)

    with E0 given by the decomposition (117). So we obtain the local equation

    y[iω1C(y){E(x)+yΦ0(x,y)}] dy=0. (136)

    The potential Φ0 may be written on the form

    Φ0(x,y)=3j=1χj(y)ejE(x)=χ(y)E(x). (137)

    From (133) and (137), we get:

    E0(x,y)=(Id+yχ(y))E(x). (138)

    Inserting E0 in (136) we obtain

    y[iω1C(y)(Id+yχ(y)]=0. (139)

    Now, we build oscillating test functions satisfying constraint (133) and use them in weak formulation (128). We define test function V(x,y)=α(x)+yβ(x,y), V(x,y)C10(Ω,C1#(Z)) and we inject in (128) test function Vε=V(x,xε), which gives:

    ΩEε(x)(x×x×V(x,xε)+2εx×y×V(x,xε)+1ε2y×y×V(x,xε))+[ω2ε5k(ϵ)+iω((1εC(xε)+ε41εR(xε))1{y<0}+ε1{y>0})]Eε(x)V(x,xε) dx=iωΓdHd×e2(e2×V(x,1,z,ξ,ζ))×e2 dσ, (140)

    with V(x,1,z,ξ,ν,ζ)=V(x,1,z,ξ,ζ) the restriction on V which does not depend on ν. The term containing the constraint, the third one, disappears. Passing to the limit ε0 and replacing the expression of V by the term α(x)+yβ(x,y), we have

    x×y×V(x,y)=x×y×[α(x)+yβ(x,y)]                        =x×y×(α(x))+x×y×(yβ(x,y))                        =x×y×(yβ(x,y)). (141)

    Since y×(y)=0, the term 2εx×y×yβ(x,y)) vanishes. Therefore, (140) becomes:

    ΩZE0(x,y)x×x×(α(x)+yβ(x,y))    +iω1C(y)E0(x,y)(α(x)+yβ(x,y) dydx    =iωΓdHd×e2(e2×(α(x,1,z)+yβ(x,1,z,ξ,ζ)))×e2 dσ. (142)

    Now in (142) we replace expression E0 giving by (138). We obtain

    ΩZ(Id+yχ(y))E(x)(x×x×(α(x)+yβ(x,y))    +iω1C(y)(Id+yχ(y))E(x))(α(x)+yβ(x,y)) dydx    =iωΓdHd×e2(e2×(α(x,1,z)+yβ(x,1,z,ξ,ζ)))×e2 dσ. (143)

    Taking α(x)=0 in (143), we obtain

    ΩZ(Id+yχ(y))x×x×E(x)yβ(x,y)               +iω1C(y)(Id+yχ(y))E(x)yβ(x,y)dydx=0. (144)

    Integrating by parts

    ΩZy{(Idyχ(y))x×x×E(x)}β(x,y)                     iωy{1C(y)(Idyχ(y))E(x)}β(x,y) dydx=0. (145)

    And since y{1C(y)(Id+yχ(y))E(x)}=0 we obtain

    ΩZy{(Id+yχ(y))x×x×E(x)}β(x,y) dydx=0, (146)

    which gives the cell problem

    y[Id+yχ(y)]=0. (147)

    From (139) and (147), the scalar function χ is the unique solution, thanks to Lax-Milgram Lemma, up to an additive constant in the Hilbert space of Z periodic function H1#(Z) of the following boundary value problem

    {y(χ(y))=0  in ZΩC,[χn]=nj  on ΩC,[χ]=0  on ΩC, (148)

    where [f] is the jump across the surface of ΩC and nj, j={1,2,3} is the projection on the axis ej of the normal of ΩC.

    Remark 4.11. (148) can be seen as an electrostatic problem. Solving (139) and (147) reduces to look for a potential induced by surface density of charges. Then χ is this potential induced by the charges on the interface of carbon fiber.

    Setting β(x,y)=0 in (143) and integrating by parts, we get

    ΩZ(Id+yχ(y))x×x×E(x)α(x)                           +iω1C(y)(Id+yχ(y))E(x)α(x) dydx                           =iωΓdHd×e2(e2×α(x,1,z))×e2 dσ, (149)

    which gives the following well posed problem for E(x)

    {θ1x×x×E(x)+iωθ2E(x)=0  in Ω,θ1x×E(x)×e2=iωHd×e2  on Γd,x×E(x)×e2=0  on ΓL, (150)

    with θ1=ZId+yχ(y) dy and θ2=Z1C(y)(Id+yχ(y)) dy.

    This concludes the proof of Theorem (124).


    5. Conclusion

    We presented in this paper the homogenization of time harmonic Maxwell equation by the method of two-scale convergence. We started by studying the time harmonic Maxwell equations with coefficients depending of ε. We remind that λ is the wave length, δ is the skin length, L is thickness of the medium and e the size of the basic cell and then ε=eL is the small parameter. We find for low frequencies the macroscopic homogenized Maxwell equations depending on the volume fraction of the carbon fibers and we find also the microscopic equation.


    6. Annexes

    A. Presentation of all cases of tables 1, 2, 3 and 4

    Table 1. For ¯ω=100 rad.s1.
    case 1 case 2 case 3 case 4 case 5 case 6 case 7 case 8
    L(m) 10-3 10-3 10-3 10-3 10-3 10-3 10-3 10-3
    e(m) 10-5 10-5 10-5 10-5 10-5 10-5 10-5 10-5
    λ(m) 106 106 106 106 106 106 106 106
    σ(S·m—1) 40000 40000 10-10 10-3 40000 104 10-10 10-3
    σ(m) 0, 1 0, 1 107 103 0, 1 0, 1 107 103
    σc(S·m-1) σ σ ˉσε7 ˉσε4 σ σ ˉσε7 ˉσε4
    σr(S·m-1) ε7ˉσ ε4ˉσ σ σ ε7ˉσ ε4ˉσ σ σ
    σa(S·m-1) ε9ˉσ ε9ˉσ ε2ˉσ ε6ˉσ εˉσ εˉσ ˉσε5 ˉσε2
    4πˉL2¯λ2 ε9 ε9 ε9 ε9 ε9 ε9 ε9 ε9
    ˉL2¯δ2 ε2 ε2 ε10 ε7 ε2 ε2 ε10 ε7
     | Show Table
    DownLoad: CSV
    Table 2. For ¯ω=106 rad.s1.
    case 1 case 2 case 3 case 4 case 5 case 6 case 1 case 8
    L(m) 10-3 10-3 10-3 10-3 10-3 10-3 10-3 10-3
    e(m) 10-5 10-5 10-5 10-5 10-5 10-5 10-5 10-5
    λ(m) 103 103 103 103 103 103 103 103
    σ(S.m-1) 40000 40000 10-10 10-3 40000 40000 10-10 10-3
    δ(m) 10-3 10-3 105 10 10-3 10-3 105 10
    σc(S·m-1) σ σ ˉσε7 ˉσε5 σ σ ˉσε7 ˉσε5
    σr(S·m-1) ε7ˉσ ε4ˉσ σ σ ε7ˉσ ε4ˉσ σ σ
    σa(S·m-1) ε9ˉσ ε9ˉσ ε2ˉσ ε6ˉσ εˉσ εˉσ ˉσε5 ˉσε2
    4πˉL2ˉλ2 ε5 ε5 ε5 ε5 ε5 ε5 ε5 ε5
    ˉL2ˉδ2 1 1 ε8 ε5 1 1 ε8 ε5
     | Show Table
    DownLoad: CSV
    Table 3. For ¯ω=1010 rad.s1.
    case 1 case 2 case 3 case 4 case 5 case 6 case 1 case 8
    L(m) 10-3 10-3 10-3 10-3 10-3 10-3 10-3 10-3
    e(m) 10-5 10-5 10-5 10-5 10-5 10-5 10-5 10-5
    λ(m) 10-1 10-1 10-1 10-1 10-1 10-1 10-1 10-1
    σ(S.m-1) 40000 40000 10-10 10-3 40000 40000 10-10 10-3
    δ(m) 10-5 10-5 103 10-1/2 10-5 10-5 103 10-1/2
    σc(S·m-1) σ σ ˉσε7 ˉσε4 σ σ ˉσε7 ˉσε4
    σr(S·m-1) ε7ˉσ ε4ˉσ σ σ ε7ˉσ ε4ˉσ σ σ
    σa(S·m-1) ε9ˉσ ε9ˉσ ε2ˉσ ε6ˉσ εˉσ εˉσ ˉσε5 ˉσε2
    4πˉL2ˉλ2 ε ε ε ε ε ε ε ε
    ˉL2ˉδ2 1ε2 1ε2 ε6 ε3 1ε2 1ε2 ε6 ε3
     | Show Table
    DownLoad: CSV
    Table 4. For ¯ω=1012 rad.s1.
    case 1 case 2 case 3 case 4 case 5 case 6 case 1 case 8
    L(m) 10-3 10-3 10-3 10-3 10-3 10-3 10-3 10-3
    e(m) 10-5 10-5 10-5 10-5 10-5 10-5 10-5 10-5
    λ(m) 10-3 10-3 10-3 10-3 10-3 10-3 10-3 10-3
    σ(S.m-1) 40000 40000 10-10 10-3 40000 40000 10-10 10-3
    δ(m) 10-6 10-6 1 10-3/2 10-6 10-6 1 10-3/2
    σc(S·m-1) σ σ ˉσε7 ˉσε4 σ σ ˉσε7 ˉσε4
    σr(S·m-1) ε7ˉσ ε4ˉσ σ σ ε7ˉσ ε4ˉσ σ σ
    σa(S·m-1) ε9ˉσ ε9ˉσ ε2ˉσ ε6ˉσ εˉσ εˉσ ˉσε5 ˉσε2
    4πˉL2ˉλ2 1 1 1 1 1 1 1 1
    ˉL2ˉδ2 1ε3 1ε3 ε3 ε 1ε3 1ε3 ϵ3 ε3
     | Show Table
    DownLoad: CSV

    -The case 1 corresponds to the air not ionized, a resin not doped and ¯σ is the effective electric conductivity in the direction of the carbon fibers. We have for the effective electric conductivity ¯σ=σc40000 S.m1, the resin conductivity is about σr1010 S.m1 and the conductivity in the air is about 1014 S.m1. So when we want to calculate the ratio in (??)-(41) depending on ε we get: σr¯σε7 and σa¯σε9.

    -In case 2, the air is not ionized, the resin is doped and ¯σ is the effective conductivity is in direction of carbon fibers. We have like the case 1 ¯σ=σc40000 S.m1. The resin conductivity is about σr103 S.m1 and the conductivity in the air is about 1014 S.m1. So σr¯σε4 and σa¯σε9.

    -In case 3, the air is not ionized, the resin is not doped and ¯σ is the effective conductivity is orthogonal to the fibers. ¯σ=σr1010 S.m1. The carbon fiber conductivity is about σc104 S.m1 and the conductivity in the air is about 1014 S.m1. σc¯σ1ε7 and σa¯σε2.

    -Case 4 corresponds to the air non ionized, the resin doped and ¯σ is the effective conductivity orthogonal to the fibers. The effective electric conductivity is ¯σ=σr103 S.m1. The carbon fiber conductivity is about σc40000 S.m1 and the conductivity in the air is about 1014 S.m1. σc¯σ1ε4 and σa¯σε6.

    -In case 5, the air is ionized, the resin is not doped and ¯σ is the effective conductivity is in the direction of the carbon fibers. This one is equal ¯σ=σc40000 S.m1, the resin conductivity is about σr1010 S.m1 and the conductivity in the air is now about 4242 S.m1. σr¯σε7 and σa¯σε.

    -Case 6 corresponds to the air ionized, the resin doped and ¯σ is the effective conductivity in direction of the carbon fibers. This one is equal ¯σ=σc40000 S.m1, the resin conductivity is about σr103 S.m1 and the conductivity in the air is now about 4242 S.m1. σr¯σε4 and σa¯σε.

    -Case 7 corresponds to the air ionized, the resin not doped and ¯σ is the effective conductivity orthogonal to the fibers. The effective conductivity is ¯σ=σr1010 S.m1, the carbon fibers conductivity is about σc40000 S.m1 and the conductivity in the air is now about 4242 S.m1. σc¯σ1ε7 and σa¯σ1ε6.

    -Case 8 corresponds to the air ionized, the resin doped and ¯σ is the effective conductivity orthogonal to the fibers. The effective conductivity is ¯σ=σr103 S.m1, the carbon fibers conductivity is about σc40000 S.m1 and the conductivity in the air is now about 4242 S.m1. σc¯σ1ε4 and σa¯σ1ε2.

    B. Structure of the equations depending of ε

    For ¯ω=100 rad.s1, we have

    Case 1

    η=9 and Σεa=ε11, Σεr=ε9,Σεc=ε2. (151)

    Case 2

    η=9 and Σεa=ε11,Σεr=ε6,Σεc=ε2. (152)

    Case 3

    η=9 and Σεa=ε12,Σεr=ε10,Σεc=ε3. (153)

    Case 4

    η=9 and Σεa=ε13,Σεr=ε7,Σεc=ε3. (154)

    Case 5

    η=9 and Σεa=ε3,Σεr=ε9,Σεc=ε2. (155)

    Case 6

    \eta =9\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{3}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{6}},\Sigma _{c}^{\varepsilon }={{\varepsilon }^{2}}. (156)

    Case 7

    \eta =9\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{5}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{10}},\Sigma _{c}^{\varepsilon }={{\varepsilon }^{3}}. (157)

    Case 8

    \eta =9\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{5}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{7}},\Sigma _{c}^{\varepsilon }={{\varepsilon }^{3}}. (158)

    For \overline{\omega} = 10^{6}~rad.s^{-1}

    Case 1

    \eta =5\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{9}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{7}},\Sigma _{c}^{\varepsilon }=1. (159)

    Case 2

    \eta =5\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{9}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{4}},\Sigma _{c}^{\varepsilon }=1. (160)

    Case 3

    \eta =5\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{10}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{8}},\Sigma _{c}^{\varepsilon }=\varepsilon . (161)

    Case 4

    \eta =5\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{11}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{5}},\Sigma _{c}^{\varepsilon }=1. (162)

    Case 5

    \eta =5\ \text{and}\ \Sigma _{a}^{\varepsilon }=\varepsilon ,\Sigma _{r}^{\varepsilon }={{\varepsilon }^{7}},\Sigma _{c}^{\varepsilon }=1. (163)

    Case 6

    \eta =5\ \text{and}\ \Sigma _{a}^{\varepsilon }=\varepsilon ,\Sigma _{r}^{\varepsilon }={{\varepsilon }^{4}},\Sigma _{c}^{\varepsilon }=1. (164)

    Case 7

    \eta =5\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{3}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{8}},\Sigma _{c}^{\varepsilon }=\varepsilon . (165)

    Case 8

    \eta =5\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{3}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{5}},\Sigma _{c}^{\varepsilon }=1. (166)

    For \overline{\omega} = 10^{10}~rad.s^{-1}

    Case 1

    \eta =1\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{7}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{5}},\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{2}}}. (167)

    Case 2

    \eta =1\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{7}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{2}},\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{2}}}. (168)

    Case 3

    \eta =1\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{8}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{6}},\Sigma _{c}^{\varepsilon }=\frac{1}{\varepsilon }. (169)

    Case 4

    \eta =1\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{9}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{3}},\Sigma _{c}^{\varepsilon }=\frac{1}{\varepsilon }. (170)

    Case 5

    \eta =1\ \text{and}\ \Sigma _{a}^{\varepsilon }=\frac{1}{\varepsilon },\Sigma _{r}^{\varepsilon }={{\varepsilon }^{5}},\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{2}}}. (171)

    Case 6

    \eta =1\ \text{and}\ \Sigma _{a}^{\varepsilon }=\frac{1}{\varepsilon },\Sigma _{r}^{\varepsilon }={{\varepsilon }^{2}},\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{2}}}. (172)

    Case 7

    \eta =1\ \text{and}\ \Sigma _{a}^{\varepsilon }=\frac{1}{\varepsilon },\Sigma _{r}^{\varepsilon }={{\varepsilon }^{6}},\Sigma _{c}^{\varepsilon }=\frac{1}{\varepsilon }. (173)

    Case 8

    \eta =1\ \text{and}\ \Sigma _{a}^{\varepsilon }=\varepsilon ,\Sigma _{r}^{\varepsilon }={{\varepsilon }^{3}},\Sigma _{c}^{\varepsilon }=\frac{1}{\varepsilon }. (174)

    For \overline{\omega} = 10^{12}~rad.s^{-1}

    Case 1

    \eta =0\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{6}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{4}},\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{3}}}. (175)

    Case 2

    \eta =0\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{6}},\Sigma _{r}^{\varepsilon }=\varepsilon ,\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{3}}}. (176)

    Case 3

    \eta =0\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{5}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{3}},\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{4}}}. (177)

    Case 4

    \eta =0\ \text{and}\ \Sigma _{a}^{\varepsilon }={{\varepsilon }^{7}},\Sigma _{r}^{\varepsilon }=\varepsilon ,\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{3}}}. (178)

    Case 5

    \eta =0\ \text{and}\ \Sigma _{a}^{\varepsilon }=\frac{1}{{{\varepsilon }^{2}}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{4}},\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{3}}}. (179)

    Case 6

    \eta =0\ \text{and}\ \Sigma _{a}^{\varepsilon }=\frac{1}{{{\varepsilon }^{2}}},\Sigma _{r}^{\varepsilon }=\varepsilon ,\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{3}}}. (180)

    Case 7

    \eta =0\ \text{and}\ \Sigma _{a}^{\varepsilon }=\frac{1}{{{\varepsilon }^{2}}},\Sigma _{r}^{\varepsilon }={{\varepsilon }^{3}},\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{4}}}. (181)

    Case 8

    \eta =0\ \text{and}\ \Sigma _{a}^{\varepsilon }=\frac{1}{\varepsilon },\Sigma _{r}^{\varepsilon }=\varepsilon ,\Sigma _{c}^{\varepsilon }=\frac{1}{{{\varepsilon }^{3}}}. (182)

    C. The transmission Maxwell problem

    Taking a test function V \in C^{1}(\Omega) with compact support in \Omega_{c}^{\varepsilon}, in weak formulation (83) associated with the problem ((63), (66), (67)). Since

    \int_{\Omega }{\nabla }\times E_{|\Omega _{c}^{\varepsilon }}^{\varepsilon }\cdot \nabla \times \overline{V}~~d\mathbf{x}={{\langle \nabla \times \nabla \times E_{|\Omega _{c}^{\varepsilon }}^{\varepsilon },\overline{V}\rangle }_{\Omega _{c}^{\varepsilon }}}, (183)

    we deduce the third equation in (84). Similarly, taking V \in C^{1}(\Omega) with compact support respectively in \Omega_{r}^{\varepsilon} and \Omega_{a}^{\varepsilon}, we obtain the first and the second equation in (84). Now, since E_{|\Omega_{a}^{\varepsilon}} \in \mathbf{H}(\mbox{curl}, \Omega_{a}^{\varepsilon}) and E_{|\Omega_{r}^{\varepsilon}} \in \mathbf{H}(\mbox{curl}, \Omega_{r}^{\varepsilon}), let V \in C_{0}^{1}(\Omega_{a}^{\varepsilon} \bigcup \Omega_{r}^{\varepsilon}) integrating by parts we get

    \begin{align} & \int_{\Omega _{a}^{\varepsilon }\bigcup{\Omega _{r}^{\varepsilon }}}{E}\cdot \nabla \times \overline{V}~~d\mathbf{x}=\int_{\Omega _{a}^{\varepsilon }}{{{E}_{|\Omega _{a}^{\varepsilon }}}}\cdot \nabla \times \overline{V}~~d\mathbf{x}+\int_{\Omega _{r}^{\varepsilon }}{{{E}_{|\Omega _{r}^{\varepsilon }}}}\cdot \nabla \times \overline{V}~~d\mathbf{x} \\ & \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ =\int_{\Omega _{a}^{\varepsilon }}{\nabla }\times {{E}_{|\Omega _{a}^{\varepsilon }}}\cdot \overline{V}~~d\mathbf{x}+\int_{\Omega _{r}^{\varepsilon }}{\nabla }\times {{E}_{|\Omega _{r}^{\varepsilon }}}\cdot \overline{V}~~d\mathbf{x} \\ & \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ +\int_{{{\Gamma }_{ra}}}{({{E}_{|\Omega _{a}^{\varepsilon }}}\times {{e}_{2}}-{{E}_{|\Omega _{r}^{\varepsilon }}}\times {{n}_{|\Omega _{r}^{\varepsilon }}})}\cdot \overline{V}~ds. \\ \end{align} (184)

    Since on every point of \Gamma_{ra} e_{2} =-n_{|\Omega_{r}^{\varepsilon}} the assumed continuity require

    {{E}_{|\Omega _{a}^{\varepsilon }}}\times {{e}_{2}}={{E}_{|\Omega _{r}^{\varepsilon }}}\times {{n}_{|\Omega _{r}^{\varepsilon }}}, (185)

    we obtain the fourth relation in (84). With the same argument on \Gamma_{cr}^{\varepsilon}, we obtain the last relation in (84). This shows that (83) implies (84). And, if {{E}^{\varepsilon }} is solution to (84) following that for any regular set \widehat{\Omega} in \Omega the Stokes's formula gives, for more details see p 57, 58 of P. Monk's book [18]:

    \forall ~~E,V\in \mathbf{H}(\text{curl},\hat{\Omega })~~\int_{{\hat{\Omega }}}{\nabla }\times E\cdot \overline{V}-E\cdot \nabla \times \overline{V}~~d\mathbf{x}={{\langle E\times {{n}_{{\hat{\Omega }}}},\overline{{{V}_{T}}}\rangle }_{\partial \hat{\Omega }}} (186)

    \mathbf{H}(\mbox{curl}, \widehat{\Omega}) has the same definition as \mathbf{H}(\mbox{curl}, \Omega) with \Omega replaced by \widehat{\Omega} and where V_{T} = (n \times V) \times n, and n_{\widehat{\Omega}} is the unit outward normal of \partial \widehat{\Omega}. For all V \in \mathbf{H}(\mbox{curl}, \Omega), V_{|\Omega_{r}^{\varepsilon}} \in \mathbf{H}(\mbox{curl}, \Omega_{r}^{\varepsilon}), V_{|\Omega_{a}^{\varepsilon}} \in \mathbf{H}(\mbox{curl}, \Omega_{a}^{\varepsilon}) and V_{|\Omega_{c}^{\varepsilon}} \in \mathbf{H}(\mbox{curl}, \Omega_{c}^{\varepsilon}). Hence, fixing any E' \in \mathbf{H}(\mbox{curl}, \Omega) according to the second equation in (84), we have {\nabla }E_{|\Omega_{r}^{\varepsilon}}^{\varepsilon} \in \mathbf{H}(\mbox{curl}, \Omega_{r}^{\varepsilon}) then applying (186) in \Omega_{r}^{\varepsilon} with E = {\nabla } E_{|\Omega_{r}^{\varepsilon}}^{\varepsilon} and V we get

    \begin{align} & \int_{\Omega _{r}^{\varepsilon }}{\nabla }\times E_{|\Omega _{r}^{\varepsilon }}^{\varepsilon }\cdot \nabla \times \overline{V}~~d\mathbf{x}=\int_{\Omega _{r}^{\varepsilon }}{\nabla }\times \nabla \times E_{|\Omega _{r}^{\varepsilon }}^{\varepsilon }\cdot \overline{V}~d\mathbf{x}+{{\langle \nabla \times E_{|\Omega _{r}^{\varepsilon }}^{\varepsilon }\times {{n}_{|\Omega _{r}^{\varepsilon }}},\overline{{{V}_{T}}}\rangle }_{{{\Gamma }_{ra}}}} \\ & \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ +{{\langle \nabla \times E_{|\Omega _{c}^{\varepsilon }}^{\varepsilon }\times {{n}_{|\Omega _{c}^{\varepsilon }}},\overline{{{V}_{T}}}\rangle }_{\Gamma _{cr}^{\varepsilon }}}. \\ \end{align} (187)

    Doing the same for \Omega_{c}^{\varepsilon}, we have

    \int_{\Omega _{c}^{\varepsilon }}{\nabla }\times E_{|\Omega _{c}^{\varepsilon }}^{\varepsilon }\cdot \nabla \times \overline{V}~~d\mathbf{x}=\int_{\Omega _{c}^{\varepsilon }}{\nabla }\times \nabla \times E_{|\Omega _{c}^{\varepsilon }}^{\varepsilon }\cdot \overline{V}~d\mathbf{x}+{{\langle \nabla \times E_{|\Omega _{c}^{\varepsilon }}^{\varepsilon }\times n_{|\Omega _{c}^{\varepsilon }}^{\varepsilon },\overline{{{V}_{T}}}\rangle }_{\Gamma _{cr}^{\varepsilon }}}. (188)

    Finally for \Omega_{a}^{\varepsilon}, we have

    \begin{align} & \int_{\Omega _{a}^{\varepsilon }}{\nabla }\times E_{|\Omega _{a}^{\varepsilon }}^{\varepsilon }\cdot \nabla \times \overline{V}~~d\mathbf{x}=\int_{\Omega _{a}^{\varepsilon }}{\nabla }\times \nabla \times E_{|\Omega _{a}^{\varepsilon }}^{\varepsilon }\cdot \overline{V}~d\mathbf{x}+{{\langle \nabla \times E_{|\Omega _{a}^{\varepsilon }}^{\varepsilon }\times {{e}_{2}},\overline{{{V}_{T}}}\rangle }_{{{\Gamma }_{d}}}} \\ & \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ -{{\langle \nabla \times E_{|\Omega _{a}^{\varepsilon }}^{\varepsilon }\times {{e}_{2}},\overline{{{V}_{T}}}\rangle }_{\Gamma _{ra}^{\varepsilon }}}. \\ \end{align} (189)

    Summing the relations above since in every point of \Gamma_{ra} n_{| \Omega_{r}^{\varepsilon}} =-e_{2} and in every point of \Gamma_{cr}^{\varepsilon} n_{| \Omega_{c}^{\varepsilon}} =-n_{| \Omega_{r}^{\varepsilon}}, it comes

    \begin{align} & \int_{\Omega }{\nabla }\times {{E}^{\varepsilon }}\cdot \nabla \times \overline{V}~~d\mathbf{x}=\int_{\Omega }{\nabla }\times \nabla \times {{E}^{\varepsilon }}\cdot \overline{V}~~d\mathbf{x}+<[\nabla \times {{E}^{\varepsilon }}\times n],\overline{{{V}_{T}}}{{>}_{{{\Gamma }_{ra}}}} \\ & \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ +{{\langle [\nabla \times {{E}^{\varepsilon }}\times n],\overline{{{V}_{T}}}\rangle }_{{{\Gamma }_{cr}}}}-i\omega \int_{{{\Gamma }_{d}}}{{{H}_{d}}}\times {{n}^{\varepsilon }}\cdot {{{\overline{V}}}_{T}}~d\sigma . \\ \end{align} (190)

    According to (83) and the first, second and third equations in (84) we have

    \begin{align} & {{\langle \nabla \times E_{|\Omega _{c}^{\varepsilon }}^{\varepsilon }\times {{n}_{|\Omega _{c}^{\varepsilon }}},\overline{{{V}_{T}}}\rangle }_{\Gamma _{cr}^{\varepsilon }}}-{{\langle \nabla \times E_{|\Omega _{r}^{\varepsilon }}^{\varepsilon }\times {{n}_{|\Omega _{r}^{\varepsilon }}},\overline{{{V}_{T}}}\rangle }_{\Gamma _{cr}^{\varepsilon }}} \\ & +{{\langle \nabla \times E_{|\Omega _{r}^{\varepsilon }}^{\varepsilon }\times {{n}_{|\Omega _{r}^{\varepsilon }}},\overline{{{V}_{T}}}\rangle }_{{{\Gamma }_{ra}}}}+{{\langle \nabla \times E_{|\Omega _{a}^{\varepsilon }}^{\varepsilon }\times {{e}_{2}},\overline{{{V}_{T}}}\rangle }_{{{\Gamma }_{ra}}}}=0, \\ \end{align} (191)

    for all V \in \mathbf{H}(\mbox{curl}, \Omega) which causes the last two equalities in (84) and concludes the first part of the proof.

    Reciprocally, integrating by parts (84) we have:

    \begin{align} & \forall ~~V\in {{\mathbf{X}}^{\varepsilon }}(\Omega ),~~\int_{\Omega }{\nabla }\times {{E}^{\varepsilon }}\cdot \nabla \times \overline{V}~d\mathbf{x}+\int_{\Omega _{a}^{\varepsilon }}{(-{{\omega }^{2}}{{\varepsilon }^{\eta }}+i\omega \Sigma _{a}^{\varepsilon })}{{E}^{\varepsilon }}\cdot \overline{V}~d\mathbf{x} \\ & \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ =-i\omega \int_{{{\Gamma }_{d}}}{{{H}_{d}}}\times {{n}^{\varepsilon }}\cdot {{{\overline{V}}}_{T}}~d\sigma , \\ \end{align} (192)

    and

    \forall ~~V\in {{\mathbf{X}}^{\varepsilon }}(\Omega ),~~\int_{\Omega }{\nabla }\times {{E}^{\varepsilon }}\cdot \nabla \times \overline{V}~d\mathbf{x}+\int_{\Omega _{r}^{\varepsilon }}{(-{{\omega }^{2}}{{\varepsilon }^{\eta }}{{\epsilon }_{r}}+i~\omega \Sigma _{r}^{\varepsilon })}{{E}^{\varepsilon }}\cdot \overline{V}~d\mathbf{x}=0, (193)

    and

    \forall ~~V\in {{\mathbf{X}}^{\varepsilon }}(\Omega ),~~\int_{\Omega }{\nabla }\times {{E}^{\varepsilon }}\cdot \nabla \times \overline{V}~d\mathbf{x}+\int_{\Omega _{c}^{\varepsilon }}{(-{{\omega }^{2}}{{\varepsilon }^{\eta }}{{\epsilon }_{c}}+i~\omega \Sigma _{c}^{\varepsilon })}{{E}^{\varepsilon }}\cdot \overline{V}~d\mathbf{x}=0. (194)

    By adding these three integrals, we get the variational formulation (83) associated with the problem ((63), (66), (67)).

    Taking the divergence of the first three equations of (84) we get (69).


    Conflicts of Interest

    All authors declare no conflicts of interest in this paper.


    [1] T. Abboud and I. Terrasse, Modélisation des phénomènes de propagation d'ondes, Centre Poly-Média de l'école Polytechnique, 2007.
    [2] Y. Amirat, K. Hamdache and A. Ziani, Homogénéisation d'équations hyperboliques du premier ordre et application aux écoulements missibles en milieux poreux, Ann. Inst. H. Poincaré, 6 (1989), 397-417.
    [3] G. Allaire, Homogenization and Two-scale Convergence, SIAM Journal on Mathematical Analysis, 23 (1992), 1482-1518.
    [4] G. Allaire and M. Briand, Multiscale convergence and reiterated homogenization, Roy.Soc.Edinburgh, 126 (1996), 297-342.
    [5] Y. Amirat and V. Shelukhin, Homogenization of time-harmonic Maxwell equations and the frequency dispersion effect, J.Maths.Pures.Appl., 95 (2011), 420-443.
    [6] A. Back and E. Frenod, Geometric Two-Scale Convergence on Manifold and Applications to the Vlasov Equation Discrete and Continuous Dynamical Systems-Serie S. Special Issue on Numerical Methods based on Homogenization and Two-Scale Convergence, 8 (2015), 223-241.
    [7] S. Berthier, Optique des milieux composites, Ed. Polytechnicia, 1993.
    [8] D. Cionarescu and P. Donato, An introduction to homogenization, Oxford University Press. , 1999.
    [9] M. Costabel, M. Dauge and S. Nicaise, Corner Singularities of Maxwell interface and Eddy current problems, Advances and Applications, 147 (2004), 241-256.
    [10] N. Crouseilles, E. Frenod, S. Hirstoaga and A. Mouton, Two-Scale Macro-Micro decomposition of the Vlasov equation with a strong magnetic field, Mathematical Models and Methods in Applied Sciences, 23 (2012), 1527-1559.
    [11] E. Frénod, P. A. Raviart and E. Sonnendrücker, Asymptotic Expansion of the Vlasov Equation in a Large External Magnetic Field, J. Math. Pures et Appl. 80, (2001), 815-843.
    [12] S. Guenneau, F. Zolla and A. Nicolet, Homogenization of 3D finite photonic crystals with heterogeneous permittivity and permeability, Waves in Random and Complex Media, 17 (2007), 653-697.
    [13] P.R.P. Hoole and S.R.H. Hoole, Guided waves along an unmagnetized lightning plasma channel, IEEE Transactions on Magnetics, 24 (1998), 3165-3167.
    [14] P.R.P. Hoole, S.R.H. Hoole, S. Thirukumaran, R. Harikrishnan, K. Jeevan and K. Pirapaharan, Aircraft-lightning electrodynamics using the transmission line model part Ⅰ: review of the transmission line model, Progress In Electromagnetics Research M, 31, (2013), 85-101.
    [15] P. Laroche, P. Blanchet, A. Delannoy, and F. Issac, Experimental Studies of Lightning Strikes to Aircraft, JOURNAL AEROSPACELAB, 112 (2012).
    [16] M. Leboulch, Analyse spectrale VHF, UHF du rayonnement deséclairs, Hamelin, CENT.
    [17] J. C. Maxwell, A dynamical theory of the Electromagnetic Field, Phisophical transacting of the Royal Society of London, (1885), 459-512.
    [18] P. Monk, Finite Element Methods for Maxwell's Equations, Oxford Science publication, Numerical Mathematics and scientific computation, Clarendon Press-Oxford, 2003.
    [19] J. C. Nédélec, Acoustic and electromagnetic equations; integral representations for harmonic problems, Springer-Verlag, Berlin, 2001.
    [20] M. Neuss-Radu, Some extensions of two-scale convergence, Comptes rendus de l'Academie des sciences, 322 (1996), 899-904.
    [21] G. Nguetseng. A General Convergence Result for a Functional Related to the Theory of Homogenization, 20 (1989), 608-623.
    [22] G. Nguetseng, Asymptotic Analysis for a Stiff Variational Problem Arising in Mechanics, SIAM Journal on Mathematical Analysis, 21 (1990), 1394-1414.
    [23] S. Nicaise, S. Hassani and A. Maghnouji, Limit behaviors of some boundary value problems with high and/or low valued parameters, Advances in differential equations, 14 (2009), 875-910.
    [24] O. Ouchetto, S. Zouhdi and A. Bossavit et al. , Effective constitutive parameters of periodic composites, Microwave conference, European, 2 (2005).
    [25] H.E. Pak, Geometric two-scale convergence on forms and its applications to Maxwell's equations, Proceedings of the Royal Society of Edinburgh, European, 135A (2005), 133-147.
    [26] N. Wellander, Homogenization of the Maxwell equations: Case Ⅰ. Linear theory, Appl Math, 46 (2001), 29-51.
    [27] N. Wellander, Homogenization of the Maxwell equations: Case Ⅱ. Nonlinear conductivity, Appl Math, 47 (2002), 255-283.
    [28] N. Wellander and B. Kristensson, Homogenization of the Maxwell equations at fixed frequency, Technical Report, (2002), 1-37.
    [29] Pr. Welter, Cours : Matériaux diélectriques, Master Matériaux, Institut Le Bel.
  • Reader Comments
  • © 2017 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(5601) PDF downloads(1123) Cited by(0)

Article outline

Figures and Tables

Figures(2)  /  Tables(4)

Other Articles By Authors

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog