R Torchio

D Voltolina

F Moro

P Alotto

P Bettini

G Meunier

J.-M Guichon

O Chadebec

Volume Integral Equation Methods for Axisymmetric Problems with Conductive and Magnetic Media

Keywords: Axisymmetric, magnetic media, low-rank approximation, hierarchical matrix, induction heating

Two new different integral equation methods for the study of axisymmetric harmonic electromagnetic problems are presented. Both the proposed formulations allow for the modeling of conductive and magnetic media and inhomogeneous media can be considered as well. Thanks to the use of low-rank approximation techniques a drastic reduction of the required memory and computational cost is achieved with a truly negligible loss of accuracy. The numerical features of the two proposed integral methods are widely discussed and compared. Moreover, most of the considerations hold for the case of 3-D integral methods. Sample codes of the two formulations are made available at https://github.com/UniPD-DII-ETCOMP/DenseMatrixMarket.

I. INTRODUCTION

When complex electromagnetic (EM) devices are embedded in a large air domain or when small air gaps are present integral equations (IE) formulations can be more convenient than differential ones [START_REF] Alotto | A Boundary Integral Formulation on Unstructured Dual Grids for Eddy-Current Analysis in Thin Shields[END_REF]. Several electromagnetic applications, such as induction hardening [START_REF] Spezzapria | Numerical Simulation of Solid-Solid Phase Transformations During Induction Hardening Process[END_REF], induction cooking [START_REF] Acero | Quantitative Evaluation of Induction Efficiency in Domestic Induction Heating Applications[END_REF], EM levitation for melting devices [START_REF] Sergeant | Magnetic Shielding of Levitation Melting Devices[END_REF], high temperature superconducting bearings [START_REF] Yan | Numerical Prediction of Levitation Properties of HTS Bulk in High Magnetic Fields[END_REF], [START_REF] Naseh | Levitation Force Maximization in HTS Magnetic Bearings Formulated by a Semianalytical Approach[END_REF], and thermonuclear fusion devices [START_REF] Freidberg | Plasma Physics and Fusion Energy[END_REF], allow for an axisymmetric solution of the EM problem. Therefore, when the axisymmetric assumption actually holds, this property can be exploited to reduce the computational cost required by the numerical simulations of 3-D EM devices. Furthermore, on the contrary of differential methods, when devices with moving source coils are studied (e.g. in induction heating hardening processes) integral formulations do not require the re-mesh of the domains, which generally requires a non-negligible time, and a small quantity of mutual coefficients must be re-computed instead.

In this paper, two different IE methods for the study of axisymmetric EM problems involving conductive and magnetic media are presented. The first relies on the amperian interpretation of the magnetization phenomena while the second relies on the coulombian interpretation [START_REF] Fano | Electromagnetic Fields, Energy, and Forces[END_REF]. These two methods are derived from the 3-D Partial Element Equivalent Circuit (PEEC) formulations proposed in [START_REF] Torchio | A 3-D PEEC Formulation Based on the Cell Method for Full-Wave Analyses With Conductive, Dielectric, and Magnetic Media[END_REF] and [START_REF] Torchio | An Extension of Unstructured-PEEC Method to Magnetic Media[END_REF], respectively, and compared in [START_REF] Torchio | A Volume PEEC Formulation Based on the Cell Method for Electromagnetic Problems from Low to High Frequency[END_REF]. Such as 3-D formulations, the two proposed methods show different numerical features which are compared and discussed. To the best of the author's knowledge, no axisymmetric integral formulation has been proposed in literature for the study of both conductive and magnetic media.

