Ludovic Claudepierre

Rémi Douvenot
email: remi.douvenot@recherche.enac.fr

Alexandre Chabory

Christophe Morlaas

A Method for Computing the VOR Multipath Error -Comparisons with In-flight Measurements

Keywords: VOR error, parabolic equation, physical optic, hybridisation, multipath, electromagnetic modelling, wind turbines

This work presents the validation of a modelling method for the characterisation of the VOR error in the presence of large scatterers. More specifically, wind turbines are considered because their increasing number makes their modelling a topical problem. Here, the electromagnetic perturbation due to the wind turbines is addressed by a method hybridising parabolic wave equation and physical optics. The VOR receiver is considered such that the overall method gives an estimation of the VOR error. This method is validated by comparisons with in-flight measurements performed in the same location before and after the erection of a wind farm.

I. INTRODUCTION

For air transport navigation, VHF Omnidirectional Range (VOR) systems are crucial: they give to aircraft their direction with respect to North. However, this information can be deteriorated by nearby obstacles, which yield multipaths. Signals coming from unwanted directions are processed by the onboard VOR receiver, which modify the on-board displayed bearing.

Today, the air transport authority delivers approvals for the construction of buildings in the vicinity of VORs [START_REF] Icao | European guidance material on managing building restricted areas[END_REF]. Their rules are very restrictive, and wind turbine industrialists request looser regulations. The method presented here could be used for delivering case-by-case approvals, for wind turbines and other buildings.

The simulation method presented here follows three objectives: 1/ the simulation can be performed on a standard computer in reasonable time, 2/ the output of the model is the VOR error, which renders the method easily interpretable for the air transport community, and 3/ the model gives a prediction of the VOR error with a good accuracy. The VOR error expression requires both the amplitude and the phase of the direct and multipath signals [START_REF] Odunaiya | Calculations and analysis of signal processing by various navigation receivers architectures[END_REF]. Therefore, a good balance between computation time and accuracy is required with both amplitude and phase of the multipaths.

Many types of simulation methods have been proposed during the last decade. The uniform theory of diffraction has been proposed [START_REF] Fontan | Assessing multipath induced errors in VOR systems using ray-tracing techniques in combination with detailed terrain databases[END_REF], [START_REF] Gonzalez | Analysis of the scattering field of wind turbine using rigorous and asymptotic techniques[END_REF]. This method is accurate but is memory intensive for complex scenes. More recently, models involving parabolic equation (PE) and physical optics (PO) have been developed [START_REF] Morlaas | Propagation model for estimating VOR bearing error in the presence of windturbines -hybridation of parabolic equation with physical optics[END_REF], [START_REF] Casanova | Modeling windmill interference via hybrid parabolic equation and equivalent current approximation methods[END_REF]. The work presented here is based on these latter. Models based on the RCS [START_REF] Casanova | Wind farming interference effects[END_REF], [START_REF] Vega | Software tool for the analysis of potential impact of wind farms on radiocommunication services[END_REF] or using the knifeedge model [START_REF] Fernandes | Assessment of wind turbines generators influence in VOR aeronautical navigation systems[END_REF] have also been developed. They propose fast computations of the field scattered by wind turbines. However, the simplifications underlying these models are important, and the VOR error cannot be obtained since the phase is not available. These methods are interesting for rough calculations, but cannot give precise results.

Finally, the method preferred here is the hybridisation between parabolic equation (PE) [START_REF] Dockery | An improved impedance-boundary algorithm for Fourier split-step solutions of the parabolic wave equation[END_REF] and physical optics (PO) [START_REF] Lee | Fourier transform of a polygonal shape function and its application in electromagnetics[END_REF] to calculate the multipaths. The hybridisation method is detailed in Claudepierre et al. [START_REF] Claudepierre | A deterministic VOR error simulation method -Application to windturbines[END_REF]. The main advantages are to considered the antenna radiation pattern and the relief between the VOR antenna and the wind turbines. Moreover, the ground composition is taken into account. The major drawback is that 3D effects are neglected. The relevancy to use PO for the simulation of the field scattered by the wind turbines has been demonstrated by Claudepierre et al. [START_REF]Influence of the lightning protection of blades on the field scattered by a windturbine[END_REF].

The simulation method is briefly presented in section II. In section III, the measurement campaign is introduced, and the choice of the relevant data to validate the simulator is presented. Finally, in section IV, the simulator is validated by comparisons to in-flight measurements and used for a deeper analysis of the scenario.

II. SIMULATION METHOD

The modelling of the computation follows four steps, as illustrated in Figure 1: first, the field is calculated on the scatterers using a PE method. One PE simulation is performed for each scatterers (except if they are aligned). Second, the currents on the scatterers are calculated from the PE incident field with an hybridisation process detailed in [START_REF] Claudepierre | A deterministic VOR error simulation method -Application to windturbines[END_REF] that takes into account the dielectric parts of the wind turbines. Third, the field radiated by the currents is computed on the receiver by PO. Finally, a VOR error given below is applied by considering the multipaths and the direct field analytically computed.

Note that the incident field on each scatterer is calculated by neglecting their mutual interactions. Indeed, in the case of wind turbines, their mutual influence has been demonstrated to be negligible at the VOR frequency [START_REF]Assessment of the shadowing effect between windturbines[END_REF].

The block diagram of the method is given in Figure 2. From the geometry of the scenery, the antenna properties, and the composition of the scatterers, the VOR error can be calculated at the required observation points. One PE calculation and one PO simulation is performed per scatterer. Finally, following [START_REF] Odunaiya | Calculations and analysis of signal processing by various navigation receivers architectures[END_REF], the VOR error is calculated by combining the scattered

ε c = tan -1     - N n=1 a n cos(θ n ) sin(ϕ n ) 1 + N n=1 a n cos(θ n ) cos(ϕ n )     , (1)
where a n , θ n , and ϕ n are the amplitude, the phase, and the azimuth of the n th multipath, respectively. These values are normalised by the amplitude, the phase, and the azimuth of the direct signal. This equation is valid for the conventional VOR. Another expression depending on the same parameters is available for the Doppler VOR [START_REF] Odunaiya | Calculations and analysis of signal processing by various navigation receivers architectures[END_REF].

III. MEASUREMENT CAMPAIGN

A. Description of the measurement campaign

The case of nine wind turbines (ENERCON E-70) erected at about 5 km from the conventional VOR of Boulogne-sur-Mer (France) is studied (Figure 3). The frequency is 113.8 MHz. The relief between the VOR and wind turbines is accounted and a dielectric ground (ǫ r = 25, σ = 0.02) is considered.

The wind turbines, represented by yellow circles in Figure 3 wind turbine model has been validated by Claudepierre et al. [START_REF]Influence of the lightning protection of blades on the field scattered by a windturbine[END_REF].

The exact aircraft positions during the flight measurements have been given by an on-board GPS receiver. Moreover, the exact orientations of the hubs during the measurement are also available.

The flight inspections have been performed along six radials (blue lines in Figure 3) at three different times: without wind turbines in May 2009, with only the masts erected in July 2012, and with the complete wind turbines in November 2012.

B. Selection of the relevant measurements

A first work before validating the simulation method with measurements is to assure the relevancy of the comparison. In this paper, this process is illustrated with two radials ϕ = -10 • and ϕ = 78 • .

The VOR error for the three measurements campaigns on the two radials ϕ = -10 • , and ϕ = 78 • are plotted with respect to distance in Figures 4 and5, respectively.

An error already exists in the first measurement campaign without the wind turbines. This error due to the surrounding environment is not negligible: up to 0.6 • of error on ϕ = 78 • and up to 1 • of error on ϕ = -10 • . The mean and variance of the absolute error are given in table I for the three measurement campaigns corresponding to no wind turbines, only the masts, and the complete wind turbines, respectively.

On the radial ϕ = 78 • , no significantly change in the mean and variance of the error is observed. It implies that the VOR error induced by the wind turbines is hidden by the surrounding disturbances. Comparing simulations with these data is irrelevant. By contrast, on the radial ϕ = -10 • , the mean and the variance of the error significantly increase with the erection of the wind turbines. Therefore, only the radial ϕ = -10 • can be used for the validation. This result can be geometrically interpreted from Fig. 3. On the radial ϕ = -10 • , the wind turbine impact is more significant because the aircraft overflies the wind farm, which is not the case on the radial ϕ = 78 • .

IV. RESULTS

A. Validation of the simulation method

The simulation method is validated by comparison with the VOR error measured along the radials ϕ = -10 turbines are plotted in Figure 7. For a distance greater than 20 NM, there is a good agreement between the measurements and the simulations. Below 20 NM, the disturbance due to the surroundings prevails, see Figure 4. Some inaccuracies still remain. Some fast variations are not correctly modelled, and a sign difference sometimes appear. It may come from the VOR receiver which is in reality more complicated than equation (1) used for the simulation. However, the simulation gives an accurate approximation of the error due to the wind turbines. Thus, the simulator is considered validated. More validations are presented in [START_REF] Claudepierre | A deterministic VOR error simulation method -Application to windturbines[END_REF].

B. Quantification of the contribution of each part of the wind turbines

To obtain more information about the influence of wind turbines on the VOR error, the simulation method is able to isolate the contributions of every part of the wind turbines in a post-processing step. For instance, the simulation result presented in Figure 7 is detailed in this section to give a better understanding of the scattering phenomena involved.

The respective contributions of the wind turbines blades, hubs, and masts to the VOR error are plotted in Figure 8 with respect to distance. The blades and the hubs add less than 0.1 • to 0.3 • error each. This is very low compared to the total error. So the difference between the two measurement campaigns (masts only and complete wind turbines) is mainly due the differences between the two flight trajectories during the two measurement campaigns.

This also shows that to study the effect of additional wind turbines in the same zone, modelling the masts only would give a good approximation of the total error. Note that this conclusion is only valid for this type of wind turbine and on that trajectory. More simulations would be necessary to extend the validity of this statement.

In conclusion, the simulation method presented in this paper allows to discriminate the contributions to the VOR error of the different scatterers and their inner parts. For instance, it is shown here that far from the VOR, the masts are the main contributors to the error.

V. CONCLUSION

The simulation method presented here to predict the VOR error due to large scatterers is based on an hybridisation between parabolic equation and physical optics. A processing step to take into account the VOR receiver is also used. In this paper, this method has been applied to wind turbines. A validation with in-flight measurements around a conventional VOR in Boulogne-sur-Mer and nine wind turbines nearby has been performed.

The simulator ability to compute the VOR error due to complex obstacles such as wind farms with a great accuracy has been shown. This simulator can also discriminate the contributions of the different obstacles and their different parts to the VOR error.

The simulator can be used before flight inspections. A simulation of the VOR error all around the VOR indicates the radials to consider if one wants to measure the maximum error during a flight inspection. It can also be used beforehand to approve or turn down a wind farm building project. Note that the reasonable simulation time allows to study several scenarii for a single project.

Figure 2 :Figure 1 :

 21 Figure 2: Synoptic of the simulator.

  , are composed by a metallic conic mast (height = 98 m, top diameter = 2 m, and bottom diameter = 7.5 m), a hub modelled by a metallic rectangular parallelepiped (4 × 11 × 4 m), and blades modelled by 2 parallel dielectric slabs (length = 35 m, maximum width = 1 m, thickness = 44 mm, and ǫ r = 5) with 0.65 m spacing. A metallic rod representing the lightning protection (diameter = 10 mm) is placed between them. This