The case of dielectric media is not considered here. Indeed, when the axisymmetric assumption holds, the polarization currents are forced to be divergence free. Thus, the main effect of dielectric media (i.e. the enhancement of capacitive effects due to the bound electric charges) is not present. However, Manuscript received December 1, 2012; revised August 26, 2015. Corresponding author: M. Shell (email: http://www.michaelshell.org/contact.html).

following [START_REF] Torchio | A 3-D PEEC Formulation Based on the Cell Method for Full-Wave Analyses With Conductive, Dielectric, and Magnetic Media[END_REF] and [START_REF] Torchio | An Extension of Unstructured-PEEC Method to Magnetic Media[END_REF], both the proposed methods can be easily extended to the case of dielectric media.

Like all the IE methods, the main drawback of the two proposed formulations is the generation of full matrices. However, thanks to the use of low-rank approximation techniques, the memory requirement and the computation time can be drastically reduced to a very small percentage of the cost required by the full matrix problem. Indeed, thanks to the smoothness of the axisymmetric kernel, a very low compression ratio can be obtained, e.g. less then 3% (lower than the ones usually attained in 3-D problems [START_REF] Torchio | PEEC-Based Analysis of Complex Fusion Magnets During Fast Voltage Transients With H-Matrix Compression[END_REF], [START_REF] Voltolina | High-Performance PEEC Analysis of Electromagnetic Scatterers[END_REF]). Thus, the use of low-rank approximation techniques widely increases the applicability of the proposed formulations.

The rest of the paper is organized as follows. First, in section II, the amperian method is formulated and then discretized. Then, in section III, a different method based on the coulombian interpretation is proposed. In section IV, the numerical features of the two methods are discussed and compared. Most of the discussions still hold for the related 3-D methods proposed in [START_REF] Torchio | A 3-D PEEC Formulation Based on the Cell Method for Full-Wave Analyses With Conductive, Dielectric, and Magnetic Media[END_REF] and [START_REF] Torchio | An Extension of Unstructured-PEEC Method to Magnetic Media[END_REF]. Then, in section V, the case of a realistic induction cookware, presented in [START_REF] Moro | Impedance design of cooking appliances with multilayer induction-efficient cookware[END_REF], is considered.

II. AMPERIAN METHOD

The 3-D conductive and magnetic domains are defined as Ω c and Ω m , respectively. These domains can intersect and their union is defined as the active domain Ω a = Ω c ∪ Ω m . These domains are obtained from the axial revolution of 2-D domains Ωc , Ωm , and Ωa , i.e.

Ω c = Ωc × [0, 2π], Ω m = Ωm ×[0, 2π],

A. Formulation

As in [START_REF] Torchio | A 3-D PEEC Formulation Based on the Cell Method for Full-Wave Analyses With Conductive, Dielectric, and Magnetic Media[END_REF], the vector potential A and the scalar electric potential ϕ e are introduced together with

E(r) = -iωA(r) -∇ϕ e (r) + E ext (r), (1)
B(r) = ∇ × A(r) + B ext (r), ( 2
)
where E is the electric field, B is the magnetic flux density field, r is the field point, i is the imaginary unit, and ω is the angular frequency. E ext and B ext are the external electric and magnetic flux density fields. Equations ( 1) and ( 2) are complemented by the following constitutive relationships:

E(r) =ρ c (r)J c (r), r ∈ Ω c (3) B(r) =α m (r)M(r), r ∈ Ω m (4)
where J c is the conduction current density, M is the magnetization, ρ c is the electric resistivity and

α m (r) = µ 0 µ r (r)(µ r (r) -1) -1 , (5)
in which µ 0 is the permeability of vacuum and µ r is the relative permeability.

When the axisymmetric assumption holds, the EM quantities can be written in terms of cylindrical coordinates ( r, φ, z). The vector EM quantities only exhibits some (cylindrical) coordinates, i.e.

J c = J c φ φ, J a = J a φ φ, K a = K a φ φ, E = E φ φ, A = A φ φ, B = B r r+B z z, M = M r r + M z z. ( 6
)
where J a = ∇ × M and K a = M × n are the volume and surface amperian currents (in which n is the outgoing unit normal of Γ m ). It is worth noting that the axisymmetric condition ( 6) assumes the magneto-quasistatic hypothesis, i.e. capacitive and time-delay effects are neglected. By imposing the Lorenz gauge condition ∇ • A = -iωε 0 µ 0 ϕ e , where ε 0 is the permittivity of vacuum, an integral expression of A can be derived by combining ( 1) and (2) with the axisymmetric Maxwell's equations [START_REF] Datta | Electric and Magnetic Fields from a Circular Coil Using Elliptic Integrals[END_REF], i.e.

A φ (r) = µ 0 Ωc g(r, r )J c φ (r )d Ω + Ωm g(r, r )J a φ (r )d Ω + Γm g(r, r )K a φ (r )d Γ , (7)
where r is the integration point and g(r, r ) is the axisymmetric Green's function [START_REF] Alexander | Generalized axially symmetric potential theory[END_REF]. The field and the integration points can be written in cylindrical coordinates, i.e. r = (r, z, φ) and r = (r , φ , z ). The axisymmetric Green's function (which is the solution of the Poisson equation for a ring carrying a uniform current [START_REF] Simpson | Simple Analytic Expressions for the Magnetic Field of a Circular Current Loop[END_REF]) is defined as

g(r, r ) = 4r K (2 -k 2 )K(k 2 ) -2E(k 2 ) 4πk 2 , (8)
where

K = (r + r) 2 + (z -z) 2 , k 2 = 4r rK -2 , ( 9
)
and K and E are the elliptic integrals of the first and second kind, respectively (defined with the same convention of [START_REF] Datta | Electric and Magnetic Fields from a Circular Coil Using Elliptic Integrals[END_REF]). The axisymmetry assumption (6) enforces the solenoidality of J c . Thus, since the volume and surface electric charges vanish (which are the sources of ϕ e ), the scalar electric potential is null, i.e. ϕ e (r) = 0.

Finally, under the axisymmetry assumption, combining (1) and ( 2) with ( 3) and (4), respectively, the final continuum equations for an axisymmetric EM problem are where A is given by the integral expression [START_REF] Freidberg | Plasma Physics and Fusion Energy[END_REF].

ρ c (r)J c (r) = -iωA(r) + E ext (r), r ∈ Ω c , ( 10
) α m (r)M(r) =∇ × A(r) + B ext (r), r ∈ Ω m , (11)

B. Discretization

The 2-D domains Ωc and Ωm in the (r, z) plane (i.e. φ = 0) are discretized into two quadrilateral grids defined as: GΩ c for Ωc , and GΩ m for Ωm . GΩ c consists of n c nodes, e c edges, and f c faces. Similarly, GΩ m consists of n m nodes, e m edges, and f m faces. The two grids overlap when a conductive magnetic media is considered.

The conduction current density vector is then expanded as

J c (r) = fc k w φ k (r)j c k , (12)
where j c k is the flux of J c through the kth quadrilateral face of GΩ c and w φ k is the azimuthal vector shape function related to f c k . The corresponding arrays of DoFs j c = (j c k ) is then introduced. The shape function w φ k is defined as

w φ k (r) = A -1 k φ, r ∈ f k , (13)
where A k is the area of f c k .

The magnetization is instead discretized as

M(r) = em k w e k (r)m k , (14)
where m k = em M • dl is the line integral of M along the kth edge of the mesh e m k and w e k is the axisymmetric edge shape function related to e m k . The array of DoFs m = (m k ) is introduced. The vector shape function w e k is obtained as the axial revolution of the (traditional) edge vector shape function. With respect to Fig. 1,w e k is given by

w e k (r) = w e k (r, z, φ) = d ZR e k , r ∈ ∆ k 0, elsewhere , (15)
where e k is the unit vector in the direction of the e m k and ∆ k is the support of w e k , i.e. the quadrilateral faces having e m k . It is useful to introduce also the expansions of J a and K a . J a is expanded similarly to [START_REF] Torchio | PEEC-Based Analysis of Complex Fusion Magnets During Fast Voltage Transients With H-Matrix Compression[END_REF] and the vector array j a = (j a k ) is introduced where j a k is the flux of J a through the kth quadrilateral face of GΩ m . Analogously, the surface amperian currents K a are expanded as

K a (r) = e b m k w φs k (r)k a k , (16)
where k a k is the flux of K a through the kth boundary edge e b m k of GΩ m . Then, the array of DoFs m k and defined as

k a = (k a k ) is
w φs k (r) = L -1 k φ, r ∈ e b k , (17)
where L k is the length of e b m k . The dual grids GΩ c and GΩ m obtained from the barycentric subdivision of GΩ c and GΩ m are now defined. These grids consist of dual geometric entities: GΩ c consists of dual ñc nodes, ẽc edges, and fc faces, while GΩ m consists of dual ñm nodes, ẽm edges, and fm faces. The geometric entities of GΩ c and the dual grid GΩ c are exemplified in Fig. 4 for a rectified geometry. Thus, the following arrays of DoFs associated to the geometric entities of GΩ c can be defined

• e = (e k ) on dual nodes ñc , where

e k = 2πr ñc E φ (r ñc ), (18)
• a c = (a c k ) on dual nodes ñc , where

a c k = 2πr ñc A φ (r ñc ), (19)
in which r ñc = (r ñc , φ ñc , z ñc ) is the position of the kth dual node ñc k . The array of DoFs e 0 related to E ext is defined likewise e. It is worth noting that, due to the axisymmetric assumption, the dual nodes ñc are in a one-to-one relationship with the 3-D (dual) edges ẽ3D c obtained as the axial revolution of ñc , as represented in Fig. 4. Indeed, ( 18) and ( 19) are equal to the line integrals of E and A along the kth 3-D edge ẽ3D c k . Other arrays of DoFs related to B and A can also be introduced for the dual magnetic grid GΩ m :

• a m = (a m k ) on dual nodes ñm , where a m k = 2πr ñm A φ (r ñm ); (20)
a m can be subdivided into a m = [a mv , a ms ], where a mv is related to the internal dual nodes and a ms is related to the boundary dual nodes;

• b = (b k ) on dual edges ẽm , where b k = em k 2πr( φ × B) • dl, (21)
i.e. b k is the flux of B trough the 3-D (dual) face f 3D m k obtained as the axial revolution of ẽm k (see Fig. 4).

The array of DoFs b 0 related to B ext is defined likewise b. Expansions ( 12), ( 14) and ( 16) can now be inserted in [START_REF] Torchio | An Extension of Unstructured-PEEC Method to Magnetic Media[END_REF] and [START_REF] Torchio | A Volume PEEC Formulation Based on the Cell Method for Electromagnetic Problems from Low to High Frequency[END_REF]. Then, the resulting equations are tested following the Galerkin scheme, i.e. ( 10) is tested with w φ k , for k = 1, • • • , f c , and ( 11) is tested with w e k , for for k = 1, • • • , e m . Thus, equations which represent the discrete form of ( 1) and ( 2) are

e + iωa c = e 0 , (22)
b

+ GΩ m a m = b 0 , (23)
where GΩ m is the dual edge to dual node (dual) incidence matrix of dimension (ẽ m × ñm ).

The discrete forms of ( 3) and ( 4) result in

e = R c j c , b = Fm, (24)
where R and F are sparse constitutive matrices:

R c,kh = Ωc 2πrw φ k (r) • w φ h (r)ρ c (r)d Ω, (25)
F kh = Ωm 2πrw e k (r) • w e h (r)α m (r)d Ω. ( 26
)
The arrays a c and a m = [a mv , a ms ], are instead obtained as

a c =L c c j c + L c m j a + L c ms k a , a mv =L m c j c + L m m j a + L m ms k a , a ms =L ms c j c + L ms m j a + L ms ms k a , (27)
in which the volume, surface, and volume-surface inductance matrices are defined as

L α β ku =µ 0 Ωβ 2πrw φ u (r) • Ωα w φ k (r )g(r, r )d Ω d Ω, L ms α ku =µ 0 Ωα 2πrw φ u (r) • Γm w φs k (r )g(r, r )d Γ d Ω, L ms ms uk =µ 0 Γm 2πrw φs u (r) • Γm w φs k (r )g(r, r )d Γ d Γ,
where α = c, m, β = c, m, L α β ku = L β α uk , and L ms α ku = L α ms uk . Moreover, j a , k a , and m satisfy:

j a = CΩ m m, k a = CΩ m m, j a k a = C a Ωm m, (28)
where CΩ m is the (f m × e m ) face to edge incidence matrix, CΓ m is the (e b m × e m ) selection matrix which selects the boundary edges e b m of GΩ m , and

C a Ωm = CΩ m CΓ m = GT Ωm .
Finally, the discrete EM problem can be represented as the following algebraic system of equations

  Rc + iωL c c iω[L c m , L c ms ]C a Ωm -GΩ m L m c L ms c F -GΩ m L m m L m ms L ms m L ms ms C a Ωm   jc m = e0 b0 , (29)
which can be symmetrized by multiplying the second row by -iω. System (29) naturally enforce the EM properties of the magnetization, but not in a numerically strong sense. For instance, in steady state condition and for homogeneous non-conductive magnetic media, only the surface amperian currents exist while the volume amperian currents are zero (∇×M = J a ≈ 0). This condition is not strongly enforced by (29) and a distribution of J a close to the boundary is obtained also when a zero frequency problem is solved. Anyhow, this usually does not particularly affect the global quality of the solution but it can be an issue when high accuracy is required. However, this condition can be enforced by imposing:

m = GΩ m ψ m , (30)
where GΩ m is the edge to node (e m × n m ) incidence matrix and ψ m is the array of DoFs related to a discrete scalar magnetic potential associated to the primal nodes n m of GΩ m (one node for each magnetic sub-domain must be removed and considered as reference). Thus, when ∇ × M holds, (30) can be applied to (29) and the second row of (29) can be projected into a new set of equations by using G T Ωm (without one row for each magnetic sub-domain). The resulting system has a reduced dimension and, since j a = 0, the evaluation of L m c , L c m , and L m m is not needed anymore.

III. COULOMBIAN METHOD

In the following an alternative method for the study of axisymmetric EM problems is proposed. In section II, the magnetization phenomena have been considered following the amperian interpretation. Here, instead, the magnetic currents and charges are introduced and different electric and magnetic potentials are adopted.

A. Formulation

Equations ( 1) and ( 4) can be replaced by the following

E(r) = -iωA e (r) -∇ϕ e (r) -ε -1 0 ∇ × A m (r) + E ext (r), (31)
H(r) = -iωA m (r) -∇ϕ m (r) + µ -1 0 ∇ × A e (r) + H ext (r), (32)
where A e and A m are the electric and magnetic vector potentials [START_REF] Torchio | An Extension of Unstructured-PEEC Method to Magnetic Media[END_REF], ϕ m is the scalar magnetic potential, H(r) is the magnetic field, and H ext (r) is the external magnetic field. The vector potentials A e and A m are related to A by

A = A e + (iω) -1 ∇ × A m , (33)
and they are defined as in [START_REF] Torchio | An Extension of Unstructured-PEEC Method to Magnetic Media[END_REF] together with ϕ m . The Lorenz gauge conditions are also imposed:

∇ • A e = -iωε 0 µ 0 ϕ e and ∇ • A m = -iωε 0 µ 0 ϕ m .
Relations [START_REF] Naseh | Levitation Force Maximization in HTS Magnetic Bearings Formulated by a Semianalytical Approach[END_REF] still hold and they are complemented with

A e = A e φ φ, H = H r r + H z z, J m = J mr r + J mr z, A m = A mr r + A mr z. ( 34
)
where J m = iωµ 0 M is the magnetic current density.

As noticed in section II, when the axisymmetric assumption holds, J c is divergence free and ϕ e (r) = 0.

The integral expression of A e can be formulated in cylindrical coordinates as shown in [START_REF] Datta | Electric and Magnetic Fields from a Circular Coil Using Elliptic Integrals[END_REF]:

A e φ (r) = µ 0 Ωc g(r, r )J c φ (r )dΩ, (35)
with the same definition given for [START_REF] Freidberg | Plasma Physics and Fusion Energy[END_REF].

The definition of ϕ m can be also derived from [START_REF] Datta | Electric and Magnetic Fields from a Circular Coil Using Elliptic Integrals[END_REF]:

ϕ m (r) = Ωm ḡ(r, r ) µ 0 m (r )d Ω + Γm ḡ(r, r ) µ 0 ς m (r )d Γ, (36)
where, m and ς m are the volume and surface bound magnetic charges related to J m by means of the continuity relations:

∇ •J m (r ) = -iω m (r ), J m (r )•n = -iως m (r ). (37)
In (36), ḡ is the axisymmetric Green's function of the scalar potential. The axisymmetric Green's function ḡ is the solution of the Poisson equation for a uniformly charged ring [START_REF] Ciftja | The electrostatic potential of a uniformly charged disk as the source of novel mathematical identities[END_REF], [START_REF] Huré | Solutions of the axi-symmetric Poisson equation from elliptic integrals[END_REF], that is

ḡ(r, r ) = 4r 4πK K(k 2 ), (38)
where the quantities are defined as in [START_REF] Fano | Electromagnetic Fields, Energy, and Forces[END_REF].

The axisymmetric integral expression of A m can be extrapolated from the one of the scalar magnetic potential, i.e.

A mα (r) = ε 0 Ωm ḡ(r, r )J mα (r )dΓ, (39)
where α = r, z. Equations ( 31) and (32) are complemented by (3) and the following constitutive equation:

H(r) =ρ m (r)J m (r), r ∈ Ω m , (40)
where ρ m = (iωµ 0 (µ r -1)) -1 is the magnetic resistivity. Finally, combining (31) with (3), and (32) with (40), the following final continuum equations for an EM axisymmetric problem are obtained

ρ c (r)J c (r) = -iωA e (r, J c (r )) -ε -1 0 ∇ × A m (r, J m (r )) + E ext (r), r ∈ Ω c , (41)
ρ m (r)J m (r) = -iωA m (r, J m (r )) -∇ϕ m (r, J m (r )) + µ -1 0 ∇ × A e (r, J c (r )) + H ext (r), r ∈ Ω m , (42)
where A e , ϕ m , and A m are given by ( 35), (36), and (39), respectively.

B. Discretization

Conductive and magnetic domains are represented by means of 2-D domains Ωc and Ωm in the (r, z) plane (i.e. φ = 0) as in section II. Again, the two quadrilateral grids GΩ c and GΩ m are introduced with the same definitions given above. Moreover, the dual grids GΩ c and GΩ m are also introduced and the conduction current J c is expanded as in [START_REF] Torchio | PEEC-Based Analysis of Complex Fusion Magnets During Fast Voltage Transients With H-Matrix Compression[END_REF].

The only difference concerning conductive media between the two formulations is in the definition of the arrays of DoFs related to the vector potentials. In this formulation A is replaced by A e and A m and the corresponding arrays of DoFs are introduced:

• a c m = (a c m,k ) on primal edges e c , where

a c m,k = e m,k A m (r) • dl. (44)
The magnetic current density J m is expanded as where

J m (r) = em k w f k (r)j m,k , (45)
j m,k = em 2πr( φ × J m (r)) • dl,
is the flux of J m through the kth 3-D face f 3D m k represented in Fig. 3 and related to the kth dual edge e m k (i.e. f 3D m k is the axial-revolution of e m k , as shown in Fig. 3). Thus, the corresponding array of DoFs j m = (j m k ) is introduced. w f k is the (scaled) face vector shape function related to e m k and exemplified in Fig. 1. With respect to Fig. 1,w f k is given by

w f k (r) = w f k (r, z, φ) = d A 3D k R ( φ × e k ), r ∈ ∆ k 0, elsewhere , (46)
where A 3D k is the area obtained from the axial revolution of e m k (i.e. the area of f 3D m k and A 3D k = 2πZ). The expansion (45) combined with (37) induces an expansion for m and ς m . Thus, the continuity equations (37) in discrete form become

CΩ m j m = -iωq v m , CΓ m j m = -iωq s m , C a Ωm j m = -iω q v m , q s m T , (47)
where

q v m = (q v m,k ), for k = 1, • • • , f m , and q s m = (q s m,k ), for k = 1, • • • , e b
m , store the (net) volume and surface bound magnetic charges, respectively.

Additional arrays of DoFs living on the geometric entities of GΩ m and GΩ m can be introduced for H, A e , A m , and ϕ m :

• h = (h k ) on dual edges ẽm , where

h k = ẽm,k H(r) • dl, (48)
• a m m = (a m m,k ) on dual edges ẽm , where

a m m,k = ẽm,k A m (r) • dl, (49)
• φ m = (φ m,k ) on dual nodes ñm , where

φ v m,k = ϕ m (r ñm ), (50)
φ m can be subdivided into φ m = [φ v m , φ s m ],
where φ v m is related to the internal dual nodes and φ s m is related to the boundary dual nodes,

• a m e = (a m e,k ) on nodes n m , where a m e,k = 2πr nm A e φ (r nm ),

The array of DoFs h 0 related to H 0 is defined likewise h. Expansions ( 12) and ( 45) are now inserted in (41) and (42). Then, the resulting equations are tested following the Galerkin scheme, i.e. ( 41) is tested with w φ k , for k = 1, • • • , f c , and (42) is tested with w f k , for for k = 1, • • • , e m . The following equations represent the discrete form of (31) and ( 32):

e + iωa c e + M 1/ε0 CΩ c a c m = e 0 , (52)
h + iωa m m + GΩ m φ m -M 1/µ0 GΩ m a m e = h 0 , (53)
where CΩ c is the face to edge conductive incidence matrix of dimension (f c × e c ), while M 1/ε0 and M 1/µ0 are projection matrices defined as

M 1/ε0,kh = 1 ε 0 Ωc 2πrw φ k (r) • w φ h (r)d Ω, (54)
M 1/µ0,ku = 1 µ 0 Ωm 2πrw f u (r) • w f k (r)d Ω. (55)
The discrete form of ( 3) is again given by the first equation in (24). Instead, the one of ( 40) is given by

h = R m j m , (56)
where R m is a sparse constitutive matrix:

R m ku = Ωm 2πrw f u (r) • w f k (r)ρ m (r)d Ω. (57)
The others arrays a c e , a c m , a m m , a m e , and φ m are instead obtained as

a c e = L c c j c , a m m = L m m j m , a c m = N c m j m , a m e = N m c j c , [φ v m , φ s m ] = P m q v m q s m = P vv m P vs m P sv m P ss m q v m q s m , (58)
in which

P vv m ku = 1 A k V 3D u µ 0 fm k fm u ḡ(r, r )d Ω d Ω, (59)
where V 3D u is the volume of the uth 3-D volume v 3D mu (i.e. the axial revolution of f mu ) and A k is the area of f m k ;

P sv m ku = 1 L k V 3D u µ 0 e b m k fm u ḡ(r, r )d Ω d Γ, ( 60
)
where L k is the length of the boundary edge e b m k ;

P vs m ku = 1 A k A 3D u µ 0 fm k e b mu ḡ(r, r )d Γ d Ω, (61)
where A 3D u is the area of the uth 3-D face f b,3D mu (i.e. the axial revolution of e b m,u );

P ss m ku = 1 L k A 3D u µ 0 e b m,k e b m,u ḡ(r, r )d Γ d Γ. ( 62
)
The magnetic inductance matrix L m m is instead given by

L m m ku = ε 0 Ωm 2πrw f u (r) • Ωm w f k (r )ḡ(r, r )d Ω d Ω, ( 63
)
and matrices N m c and N c m (which represent the coupling between the electric and magnetic domains) are defined by

N m c ku = µ 0 2πr n m,k φ • Ωc w φ u (r )g(r, r )d Ω , (64)
where r nm is the coordinate along r of n m,k ;

N c m ku = ε 0 ec e • Ωm w f k (r )ḡ(r, r )d Ω dl, ( 65
)
where e is the unit vector in the direction of e c . Finally, the discrete EM problem derived from ( 3), ( 31), (32), and (40) can be represented as the following algebraic system of equation

R c + iωL c c M 1/ε0 CΩ c N c m -M 1/µ0 GΩ m N m c R m + iωL m m + 1 iω GΩ m P m C a Ωm × j c j m = e 0 h 0 .
(66) Due to the symmetry of ( 31) and (32), matrix blocks (1, 2) and (2, 1) should satisfy

M 1/ε0 CΩ c N c m = -(-M 1/µ0 GΩ m N m c ) T . ( 67
)
However, since N c m and N m c are only numerical approximations of integral discrete operators, this property is numerically lost. However, (67) can be enforced resulting in a more accurate method which strongly satisfies the physical properties of the EM fields. Moreover, by substituting the matrix block (1, 2) with -(2, 1) T , matrices M 1/µ0 and N m c are not needed anymore (or vice-versa). Thus, multiplying the second row of (66) by -1, the system can be symmetrized.

Finally, when homogeneous magnetic media are considered ∇ • J m = 0 (i.e. m = 0). This property is enforced by (66) but not in a numerically strong sense. However, the vanishing of m can be strongly enforced by imposing

j m = GΩ m t m , ( 68
)
where t m is the array of DoFs related to a discrete vector potential associated to the nodes n m (3-D edges n 3D m ) of GΩ m , (except one for each magnetic sub-domain considered as reference). Then, (68) is applied to (66) and GΩ m (without one column for each magnetic sub-domain) is used to project the second matrix row of (66) into a new set of equations.

The resulting system has a reduced dimension and it strongly enforces CΩ m j m = 0, i.e. q v m = 0 and ∇•J m = 0. Moreover, matrices P vv m , P vs m , and P sv m not needed anymore.

IV. COMPARISON BETWEEN THE METHODS

The method proposed in section II is based on the amperian interpretation of the magnetization phenomena. When this approach is adopted, the inductance matrices are the only full integral matrices to be evaluated, resulting in a relatively easy implementation. Moreover, several analytical expressions for the self and mutual inductances of loops with rectangular cross section exist [START_REF] Chitarin | An integral formulation for eddy current analyses in axisymmetric configurations[END_REF]. Thus, all the inductance matrices can be efficiently evaluated using analytical formulas.

The main drawback of this method consists in a poor numerical performance when magnetic media with relative high permeability are involved. Especially when they have a non-zero conductivity or when they are very close to conductive media. Indeed, the constitutive relation (4) depends on α m that, with the increase of the relative permeability, tends to the value µ 0 (i.e. α m → µ 0 with the increasing of µ r ). For instance, from a numerical point of view, the difference between a media with µ r = 1,000 and µ r = 10,000 is left to the last significant digits of α m and, due to round-off errors [START_REF] Turing | Rounding-Off Errors in Matrix Process[END_REF], this numerical information can be completely lost when floating-point arithmetic is used. When F is summed to other integral matrices (which are affected by numerical approximations due to the integration of the Green's function) this numerical issue is enhanced. Indeed, the numerical approximations which affect the integral matrices are orders of magnitude greater than the last significant digits of α m . This problem is alleviated when (30) is enforced. However this is only possible when non-conductive magnetic media are considered.

The effort required by the implementation of the coulombian method described in section III is slightly more expensive since different integral matrices must be evaluated and matrix L m m does not allow for an easy analytical evaluation. Moreover, when properties (30) and (68) are enforced, the amperian approach leads to a magnetic matrix equation with a lower sparsity ratio (∼ n b m /n m ) while the magnetic equation of (66) remains full. However, system (66) does not suffer from the numerical issue described above for conductive magnetic media with high permeability. Moreover, when the frequency is sufficiently low, L m m ≈ 0. Thus, system (66) can be further simplified by neglecting the presence of the term iωL m m in the matrix block (2, 2), which leads to a relatively sparse magnetic equation when property (68) is enforced (sparsity ratio ∼ n b m /n m ). Both systems (29) and (66) are suitable for the application of low-rank approximation techniques [START_REF] Torchio | PEEC-Based Analysis of Complex Fusion Magnets During Fast Voltage Transients With H-Matrix Compression[END_REF], allowing a drastic reduction of the storage and the computational cost arising from the handling of full matrices. In the numerical studies, the low-rank approximation techniques described in [START_REF] Voltolina | High-Performance PEEC Analysis of Electromagnetic Scatterers[END_REF] and implemented in HLIBpro library [START_REF] Kriemann | .7.1)[END_REF] has been applied to both the proposed methods. HLIBpro relies on the hierarchical (H)-matrix representation coupled with adaptive cross approximation (ACA). Such library also provides a partitioning algorithm, based on geometric bisection, which splits the DoFs into disjoint clusters and builds a binary tree which terminates as soon as a prescribed condition is reached. However, in order to obtain matrix blocks which are actually low-rank, it is important to force the first clustering partition of the unknowns (which is viable thanks to the functionality of HLIBpro). Thus, when (29) (or (66)) is solved, j c and m (or j c and j m ) must be subdivided into two different clusters, and m (or j m ) should be further subdivided into boundary and internal unknowns. Then, as described in [START_REF] Voltolina | High-Performance PEEC Analysis of Electromagnetic Scatterers[END_REF], the partitioning algorithm proceeds with the construction of the cluster tree starting from these three master cluster nodes. A completely equivalent discussion holds for ψ m (or t m ) when (30) (or (68)) is enforced. When these measures are embraced, thanks to the smoothness of the Green's kernel, very low compression ratios can be reached (e.g. less than 3%, as shown in the numerical studies). Finally, it is worth noting that both (29) and (66) can be adopted for the solution of steady state problems: (29) can be directly solved when ω = 0 while (66) only requires a scaling of the unknowns j m by iω.

Sample MATLAB R codes of the two formulation are available at https://github.com/UniPD-DII-ETCOMP/ DenseMatrixMarket.

V. NUMERICAL RESULTS

Two algorithms based on the two proposed integral axisymmetric formulations have been implemented in MATLAB R and parallel MEX-FORTRAN functions based on OpenMP libraries. Moreover, the code has been also combined with HLIBpro library in a matrix-free logic. Thus, when HLIBpro is adopted, systems (29) and (66) are never actually assembled and stored. In the following, the code based on the formulation of section II is defined as A-code (Amperian code) while the one based on section III is defined as C-code (Coulombian code). Self inductance coefficientsare evaluated following [START_REF] Chitarin | An integral formulation for eddy current analyses in axisymmetric configurations[END_REF], whereas mutual coefficients are numerically evaluated with 8 and 9 quadrature points. Simulations were run on a Linux machine equipped with a Xeon E5-2643 v4 processor (dual 6-core, 24 thread, @3.40 GHz) and with 512 GB of RAM.

The case of the induction cookware presented in [START_REF] Moro | Impedance design of cooking appliances with multilayer induction-efficient cookware[END_REF] is considered. The geometry of the problem is completely described by Fig. 4. With respect of Fig. 4, layers l 1 , l 2 , l 3 , and l 4 constitute the bottom part of the pot, layer l 5 is the external current-driven coil (the problem excitation), l 6 is the magnetic flux concentrator, and l 7 is the shield. The material parameters are the ones in Table I. With the aim of investigating the numerical behavior of the A and C-codes under different conditions, the magnetic steel (layer l 3 in Fig. 4) is modeled in two different ways: a) by non conductive magnetic media, and b) by conductive magnetic media with a conductivity of 2 MS/m. The frequency of the problem, f , is swept from 20 kHz to 100 kHz. The equivalent impedance of the device, with respect to a uniform (external) current flowing in the external coil, is modeled as Z eq = R + iωL and evaluated with the two proposed methods. Due to the strong skin effects (skin depth of δ ≈ 0.06 mm for the worst case), a great amount of mesh elements is required for the discretization of the conductive layers. Thus, the mesh of the induction cookware model consists of 54,120 conductive quadrilateral faces and 16,680 magnetic quadrilateral faces. With the aim of allowing a fair comparison, the same mesh is adopted for the A and C-codes. In Fig. 5, the results obtained from A-code and C-code are compared with the ones from the commercial Finite Element Method (FEM) software COMSOL R , which allows for a very efficient solution of this kind of device and therefore is considered as a reference. The axisymmetric FEM model consists of 64,847 mesh elements and third order basis function.

Fig. 6 shows the relative error of the magnitude of Z eq obtained from A-code and C-code with respect to FEM for both cases (a) and (b).

Results in Fig. 5 and Fig. 6 show an excellent agreement when case (a) is considered. For case (a), the simulations have been carried out with and without the imposition of (30) and (68). However, as expected, even if a little distortion on the distribution of M is obtained when (30) and (68) are not enforced, the results in terms of Z eq differ less than 0.05%.

For case (b), C-code still shows an excellent agreement with FEM. A-code instead shows a lower accuracy, especially for the real part of Z eq which is more affected by the actual distribution of M, which is not enforced in a numerically strong sense. Indeed, it is worth noting that for case (b), since a conductive magnetic media is considered, (30) cannot be enforced since M is not curl free. Moreover, due to roundoff errors, the A-code is very sensitive to the accuracy of the numerical integrations. However, due to the numerical nature of round-off errors, a higher accuracy in the computation of the integral matrices does not ensure to avoid this numerical issue. For the sake of comparison, the same number of Gauss points has been chosen for the A and C-codes.

In Table II the computational details of the two methods are reported for case (b) and f = 100 kHz. For the A and the C-codes, the computational times and the Peak Memory Usage (PMU) are reported with and without the adoption of low-rank approximation techniques. The results of Table II have been obtained without the imposition of (30) and (68), the compression tolerance described in [START_REF] Voltolina | High-Performance PEEC Analysis of Electromagnetic Scatterers[END_REF] is set to ε ACA = 10 -3 (much lower than the ones adopted in [START_REF] Voltolina | High-Performance PEEC Analysis of Electromagnetic Scatterers[END_REF] for a 3-D case) and the admissibility parameter is η = 2. With this choice of ε ACA , the relative error introduced by the low-rank approximation on the solution is less than 1%. Thus, it is safe to say that the errors introduced by the low-rank approximation do not noticeably affect the accuracy (results shown in Fig. 5 are the ones obtained from HLIBpro).