Figure 3 :

 3 Figure 3: Wind farm at 5 km from the VORC of Boulognesur-Mer (France) and measurement radials (blue lines).

Figure 4 :Figure 5 :

 45 Figure 4: VOR error measured during the flight inspection campaigns without the wind turbines (blue), with the masts only (red), and with the complete wind turbines (green). The radial is ϕ = -10 • .

Figure 6 :

 6 Figure 6: VOR error due to the wind turbine masts measured (red) and simulated (black) on the radial ϕ = -10 • .

Figure 7 :

 7 Figure 7: VOR error due to the wind turbines measured (red) and simulated (black) on the radial ϕ = -10 • .

Figure 8 :

 8 Figure 8: Simulated VOR error ( • ) due to the wind turbine hubs (blue), blades (black), and masts (green) on the radial ϕ = -10 • .

Table I :

 I • . The simulated and measured VOR errors due to the masts only are plotted in Figures 6 with respect to distance. In the same way, the simulated and measured VOR error for the complete wind Mean and variance of the absolute value of the VOR error for the three flight inspection measurements.

	Windturbines	None	Masts only	Complete
	Rad. 78 •	Mean Variance	0.26 0.028	0.37 0.069	0.27 0.056
	Rad. -10 •	Mean Variance	0.31 0.056	0.40 0.09	0.53 0.25

ACKNOWLEDGMENT

The authors would like to thank the French Air Navigation Service Provider (DSNA/DTI) for the funding.


Engang Wang
email: egwang@mail.neu.edu.cn

Progress of some techniques on electromagnetic metallurgy

Keywords: Fluid flow, V-EMBr, EMS, Superalloy, Solidification, Continuous casting

Some of the electromagnetic metallurgy techniques are introduced in this paper. One is the vertical electromagnetic brake (V-EMBr) which have some different effects with Ruler-EMBr and FC-Mold on the control of fluid flow in continuous casting mould. Another is the vertical continuous casting of nickel-based superalloy with high frequency electromagnetic field in the mould and lower frequency electromagnetic stirring (EMS) in second cooling zone of pilot continuous casting machine, that have apparently improved the surface and internal quality of superalloy billets.

Introduction

During the past decades, the Division of Electromagnetic Metallurgy on Steel, Key Laboratory of Electromagnetic Processing of Materials (Ministry of Education), Northeastern University, P. R. China, has focused on the research of EPM field. The different electromagnetic field, such as high frequency, lower frequency electromagnetic fields and DC field, were respectively applied to investigate their influence on the fabrication processing, qualities and properties of materials. In recent years, we also have some progress about the high magnetic field on the fabrication of Fe-based magnetic nanocomposites [START_REF] Lin | Proceeding of 8th Pacific Rim International Congress on Advanced Materials and Processing (PRICM8)[END_REF] , the solidification, precipitation and properties of Cu based composites [START_REF] Lei | [END_REF][3][4] , and the amorphous state of bulk metallic glass during annealing process [5] . In this paper, a new pattern of electromagnetic brake (EMBr) device is introduced, which was firstly proposed in EPM 2009 conference [START_REF] Wang | Proceeding of the 6th Int. Conference on Electromagnetic Processing of Materials[END_REF] and it is called the vertical electromagnetic brake (V-EMBr). The metallurgical effects of V-EMBr to control the fluid flow in continuous casting mould are presented by numerical simulation methods. On the other hand, the metallurgical effects of continuous casting of nickel-based superalloy, which is generally fabricated by ingot casting, are presented with high frequency electromagnetic field in the mould and lower frequency electromagnetic stirring (EMS) in second cooling zone of pilot continuous casting machines.

Effect of V-EMBr in continuous casting mould

In continuous casting of steel, some patterns of EMBr, such as EMBr Ruler and FC-Mould, are used to control the fluid flow in the mould, which generally have the level patterns of magnetic poles in the wide side of the mould (Fig. 1). However, these magnetic poles are fixed in a certain height of the mould, which generally couldn't fit the intermittently change of the depth and angle of submerged entry nozzle (SEN) during the casting process, so it will reduce the effect of EMBr and even bring some negative effect. The proposed new pattern of V-EMBr, which has two pairs of magnetic poles vertically installed on the narrow sides of the mould from the meniscus to the impinge point of melt from the outlet of SEN (Fig. 2). The important characteristic of V-EMBr is that it covered the meniscus and the initial solidifying shell near the narrow sides of mould which generally brings the surface fluctuation, powder flux engulfment, and subsurface defects in slabs owing to the bubbles and inclusions are captured by initial solidifying shell, moreover their control effect to fluid flow is not affected by the change of melt surface level, the outlet depth and angle of SEN. The numerical simulation results show that the flow of molten steel in the mold is enhanced as the casting speed increasing, it presents that the vortex center in lower recirculation zone is moved down (Fig. 3a), and the flow velocity of free surface is gradually increased (Fig. 4a), and the velocity and turbulent kinetic energy near the impact points of molten steel jet from the SEN in the narrow of mould is also gradually increased (Fig. 5a). That mean the chance is increase to have the powder flux engulfment in the mold surface and involve inclusions and bubbles to the slabs. After applying V-EMBr under different casting speed in continuous casting process, the fluid flow in upper recirculation zone is remarkable reduced and the vortex center in lower recirculation zone is obviously moved upward (Fig. 3b). The flow velocity of molten steel in the free surface of mould is significantly reduced (Fig. 4b). The flow field near the impact points in the narrow side of mould forms a regular flow field moved upward and downward (Fig. 5b) contrast to irregular divergence flow without EMBr (Fig. 5a), it proves that proposed V-EMBr can significantly reduce the impact of molten steel jet to the narrow side and the meniscus of mould, and further prevent the surface fluctuations and powder flux engulfment in the mold. The numerical simulation results proved that the V-EMBr device can play a effect control to the molten steel in the continuous casting mould. The numerical simulation results show that when the submergence depth of SEN changes from 200mm to 170mm, which is the conventional operation in the continuous casting, the flow velocity of molten steel in the free surface is increased, especially the region which velocity is over 0.4m/s is also increased (Fig. 6a). That mean the chance to bring the powder flux engulfment near the meniscus of mold is increased. After applying the V-EMBr (B ave is 0.3T), the flow velocity of molten steel in the free surface is not over 0.4m/s as the SEN depth is 170mm and 200mm (Fig. 6b), that is helpful to prevent the surface fluctuations and powder flux engulfment in the mold, and it also proved that the control effect to fluid flow of V-EMBr is a few affected by the depth and outlet angle of SEN.

Fabrication of superalloy strands in continuous casting with electromagnetic fields

Incoloy 800 is a kind of Ni-based superalloy, which is widely used in steam generator, nuclear reactor, petrochemical industrials and solar energy systems due to its high strength at high temperatures and the ability to resist oxidation and high-temperature corrosion. However, Incoloy 800H strands is traditionally fabricated by ingot casting, which not only brings a higher production cost and also has a series of quality problems such as serious segregation, coarse grain and serious cracks. In this paper, the fabrication of rectangular strands and billets of Incoloy800 superalloy were carried out in a vertical pilot continuous caster in Northeastern University, China. To improve the quality of Incoloy800 strands, the electromagnetic fields are applied, one is the high frequency electromagnetic field in mould to improve its surface quality and another is electromagnetic stirring (EMS) in the secondary cooling zone to control the solidifying structure and segregation of elements in the strands. Fig. 7 shows the solidification macrostructure of Incoloy800 slabs (225mm×100mm) casted in the vertical pilot continuous caster without and with EMS. It is shown that the solidification structure of Incoloy800 transverse section without EMS is mainly consists of the unidirectional column grains and a few fine equiaxed grains in the surface of strands, and the serious transgranular phenomenon is appeared in the center of strands (Fig. 7a). When EMS is applied, the transgranular phenomenon is disappeared and the equiaxed ratio in the transverse section of Incoloy800 strands is increased to 31% (Fig. 7c) and the average equiaxed grain size in the central area of slabs is apparently refined. Fig. 8 shows the macrostructure of Incoloy800HT billets (100mm×100mm) in longitudinal section casted in the vertical pilot continuous caster without and with EMS. It is also shown the transgranular phenomenon no EMS is applied in the center of strands and the small shrinkage or porosity and internal crack between the column grains (Fig. 8a). When the EMS is applied, the growth of columnar dendrites were apparently restrained and many of fine equiaxed grains are formed in the center of billets, the transgranular phenomenon, shrinkage and internal crack are eliminated, the internal quality of the Incoloy800HT billet is remarkably improved. In order to improve the serious surface defects of Incoloy800HT billets (100mm×100mm) such as oscillation mark (Fig. 9a), etc, a high frequency electromagnetic field (20Hz) was applied near the meniscus zone of the mould during continuous casting process. The experimental results show that appearance of Incoloy800HT billets is significantly improved with high frequency electromagnetic field, the oscillation marks are suppressed and surface defects can be eliminated.

Fig. 1

 1 Fig.1 EMBr-Ruler and FC-Mold Fig.2 Vertical Electromagnetic Brake (V-EMBr)

Fig. 5

 5 Distribution of velocity and turbulent kinetic energy in narrow side of mold with and without V-EMBr (B ave =0.3T) Depth of SEN = 170mm Depth of SEN = 200mm a) Without EMBr b) With V-EMBr Fig.6.The velocity distribution of molten steel on the free surface of mold with and without V-EMBr (B ave = 0.3T)

  EMS Fig. 7 Macrostructure of continuous casting Incoloy800 slabs (225mm×100mm)

Conclusions 1 )

 1 A new pattern of electromagnetic brake (so-called V-EMBr) is proposed and the numerical simulation results show that the V-EMBr can significantly reduce the impact of molten steel jet to the narrow side and the meniscus of mould and the velocity of molten steel in the free surface of mould, that further prevent the surface fluctuations and powder flux engulfment in the mold. Moreover the control effect of V-EMBr to fluid flow is few affected by the depth and outlet angle of SEN. 2) By applying of high frequency electromagnetic field and electromagnetic stirring, the internal and external quality of Incoloy800HT billets could be improved greatly in the continuous casting. a) No EMC b) EMC Fig. 9 Macrostructure of billet a) No EMS b) EMS Fig. 8 Macrostructure of Incoloy800HT billets the longitudinal section

Acknowledgment

The author gratefully acknowledge all of colleagues and students in the Division of Electromagnetic Metallurgy on Steel, Key Laboratory of EPM, Northeastern University, P. R. China. This work was financially supported by the National Nature Science Foundation of China (No.51574083), the Key Project of the Ministry of Education of China (No. 311014), the Program of Introducing Talents of Discipline to Universities (the 111 Project of China, No. B07015) and the National High-Tech Research and Development Programme of China (No. 2001AA337040)


Philippe Brocard
email: brocard@recherche.enac.fr

Daniel Salos

Olivier Julien

Mikael Mabilleau

Performance Evaluation of Multipath Mitigation Techniques for Critical Urban Applications Based on a Land Mobile Satellite Channel Model

Keywords: multipath mitigation, urban environment, error model, rail users I