With respect Table II, for the uncompressed case the problem is solved with an LU decomposition of the matrix system. Instead, the compressed case is solved by means of GMRES solver with a H-LU preconditioner. Thus, the solution time includes both the H-LU and the GMRES solution. For both the uncompressed and compressed cases, the PMU almost coincides with the memory required for the storage of the system matrices and LU factorizations. The H-LU is generated by imposing the same tolerance ε ACA and its storage requires almost the same memory of the compressed system. As can been seen from Table II, due to the big size of the problem, the use of low-rank approximation techniques is mandatory. Indeed, the adoption of HLIBpro library allows for a drastic reduction of the computational cost required by the IE methods which become comparable with the highly optimized FEM commercial software for this kind of application.

VI. CONCLUSION

Two axisymmetric integral equations methods have been proposed. Both methods allow for conductive and magnetic media and inhomogeneous properties of the materials can be considered as well. The numerical features of the two methods have been widely discussed and compared. Most of the considerations of section IV still hold for the 3-D methods presented in [START_REF] Torchio | A 3-D PEEC Formulation Based on the Cell Method for Full-Wave Analyses With Conductive, Dielectric, and Magnetic Media[END_REF] and [START_REF] Torchio | An Extension of Unstructured-PEEC Method to Magnetic Media[END_REF] from which the axisymmetric methods are derived. Low-rank approximation techniques based on HLIBpro library have been successfully combined with both the proposed methods obtaining a drastic reduction of the computational cost with a truly negligible loss of accuracy. Thanks to the regularity of the axisymmetric Green's kernel, the compression ratio attained in the numerical studies is lower than the ones obtained for integral 3-D methods as in [START_REF] Torchio | PEEC-Based Analysis of Complex Fusion Magnets During Fast Voltage Transients With H-Matrix Compression[END_REF] and [START_REF] Voltolina | High-Performance PEEC Analysis of Electromagnetic Scatterers[END_REF]. As demonstrated by the numerical results, the use of HLIBpro library widely increases the applicability of the proposed formulations that become competitive with finite element method approaches and which may become preferable when large air domains or small air gaps are considered.

  and Ω a = Ωa ×[0, 2π]. The boundaries of Ω c and Ω m are defined as Γ c = ∂Ω c and Γ m = ∂Ω m , respectively, and analogously Γ c = Γc × [0, 2π], Γ m = Γm × [0, 2π], where Γc = ∂ Ωc and Γm = ∂ Ωm .