The use of GNSS for critical terrestrial applications requires adapted integrity monitoring algorithms for either safety or liability issues. Prior to the design of integrity monitoring algorithms, it is necessary to characterize the nominal error model of every error source affecting the pseudorange measurements. This paper focuses on describing a way to calculate the contribution of the multipath to the nominal error model by a set of simulations based on an urban channel model and a realistic GNSS tracking loop simulator. To improve the accuracy of the positioning, multipath mitigation techniques are usually implemented in the receivers. This paper assesses the performances of the Narrow Correlator, the Double-Delta and the A Posteriori Multipath Estimation techniques in an urban environment and provides their multipath residual nominal error model as a function of the satellite elevation angle. A relationship between quality monitoring indicators and the multipath estimation technique (APME) is also discussed. The results of this paper are obtained by simulating an actual urban railway environment and therefore can be taken as inputs for the design of GNSS integrity monitoring algorithms for rail users. The same approach can be used for the determination of the nominal error model for other terrestrial application subject to simulating an adapted channel.

INTRODUCTION

Due to the modernization of current GNSS, the development of new GNSS like Galileo and Beidou, and the improvement of Satellite Based Augmentations Systems (SBAS), GNSS is a promising alternative for terrestrial rail and road position monitoring systems that require costly ground infrastructures. In the rail domain, the European Train Control System (ETCS) [1] is the automatic signaling, control and train protection system that is currently being deployed for an improved interoperability in Europe. In ETCS level 2 and 3 the vehicles have to self-monitor their position and speed, based on a combination of radiobeacons (Eurobalises) installed along the railway that provide reference positions, and odometry. Several European projects have investigated the possible introduction of GNSS in ETCS to reduce the number of radiobeacons that currently has to be at least one per 2.5 km [1]. Among these projects, it can be mentioned INTEGRAIL (ESA, 2001), GRAIL (6 th FP, GSA, 2005), GRAIL-2 (7 th FP, GSA, 2010), 3INSAT (ESA, 2012) and SATLOC (7 th FP, GSA, 2012). The memorandum of understanding between the European Commission and the European Rail sector Association (ERA) of 2012 states that "GNSS can play a major role in the rail sector for fleet management and rail safety". To speed up the standardization process for the use of GNSS in ETCS, the Union Industry of Signalling (UNISIG) created a working group called Satellite Positioning Working Group in cooperation with the European GNSS Agency (GSA). The main challenge for the introduction of the GNSS in train control is the very low Tolerable Hazard Risk (THR), which, for the whole signaling system, shall not be over 2.10 - 9 /h to fulfill the SIL4 requirement [1]. Therefore, the tolerable uncertainty on the GNSS position solution shall be even lower, as referred in [START_REF] Mocek | Galileo Safety-of-Life Service Utilization for Railway Non-Safety and Safety Critical Applications, AA(Railway Infrastructure Administration, TUDC -Laboratory of Intelligent Systems)[END_REF] which mentions a THR between 1.10 -10 /h and 1.10 -11 /h, [START_REF] Mertens | LOCPROL: A low cost Train Location and Signalling system for[END_REF] that mentions a THR between 1.10 -11 and 1.10 - 12 /h, [START_REF] Zheng | The Effects of Railway Track Database Quality on the Performance of Tightly Coupled GNSS/Track Database Train Positioning System[END_REF] that sets the Integrity risk at 4.10 -12 /h, and [START_REF] Rispoli | Recent progress in application of GNSS and advanced communications for railway signaling[END_REF] THR = 1.10 -11 /h. These operational requirements are far more stringent than the performances of the existing standalone GNSS, therefore the use of an adapted integrity monitoring concept is required. In the road domain, GNSS is one of the recommended technologies in the EU directive for Electronic Toll Collection (ETC) [START_REF] Parliament | of 29 April 2004 on the interoperability of electronic road toll systems in the Community[END_REF] and GNSS based ETC systems already exists for heavy good transportation in Germany (Toll Collect) and Slovakia (MYTO). The norm CEN ISO TS 17444-1 provides different metrics that allow to define and assess the performance of ETC systems in terms of overcharging and undercharging. GNSS integrity monitoring systems can be used in ETC in order to enhance the liability of the charging system. GNSS based ETC and Train Control are two examples of terrestrial applications that require their own integrity monitoring concepts to augment the GNSS and to fulfill their operational requirements. The design of such algorithm shall be adapted to the system and take into account the available sensors, augmentations systems, and the integrity requirements to real-time monitor the protection levels. Prior to the design of any integrity monitoring algorithm, it is necessary to characterize the nominal error model of every error sources affecting the pseudorange measurements. For both applications, the vehicles are likely to operate in urban area where the multipath is a significant contributor to the error measurement. This paper aims at characterizing the nominal code pseudorange multipath error standard deviation in urban environments for rail users with dual constellation GPS/Galileo single-frequency L1/E1 receivers. To do so, the characterization is based on the use of the wideband Land Mobile Satellite (LMS) channel model developed by the German Aerospace Center (DLR). This model is the reference wideband model for the ITU (International Telecommunication Union). Based on an urban railway channel scenario, the code measurements' statistics and residual errors are investigated using standard tracking loops as well as several multipath mitigation techniques. The channel simulated in this paper represents a typical urban railway environment, therefore the provided results apply to rail users. The same approach could be used to determine the nominal pseudorange error model for ETC by modifying the scenario.

II. DESCRIPTION OF SIMULATION TOOLS

The approach consists in connecting a simulator of urban channel that generates time series of signal with a realistic GNSS receiver simulator.

A. LMS channel model

The Land Mobile Channel (LMS) model [START_REF] Lehner | A Novel Channel Model for Land Mobile Satellite Navigation[END_REF] which can be downloaded for free, is a generative wideband urban channel model based on a measurement campaign conducted in Munich in 2002. Compared to a narrowband model such as the one developed in [START_REF] Prieto-Cerdeira | Versatile two-state land mobile satellite channel model with first application to DVB-SH analysis[END_REF], a wide band model assumes that the coherence time of the channel is negligible compared to the chip duration of the Pseudo Random (PRN) code. A wideband model is preferred in this paper in order to realistically simulate the impact of multipaths on the receivers tracking loops and therefore on the pseudorange measurement. The model simulates the received signal for a fixed satellite whose relative position to the user is characterized by its azimuth and its elevation. The used LMS model simulates the channel for a vehicle moving along a street with predefined parameters. The model generates differently the direct satellite to receiver Line-of-Sight (LOS) and the multipaths. The properties of the LOS are deterministically generated by simulating realistic physical phenomena. The diffraction of the LOS over the edges of the buildings, the shadowing induced by trees and the diffraction over lamp posts are considered. The echoes are statistically generated based on several distributions deduced from the measurement campaign. The number of reflections, the geometric distribution of the reflectors (from which the phase and delay are deduced), the power distribution of the echoes, their lifespan (space life duration) and the distribution of the time variation of their amplitude (fading) are used to generate the echoes.

At every sample ( 𝑘 ), each received LOS and echo is characterized by its amplitude 𝛼 𝑛 (𝑘), code delay 𝜏 𝑛 (𝑘) and carrier phase 𝜑 𝑛 (𝑘). The model of the received signal is:

𝑠(𝑘) = ∑ 𝛼 𝑛 (𝑘)𝑑(𝑘𝑇 𝑠 -𝜏 𝑛 (𝑘))𝑐(𝑘𝑇 𝑠 -𝜏 𝑛 (𝑘))cos (𝜑 𝑛 (𝑘)) 𝑁 𝑛=1 (1)
Where c is the spreading code, d is the GNSS data and 𝑇 𝑠 is the sampling period. A slight modification was added to the original LMS model in order to consider the dynamic of the mobile in the generation of the code delay. In the original model, the delay of the LOS remains null even when considering a fixed satellite and a moving vehicle. To model the impact of the dynamic on the delay, a delay variation was added to both the LOS and echo to LOS relative delay. This phenomenon is referred to as code Doppler which has the same origin as the carrier Doppler. The additional delay variation for two consecutive samples k and k+1 is expressed by:

𝜏(𝑘 + 1) -𝜏(𝑘) = 𝑣 𝑟 𝑐 . 𝑇 𝑠 (2)
Where 𝑣 𝑟 is the derivative of the distance between the receiver and the satellite.

B. GNSS tracking simulator

A realistic GNSS receiver simulator, developed by ENAC and referred to as GeneIQ, uses the outputs of the LMS channel model, consisting of the times series containing the amplitude, code delay, carrier phase and Doppler of each LOS and reflections. This simulator is based on the modeling of the receiver correlator outputs which expression is given in Eq.

(3). The Doppler is not an explicit output of the channel model, even if it is used by the model for the phase generation. As a consequence, the echoes' Doppler has been calculated differentiating their phase measurements on two consecutive epochs. The correlation between each component of the signal and the local replica are cumulated to form a composite correlator output. The phase is assumed to vary linearly over the integration interval. The expression of the In-phase (I) and Quadra phase (Q) correlator outputs are:

𝐼 𝑃 (𝑘) = ∑ 𝛼 𝑛 (𝑘) 2
Where  K cc is the autocorrelation function of the spreading code.  𝜏(𝑘) is the delay of the local spreading code  𝑓 ̂(𝑘) is the frequency of the local carrier  𝜑 ̂(𝑘) is the phase of the local carrier

C. Validation of the LMS/Receiver coupling

The simulation platform, i.e. the combination of the LMS channel model followed by the GeneIQ receiver simulator, has been validated. In a first validation step, the LOS component of the generated received signal has been isolated from the echoes and has been processed by the simulator. This allows to check the right generation of the dynamic of the LOS signal by comparing the behavior of the tracking loops with their theoretical behavior. In this case, the received signal has been simulated for a vehicle with a constant velocity and a fixed satellite position. When the satellite to receiver distance is linearly varying, it generates a constant Doppler shift on the carrier phase and code delay, which leads to a biased estimation of the delay and phase by first orders DLL and PLL. The code delay and carrier phase estimated by higher orders tracking loops are not biased due to a constant Doppler shift. This paper only describes the validation process of the DLL. The theoretical expression of the steady state error of a first order DLL is the following (in chips) [START_REF] Van Dierendonck | Global Positioning System: Theory & Applications[END_REF]:

E ∞ = ∆𝑓 𝑐 4𝐵 𝐿 = 𝑣 0 𝑇 𝑐 4𝐵 𝐿 𝑐 (4)
Where  v 0 is the radial speed of the receiver in m/s  B L is the one-sided loop bandwidth of the loop filter in Hz

 ∆f c = v 0 c.T c
is the code Doppler in chips/s,  T c is the chip duration in second In the chosen scenario, the velocity of the vehicle is set to 20 km/h. The satellite has an azimuth of 40° and an elevation of 50°. The radial velocity is v 0 = 20/3.6 * cos(50) * cos (45) = 2.525 ms -1 . The loop bandwidth is set to 1 Hz which is a typical value for a DLL in a GNSS receiver. Using Eq. ( 4), we find a steady state error of 2.16*10 -4 chips, which corresponds to the simulation results shown in Fig. 1. The validation of the tracking loops can also be done by means of testing their behavior in presence of additive white Gaussian noise. In this set of simulations, the receiver simulator only is tested, and therefore the LMS channel is not connected to it. An additive white Gaussian noise is added to a simulated GPS L1 C/A signal. The simulations are conducted over 20 seconds of signal and averaged 10 times. These simulations allow to test the right consideration of the predetection bandwidth of the loops which has no impact on the steady state error due to the dynamic, and therefore was not tested. The theoretical expression of the DLL tracking error due to the thermal noise for an Early-Minus-Late Power (EMLP) discriminator which expression is given in [START_REF] Lehner | A Novel Channel Model for Land Mobile Satellite Navigation[END_REF], is [START_REF] Van Dierendonck | Global Positioning System: Theory & Applications[END_REF], [START_REF] Van Dierendonck | Theory and Performance of Narrow Correlator Technology in GPS Receiver[END_REF]:

σ τ 2 = 𝐵 𝐿 𝑑 2 𝐶 𝑁 0 (1 + 2 (2 -𝑑) 𝐶 𝑁 0 𝑇 𝑖𝑛𝑡 ) (5)
Where  v 0 is the radial speed of the receiver  𝑇 𝑖𝑛𝑡 is the accumulation time

 C N 0
is the ratio of the carrier power to the noise power spectral density in dB-Hz  d is the chip spacing, was set to 0.5 chips B L,DLL = 1Hz in the simulations. d is set to 0.5 chips Fig. 2 compares the theoretical expression of the standard deviation of a conventional DLL due to the additive white Gaussian noise with its simulated value. The simulated behavior of the loop is well characterized by the theoretical expression. Both PLL and FLL have been tested using the same approach but are not detailed in the paper. The generation of the LOS and the right behavior of the tracking loops have been verified. The right consideration of the multipaths also has to be tested. To check that the reflections are well treated in the receiver, the LOS was set to zero and the LOS characteristics were injected in the multipath generation module. The same validations process has been conducted again on the LOS treated as an echo.

III. PARAMETERS USED FOR THE SIMULATIONS

A. Parameters for the channel generation

The first step consists in setting the parameters for the generation of the virtual street. A typical railway in downtown Toulouse +43° 35' 58.04", +1° 27' 29.82", as illustrated in Fig. 3, was taken as a reference and reproduced in the simulator. In urban area, the railway is several meters wide because of the proximity with train station. The rows of trees and housefronts are respectively placed at a distance of 12 and 22 meters from the center of the railway. The average height of the buildings was set to 15 m and the height of the antenna was set to 4.3 m which corresponds to the standard loading gauge for trains that carry passengers (GA, GB, GB1 as defined by UIC in [START_REF]Rules governing application of the enlarged GA, GB, GB1[END_REF]). The typical mask angle of the buildings for this scenario is 26° for an azimuth of 90°.

According to [START_REF] Perez-Fontan | Statistical modeling of the LMS channel[END_REF], to represent the urban environment with a sufficient time resolution, the sampling frequency shall respect the relation:

f 𝑠 ≥ 8. v 𝜆 (6)
The velocity tested for the vehicle in urban area has been set to 20km/h, and therefore a 1 kHz sampling frequency that verifies the relation was chosen for the simulations. A low velocity is preferred to avoid the multipaths to be filtered by the tracking loops due to the difference between their Doppler shift and the Doppler shift of the LOS. The satellite position relative to the user is set in terms of azimuth and elevation, and will vary in order to cover most of the possible scenarios and to characterize the pseudorange error as a function of the satellite elevation.

B. GNSS signals

The characteristics of the GNSS signals that are investigated in this paper are described in this section.

1) GPS L1 C/A

The GPS L1 C/A has the following characteristics [START_REF]-D: Navstar GPS Space Segment/Navigation User Interfaces[END_REF]: Its carrier is located at 1575.42 MHz. The carrier is modulated by a pseudorandom code (PRN) and a navigation message by Binary Phase Shift Keying (BPSK) with a Non-return to zero (NRZ) waveform. The spreading code is 1023 chips long and is transmitted at a rate of 1.023 Mchips/s. The navigation message has a bit rate of 50 bit/s.

2) Galileo E1 OS The Galileo E1 OS has the following characteristics [START_REF]European GNSS (Galileo) Open Service -Signal-in -Space[END_REF]: Its carrier is located at 1575.42 MHz. The signal is the sum of a data channel which waveform is CBOC(6,1,1/11,+) and a pilot channel modulated with CBOC(6,1,1/11,-). The spreading code is 4092 chips long and is transmitted at a rate of 1.023 Mchips/s. The pilot component is modulated by a secondary code. The data component has a bit rate of 250 bit/s.

C. Receiver model

The aim is to model the standard deviation of the code pseudorange error due to the multipaths only. In this case, no thermal noise due to the electronic component will be added to the correlator outputs. Our attention will focus on the DLL error because a basic GNSS receiver used the code measurement as the input for the navigation filter. The normalized EMLP discriminator is preferred as it is noncoherent and therefore less affected by an abnormally large phase tracking error. It is thus more suitable to be used in urban environment.

𝐷 𝐸𝑀𝐿𝑃 = 𝐼 𝐸 2 + 𝑄 𝐸 2 -𝐼 𝐿 2 -𝑄 𝐿 2 𝐼 𝐸 2 + 𝑄 𝐸 2 (7)
I E , Q E , I L , Q L are respectively the Early In-phase, Early Quadra-phase, Late In-phase, Late Quadra-phase correlator outputs of a conventional DLL.

The loop bandwidth of the DLL was set to 1 Hz which is a typical value for this loop in a GNSS receiver.

The accumulation duration for the DLL is 20 ms which correspond to the duration of one data bit for GPS L1 C/A and is typically selected when bit-synchronization has been achieved. For a fair comparison, the same accumulation duration is used for the processing of the pilot channel of Galileo E1 OS. To integrate over 20 ms, the synchronization with the secondary code is assumed.

In order to avoid the estimation of the delay to be biased due to the Doppler shift, a second order DLL is used. An FLL is preferred over a PLL for an improved robustness of the carrier tracking at the cost of a reduced precision. The bandwidth of the second order FLL is set to 5 Hz, a typical value for such a loop in actual GNSS receivers and a non-coherent arctangent discriminator is used to estimate the frequency error based on the correlator outputs. This discriminator is defined as [START_REF] Curran | Improving the Design of Frequency Lock Loops for GNSS Receivers[END_REF]:

𝐷 𝐴𝑇𝐴𝑁 = 1 𝑇 𝐼 𝑈𝑊 [𝑎𝑟𝑐𝑡𝑎𝑛 ( 𝑄 𝑃 (𝑘) 𝐼 𝑃 (𝑘) ) -𝑎𝑟𝑐𝑡𝑎𝑛 ( 𝑄 𝑃 (𝑘 -1) 𝐼 𝑃 (𝑘 -1) )] (8)
Where