Figure 1 .

 1 Figure 1. Axial-symmetric vector edge (red) and face (blue) shape functions related to edge em k . An equivalent representation holds for edges along the r direction.

Figure 2 .

 2 Figure 2. Primal and dual geometric entities for the conductive grid. Rectified axisymmetric geometry introduced. w φs k is the azimuthal vector shape function related to e b m k and defined as

  ) on dual nodes ñc , where a c e,k = 2πr ñc A e φ (r ñc ),

Figure 3 .

 3 Figure 3. Primal and dual geometric entities of the magnetic grid for the coulombian formulation. Rectified axisymmetric geometry.

Figure 4 .

 4 Figure 4. Axial-symmetric model of the induction cookware. Dimensions are in mm. l 1 ) Stainless Steel, l 2 ) Aluminum, l 3 ) Magnetic Steel, l 4 ) Stainless Steel, l 5 ) External coil, l 6 ) Ferrite (flux concentrator), l 7 ) Aluminum (shield).

Figure 5 .

 5 Figure 5. Zeq of the induction cookware. Results obtained from A-code, C-code, and FEM for cases (a) and (b).

Figure 6 .

 6 Figure 6. Relative error of |Zeq| for A-code and C-code for cases (a) and (b).


Coralie Valle

Baptiste Martin

Françoise Debart

Jean-Jacques Vasseur

Isabelle Imbert

Bruno Canard

Bruno Coutard

Etienne Decroly
email: etienne.decroly@afmb.univ-mrs.fr

The C-Terminal Domain of the Sudan Ebolavirus L Protein Is Essential for RNA Binding and Methylation

d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

INTRODUCTION

Ebola virus is an emerging virus causing severe epidemics such as the unprecedented devastating outbreak of 2015 in West Africa, which killed more than 11,000 people [START_REF] Shultz | Distinguishing epidemiological features of the 2013-2016 West Africa Ebola virus disease outbreak[END_REF]. Although Ebola faded from public attention, outbreaks continue to occur especially in Democratic Republic of the Congo (DRC) where the current epidemic, which has been declared in August 2018, has already killed more than 1,000 citizens [START_REF]Ebola | Ebola situation reports: Democratic Republic of the Congo[END_REF]. This virus infects human and non-human primates to provoke a massive viremia leading to haemorrhagic fever, which is fatale in most cases. Despite a vaccine is currently in development [START_REF] Henao-Restrepo | Efficacy and effectiveness of an rVSV-vectored vaccine expressing Ebola surface glycoprotein: interim results from the Guinea ring vaccination cluster-randomised trial[END_REF] [START_REF] Wells | Ebola vaccination in the Democratic Republic of the Congo[END_REF] and some therapeutic antibodies have shown antiviral activities [START_REF] Saphire | Antibody-mediated protection against Ebola virus[END_REF], we urgently need to develop new therapeutic compounds to treat infected people and contact patients.

Ebolavirus genus contains five strains: Zaire ebolavirus (Ebola virus, EBOV), Sudan ebolavirus (Sudan virus, SUDV), Tai Forest ebolavirus (Tai Forest virus, TAFV), Bundibugyo ebolavirus (Bundibugyo virus, BDBV), and Reston ebolavirus (Reston virus, RESTV). These filamentous viruses together with Marburg marburgvirus (Marburg virus, MARV), also responsible of acute haemorrhagic fever, and Lloviu cuevavirus form the Filoviridae family [START_REF] Kuhn | Proposal for a revised taxonomy of the family Filoviridae: classification, names of taxa and viruses, and virus abbreviations[END_REF]. They are enveloped virus with nonsegmented negative strand RNA genome (NNS) belonging to the Mononegavirales order. These mononegaviruses, which contain several important human pathogens such as measles, RSV and rabies viruses, display similar genetic organization and share common replication strategies.