𝑈𝑊(𝑥) = { 𝑥 -𝜋 for 𝑥 ≥ 𝜋/2 𝑥 for 𝜋/2 > 𝑥 ≥ -𝜋/2 𝑥 + 𝜋 for -𝜋/2 > 𝑥 (9)
The accumulation duration of the FLL is set to 10 ms because the FLL requires two consecutives correlator outputs to update the frequency of the NCO. Therefore, the corrected frequency can be used for the generation of the correlator output for the DLL.

Two GPS/Galileo dual constellation receivers are considered in this paper. Both processes L1 C/A and GALILEO E1 OS but differ from each other by their front-end filters and the delays between the early and late channels, referred to as the early-late spacing, and noted d. The rule-of-thumb for the selection of d is d=1/BW, where BW is the one-sided bandwidth of the front-end filter in MHz.

 The first receiver is a Narrowband BPSK(1)/BOC(1,1) receiver with a front-end bandwidth of 4 MHz and a chip spacing d = 0.5

 The second receiver is a Wideband BPSK(1) /CBOC(6,1,1/11) receiver, with a front-end bandwidth of 16 MHz and a chip spacing d = 0.125. The chosen bandwidth is large enough to take into account the contribution of the BOC(6,1) in the power spectral density of the CBOC(6,1,1/11).

The aim of this study is to determine the nominal standard deviation in a representative urban scenario due to the multipaths. The nominal scenario for the multipath error only occurs when the receiver is actually tracking the signal, and when the pseudorange errors are only due to the multipath interference phenomenon. The samples corresponding to the transient state of the loops, when the local replica is not locked on the LOS dynamic shall be excluded. Moreover the subsequent samples corresponding to a loss of lock shall also be removed not to impact the standard deviation. The receiver can self-monitor the effectiveness of its tracking thanks to detectors discussed in [START_REF] Van Dierendonck | Global Positioning System: Theory & Applications[END_REF]. In this paper, as a FLL is used to track the carrier, a 𝐶/𝑁 0 estimator is preferred over a phase lock detector:

𝐶 𝑁 0 ⁄ ̂= 𝐴 -1 + √𝐴. (𝐴 -1)
Where

𝐴 = [𝑚𝑒𝑎𝑛(𝐼 𝑃 2 + 𝑄 𝑃 2 )] 2 𝑣𝑎𝑟(𝐼 𝑃 2 + 𝑄 𝑃 2 ) (10)
In [START_REF]Rules governing application of the enlarged GA, GB, GB1[END_REF], the mean and variance are computed over 1 second. When the estimator is over a set threshold that has to be determined, it is decided that the signal is tracked and the samples shall be taken into account for the nominal standard deviation evaluation. The updating rate of the 𝐶/𝑁 0 estimation shall be sufficiently high to quickly detect any loss of lock but shall not be too low for computation issues. 0.1 second appears to be a good trade-off between sensitivity and complexity. If the 𝐶/𝑁 0 is over a set threshold, the latest 0.1 seconds of samples will be added to the "nominal case" vector by simple concatenation.

IV. MULTIPATH MITIGATION TECHNIQUES

Multipaths are a major challenge for the use of GNSS in urban environment, and particularly for critical and highaccuracy applications. Several reflections, that are delayed, phase shifted and attenuated compared to the LOS are processed by the receiver. These echoes distort the shape of the correlation function leading to an error on the estimation of the code delay by the DLL. The mitigation of multipath interference can be done at different level in the GNSS receiver: antenna, signal processing or position computation. The use of choke ring antennas, which filter out most of the multipath signal that are Left Hand Circular Polarized (LHCP), while allowing the direct signal which is Right Hand Circular Polarized (RHCP) to pass through, attenuates the multipath interferences. The mitigation can also be done in the GNSS receiver by modifying the design of the DLL or by adding a multipath estimation module. This paper focuses on this approach. Several techniques have been proposed in order to mitigate at the signal processing level the code pseudorange multipath error magnitude. The paper focuses on three methods: the Narrow Correlator, the Double-delta and the A Posteriori Multipath Estimation (APME). It has to be mentioned that none of these techniques mitigate the problem of NLOS reception or improve the quality of the tracking of the carrier phase.

A. Narrow Correlator (NC)

The Narrow correlator (NC) [START_REF] Van Dierendonck | Theory and Performance of Narrow Correlator Technology in GPS Receiver[END_REF], patented by NovAtel in 1991, is based on the principle of narrowing the chip spacing (d) of a conventional DLL. For a conventional DLL with an infinite bandwidth receiving an unique echo attenuated by 1/α, the maximum magnitude of the theoretical error envelope is αd/2. Therefore the reduction of d improves the multipath resilience of the DLL. The NC also reduces the standard deviation of the DLL due to the thermal noise according to [START_REF] Rispoli | Recent progress in application of GNSS and advanced communications for railway signaling[END_REF]. The limiting factor for the choice of d is the bandwidth available for the front-end filter [START_REF] Betz | Extended theory of early-late code tracking for a bandlimited GPS receiver[END_REF]. The estimation of the code delay is done by mean of the EMLP discriminator (D EMLP ).

B. Double-Delta (∆∆)

The Double-Delta (∆∆) technique uses more correlators (5 instead of 3) than a conventional early-late DLL [START_REF] Mcgraw | GNSS Multipath Mitigation Using Gated and High Resolution Correlator Concepts[END_REF]. A very early E2(τ) and a very late L2(τ) additional correlators are used, delayed by respectively -d and d compared to the prompt channel. Several different implementation of the ∆∆ technique are patented, using different discriminators, among them we can cite:  Ashtech's Strobe Correlator [START_REF] Irsigler | Comparison of multipath mitigation techniques with consideration of future signal structures[END_REF] which uses the discriminator:

𝐷 = 2(𝐸(𝜏) -𝐿(𝜏)) -[𝐸2(𝜏) -𝐿2(𝜏)] (11)
Where 𝐸(𝜏) and 𝐿(𝜏) are the Early and Late correlator of the conventional DLL.

 Rockwell Collin's High Resolution Correlator (HRC) [START_REF] Mcgraw | GNSS Multipath Mitigation Using Gated and High Resolution Correlator Concepts[END_REF], NovAtel's Pulse Aperture Correlation (PAC) [START_REF] Jones | Theory and Performance of the Pulse Aperture Correlator[END_REF] which use the discriminator:

𝐷 = 𝐸(𝜏) -𝐿(𝜏) -0.5 * [𝐸2(𝜏) -𝐿2(𝜏)] (12)
Both discriminators are equivalent to an early minus late discriminator, with:

𝐾 𝑐𝑐 ∆∆ (𝜏) = 𝐾 𝑐𝑐 (𝜏) -[𝐾 𝑐𝑐 (𝜏 + 𝑑/2) + 𝐾 𝑐𝑐 (𝜏 -𝑑/2)] (13)
The correlator function 𝐾 𝑐𝑐 ∆∆ (𝜏) is much narrower than 𝐾 𝑐𝑐 (𝜏), and thus, performs better in presence of multipath than a conventional DLL. The discriminator that is used in the paper is a non-coherent ∆∆ discriminator:

D = I E 2 (τ) + Q E 2 (τ) -I L 2 (τ) -Q L 2 (τ) -0.5 * [I E2 2 (τ) + Q E2 2 (τ) -I L2 2 (τ) -Q L2 2 (τ)] (14)
The double delta is implemented on the wideband GNSS receiver because its implementation increases the complexity of the receiver and the narrowband receiver is assumed to be a mass market low cost receiver. The NC does not increase the computation burden of the receiver and may be suitable for low cost receivers even if it requires the use of a wider frontend filter, which increases the power of the incoming signal processed.

C. A Posteriori Multipath Estimation (APME)

The A Posteriori Multipath Estimation (APME), patented by Septentrio [START_REF] Sleewaegen | Method and apparatus for processing signals for ranging applications[END_REF], is a multipath estimation technique described in [START_REF] Sleewaegen | Mitigating Short-Delay Multipath: a Promising New Technique[END_REF] and designed for GPS L1C/A. The APME does not require any modification on the design of the DLL, but relies on a multipath estimation module instead. The estimation module requires an extra very late correlator output, delayed by d with respect to the prompt channel. The amplitude of the multipaths affecting the narrow correlator is estimated by the following formula:

𝑀𝑃 = -0.42. 𝜆 𝑐 (1 - 𝛾 +2 𝐼 𝐿2 𝛾 0 𝐼 𝑃 1 1 -𝑑 ) (15)
γ +2 and γ 0 are corrections terms that accounts for the rounding of the correlation peak due to the filtering of the secondary lobes of the incoming signal by the front-end. For a front-end with a 16 MHz bandwidth, γ 0 = 1.013.

The estimation is filtered in a low-pass filter, with a noise equivalent bandwidth of 0.1 to 1 Hz to filter out the noise [START_REF] Sleewaegen | Method and apparatus for processing signals for ranging applications[END_REF].

The obtained correction is substracted to the code pseudorange measurements of the conventional narrow correlator DLL to reduce the multipath error.

The design of the low pass filter implemented to filter-out the noise on the estimate has a major impact on the quality of the multipath estimation. The filter output shall have a dynamic which is close to the dynamic of the multipath tracking error that results of a combination of both the multipath coherence time and the dynamic of the tracking loop. The optimal bandwidth of the low pass filter has been investigated by means of Monte-Carlo simulations on the LMS. For a vehicle driving at 20 km/h over 1 km, for several satellite positions, the low pass filter that minimizes the residual pseudorange error in the least square sense was obtained as shown in Table I. The filter implemented is a second order Butterworth filter, because it matches well with the dynamic of the second order DLL. A Bandwidth of 1.28 Hz and a linear Gain of 3.2 performs the best overall. Fig. 4 illustrates the multipath error estimation process and the reduction of the multipath interference by the APME. The APME residual error has been translated around -1 meter for a better visibility. From 20 seconds to 30 seconds the variance reduction obtained by the APME compared to a NC can be observed. Basic simulations showed that the FLL has to be replaced by a PLL when implementing the APME. A phase offset can result in a low power for the In-phase correlators compared to the noise level, which can result in an unbounded multipath estimation, derived from the quotient of two terms of noise. Nevertheless this technique could be extended to the non-coherent case by squaring the in-phase and quadra-phase correlators of the DLL. This study focuses on the implementation described in [START_REF] Sleewaegen | Mitigating Short-Delay Multipath: a Promising New Technique[END_REF].

D. Theoretical performances

The study of the error envelope is a classical approach to evaluate the performances of the existing multipath mitigation techniques. The envelope is obtained by assuming that the signal is affected by one reflection, with a fixed attenuation, in phase or antiphase compared to the LOS. The composite discriminator is null for a delay that is neither zero, nor superior to 1+d/2 chip, leading to an error on the code pseudorange. The multipath error envelope of the narrowband (BW = 4 MHz) receiver model, that processes BPSK(1) and BOC(1,1) with a conventional early-late DLL is plotted in Fig. 5. The Signal to Multipath Ratio (SMR) is set to 20 dB in these simulations. According to Fig. 5 BOC(1,1) outperforms BPSK(1) for long to medium multipath delays. The multipath error envelope of the dual constellation wideband (BW = 16 MHz) receiver model with multipath mitigation techniques is plotted in Fig. 6 for BPSK(1) and Fig. 7 for CBOC(6,1,1/11). The benefit of reducing the chip spacing from 0.5 to 0.125 reduces the maximum magnitude of the theoretical error by 75% for BPSK. For BPSK, the ∆∆ and APME outperforms the Narrow correlator for short, medium, and long delays multipaths. The APME does not mitigate medium to long delay multipaths, but has the envelope with the lowest magnitude for the very short delays multipaths (less than 0.1 chips). This characteristic makes this technique promising because, as stated in [START_REF] Sleewaegen | Mitigating Short-Delay Multipath: a Promising New Technique[END_REF], the distribution of the reflections is higher for short and very short delays. As for the CBOC(6,1,1/11), which envelope is plotted in Fig. 7, the ∆∆ improve the multipath envelope for reflections with short to medium delays. The estimation step of the APME that is described in [START_REF] Sleewaegen | Mitigating Short-Delay Multipath: a Promising New Technique[END_REF] is designed for the BPSK signal, and does not behave well for CBOC(6,1,1/11). Therefore the performances of APME are not evaluated on this signal. The error envelope approach has a limited scope because it does not consider the multiplicity of the reflections, the distribution of the code delays which is not uniform, and the filtering by the carrier loop due to the Doppler shift of the echoes. Among the techniques which are not investigated in this paper are the Early1/Early2 [START_REF] Irsigler | Comparison of multipath mitigation techniques with consideration of future signal structures[END_REF] tracker that processes the tracking in the area where the correlation function is less distorted, and the NovAtel's MEDLL [START_REF] Van Nee | The multipath estimating delay lock loop: approaching theoretical accuracy limits[END_REF] that uses several correlator to determine the amplitude, delay and phase of all the received signals, including the multipaths. The former is sensitive to the thermal noise because of the low magnitude of the correlator at the base of the correlation function and the latter increases the complexity of the receiver whereas achieving comparable level of performances as the studied techniques.

V. SIMULATION RESULTS

The simulations were conducted for a train travelling in the virtual street generated by the LMS at a velocity of 20 km/h. The satellite is positioned at eight different elevations (10, 20, 30, 40, 50, 60, ). As the topography of the street is symmetric, it is possible to limit the study at these three azimuths. The channel is generated and processed over 1km for each satellite position.

A. Threshold determination for the 𝐶/𝑁 0 estimator

The nominal multipath error model has been calculated using the samples for which the GNSS receiver is tracking the signal. In this paper, the 𝐶/𝑁 0 estimator given in [START_REF]Rules governing application of the enlarged GA, GB, GB1[END_REF] has been used as tracking detector. This estimator can easily be implemented in an actual GNSS receiver. First it is necessary to set the 𝐶/𝑁 0 threshold at which the receiver is considered to be tracking, in presence of multipath only. Considering four thresholds [START_REF] Sleewaegen | Method and apparatus for processing signals for ranging applications[END_REF][START_REF] Meguro | GPS Multipath Mitigation for Urban Area Using Omnidirectional Infrared Camera[END_REF]30 and 35 dB-Hz), the distribution of the residual error obtained after the sample selection is studied for the specific case of the narrowband BPSK(1) receiver. For a satellite which azimuth is firstly 45° and then 90°, for each elevation, and each threshold, the residual samples have been concatenated. The mean pseudorange error is plotted in Fig. 8. For elevations higher than 20°, the mean pseudorange error is sub-decimetric and does not vary with the 𝐶/𝑁 0 threshold. For very low elevations (i.e. 20° and 10°), the distribution of the pseudorange error is not centered and its 95th percentile confidence interval is abnormally important as illustrated in Fig. 9.

The detector with a thresholds set to 30 and 35 dB-Hz effectively detect the loss of lock and exclude all the samples as illustrated in Fig. 10. The detector with a threshold set to 20 or 25 dB-Hz detects most of the biased samples, but 3% of the samples still remain not detected and are considered relevant. The magnitude of the bias is sufficiently important (several meters) to lead to a positioning failure and jeopardize the safety of the signaling system. Fig. 10 shows the percentage of samples over the threshold, and therefore the availability of the satellite. It can be inferred from Fig. 10 that setting the threshold to 35 dB-Hz is not advisable because it reduces the availability of the satellite to less than 80 % even when the signal is tracked and that the pseudorange error is centered and has 95 th percentile confidence interval lower than 2 m. The 𝐶/𝑁 0 mask is set to 30 dB-Hz as it represents the best trade-off between the ability to detect a loss of lock and signal availability. The case of the 0° azimuth is not taken into account for the threshold determination as the LOS is never obstructed and the receiver is tracking for the majority of the samples. The same approach has been followed for each modulation and multipath mitigation technique but is not detailed in the paper.

The same conclusion has been done for the other receivers and modulations. In Fig. 10, it can be observed that the percentage of effective tracking for nonzero azimuths decreases with the elevation and is negligible for satellite elevation angles of 10° and 20°. To explain this phenomenon, the Signal to Multipath amplitude Ratios (SMR) histograms are plotted in Fig. 11, Fig. 12 and Fig. 13 which respectively correspond to satellite elevations of 10°, 20° and 30°. The SMR is defined as the ratio between the amplitude of the LOS and the amplitude of the strongest echo. The bimodality of the SMR distribution is characteristic of the two configurations that occur. The mode corresponding to the low SMR is the consequence of the reception of a combination of a NLOS and reflections. The mode corresponding to the high SMR is the consequence of the reception of a LOS and reflections. The NLOS reception is characterized by very low or negative SMR which lead to losses of lock of the tracking loops. It can be inferred from Fig. 11, Fig. 12 and Fig. 13 that the number of samples with very low or negative SMR increases as the elevation decreases. Therefore the amount of samples that are over the 𝐶/𝑁 0 mask increases with the elevation.

B. Evaluation of the performances of the different receiver's architectures

For the simulations, the channel is generated once and each receiver's architectures and signals are processed on the same channel. The raw standard deviations of the pseudorange error are given for every azimuth. Fig. 14 illustrates the standard deviation of the code pseudorange error due to multipath. When the satellite is orthogonal to the housefronts (Azimuth = 90°) the GNSS receiver is more affected by multipath interference because this scenario is geometrically more prone for the reflections on the housefronts. According to Fig. 14, the order of magnitude of the standard deviation of the DLL error in urban is submeter. It could be inferred from that both BPSK(1) and BOC(1,1) perform the same except for the elevations that are lower than or equal to 40°. These elevations are close to or lower than the mask angle of the buildings, which means that the received signal is most of the time a NLOS embedded in statistically generated reflections. In this scenario the LOS to multipath power ratio is sufficiently low for the DLL not to be able to track the signal of interest. In the case of BPSK(1) for elevation under than 30° and azimuth different of 0°, no samples has an estimated 𝐶/𝑁 0 over the threshold that was set to 30 dB-HZ. This phenomenon justifies the lack of points for BPSK(1) under 30° of elevation. When using the BOC(1,1) modulation, a tiny percentage of samples of the signal is tracked for low elevations. Narrowing the chips spacing of the DLL improves the robustness of the code tracking as illustrated in Fig. 15, where the channel is processed by the wideband receiver. For elevation higher than or equal to 40°, for BPSK(1), the mean reduction of the standard deviation of the multipath code error obtained by narrowing the correlators from 0.5 to 0.125 is 51.6%. This gain is more representative of the real-life performances than the theoretical gain of 75% on the multipath envelope because it takes into account the tendency of the reflection to be short delay and the multipath filtering phenomenon by the FLL due to the Doppler shift of the reflections. On the same set of data, for BPSK(1), the Double-Delta reduces the standard deviation of the code error by 29.3% on average compared to the Narrow-Correlator. The APME performs a 37.8% reduction compared to the Narrow Correlator. When the tracking is sporadic, i.e. for low elevations, the PLL that is used for the APME instead of a FLL is not able to track efficiently which explains the lack of measurement for Azimuths lower than 30° for this technique.

For every receiver's architecture and every modulation, the 95th percent confidence interval of the pseudorange error is given in annex. The APME outperforms the Double-Delta for high elevations in terms of reduction of the error magnitude as illustrated in Appendix in Table V Reducing the chips spacing also improves the variance of the DLL for CBOC(6,1,1/11). The double delta only improves the tracking quality for low elevations as illustrated in Fig. 16.

The gain observed on the multipath envelope is not observed on the simulations on the LMS channel for high elevations.

C. Nominal multipath pseudorange standard deviations 1) Overbounding

The design of integrity monitoring algorithms requires the characterization of the actual error distribution by a known distribution. Integrity monitoring systems usually assume that the error follows a zero mean Gaussian distribution. The distribution of the raw multipath pseudorange error obtained by processing the LMS channel output is neither Gaussian nor centered. Thanks to the 𝐶/𝑁 0 discrimination, with a threshold of 30 dB-Hz, the mean of the residual code error is subdecimeter as shown in Fig. 8. The centered distribution still remains not Gaussian. The well-known Cumulative Distribution Function (CDF) overbounding technique described in [START_REF] Decleene | Defining Pseudorange Integrity -Overbounding[END_REF] is used in this paper. The CDF overbounding technique requires the distribution to be unimodal, centered, and symmetric. The unimodality was checked by watching the histograms of the observed pseudorange error. All three assumptions are fulfilled after the 𝐶/𝑁 0 discrimination. In real-life, the vehicle may not travel along a straight street such as what is generated in the LMS, but may face buildings when turning, or bridges, so that the satellites located in front of it can be masked. Therefore proposing a model in which the nominal multipath standard deviation is a function of the azimuth would not be relevant. For each simulated elevations, the pseudorange error from the range of possible azimuth shall be concatenated before processing the overbounding.

2

) Simulation results

The results of the CDF overbounding on the pseudorange error for the narrowband receiver are given in Fig. 17. Using multipath mitigations techniques reduces the nominal standard deviation due to multipaths which can improve the protection levels and therefore improve the availability of the integrity monitoring algorithm. For elevations lower than the mask angle of the buildings, the only samples that are above the 𝐶/𝑁 0 threshold are those obtained for a null azimuth. For this azimuth, the multipath interference is lower which results in lower sigmas as illustrated in Fig. 18. Few samples with low elevation and nonzero azimuth are effectively tracked for the wideband CBOC receiver and the overbounding process, leads to a high nominal standard deviation (as shown in Fig. 19).

VI. FUTURE WORK

The use of GNSS in critical applications requires an augmentation system for integrity monitoring. Ground and Satellite augmentations are not sufficient for critical terrestrial applications because they just provide alerts in case of SiS fault whereas, in urban area, the main faults are due the NLOS reception and the strong multipath interferences. Receiver Autonomous Integrity Monitoring (RAIM) techniques are able to detect and exclude faulty measurements in those situations. The classic RAIM algorithms usually assume that one pseudorange measurement is biased. In urban environment, several pseudorange may be biased simultaneously due to the receiver's environment. The use of RAIM that can detect multiple failures makes it possible to overcome the problem of simultaneous faults. The trade-off of such an approach is the augmentation of the protection levels and the calculation burden in the receiver. Another approach consists in pre excluding the faulty measurements before providing them to the integrity monitoring module. The NLOS detection may be feasible by means of Shadowmatching, [START_REF] Obst | Urban Multipath Detection and mitigation with Dynamic 3D Maps for Reliable Land Vehicle Localization[END_REF], fisheye camera [START_REF] Meguro | GPS Multipath Mitigation for Urban Area Using Omnidirectional Infrared Camera[END_REF], or by applying masks on the 𝐶/𝑁 0 or the elevations. The detection of large pseudorange error due to the multipath interference (LOS + echoes) is necessary. Indeed, the multipath mitigation techniques that improve the availability of the integrity monitoring algorithms may reduce the magnitude of large multipath errors but do not necessarily prevent from potential abnormal large errors. A way to monitor the likelihood of tracking problems is the monitoring of the correlation function distortion. It can be done by means of by signal quality indicators such as the 𝐶/𝑁 0 estimators or the metrics that are used for the detection of evil-waveforms as described by [START_REF] Irsigler | Development of a Real-Time Multipath Monitor Based on Multi-Correlator Observations[END_REF]. One of the proposed metrics in [START_REF] Irsigler | Development of a Real-Time Multipath Monitor Based on Multi-Correlator Observations[END_REF] is the simple ratio test which expression is:

𝑚 = 𝐼 𝑋 𝐼 𝑌 (16)
Where X and Y are the fraction of chip by which the local replica is delayed compared to the prompt channel. An analogy can be done between the differential ratio test and the multipath estimator of the APME. This estimator can be rewritten on the following form:

𝑀𝑃 = 0.42 𝛾 2 𝛾 0 (1 -𝑑) [ 𝐼 2 𝐼 0 - 1/𝛾 2 (1 -𝑑) 1/𝛾 0 ] (17)
The APME is a particular case of Simple Ratio Test modified to be zero mean when the signal is only affected by thermal noise (by subtracting its mean), and that corrects the distortion of the correlation function due to the front-end filter thanks to the γ terms. The important ability of the APME to mitigate multipath interference can therefore be used for multipath and (potential failure) detection. Based on this assumption, further work will be focused on proposing and testing new signal quality indicators.

VII. CONCLUSIONS

A method for the estimation of the nominal pseudorange error measurement due to multipaths in urban area has been proposed. The method is based on the connection between a realistic urban channel model and a GNSS receiver simulator. This simulation platform has been validated comparing its results with its theoretical behavior in well-known scenarios. Several GNSS receiver architectures have been tested, including high end receivers that apply multipath mitigation techniques. The simulations have made it possible to compare the performances of these algorithms when processing a signal affected by reflections only. The Narrow Correlator has been demonstrated to improve the standard deviation due to multipath by 51.6% with respect to the a conventional DLL with a half chip Early-Minus late spacing. The Double-Delta and APME reduce the standard deviation due to multipath by respectively 29.3% and 37.8% with respect to the Narrow Correlator. For the different receivers and modulations, the distribution of the multipath error has been overbounded to obtain the nominal error model as a function of the satellite elevation. The simulations showed the improvement in the accuracy of the pseudorange measurement and in the nominal model when using multipath mitigation techniques for a better availability of the integrity monitoring algorithm in urban areas. The possible detection at the signal processing level of pseudorange measurement affected by non-nominal multipath interference has been discussed.

Appendix

This appendix provides the 95th percentile confidence intervals of the pseudorange error (in meters) for the various configurations of receivers and signals that are tested in the simulations.

Narrowband receiver

  𝑑(𝑘𝑇 𝑠 -𝜏 𝑛 (𝑘))𝐾 𝑐𝑐 (𝜏 𝑛 (𝑘) 𝑁 𝑛=1 -𝜏(𝑘))𝑠𝑖𝑛𝑐 (𝜋 (𝑓 𝑛 (𝑘) -𝑓 ̂(𝑘)) 𝑇 𝑖 ) cos (𝜑 𝑛 (𝑘) -𝜑 ̂(𝑘)) 𝑄 𝑃 (𝑘) = ∑ 𝛼 𝑛 (𝑘) 2 𝑑(𝑘𝑇 𝑠 -𝜏 𝑛 (𝑘))𝐾 𝑐𝑐 (𝜏 𝑛 (𝑘) 𝑁 𝑛=1 -𝜏(𝑘))𝑠𝑖𝑛𝑐 (𝜋 (𝑓 𝑛 (𝑘) -𝑓 ̂(𝑘)) 𝑇 𝑖 ) sin (𝜑 𝑛 (𝑘) -𝜑 ̂(𝑘))

Fig. 1 .

 1 Fig. 1. Steady state error of a first order DLL for a vehicle at a constant speed on the LMS

Fig. 2 .

 2 Fig. 2. DLL tracking error due to thermal noise

Fig. 3 .

 3 Fig. 3. Street model used for the channel generation