The genome of filoviruses is about 19 kb long and codes for seven major proteins: the nucleoprotein (NP), viral proteins VP35 and VP40, the glycoprotein (GP), viral proteins VP30 and VP24 and the "large" protein L [START_REF] Sanchez | Sequence analysis of the Ebola virus genome: organization, genetic elements, and comparison with the genome of Marburg virus[END_REF] [START_REF] Elliott | Ebola protein analyses for the determination of genetic organization[END_REF]. The viral life cycle is initiated by the interaction of the GP envelope protein with cell surface determinants allowing virus entry into the cell by macropinocytosis. The NPencapsidated viral genome is next released in the host cell cytoplasm where the virus replication cycle is initiated. The large protein L, associated with NP, VP35 and VP30 first transcribes viral mRNAs by using a discontinuous transcription mechanism that generates seven monocistronic capped and polyadenylated RNAs [START_REF] Mühlberger | Termini of all mRNA species of Marburg virus: sequence and secondary structure[END_REF] [START_REF] Weik | Ebola virus VP30-mediated transcription is regulated by RNA secondary structure formation[END_REF]. The polymerase complex initiates the mRNA synthesis at a conserved genestart (GS) sequence and transcribes the genes to the gene-stop sequence (GE) [START_REF] Whelan | Transcription and replication of nonsegmented negative-strand RNA viruses[END_REF]. The poly-A tail is then added by a slippage mechanism of the polymerase occurring in a conserved poly(U)6 tract in the GE sequence [START_REF] Jordan | Initiation, extension, and termination of RNA synthesis by a paramyxovirus polymerase[END_REF]. The RNA-dependent RNA polymerase (RdRp) then scans the intergenic sequences and reinitiates the RNA transcription at the next GS sequence. The viral mRNAs of mononegaviruses are co-transcriptionally capped by an unusual sequence of capping reactions. This process was elucidated for vesicular stomatitis virus (VSV) and takes place when the nascent RNA chain has reached 31 nucleotides length [START_REF] Tekes | A freeze frame view of vesicular stomatitis virus transcription defines a minimal length of RNA for 5' processing[END_REF]. The cap synthesis starts with the formation of a covalent adduct between a conserved histidine residue of the polyribonucleotidyltransferase (PRNTase) embedded into the Cap domain and the 5'-phosphate of viral mRNA. In presence of GDP, the PRNTase transfer a GMP moiety to the 5′-diphosphate of the covalently bound RNA, yielding to the formation of the cap structure (GpppN). The 5' cap structure is subsequently methylated at the 2'O position of the first nucleotide (N1) and at the N7 position of the guanosine ( m GpppNm or cap-1) [START_REF] Shatkin | Capping of eucaryotic mRNAs[END_REF]. This cap-1 structure is undistinguishable from cellular caps and plays a critical role for virus replication as it protects viral mRNAs from cellular 5' exonucleases and it allows the eIF4e-dependent initiation of mRNA translation. Additionally, the cap-2'O methylation is a self-marker that limits the detection of infecting viral RNAs by the innate immunity sensors of the retinoic acid-inducible gene-I (RIG-I) like family. Thus, viruses expressing RNAs lacking cap 2'O methylation ( m GpppN-RNA or cap-0) are early detected by RIG-I receptors in infected cells leading to the secretion type I interferon (INFα/β).

Furthermore, INFs stimulate the expression of INF-stimulated genes (ISG) such as INF Induced Protein with Tetratricopeptide Repeats (IFIT) that sequesters miss-capped RNAs (ie: cap-0) leading to restriction of virus propagation (for review see [START_REF] Decroly | Conventional and unconventional mechanisms for capping viral mRNA[END_REF]). Capping enzymes such as cap-MTases are thus considered as potent antiviral targets as the inhibition of N7 MTase is supposed to limit RNA translation into viral proteins and 2'O MTase inhibitors should unmask viral RNA to the intracellular self-sensors.