Fig. 5 .Fig. 6 .Fig. 7 .

 567 Fig. 5. Multipath error envelope for the narrowband BPSK(1)/BOC(1,1) receiver

  70 and 80°) and three different azimuths (0, 45

Fig. 8 .Fig. 9 .Fig. 10 .

 8910 Fig.8. Mean error of the code pseudorange error for the narrowband receiver processing BPSK(1)

Fig. 11 .Fig. 12 .Fig. 13 .

 111213 Fig. 11. Signal to multipath amplitude ratio for a satellite elevation of 10°

Fig. 14 .

 14 Fig. 14. Standard deviation of the code pseudorange error for the narrowband receiver

  (1), Azimuth = 0°W C-BPSK(1), Azimuth = 45°W C-BPSK(1), Azimuth = 90°W C-BOC(1,1), Azimuth = 0°W C-BOC(1,1), Azimuth = 45°W C-BOC(1,1), Azimuth = 90°

Fig. 15 .

 15 Fig. 15. Standard deviation of the code pseudorange error for the wideband receiver for the BPSK(1) modulation

Fig. 16 .

 16 Fig. 16. Standard deviation of the code pseudorange error for the wideband receiver for the CBOC(6,1,1/11) modulation

Fig. 17 .

 17 Fig. 17. Nominal multipath standard deviation for the narrowband receiver obtained after CDF overbounding

Fig. 18 .Fig. 19 .

 1819 Fig. 18. Nominal multipath standard deviation for the wideband receiver obtained after CDF overbounding for the BPSK(1) modulation

ACKNOWLEDGMENT

TABLE I .

 I OPTIMAL FILTER FOR THE MULTIPATH ERROR ESTIMATION

	Elevation	0°	Azimuth 45°	90°
	50°	BW =1.28 Gain = 3.18	BW =1.24 Gain = 2.91	BW =1.28 Gain = 2.85
	60°	BW =1.48 Gain = 2.79	BW =1.28 Gain = 3.15	BW =1.08 Gain = 3.29
	70°	BW =1.28 Gain = 3.28	BW =1.08 Gain = 3.54	BW =1.28 Gain = 3.28
	80°	BW =1.28 Gain = 3.24	BW =1.28 Gain = 3.18	BW =1.28 Gain = 3.35

Fig. 4. Principle of the multipath mitigation by the APME. The residual error of the APME has been translated of -1 meter for a better visibility.

TABLE IV .

 IV 95 TH PERCENTILE CONFIDENCE INTERVAL FOR THE WIDEBAND RECEIVER (NC) PROCESSING BPSK(1)

					Elevation	0°	Azimuth 45°	90°
					70°	0.3138	0.506	0.8578
					80°	0.4739	0.4385	0.7583
					TABLE V.	95 TH PERCENTILE CONFIDENCE INTERVAL FOR THE
					WIDEBAND RECEIVER WITH ∆∆ PROCESSING BPSK(1)
					Elevation	0°	Azimuth 45°	90°
					10°	0.1669	-	-
					20°	0.2314	0.1715	-
					30°	0.227	0.1822	0.7064
					40°	0.2864	0.3055	0.7665
					50°	0.2979	0.4393	0.631
					60°	0.2454	0.3158	0.5389
					70°	0.2581	0.3838	0.6459
					80°	0.3859	0.3554	0.6489
					TABLE VI.	95 TH PERCENTILE CONFIDENCE INTERVAL FOR THE
					WIDEBAND RECEIVER WITH APME PROCESSING BPSK(1)
					Elevation	0°	Azimuth 45°	90°
					10°	0.4423	-	-
					20°	0.5887	-	-
					30°	0.4158	0.2529	0.5538
					40°	0.3484	0.2691	0.5439
					50°	0.3039	0.3797	0.528
					60°	0.2259	0.2507	0.4578
					70°	0.2085	0.3506	0.4598
					80°	0.2898	0.2721	0.4299
					TABLE VII.	95 TH PERCENTILE CONFIDENCE INTERVAL FOR THE
					WIDEBAND RECEIVER (NC) PROCESSING CBOC(6,1,1/11)
					Elevation	0°	Azimuth 45°	90°
	Elevation	0°	Azimuth 45°	90°	10° 20°	0.1683 0.2467	3.9713 1.1279	4.4582 3.6121
	10°	0.4166	0.5763	-	30°	0.2209	0.1529	1.6784
	20°	0.6038	0.3563	-	40°	0.3254	0.2329	0.6924
	30°	0.5135	0.3426	1.0606	50°	0.2453	0.3004	0.5097
	40°	0.4989	0.4309	0.9723	60°	0.1689	0.1935	0.3729
	50°	0.3692	0.5235	0.8597	70°	0.1613	0.2589	0.4453
	60°	0.3018	0.4023	0.7164	80°	0.2485	0.2409	0.4143

TABLE VIII .

 VIII 95 TH PERCENTILE CONFIDENCE INTERVAL FOR THE WIDEBAND RECEIVER WITH ∆∆ PROCESSING CBOC(6,1,1/11)

	Elevation	0°	Azimuth 45°	90°
	10°	0.1406	0.3521	0.2685
	20°	0.1803	0.2664	0.3678
	30°	0.17	0.1516	0.4512
	40°	0.1799	0.2166	0.497
	50°	0.2077	0.307	0.5237
	60°	0.1883	0.226	0.4262
	70°	0.1879	0.3049	0.4952
	80°	0.2989	0.2743	0.5177

This study was funded by the European GNSS Agency (GSA) as a part of the GENIUS project, Ecole Nationale de l'Aviation Civile (ENAC) and Egis Avia.


B Nacke
email: nacke@etp.uni-hannover.de

V Kichigin

V Geza

I Poznyak

Continuous melting and pouring of an aluminum oxide based melt with cold crucible

Keywords: Electromagnetic melting, cold crucible melting, skull melting, oxide melting, continuous pouring

Continuous melting and pouring of an aluminum oxide based melt with cold crucible

B Nacke, V Kichigin, V Geza, I Poznyak

Introduction

At the Institute of Electrotechnology (ETP) an inductive skull melting installation is available to conduct experiments with all sorts of low conductive materials. In particular, the installation is used for production of optical glasses, melting of oxides, chemical high temperature reactions and directional solidification or pulling of crystals [START_REF] Nacke | Induction skull melting of oxides and glasses in a cold crucible[END_REF][START_REF] Kudryash | Experimental investigation of induction melting in cold crucible for high temperature processing of oxides using HF transistor power supply[END_REF][START_REF] Nacke | Simulation of 3D heat and mass transfer in glass and oxide melts using the inductive skull melting technology[END_REF]. The needed process temperature, especially in case of oxides, often exceed 2000°C. With the inductive skull melting technology it is possible to handle process temperatures of more than 3000°C by using an inductor crucible when melting glasses or a cold crucible when melting oxides. The inductor crucible can be regarded as a single turn coil combined with a crucible into a single unit, whereas a separate inductor is needed when a cold crucible is used. In both cases, massive water cooling of all compoments is required. Between the hot melt inside the crucible and the water-cooled crucible walls, a layer of solid material is formed, which is commonly referred to as skull layer. The skull layer provides protection against contamination of the melt and damage to the crucible in general. The existing installation was extended in the last years by a new cold crucible for continuous melting and pouring of high temperature oxides [START_REF] Kichigin | Investigation of hydrodynamics of alumina oxide in a cold crucible at continuous melting and discharging[END_REF]. A special design of the cold crucible with two separated working areas allows continuous melting, superheating and pouring in one device. Fig. 1 shows the principle of the cold crucible with the inductor, the water-cooled bottom, an additional separator to create a separate chamber for superheating and the discharging hole for pouring.

Design of the cold crucible for continuous pouring Fig. 2 shows the design of the induction furnace. The furnace consists of the single-turn inductor with inner diameter of 340 mm and height of 250 mm. The water-cooled crucible has inner diameter of 250 mm and height of 505 mm. Every section of the crucible is separately water-cooled. The cross section of the crucible is mostly of round shape, only the front sections are mounted in a straight line. The middle front section has a discharging hole. In the upper part the sections are fixed on the pertinax plate by a bolt connection, at the bottom they are fixed by a copper ring. A water-cooled bottom is located inside the crucible.

A barrier (separator) is located inside the crucible in front of the discharging hole. It has no contact with the crucible. The barrier is produced from copper profile and is water-cooled too. The separator provides the division of the crucible into two chambers (Fig. 3) and allows overheating the melt inside the pouring chamber.

The induction system with cold crucible is located inside the watercooled working chamber (see Fig. 2). The crucible is tilted by a certain angle to allow pouring the melt through the discharging hole.

The inductor is connected to the water-cooled current feeder of the high frequency transistor generator with an output power of 300 kW and an output voltage of 450 V which is doubled by the capacitor transformation to around 900 V. Frequency can vary in the range from 80 to 350 kHz.

Experiments

Experiments with new designed cold crucible for melting of aluminium oxide have been carried out at ETP. The melting process includes several steps: preparation of the oxide material; starting of the process; induction melting to form the melt and pouring. The oxide material is filled into the cold crucible up to the level lower than the discharging hole. The oxide in the solid state is good electrical insulator and but its electrical conductivity increases at temperatures near the melting point. Wherein coupling to the electromagnetic field is appearing. For starting of the melting process a little quantity of pure metal is embedded inside the oxide material. The small pieces of this starting metal melt down by induction heating and transfer the induced heat to the surrounding oxide. As soon as the oxide has itself reached high enough for coupling temperature it is able to absorb the electromagnetic energy directly because the electrical conductivity is increased. The melting pool is extended by charging of additional oxide powder.

Induction melting is carried out at frequency of 300 kHz with the generator output voltage of 400 V and power of 200 kW. Fig. 4 shows the top view to the crucible during the melting process before pouring. The oxide powder is permanently charged to the melting zone of the crucible through the water-cooled Fig. 6: Calculation area and boundary conditions charging tube. The barrier prevents flow of non-melted material to the superheating zone during melting. The flow of the melted oxide from the melting zone to the superheating zone takes place under the separator. The pouring stream of the melt is shown in Fig. 5. The main aim of the pouring process is to receive the continuous melt jet through the discharging hole in the cold crucible. This is mainly influenced by the melt flow and the temperature distribution in the melt. In order to receive the continuous melt jet, permanent transfer of the melted oxide from the melting zone to the superheating zone must be provided. The temperature in the melt must be high enough for melting the charged oxide and the pressure on the melt surface generated by the non-melted oxide must be constant. The melt surface temperature is kept up around 2300°C during the pouring.

Mathematical model

In order to optimize the flow from the melting to the superheating chamber and secondly to optimize also the pouring behaviour, a 3D numerical model for the cold crucible with separator has been developed. Calculation of electromagnetic field based on Maxwell equations was realized in ANSYS software [START_REF][END_REF]. The obtained electromagnetic forces and Joule heat are imported to the fluid dynamic analysis. The melt flow is described by the Navier-Stokes equation:

g T f p u p u u t u                    1 1
Here uvelocity; ppressure; Ttemperature; f -Lorentz forces density; β -thermal expansion coefficient; νkinematic viscosity. The thermal process in the melt is described by the heat transfer equation:

q T c T u t T          
Here λ -thermal conductivity, c -specific heat, q -Joule heat density.

Geometry of the melt and the boundary conditions are shown in Fig. 6. Constant temperature T=T MELT was assumed for the bottom and the walls with no slip for the flow field which corresponds to the solid-liquid interface between the melt and the skull. The top surface divided in pouring and melting zones is treated differently. In the pouring zone a radiation boundary conditions with free slip are used. In the melting zone the fixed temperature T=T MELT and flow rate are used.

Results of numerical simulation

Fig. 7 shows the temperature distribution in the vertical cross section of the crucible. In the middle part of the crucible the temperature isolines are horizontal. The arrows show the direction of melt flow which is in this zone from right side to the left side. In front of the pouring zone the flow arrows change its direction to the bottom. This effect is explained by temperature differences between the separator and the melt. In the pouring zone the flow arrows change its direction again to the top and then reach the discharging hole. Fig. 8 shows the streamlines of velocity in the melt. Along the whole perimeter of the crucible buoyancy vortices are observed which are caused by strong radial temperature gradients. These vortices are narrow in radial direction and long in vertical direction. In the region of outflow the vortex is more intensive and larger in radial direction. An additional vortex is formed directly under the surface in the superheating zone, it has velocity vectors directed downwards near the separator due to the cooling of melt.

Comparison of experimental and numerical values of power losses

The numerical simulation of the induction melting process in cold crucible allows obtaining the power losses in each induction furnace unit and the induced power in the melt. For this purpose a 3D numerical model of the induction furnace was developed so that the system geometry and dimensions are fully representing the real induction furnace.

The numerical model was used to predict the power losses in the furnace system without any material in the crucible in order to separate the losses in the crucible and the induced power in the melt. The comparison of obtained experimentally energetical parameters of the induction furnace and results of numerical modelling, presented in Table 1, show very good agreement. Experiment allows measuring the power losses in the cold crucible and in the melt only together, while numerical simulation provides the data separately.

Conclusions

A new type of cold crucible furnace has been developed for continuous melting and pouring of oxide materials. The furnace has two separate chambers for melting and superheating but it is powered by one induction coil only. The superheating zone in the second chamber is achieved by using an additional barrier which provides the separation but also the superheating of the melt before pouring. With a special numerical simulation model the melt flow inside the double-chamber crucible has been investigated and optimized. The second numerical model is able to calculate the energetical parameters of the furnace and allows optimization of energy efficiency of the melting and pouring processes. The new furnace has been experimentally investigated and meanwhile tested in practice. It confirms that continuous melting and pouring of oxide materials in one cold crucible is feasible and efficient.

Fig. 1 :Fig. 2 :

 12 Fig. 1: Principle of cold crucible for continuous pouring

Fig. 3 :

 3 Fig. 3: Cold crucible with separator

Fig. 7 :

 7 Fig. 7: Temperature field and direction of melt flow

Table 1 :

 1 Comparison of experimental and numerical data

	Parameter	Experiment	Numerical modelling
	Diameter of the melt, cm	25	25
	Height of the melt, cm	14	14
	Frequency, kHz	298	298
	Power losses in the inductor, kW	3.4	3.4
	Power losses in the cold crucible and melt, kW	76	15.8 (crucible) 63,2 (melt)
	Power losses in the bottom, kW	3.8	1.4
	Total power losses, kW	83.2	83.8


L Nastac
email: lnastac@eng.ua.edu

D Zhang

S Jia

Y Xuan

Advances on Experimental and Numerical Modeling of Al-based Alloys and Nanocomposites Fabricated via Ultrasonic and Electromagnetic Processing

Keywords: Ultrasonic and EM stirring, Cavitation, A356 alloy, Metal-Matrix-Nanocomposites, Solidification microstructure, Numerical modeling of ceramic nanoparticle dispersion

HAL is

Introduction

Aluminum matrix nanocomposites have the potential to offer outstanding properties, including low density, high specific strength, high specific stiffness, excellent wear resistance and controllable expansion coefficient, which make them attractive for numerous applications in aerospace, automobile, and military industries field [1][2][3][4][5][6][7][8][9][START_REF] Jia | Advances in the Science and Engineering of Casting Solidification: An MPMD Symp[END_REF]. As shown in literatures [START_REF] Mussert | [END_REF]12], it is beneficial to use nano-sized ceramic particles to fabricate metal matrix nanocomposites (MMNCs) since similar ductility and higher strength with the metal matrix can be achieved. Currently, there are several fabrication methods of MMNCs, including mechanical alloying with high energy milling [START_REF] Chen | Advanced Particulate Materials And Processes International Conference[END_REF], ball milling [START_REF] Filho | [END_REF], nano-sintering [15], spray deposition, electrical plating, sol-gel synthesis, laser deposition, etc. The mixing of nano-sized ceramic particles is normally lengthy, expensive, and energy consuming.