The structure of VSV L protein, recently solved by cryo-electromicroscopy (cryo-EM) at 3.8 Å, reveals its organization in 5 main topological domains [START_REF] Liang | Structure of the L Protein of Vesicular Stomatitis Virus from Electron Cryomicroscopy[END_REF]. The RdRp, that contains three conserved regions (CRI, CRII, and CRIII), is intimately associated with the PRNTase domain (CRIV) forming a "donut-like" structure observed by negative staining [START_REF] Rahmeh | Molecular architecture of the vesicular stomatitis virus RNA polymerase[END_REF]. The RdRp is organised as a right-hand "fingers-palm-thumb" structure, typical of polymerases. The Cap domain (PRNTase), which harbours an original fold, projects a loop near the catalytic site of the polymerase. This structural interplay suggests that the PRNTase domain participates to the initiation of polymerisation similarly to polymerases having a priming loop. In addition, the cryo-EM structure corroborates biochemical analysis showing the interplay between RdRp and capping activities [START_REF] Tekes | A freeze frame view of vesicular stomatitis virus transcription defines a minimal length of RNA for 5' processing[END_REF]. The EM studies shows that the "donuts-like" structure (RdRp and PRNTase domains) is followed by 3 globular structures corresponding to the connector domain (CD), the MTase domain and a small C-terminal domain (CTD), respectively. The CD forms a globular domain, separating the PRNTase domain from the MTase domain, which is composed of a bundle of eight helices. It is likely that the CD domain plays an organizational role as it interacts with the viral P protein. Conversely, the MTase domain contains a K-D-K-E catalytic tetrad, characteristic of 2'O MTases. The VSV cryo-EM structure reveals its organization in a typical Rossman fold of S-adenosyl methionine (SAM) dependent methytransferases. The 2.2 Å X-ray structure of human metapneumovirus (hMPV) MTase+CTD domain confirms the methyltransferase fold with a SAM-binding site positioning the released methyl group face to the catalytic tetrad (Paesen et al.). In contrast, the MTase+CTD lacks a canonical cap-binding site like observed in flavivirus MTases [START_REF] Egloff | An RNA cap (nucleoside-2'-O-)-methyltransferase in the flavivirus RNA polymerase NS5: crystal structure and functional characterization[END_REF][19] [START_REF] Zhou | Structure and Function of Flavivirus NS5 Methyltransferase[END_REF]. Thus, it is likely that the RNA substrate accommodates an unusually narrow RNA-binding groove formed by the MTase domain overlayed by the CTD enriched in basic amino acids.

The MTase activity of VSV L protein was demonstrated in vitro for VSV L protein that sequentially methylates the RNA cap structure at the ribose 2'O and the guanosine N7 positions [START_REF] Rahmeh | Ribose 2'-O Methylation of the Vesicular Stomatitis Virus mRNA Cap Precedes and Facilitates Subsequent Guanine-N-7 Methylation by the Large Polymerase Protein[END_REF]. The catalytic activity of the MTase domain has been confirmed in vitro using MTase+CTD constructs of hMPV and filovirus (SUDV) [START_REF] Paesen | X-ray structure and activities of an essential Mononegavirales L-protein domain[END_REF] [START_REF] Martin | The methyltransferase domain of the Sudan ebolavirus L protein specifically targets internal adenosines of RNA substrates, in addition to the cap structure[END_REF]. The MTase+CTD domains of these viruses have also been reported to methylate the 2'O position of the N1 of capped and triphosphate RNAs. The subsequent methylation of the N7 guanosine was evidenced using longer RNA substrates using SUDV MTase+CTD. Besides this, we also demonstrated recently that the SUDV MTase and hMPV MTase (in a lesser extend) harbour a cap-independent methyltransferase activity inducing internal adenosine 2'O methylations. The role of such epitranscriptomic RNA modification is not yet fully understood, but similar internal 2'O methylations observed inside the HIV RNA genome have been shown to protect viral RNAs from detection by the RIG-like receptor Melanoma Differentiation-Associated protein 5 (MDA5) [START_REF] Ringeard | FTSJ3 is an RNA 2′-O-methyltransferase recruited by HIV to avoid innate immune sensing[END_REF].

Altogether, these results suggest that the filovirus MTase evolved towards a dual activity with distinct substrate specificities. The cap-N7 MTase activity that promotes viral protein translation, the cap-

RESULTS

The CTD domain of SUDV L protein is enriched in basic amino acids and stabilizes the MTase domain.

The C-terminal part of SUDV L protein contains a conserved MTase domain including the conserved region VI (CRVI) followed by a CTD domain of unknown function (Fig. 1A). The MTase harbours a K-D-K-E catalytic tetrad, canonic of 2'O MTases, and a SAM binding GxGxG motif highlighted in the multi-alignment of filovirus C-terminal part of L proteins (Fig. S1). The MTase domain (residues 1693 to 2036) is followed by the CTD (residues 2037 to 2210), which greatly varies in length within the Mononegavirales order (from ~120 residues in Pneumoviridae to ~240 residues in Rhabdoviridae) (Fig. S2). Despite this weak conservation, CTDs are enriched in basic residues in most of Mononegavirales.

Calculated isoelectric points of each CTD vary between 8.33 and 9.94 except for rubelavirus and bornavirus (Fig. S2). The alignment of filovirus CTDs also reveals the presence of 3 strictly conserved basic amino acids (arginine and lysine, indicated by *) and 2 other conserved positions with either arginine or lysine residues (indicated by x) (Fig. S1). Two conserved aromatic residues (histidine and phenylalanine, indicated by *), which might participate to RNA recruitment by interacting with nucleobases by staking, are also present in the CTD. We also analysed the conservation of residues in the protein alignment of the MTase+CTD domains of Mononegavirales. Results are represented onto the hMPV structure (Fig. 1B&C). The conservation model indicates that the most conserved residues are localized around the catalytic pocket, the SAM binding site and in a groove enriched in basic amino acids previously proposed to accommodate the RNA substrate [START_REF] Paesen | X-ray structure and activities of an essential Mononegavirales L-protein domain[END_REF]. These observations together with the structural analysis of hMPV MTase+CTD showing that MTase domain is overlayed by the CTD to form a narrow RNA binding groove suggest that the CTD plays a key role in RNA recognition.

The CTD domain L protein is essential for the SUDV MTase activity.

We expressed in bacteria the MTase adjoining the CTD (MTase+CTD) and the SUDV MTase domain (MTase-ΔCTD). Both proteins were purified by affinity chromatography on NTA column combined with exclusion chromatography. The purity of both proteins was assessed by SDS-PAGE analysis upon Coomassie blue staining. Fig. 2A shows that MTase+CTD and SUDV MTase-ΔCTD migrate at their expected apparent molecular weight of 58 and 39 kDa, respectively. Their identity was confirmed by MALDI-TOF (not shown). Both proteins were also analysed by thermal shift assay (TSA) in order to assess their folding and compare their stability (Figure 2B). The proteins revealed typical denaturation curve of folded proteins with a Tm value of 54.3°C and 46.0°C for the MTase+CTD and the MTase-ΔCTD, respectively. These data suggest that the CTD domain increases the stability of the MTase domain. We next demonstrated the involvement of the CTD in the MTase activity by comparing the enzymatic activity of both proteins on a capped RNA (GpppN-RNA) mimicking the conserved 5' sequence of SUDV 5' transcripts. The 13-nucleotides length-capped RNA substrate was incubated with both enzymes in presence of radiolabeled methyl donor (3H-SAM). The reaction products were filtered through a DEAE membrane, and the radioactivity transferred to the RNA substrate was determined (filter binding assay). Fig. 2C shows that the MTase-ΔCTD domain was not active on the short capped-RNA substrate whereas the MTase+CTD domain induced the methylation of the capped-RNA substrate.

We also compared RNA binding properties of MTase+CTD and MTase-ΔCTD proteins by fluorescence polarization assay (FP). For this purpose, the capped RNA substrate was labelled at its 3'end by ligation of a pCp-Cy5 residue (GpppN-RNA-pCp-Cy5) and was incubated with increasing concentration of MTase+CTD or MTase-ΔCTD. The interaction of RNA with both proteins was next determined by FP measurement. Fig. 2D demonstrates that MTase-∆CTD does not interact with the capped RNA substrate. Conversely, the MTase+CTD induces an increased FP signal confirming its interaction with the synthetic capped RNA substrate with a Km value of 740 nM. Altogether, these results demonstrated that the CTD domain of SUDV plays a key role in RNA recognition and, in turn, promotes the MTase activity.

Effect of a single mutation of conserved residues in the SUDV CTD on the MTase activity.

We previously demonstrated that the SUDV MTase+CTD induces 2'O methylation of internal adenosines in addition to N7 and 2'O methylations of the cap [START_REF] Martin | The methyltransferase domain of the Sudan ebolavirus L protein specifically targets internal adenosines of RNA substrates, in addition to the cap structure[END_REF]. To further characterize the role of the CTD in the RNA recruitment and unravel the involvement of CTD in the regulation of the different MTase activities carried by MTase+CTD, we mutated conserved residues of filoviruses (Fig. S1) into alanine. The mutant proteins were expressed in bacteria and purified by affinity chromatography on NTA column. Figure 3A shows that the mutated proteins migrate at a molecular weight similar to the WT MTase+CTD. The mutated proteins were also analysed by TSA and show typical denaturation curve of folded proteins (Fig. S3) with Tm value ranging from 50.6 °C to 55.6 °C (Table 1). Table 1 indicates that most mutants display Tm values similar to the wild-type MTase+CTD protein. In contrast, we observed that the mutant R2068A has a slightly decreased melting temperature (ΔTm = -3,7°C) suggesting that this mutation affects the protein stability. We next tested the effect of mutations in the CTD on the RNA binding properties by FP assay. We performed FP measurements using the GpppG-RNA-Cy5 and mGpppGm-RNA-Cy5 SUDV-RNA (Table 1 and Fig. S4). Mutations K2043A and K2189A lead to an alteration of interactions between the substrate (GpppG-SUDV12) and the protein and a loss of recognition for the product m GpppGm-SUDV12 (Kd > 9µM). We next tested the effect of mutations on cap-dependent and cap-independent MTase activities using different RNA substrates. By incubating the WT and mutated MTase+CTD proteins with either GpppGm(Am)-SUDV12, m GpppA(Am)-SUDV12 or m GpppGm-SUDV12 in presence of radiolabelled methyl donor (3H-SAM), we determined the effect of CTD mutations on cap-N7 MTase, cap-2'O MTase and internal 2'O-A MTase activities, respectively (Fig. 3B). We identified three categories of mutations. (i) Mutations leading to the reduction (H2112A and F2113A) or the lost (R2068A, K2118A and K2189A) of all MTase activities. Interestingly, the RNA recognition feature of R2068A, K2118A and K2189A is also strongly impaired and R2068A mutation is also associated to a reduced protein stability compared to the WT (ΔTm -3.7°C). (ii) Mutations decreasing the 2'O MTase activities (cap-2'O MTase and internal 2'O-A MTase), but not the cap-N7 MTase activity (K2043A). This result could be correlated with the alteration of recognition observed in FP, where the interaction with the RNA substrate m GpppGm-SUDV12 is more impaired in comparison to the interaction between mutated protein and substrate GpppG-SUDV12.

DISCUSSION

The mononegavirus L protein plays a key role in replication/transcription and RNA capping. Cryo-EM and biochemical studies have evidenced its organization in 5 main topological domains [START_REF] Liang | Structure of the L Protein of Vesicular Stomatitis Virus from Electron Cryomicroscopy[END_REF] embedding the RdRp, the PRNTase and the MTase activities. Nevertheless, the role of the Connector Domain (CD) and the CTD, framing the MTase domain, is still elusive. In this work, we focused on the regulation role of the CTD on the SUDV MTase activities.

The CTD domain, the less conserved domain of mononegavirus L protein, is enriched in basic amino acid suggesting their involvement in RNA recruitment. This hypothesis is also supported by the X-ray structure of hMPV MTase+CTD showing that the CTD together with the MTase form a kind of conserved RNA binding groove at immediate proximity of the MTase catalytic site (Fig. 1B &C). In this work, we experimentally demonstrated that the CTD domain is essential for the RNA recognition as the recognition of capped RNA is strongly impaired in absence of the CTD (Fig. 2D). The key role of the CTD in RNA recognition was further confirmed by alanine mutagenesis of the CTD. A single mutation of conserved basic and hydrophobic residues was engineered and different mutated proteins were produced and purified. We next demonstrated by TSA that these mutations barely affect the protein folding and stability (Table 1). In addition, all these single mutations decrease the recognition of capped RNAs (Table 1), but less than observed with the MTase-∆CTD protein. Our results indicate that some mutations strongly impair the interaction with m GpppGm-SUDV12, more than with GpppG-SUDV12 (K2043A and K2189A) suggesting that these mutations might regulate in a different way the MTase activity targeting the cap structure or adenosine residues within RNA sequence. Thus, all these results demonstrate the crucial role of CTD conserved residues in MTase substrate recognition. These data agree with the X-ray structure of hMPV and the cryo-EM structure of VSV full length protein indicating that the canonical cap binding domain is virtually absent in these proteins. Conversely, it was suggested that the RNA substrate might accommodate a groove enriched in basic amino acids at the interface between the MTase and the CTD domain.

We next evaluated the effect of CTD deletion and single mutations on the MTase activities. The cap-N7 MTase and the cap-2'O MTase activities were followed by incubating the different mutants with specific RNA substrates previously described [START_REF] Martin | The methyltransferase domain of the Sudan ebolavirus L protein specifically targets internal adenosines of RNA substrates, in addition to the cap structure[END_REF]. In addition, we followed the methylation of adenosines within RNA sequence using a cap-1 RNA substrate ( m GpppGm-SUDV12). We observed that MTase activities are completely abrogated by the deletion of the CTD that hampers the RNA recognition. In contrast, single mutations in the CTD induce more contrasted effects on the SUDV MTase activities. We identified mutations decreasing or abrogating the overall MTase activities (F2113A, H2112A, R2068A, K2118A and K2189A), as well as one mutation uncoupling the N7 MTase activity from the 2'O MTase activities (K2043A), and one mutation affecting internal methylation rather than cap methylations (H2067A). These results are consistent with the RNA binding assay as mutations decreasing efficiently RNA recognition (MTase-∆CTD > K2118A > R2068A) also strongly impair MTase activities. Conversely, we identified one mutant (K2189A), which still binds partially to GpppG-SUDV12, but is barely active, suggesting that such mutation affects the RNA position in the catalytic site. Finally, we also observed that mutants keeping partially the capacity to recognize RNA (K2043A, H2112A, F2113A, H2067A) also keep cap-N7 and/or cap-2'O MTase activities.

Altogether, these results elucidate one function of the CTD domain of SUDV, which regulates RNA substrate recognition and in turns, participate to the fine-tuning of the MTase activity. Interestingly, we recently demonstrated that the SUDV MTase induces in vitro internal adenosine 2'O methylation in addition to the N7 and 2'O methylation of the cap structure [START_REF] Martin | The methyltransferase domain of the Sudan ebolavirus L protein specifically targets internal adenosines of RNA substrates, in addition to the cap structure[END_REF]. On the other hand, other viruses belonging to the Mononegavirales order such as hMPV seem to induce preferentially cap-N7 and cap-2'O methylations beside internal methylation. The different specificities between these enzymes might depend on their capacity to recognize and position the RNA substrate into the MTase catalytic site linked to the CTD recognition properties, which is barely conserved among the mononegavirus families.

We do not know yet whether SUDV induces the methylation of adenosines inside its own genomic or messenger RNAs in infected cells. In addition, the role of such methylations in viral RNA sensing, mRNA translation into viral protein, or RNA encapsidation is still elusive. In contrast to other viral systems such as HIV infection, it was recently demonstrated that internal 2'O methylations of RNAs induced by a cellular MTase recruited by the virus during infection (FTSJ3) limits RNA detection by the MDA5 pathway yielding to reduction of interferon secretion [START_REF] Ringeard | FTSJ3 is an RNA 2′-O-methyltransferase recruited by HIV to avoid innate immune sensing[END_REF]. The identification of mutations in the CTD of SUDV uncoupling the different MTase activities of the L protein opens the way to elucidate the function of internal methylations in filovirus infection and would participate to develop new antiviral strategies leading to an activation of the innate immunity system, especially through RIG-I/MDA5 pathway. μM of 3 H-SAM (Perkin Elmer) in an optimized MTase assay buffer (50 mM Tris-HCl variable pH, 10 mM arginine). Reactions were stopped by a 10-fold dilution in water after 3 h at 30°C. Samples were transferred to DEAE filtermats (Perkin Elmer) using a Filtermat Harvester (Packard Instruments).

Methyl transfer was then evaluated as described before [START_REF] Bollati | Recognition of RNA cap in the Wesselsbron virus NS5 methyltransferase domain: implications for RNA-capping mechanisms in Flavivirus[END_REF]. Briefly, the RNA-retaining mats were washed twice with 10 mM ammonium formate pH 8, twice with water and once with ethanol. They were then soaked with liquid scintillation fluid, allowing the measurement of 3 H-methyl transfer to the RNA substrates using a Wallac MicroBeta TriLux Liquid Scintillation Counter13.

Fluorescence polarization (PF)

Using T4 RNA ligase 1 (20 units; New England Biolabs), cyanine 5-cytidine-5-phosphate-3-(6aminohexyl)phosphate (12.5 µM, Jena Bioscience) was ligated to the 3'ends of the RNA substrates (10 µM) in T4 RNA ligase 1 buffer (NEB), 1 mM ATP (16 °C, overnight). Ligase was removed by RNA precipitation in 3 M sodium acetate supplemented with 1 µg/µL of glycogen (Thermo Scientific). The fluorescent RNA was incubated (5 min at room temperature) with increasing concentrations of CRVI+CTD, in 20 mM Tris pH 8 or pH 8,5 150 mM NaCl, 5% glycerol. Fluorescence polarization (FP) measurements were performed in a microplate reader (PHERAstar FS; BMG Labtech) with an optical module equipped with polarizers and using excitation and emission wavelengths of 590 and 675 nm, respectively. Dissociation constants (Kd) were determined using Hill slope curve fitting (Prism).

Thermal shift assay (TSA)

A mix of 4 µM protein, buffered by resuspension buffer (50 mM Tris pH 8, 150 mM NaCl, 5% glycerol) with 0,02% (v/v) SYPRO Orange dye (Thermo Scientific) was made up to a total volume of 20 µL.

Samples were placed in a semi skirted 96 well PCR plate (BioRad), sealed and heated in an Mx3005p qPCR machine (BioRad) from 25 to 95°C at a rate of 1 °C.min -1 . Fluorescence changes were monitored with excitation and emission wavelengths at 492 and 610 nm respectively.

dependent and independent 2 '

 2 O MTase activities inducing cap 2'O methylation and internal 2'O methylation limiting the viral RNA sensing by RIG-I and MDA5 respectively [25][26][27][28]. In this work, we addressed the functional role of the SUDV CTD domain by comparing the MTase activity of MTase+CTD and MTase-∆CTD constructs. After expression and purification of various SUDV MTase constructs, we demonstrated the critical role of the CTD domain for RNA substrate recognition and thereby controlling the MTase activity. The importance of CTD basic residues in RNA recognition and MTase activity was confirmed by directed mutagenesis and we identified single mutations uncoupling the different activities carried by SUDV MTase. These results highlight the key regulation role of the CTD domain in the RNA methylation process and raise the possibility to identify compounds inhibiting specifically these viral MTases by targeting the CTD.

  (iii) Mutation reducing more specifically the internal 2'O-A MTase activity (H2067A) than the cap-MTase activities. Altogether, these results demonstrate that mutations in the SUDV CTD of L protein regulate the RNA recognition and the MTase activities. We highlighted mutations uncoupling cap methylations (cap-2'O MTase and cap-N7 MTase) from internal 2'O-A methylation or N7 methylation from 2'O methylations (cap-2'O MTase and internal 2'O-A MTase) demonstrating that the CTD plays a key role in the fine tuning of the N7 and 2'O MTase activities of the L protein.