Ultrasonic Stirring Technology (UST) has been extensively used in purifying, degassing, and refinement of metallic melt [16][17][18][19][20], mainly because introducing the ultrasonic energy into a liquid will induce nonlinear effects such as cavitation and acoustic streaming. SiC nanoparticles are widely used as reinforcement particles due to their relatively good thermal and chemical stability. Ultrasonic vibration can improve the wettability between the reinforced nanoparticles and the metal matrix [21][22][23], which will assist to distribute the nanoparticles more uniformly into the metal matrix. Ultrasonic cavitation in liquids causes high speed liquid jets of up to 1000 km/h. Such jets press liquid at high pressure between the particles and separate them from each other [START_REF] Hielscher | Dans European Nano Systems Worshop -ENS 2005[END_REF].

In this article, the effects of the ultrasonically dispersed SiC nanoparticles on the microstructure of A356 nanocomposites are studied in detail.

Experimental and Modeling Approaches

Experimental Approach. A356 alloy is selected as the metallic matrix. The ceramic nanoparticles used here are β-SiC (spherical shape, average diameter of about 50 nm). The UST system used in this study is illustrated in Fig. 1. The main parameters of the ultrasonic equipment are: max. power, P = 2.4 kW and frequency, f = 18 kHz. An induction furnace with a capacity of 2.7 kg was used to melt the A356 alloy. After the alloy is melted, the Nb ultrasonic probe is inserted to about 50 mm beneath the melt surface to perform ultrasonic stirring at 1.75 kW and 18 kHz frequency. 1 wt. % SiC nanoparticles are injected into the cavitation area. The molten pool is protected by Argon gas atmosphere. The experiments are repeated several times for statistical interpretation of the results. Experiments applying UST treatment during solidification are also performed to better understand the effects of UST. The comparison experiment is done by melting the alloy and turning off the furnace to let the molten metal solidify. The other experiment is done by melting the alloy and turning off furnace, and then treating the molten alloy with UST during solidification. Modeling Approach. The geometry of the model is shown in Fig. 2(a). The ultrasonic probe has a diameter of 40 mm. The liquid aluminium is A356. The SiC nanoparticles are treated as inert-particles. The mass flow rate of the nanoparticles is 0.014 kg/s. Thus, 1.0 wt. % of nanoparticles can be injected at about 20 mm above the bottom of the furnace for 1.0 sec. The induction furnace is comprised of an Al2O3 crucible and 3 Cu induction coils. An electric current of 200A with frequency of 50 Hz is passed through the coils to induce eddy currents in the liquid aluminium alloys. All material properties required for EM modeling are available in ANSYS Maxwell. The multiphase CFD model is able to account for turbulent fluid flow, heat transfer, solidification and the complex interaction between the molten alloy and nanoparticles by using the ANSYS Fluent DDPM and κ-ω turbulence model. The CFD model is described in detail in [START_REF] Zhang | Proceedings of the "Advanced Metal Casting Technologies" Symposium, MS&T2014[END_REF][START_REF] Zhang | [END_REF]. The solution procedure is presented in the following paragraph. The nanoparticles are injected at every fluid flow time step with a mass flow rate of 0.014 kg/s in the first second. Particles are tracked at every time step after the fluid velocity field is solved. Because of the low volume fraction of the discrete phase, one-way coupling is employed, which neglects the effect of the discrete phase on the fluid turbulence. The boundary conditions used in the simulations are described in the following paragraph. After 3 seconds, the nanoparticles are distributed well in fluid. Then solidification process is added in the model. The ultrasonic probe surface is set as velocity inlet, and the interface between liquid aluminium and air is pressure outlet. The other boundaries are set as wall. All of the discrete phase boundary condition types are set as fully reflective.

Modeling and Experimental Analyses

Figs. 2(b) and 2(c) show the distribution of the SiC nanoparticles into the molten alloys after 30s. As it can be seen from these figures, the nanoparticles are reasonably well dispersed in the molten alloys, especially when induction (EM) stirring is used. Some agglomeration of nanoparticles can be observed as well. Fig. 3 shows the SiC distribution in A356 samples reinforced with UST + 1% SiC nanoparticles during solidification. The SEM with EDS analysis has been performed to reveal the distribution of the SiC nanoparticles. It can also be seen from Fig. 3 that the SiC nanoparticles are dispersed reasonably well into the A356 matrix. Fine globular grain structures can be observed in the samples obtained using UST treatment during solidification. Also, it can be seen that the Si eutectic structure (black area) is more modified. The mechanisms for UST induced grain refinement are related to pressure and temperature oscillations in the melt. These oscillations introduced by UST are likely to induce heterogeneous nucleation in the melt. In addition, they are probably enhancing dendrite fragmentation. However, the dominant mechanism for grain refinement has not yet been determined.

There are several ways that UST can alter grain nucleation [16]. Since the pressure oscillations exist in a melt under UST processing, the liquidus temperature for the melt is changed. Thus, some part of the melt is superheated and the other part is undercooled. This phenomenon occurs at high frequencies and causes increase in the amount of nuclei into the melt. Dendrites usually start melting at the root due to local temperature rise and segregation. In the melt, UST produces strong convection and shock waves which may promote dendrite fragmentation. Convection promotes dendrite fragmentation because it causes local temperature and composition variations and promotes diffusion of solute. Shock waves will induce the breakage of the melting root [27].

Fig. 4(c) shows the optical micrograph at X50 magnification for A356 samples reinforced with SiC via UST. Similar fine globular structures can be seen in Figs 4(b) and 4(c), which means that the UST is the main approach that promoted the refinement of the primary phase. The eutectic phase appears to have similar morphologies in the samples treated with only UST and samples treated with UST and the addition of nanoparticles.

Conclusions

It was determined from this study that A356 dendritic structures were obtained without UST processing and globular grain structures were observed in the A356 samples obtained using UST processing. UST is able to break up the dendritic structures and form globular grain structures, resulting in better mechanical properties. Also, it was found that the Si eutectic structure is modified in the A356 samples processed by UST. By using SEM and EDS analyses, it was shown that the nanoparticles were dispersed reasonably well in the A356 matrix; but some insignificant agglomeration still occurred. This agglomeration could happen during UST processing or during solidification. It was also demonstrated that induction stirring can help improve the fluid flow characteristics during UST processing, making the nanoparticle distribution even more uniform.

Fig. 1 .

 1 Fig. 1. A schematic of the UST and induction furnace equipment.

Fig. 2 .Fig. 3 .

 23 Fig. 2. Model geometry and particle distribution without and with EM stirring after 30s.

Fig. 4 .

 4 Optical micrographs at X50 magnification for A356 samples (a) without UST treatment; (b) with UST treatment during solidification and (c) reinforced with 1% SiC nanoparticles.