Acknowledgements

The research was made possible by funding from the European Union Seventh Framework Programme (FP7/2007-2013) under SILVER grant agreement no. 260644 and from an MRC grant (MR/L017709/1). The authors thank Diamond Light Source for beamtime (proposal mx8423), and the staff of the MX beamlines for assistance. This work was supported by a Wellcome Trust administrative support grant 075491/Z/04.

Methods

Cloning and expression

Codon-optimized CRVI-CTD synthetic gene (Biomers) was designed for expression in c2566 pRARE2. N-terminal truncated construct of CRVI-CTD was synthetized by cloning the sequence of interest from full-length CRVI-CTD into a pDEST14 plasmid with Gateway system for expression in c2566 bacteria. After transformation, bacteria were cultured at 30°C until a OD of 0.6, and induced overnight with 500 mM IPTG (Sigma) at 17°C. Finally, cultures were centrifuged at 8,000 xg for 10 min at 4°C using a Sorval Lynx 6000 centrifuge (Thermo), and pellets were conserved at -80°C until purification.

Purification

Pellets were thawed on ice, and lysed in a 10 times final culture OD volume of optimized lysis buffer (50 mM Tris pH 8, 150 mM NaCl, 5% glycerol, 30 mM imidazole, 1 mM PMSF, 100 µg/mL lysozyme, 1 µg/mL DNase, 0.1% Triton X100) supplemented with BugBuster (Merck Millipore) for MTase+CTD.

After clarification (18,000 xg, 30 min, 4 °C), the lysate was incubated with CoNTA resin (Thermo; 0.5 mL/L culture) for 1 h at 4 °C, with gentle shaking. The beads were transferred to a 25-mL column and washed with 2 x 20 mL of buffer W1 (50 mM Tris pH 8, 1 M NaCl, 5% glycerol, 30 mM imidazole) and 2 x 10 mL of buffer W2 (50 mM Tris pH 8, 150 mM NaCl, 5% glycerol). Protein was eluted in buffer E (50 mM Tris pH 8, 150 mM NaCl, 5% glycerol) supplemented with 150 mM imidazole in the case of MTase-ΔCTD and 1 M arginine for CRVI+CTD. Imidazole was removed by gel filtration using a Superdex S75 16/60 (GE Healthcare) for MTase-ΔCTD. Finally, proteins were concentrated on Amicon Ultra (EMD Millipore), and conserved in 50% glycerol at -20°C.

Synthesis of RNA substrates

RNA sequences were chemically synthesized on solid support using an ABI 394 oligonucleotidessynthesizer. After RNA elongation with 2'-O-pivaloyloxymethyl phosphoramidite ribonucleotides, 10-11 and 2'-O-methyl phosphoramidite ribonucleotides (Chemgenes, USA), the 5'hydroxyl group was phosphorylated and the resulting H-phosphonate derivative 12 was oxidized and activated into a phosphoroimidazolidate derivative to react with guanosine diphosphate (GpppRNA) 13- 14 . After deprotection and release from the solid support, GpppRNAswere purified by IEX-HPLC and validated to be >95 % pure by MALDI-TOF spectrometry. N7-methylation of the purified GpppRNA was performed enzymatically using N7-hMTase 13-14 .

MTase activity assays

To evaluate methyltransferase activities, a radioactive test was set up by mixing 4 μM of protein with 1 μM of purified and validated synthetic RNAs (see all RNAs in Table S1), 10 μM of SAM and 0.5