domain: math spi

text:
Figure 1. Base maps of blast (black) and drill (red) measurements

Over this sub domain, the averaged copper grades of the blast and the drill holes are almost identical (around 0.6%). The variograms of blast and drill holes have similar behaviours (Figure 2), a high percentage of nugget effect (around 50%) and they differ mainly by their sills (0.12 for drill holes, 0.8 for blast holes), a comprehensible property as the blast support is larger.

Methodology

The geostatistical comparison between the two types of measurements is decomposed into two steps:

1 Deconvolution & Convolution:  Starting from the drill variogram, identifying the basic structures that model its behavior and deducing the underlying "point" variogram by deconvolution;  Making the theoretical convolution of the point variogram on 15-meter long supports and checking that it correctly fits the vertical and horizontal blast variograms, except for an additional nugget effect of 0.2.

Migration & Cross variogram 

As there is no point where both drill and blast measurements are known, we make some blast holes migrate to drill hole locations and calculate the cross variogram; 

The objective is to measure the nugget effect shared by the two types of measurements.

There are not enough drill samples to distinguish between horizontal and vertical drill variograms (they are drilled along many different directions). This is the first reason why an omnidirectional variogram will be considered for the drill samples, the second one is that all the formulas at our disposal require isotropy.

Consequently, we make two comparisons between:  An omnidirectional drill variogram and a vertical blast one;  An omnidirectional drill variogram and a horizontal blast one.

The distinction is important because the formulas differ between the two cases.

General formulas

All the formulas have been known for a long time in the literature, but in different places, and some are not even published. For the convolution charts, the most useful reference is probably 1 ; for the complete fundamental formulas, refer to 2 . Concerning the extension formulas, refer to [START_REF] Matheron | Les Variables Régionalisées et leur estimation[END_REF] .

In the following we apply a procedure illustrated in [START_REF] Séguret | Spatial sampling effect of laboratory practices in a porphyry copper deposit[END_REF] where we use the following approximation of a variogram regularized over a support "l" (the distance "h" being large in comparison with the dimension of the support):

γ ( ) γ( ) γ( , ) l h h l l  (1) with 2 00 1 γ( , ) ( ) ll l l u v dudv l     (2)
is the average of the point variogram when both extremities of vector h describe the support independently. In (2), 1D integrals are used because the core diameters are small compared to the lengths.

The way this formula is applied depends on the structure of the point variogram (spherical, exponential, linear, etc…) but also on the calculation direction compared with the regularization direction. In the following, we consider two situations:



The calculation direction is parallel to the regularization direction, notation // () l h  ;



The calculation direction is perpendicular to the regularization direction, notation

() l h   .
For the structures with a range, whether asymptotically (Exponential, Gaussian) or real (Spherical), we have:

range of γ (h) l = range of γ(h) + l (3)
Note that (3) is not compatible with approximation (1) which amounts to assigning to the regularized model the same range as that of the point model. So (1) is essentially useful for comparing the sills of regularized structures. The attenuation of the nugget effect, whether "pure" or associated with a microstructure which reaches its sill long before the first variogram lags, is proportional to the ratio of the supports. In the present case study, the diameters of drill holes and blast holes are considered to be equal and we ratio the lengths but generally speaking, one has to consider the ratio of the volumes:

Step 1: deconvolution & convolution

22 drill blast drill blast l l    (5)
With l drill =3, l blast =15, the nugget effect of the blasts must be five times smaller than that of the drills. Conclusion: the support cannot explain the differences between the nuggets of the blasts and of the drills. The blast nugget is too large.

Vertical variograms -Deconvolution & convolution

The calculation direction is parallel to the blast regularization direction ( i.e. vertical).

Exponential structure

If the practical drill range is 35m, the parameter associated with the underlying point exponential structure is expressed by (3):

0 0 3 35 3 10.7 a a     If ()
 denotes a variogram normalized by its sill, the underlying point sill C 0 of the exponential structure is produced by (1):

0 (1 (3,3)) drill C C   
For the exponential structure, the charts in 1 yield:

0 (3,3) 0.087 0.055 C    
For l=15m, we deduce:

15 0 (1 (15,15)) C C   
and we obtain:

15 (15,15) 0.34 0.036 C    
We will see later if these results correspond to the experimental blast variogram, but we must first look at the linear structure which completes the model (4).

Linear structure

For h>l we have, where b is the slope of the structure 1 :

( , ) 3
l l l b   (6)
The slope b, which does not change with the support, is given by the drill samples and the difference between two supports l and l' equals: 4), the attenuation is 0.00015, a negligible quantity. When l=3m and l'=15m, the attenuation is 0.0006, still negligible. In any case, the effect of the regularization on the linear structure is negligible. This is due to the weak slope of the linear structure.

' ( , ) ( ', ') 3 ll l l l l b      (
The combination of all the regularizations is shown in figure 4a where the dotted line represents the actual model and the red line the model we should obtain with a more realistic nugget effect. One can see that apart from the problem of the nugget effect, the variation range is acceptable, even if the linear part of the theoretical structure does not appear in the vertical experimental blast variogram.

Horizontal variograms -Deconvolution & convolution

The calculation direction is perpendicular to the blast regularization direction (i.e. horizontal).

The same procedure is followed, the only difference is that approximation (1) is not acceptable and we have to use charts that produce the exact calculation (see [START_REF] Journel | Mining Geostatistics[END_REF] , chart number 11). We obtain figure 4b where the dotted line represents the actual model and the red line the model we should obtain with a more realistic nugget effect. The fit is good.

First conclusions

If we omit the problem of the nugget effect, we see that both blast and drill holes can be considered as a regularization of the same reality according to their respective supports. This result, which we did not dare to hope, surprised us pleasantly and shows that the measurements from the blast holes are not as bad as people often think, anyway the case for this company. But the approach followed up to now suffers from two uncertainties:  The analyses are done independently. Imagine that all the blast locations have been shifted from a constant equal to the range (around 100m). In that case, the correlation between blast and drills will be zero while the same coherence properties are maintained when making individual regularizations as previously; 

The analyses refer to the drill nugget assumed to be a "natural" micro structure; is this true?

To answer these questions, cross variograms must be calculated but we do not have any location with both measurements, so a migration is necessary.

Step

Migration

In order to obtain a significant number of measurements at the same location, around 1,000 blasts samples were migrated to drill locations when the migration distance did not exceed 10 meters. Figure 5a presents the scatter diagram between the migrated values and the drill ones. The correlation coefficient is low (0.4) because the nugget effects are large.

Figure 5. (a) Scatter diagram between migrated blasts and drills; (b) Direct variogram of migrated blasts (black triangles), corresponding drills (red points) and cross-variogram of both (blue stars). The cross variogram reveals a tiny negative nugget effect with no comparison with the drill or blast ones

On Figure 5b, points (resp. triangles) present the migrated blast (resp. drill) variograms. They differ slightly from the previous ones because the number of samples is smaller and the migration affects the results. In the same figure, the stars represent the cross variogram which does not show a significant nugget effect, possibly a small negative one without any magnitude in common with the effects encountered on the individual variograms.

Conclusions

It seems that the drill-holes have their own errors too, independent of the blast ones, and the two measurements share only the structured parts of the variogram: the exponential and linear structures.

Analyse of the blast error

Description of the blast sampling

Up to now the theoretical blast support has been set to 15m but in fact the blast drilling length is approximately 17m, producing a large cone from the floor of which around 5cm of material is removed by hand across the entire surface, the idea being to restore an overall volume of 15m. Without any consideration of the numerous sampling procedures, we stay at this stage and ask the question: could the error specific to the blasts be due to the arbitrary removal of material and the blast length variability?

Randomization of the blast support

Let l and l' be two different supports. One finds in [START_REF] Matheron | Les Variables Régionalisées et leur estimation[END_REF] the formula which expresses the variance of the difference between the two grades Y over l and l', called "extension variance from l to l', also equal to twice the variogram between the grades averaged over the two supports: ( ',') ll  represent the averaged variogram when two points move independently along both the supports involved. Suppose that l and l' are randomly and independently selected uniformly in an interval, for example equal to [12.5m, 17.5m]. Then one has to calculate the mathematical expectation of (8) to obtain the resulting variogram. We have:

2 2 ' ' ' ' ( ( ) ( )) [( ( ) ( )) ]) 2 ( ) 2 ( , ) ( , ) ( ', ') l l l l ll h D Y x Y x h E Y x Y x h h l l l l l l              (
[ ( , )] [ ( ', ')] E l l E l l    (9) ' ' [ ( )] [ ( , )] [ ( , )] ll h E h E l l E l      (10)
is the theoretical variogram that we want to compare to the actual experimental variogram in order to verify if the blast nugget be associated with some support-length uncertainty.

' [ ( , )]

h E l l  is a continuous function, complex to calculate as it depends on the mutual configuration of l and l', but about which we know that for h greater than the range plus l, it reaches and stays at the sill of the underlying point variogram. In practice, the only structure that we consider is the exponential; its point sill is 0.055. For the interval of support-length uncertainty [12.5m, 17.5m],we deduce from (10) that the sill is reduced by a quantity obtained by:

17.5 12.5 1 [ ( , )] ( , ) E l l l l dl l     (11)
To evaluate the range of variations, the integral (11) is approximated by a finite sum:

13,14,

E l l l l      15,16,17 1 [ ( , )] ( , ) 5 l
We use the same charts as previously to calculate the values of ( , ) ll  involved and finally (12) yields: 1 [ ( , )] 0.055 (0.295 0.305 0.325 0.337 0.352) 0.055*0.323 5

E l l        (13)
Notice that even if we randomize the blast support over a larger interval still centered around 15m, the variance reduction does not change and stays approximately equal to the sill multiplied by 0.325. If we suppose that the support fluctuation is not symmetric around 15m, but around 13m for example, the multiplicative factor for the sill reduction decreases to 0.295. In any case, we conclude that:



The uncertainty on the support length does not produce a nugget effect but a variance reduction;  This variance reduction represents approximately 30% of the underlying variogram sill ; 

The arbitrary removal of the material, as well as the uncertainty on the blast length, cannot explain an error specific to the blasts and necessarily linked to the subsequent sampling procedures.

Summary: a formal link between blast and drill holes The model supposes that the blasts and the drills have the same average because the independent residuals are of zero mean. It must be verified when using this model. It is approximately the case here (0.63 for the drills, 0.69 for the blasts).

Removing the blast error by kriging

Model ( 14) can be used to remove the blast error by "Factorial Kriging" estimation [START_REF] Chilès | Geostatistics[END_REF] . One can easily build a linear system applicable to each blast measurement, choosing a local neighborhood of surrounding blast samples. The system is presented symbolically by using matrix formalism:

In this system, R  disappears from the second member of the linear system whereby we remove, from the estimation, the part associated with the measurement error. It does not mean that in the remaining part there is no nugget effect; it means that only the "natural" part remains. In our case, the complete nugget effect has to be removed because blasts and drills do not share any micro-structure.

Deconvolution by kriging

It may be interesting to remove the effect of regularization on the blast using a kriging system which estimates, for each blast measurement, a "point" value while simultaneously removing the part of the nugget effect associated with blast errors:

The difference with the previous system is that in the second member, (capital " ") is replaced by: 15 15 ( , , ) ( , , ) ( )

m x y z x y m z p h h h h h u p h u du         (small "p").

Block estimate by cokriging drill and blast measures

Finally, one can imagine locally renewing the mine planning block model by using blasts and drills together through a cokriging system with a linked mean (same average for both measurements): These systems were tested on a realistic simulation where the truth is known; they produce good results which will be published in the near future.

S. A.

Conclusion

In this depositand more generally, in this company (other test have been done), diamond drill hole grades and blast hole grades are consistent in the sense that, apart from the nugget effect, the structured part of their respective variograms follow the theoretical laws of regularization.

Concerning the nugget effects, we discover, by cross-analyses, that there is no natural micro-structure in the underlying point grade and the large nugget effects encountered on the variograms (approximately 50% of the variance for blasts and drills) are due to blast and drill measurement errors, independent of either measurement type.

The analysis of the blast error leads to the conclusion that the error is not due to the first step of the sampling procedure, it has to be found later in the process.

As a conclusion, some linear systems are proposed for removing the nugget effects from the data, reducing the effect of convolution and, more importantly, using blasts and drills together for the short-term mine planning. These systems, among numerous different potential ones, easy to demonstrate, result directly from the formal link established here between blast and drill holes. Before using these systems, the link must be verified by adhering to the methodology presented here.

Figure 2 .

 2 Figure 2. (a) Drill hole copper grade variogram; (b) Blast copper grade variogram. Three directions are represented, 45° North (N45), 135° North (N135), and vertical (D-90). Black continuous line is the isotropic variogram

FittingFigure 3 .

 3 Figure 3. Drill hole variogram fitting. Dotted line, the experimental curve; continuous line, the model

  three times smaller than the 0.03 value deduced from the blast variogram (Figure2b). If one takes the blast nugget effect as a reference, the drill nugget effect should be 0.15, a quantity above the local sill of the variogram and not realistic.

Figure 4 .

 4 Figure 4. In blue, the points of the experimental blast variogram; dotted line, the theoretical model for the blasts deduced from the drills; in red, the theoretical model with a more realistic nugget effect (a) Theoretical regularization parallel to the vertical blast variogram (b) Theoretical regularization perpendicular to the horizontal blast variogram.

  z the indicator function equal to 0 outside the interval and 1 inside it; R(x,y,z), a "white noise" residual statistically and spatially independent from Y(x,y,z) and representing the blast error The variogram of Y blast (x,y,z) becomes:

		15 p z ( ) m		1 15	15 [0, ] 2 1 (| |) z	;
		2 [0, ] 15 1 (| |)
					( ) h	 	( ) h	 	( ) h	(15)
			blast			15	m	R
	with
			() R h	R , the nugget effect due to the blast error, having the variance 2  ;
					( ) ( h   	P	)( ) ( h   	P	)(0)	;
		15	m			15	m	15	m
		() h  , the point variogram, assumed to be isotropic;
		( ) P h	( ) p p h  		1	( | | 15)1 (| |) h h  
		15	m			2	[0,15]
								15
	Formal link
	Finally, we have:
	with	15 Y x y z Y x y z p z R x y z ( , , ) ( , , ) ( ) ( , , ) blast m   	(14)
		Y(x,y,z), the point grade assumed to be isotropic and devoid of any measurement error;
		"*" denotes a convolution product;
								
		15 Y x y z p z ( , , ) ( ) m 		 	15 ( , , ) ( m Y x y u p z u du ) 	;
	S. A. Séguret "Blast & Drill holes..." Draft author. Final paper in: Proceedings of the 7th International Con-
	ference on Sampling and Blending , Kim H. Esbensen Editor, IMP publications, DOI 10.1255/t0sf.51

Acknowledgements

The author would like to thank CODELCO Chile, for its strong support in the implementation of good geostatistical practices along the copper business value chain. I warmly thank the Editor for his excellent recommendations which improved the manuscript. Many thanks to Sebastian de La Fuentes and the complete team of R&T. Without them, nothing would have been possible.


Rachid Deriche
email: rachid.deriche@inria.fr

Computational Brain Connectivity Mapping: A Core Health and Scientific Challenge

Keywords: Diffusion MRI, MEG, EEG, Structural and Functional Brain Connectivity, Brain Mapping, Networks

One third of the burden of all the diseases in Europe is due to problems caused by diseases affecting brain. Although exceptional progress have been obtained for exploring the brain during the past decades, it is still terra-incognita and calls for specific efforts in research to better understand its architecture and functioning. To take up this great challenge of modern science and to solve the limited view of the brain provided just by one imaging modality, this article advocates the idea developed in my research group of a global approach involving new generation of models for brain connectivity mapping and strong interactions between structural and functional connectivities. Capitalizing on the strengths of integrated and complementary non invasive imaging modalities such as diffusion Magnetic Resonance Imaging (dMRI) and Electro & Magneto-Encephalography (EEG & MEG) will contribute to achieve new frontiers for identifying and characterizing structural and functional brain connectivities and to provide a detailed mapping of the brain connectivity, both in space and time. Thus leading to an added clinical value for high impact diseases with new perspectives in computational neuro-imaging and cognitive neuroscience.

Introduction

In its reports published in [START_REF] Andlin-Sobocki | Costs of disorders of the brain in europe[END_REF][START_REF] Wittchen | Ecnp/ebc report 2011 -the size and burden of mental disorders and other disorders of the brain in europe 2010[END_REF], the European Brain Council (EBC) investigated the socio-economic impact of brain diseases on European society and concluded that about one third of the burden of all diseases in Europe is due to brain diseases. Every year, over a third of the total EU population suffers from mental disorders and the burden of brain diseases (neurological, neurosurgical and psychiatric diseases together) in Europe is estimated to cost e798 billion per year. Both EBC reports project that this burden will further increase in the coming decades due to the ageing of the European society and conclude that for now and in the near future there is a very strong societal and economical need for improvement in diagnosis and therapy of brain diseases. They also call for increased efforts in basic and clinical research and strongly advocate the need for new models and computational approaches in neuroscience and neuroimaging to better understand the brain architecture and functioning.

Mapping the human brain : a core health ambition of modern science

Mapping the brain architecture and functioning is a core health ambition of the 21st century and one of the greatest challenges of modern science. The human brain is one of the most impressive and complex structures in the known universe. With approximately 85 to 100 billion neurons, the human brain is definitely an object of prodigious complexity. However, more important than the number of neurons is the extraordinary complex brain circuitry formed by the 1.000 trillion possible synaptic connections between the neurons. Since the work of pioneers such as Wernicke and Dejerine, and following Geschwind's revolutionary disconnectionist framework, many important neurological diseases and disorders including Alzheimer's, Schizophrenia, Autism, Brain trauma, Epilepsy and others have been shown to be related to pathological alterations in the connectivity of the brain (Johansen-Berg and Behrens, 2009). This is why ambitious worldwide scientific efforts to understand the hierarchical, complex, architectural and functional network organization of the brain, are combined to conduct huge data generation projects in the USA (Human Connectome Project 1 and Brain Research through Advancing Innovative Neurotechnologies (BRAIN) Initiative 2 and models integration project in Europe (FET Flagship Human Brain Project 3 . Slightly smaller but equally ambitious complementary initiatives have been launched around the globe such as China's Brain Database to identify clues to tackling brain diseases and disorders, Japan's Brain/MINDS (Brain Mapping by Integrated Neurotechnologies for Disease Studies) project, Israel Brain Technologies, Brain Canada, Aus-Brain (Australia) to name just a few. These large initiatives and brain research projects differ in scope and detail but all have been developed to federate efforts to improve our understanding of the brain.

Due to the vast amount of complex data to acquire, broad range of scales at which they have to be considered, huge variety of problems to tackle and models to integrate, projects as HCP, HBP, BRAIN etc are only a first step and there is clearly a long way to go to achieve the aims of precisely understanding the structure and functioning of the brain. Although exceptional progress has been obtained for exploring the brain during the past decades, the contribution from the exact sciences in this endeavor is still relatively less than what it should be, not to say insufficient and more appropriate models and robust methods for identifying and characterizing structural and functional brain connectivity are still lacking while very much needed. By providing the right mathematical and computational tools as well as the rigorous methodologies, exact sciences may indeed offer substantial added value for exploring the brain, still terra incognita and last continent remaining to discover. In addition, although various new insights from individual imaging modalities have been obtained for brain diseases, it is very clear that an effective research approach is to use various noninvasive imaging modalities and capitalize on the strengths of each of them to solve the limited view of the brain provided separately by any imaging modality.

BCM : Brain Connectivity Mapping

During the last decade, huge progress has been made with non-invasive and in vivo medical imaging technologies, such as diffusion and functional magnetic resonance image (dMRI, fMRI) [START_REF] Buxton | Introduction to Functional Magnetic Resonance Imaging[END_REF] as well as electro and magneto-encephalography (EEG & MEG jointly denoted as M/EEG) to reconstruct the hierarchical complex structural and functional network organization of the brain. Multi-modal imaging, integrating functional (fMRI & M/EEG) with structural (dMRI) descriptions has the potential to yield a detailed picture of brain architecture and dynamics and dramatically improve our understanding of brain connectivities. A broad range of studies have shown links between anatomical and functional whole-brain connectivity. However, there has been limited work on a systematic framework to investigate whole-brain interactions between structure and function.

A large amount of work has already been devoted to joint analysis methods for dMRI and fMRI data, Hence, a significant review study of dMRI and fMRI data fusion methodologies and their applications in cognitive and clinical neurosciences can be found in [START_REF] Zhu | Fusing dti and fmri data: A survey of methods and applications[END_REF]. The large amount of relevant literature published on dMRI and fMRI has been estimated to include hundreds of papers (at the date of March 2013) with methods categorized into three classes: fMRI assists dMRI, dMRI assists fMRI, and joint dMRI and fMRI fusion. This interest in combining dMRI and fMRI will certainly continue to grow and I envision that effective multi-modalities fusion will play increasingly important roles in neuroimaging and brain sciences in the years to come. See for instance [START_REF] Castellanos | Clinical applications of the functional connectome[END_REF] where clinical applications of the functional connectome are considered with a focus on how resting-state fMRI methods can lead to biomarker identification for brain disorders and to [START_REF] Griffa | Structural connectomics in brain diseases[END_REF] where interesting highlights and outlooks of structural connectomics for clinical applications in brain diseases are discussed.

However, if a large amount of work has been devoted to joint analysis methods for dMRI and fMRI data, much less work, not to say none, has been devoted for dMRI and M/EEG. Therefore, to advance our current fragmentary and limited understanding of the brain, and solve the limited view of the brain provided just by one imaging modality, I push forward the idea to map the brain connectivity with new generation of models and well appropriate noninvasive neuro-imaging modalities such as diffusion Magnetic Resonance Imaging (dMRI) & Electro and Magneto-Encephalography (EEG & MEG jointly as M/EEG) methods. At Inria, my research group Athena actively contributes in the field of non-invasive and multimodal brain imaging and has acquired a remarkable joint expertise in the modeling and analysis of dMRI & M/EEG data, which is clearly fundamental to non-invasively examine connections through combinations of functional and anatomical imaging techniques.

During the last decade, my research group has pushed far forward the stateof-the-art in both dMRI & M/EEG, developing new computational models and tools, tackling and solving a set of important and challenging problems to unleash the full power and multivariate information content of dMRI & M/EEG data. Combining the measurements obtained by these two types imaging modalities has the potential of providing a detailed view both in space and time of the functioning brain at a macroscopic level. Capitalizing on the strengths of dMRI & M/EEG and building on the bio-physical and mathematical foundations of our models will bring a remarkable and significant added clinical value to identify and characterize brain connectivity and positively impact brain diseases. In that spirit, I summarize in the sequel some challenges, with a particular emphasis on those related to dMRI, and target areas for the upcoming years and the rationale behind them.

Advanced dMRI for structural connectivity

To better understand the brain, it is important first to focus on its structural connectivity and to rely on diffusion MRI (dMRI) which is the unique in vivo and non-invasive technique capable of providing the structural connectivity information and investigating the complex microstructure of the cerebral white matter. Introduced in the mid 80's, dMRI is well adapted to detect, characterize and quantify possible white matter abnormalities of brain tissues that cannot be revealed by standard imaging techniques [START_REF] Jones | Diffusion MRI: Theory, Methods, and Applications[END_REF]Johansen-Berg and Behrens, 2009).

To take up the challenge to precisely and accurately recover and characterize the structural brain connectivity and map it with ever-increasing detail and richness, it is important to develop generative and ground-breaking models for the advanced acquisition, processing and analysis of dMRI data. Overall, these last twenty years have seen an explosion of intensive scientific research which has vastly improved and literally changed the face of dMRI. These ground breaking changes range from improved and more powerful gradients, scanner-technology to acquisition, post-processing and reconstruction schemes, resulting in shorter acquisitions and more accurate modeling of the brain-tissues microstructure and connectivity. However, although great improvements have been made in the last twenty years [START_REF] Assaf | Axcaliber: a method for measuring axon diameter distribution from diffusion mri[END_REF][START_REF] Zhang | NODDI: practical in vivo neurite orientation dispersion and density imaging of the human brain[END_REF][START_REF] Descoteaux | Multiple q-shell diffusion propagator imaging[END_REF][START_REF] Ozarslan | Mean apparent propagator (MAP) MRI: a novel diffusion imaging method for mapping tissue microstructure[END_REF], major improvements are still required primarily to optimally acquire dMRI data, better understand the biophysics of the signal formation, recover invariant and intrinsic microstructure features, identify biophysically important bio-markers and improve tractography. For short, there is still considerable room for improvement to take dMRI from the benchside to the bedside.

Some recents results and debates in the dMRI community perfectly illustrate this importance and the vital need to develop ground-breaking models for recovering a precise and accurate cerebral connectivity. Indeed, recently, Van Wedeen & al. (Wedeen et al., 2012) have analyzed relationships of adjacency and crossing between cerebral fiber pathways in four nonhuman primate species and in humans by using dMRI. Their results suggest that the cerebral fiber pathways adhere to a simple and well organized three-dimensional geometric structure consisting of sheets of fibers that intersect at right angles like Manhattan-like grids. However, this point of view is not totally shared within the community, as pointed out by [START_REF] Catani | Comment on "the geometric structure of the brain fiber pathways[END_REF] who argue that Wedeen's method is likely to miss fibers crossing at angles less than about 70 degrees. They claim also that this view is clearly biased by the limits of their techniques, and conclude that our brain circuitry is wired more like London's chaotic tangle than Manhattan's organized grid. These results and debates clearly illustrate the challenge and complexity of recovering a precise and accurate cerebral connectivity, and the importance and the vital need to develop accurate and precise computational models and tools to help correctly answer such kind of concerns.

Today, it is well understood that more accurate models and better signal descriptions are necessary to overcome the limitations of these state-of-the-art signal reconstruction and tractography methods. Due to the ambiguities of low dMRI spatial resolution, complexities of the underlying tissue and uncertainties of signal noise, newer and more accurate dMRI models have to be found and developed. In my research group, we have recently contributed to design optimal acquisition schemes for single and multiple q-shell in diffusion MRI [START_REF] Deriche | Optimal real-time q-ball imaging using regularized kalman filtering with incremental orientation sets[END_REF][START_REF] Caruyer | Design of multishell sampling schemes with uniform coverage in diffusion mri[END_REF]. In particular, it is our optimal sampling scheme [START_REF] Caruyer | Design of multishell sampling schemes with uniform coverage in diffusion mri[END_REF] freely available for download4 [START_REF] Caruyer | Design of multishell sampling schemes with uniform coverage in diffusion mri[END_REF], which has been used in the HCP (Human Connectome Project). We have also exploited the ability of Compressed Sensing to recover the whole 3D dMRI signal and some of its important features from a limited number of samples [START_REF] Merlet | Continuous diffusion signal, eap and odf estimation via compressive sensing in diffusion mri[END_REF]. However, and up to now, in all our acquisitions schemes as well in those developed and used in the community, the gradient directions are optimized by considering only geometrical constraints i.e without taking into account any possible a-priori information related to the data and/or to the scanner. We believe it is important to revisite the important problem of the acquisition and tackle the challenging issue to develop more appropriate and efficient sampling schemes, in particular in some specific situations where we could have some prior knowledge on data and/or the used scanner.

The dMRI signal is highly complex, hence, the mathematical tools required for processing it have to be commensurate in their complexity. The dMRI signal is often analyzed partially (single q-shell or QBI) or more globally (DSI or multi-q-shell) and in various representations such as spherical coordinates or Cartesian coordinates. Hence, it is important to explore a number of new and sophisticated mathematical approaches to improve the quality, accuracy and efficiency of modeling the dMRI signal [START_REF] Ghosh | A survey of current trends in diffusion mri for structural brain connectivity[END_REF]. Historically, the first model in dMRI is the 2nd order diffusion tensor (DTI) [START_REF] Basser | Mr diffusion tensor spectroscopy and imaging[END_REF] which assumes the EAP to be Gaussian centered at the origin. While practical, this assumption limited the ability of the DTI to describe complex, singular and intricate fiber configurations (U-shape, kissing or crossing fibers). To overcome this limitation, so-called Diffusion Spectrum Imaging (DSI) and High Angular Resolution Diffusion Imaging (HARDI) methods such as Q-ball imaging and other multi-tensors, compartment models and high order tensor models were developed to resolve the orientationnality of more complicated fiber bundle configurations [START_REF] Schultz | Higher-order tensors in diffusion imaging[END_REF][START_REF] Jones | Diffusion MRI: Theory, Methods, and Applications[END_REF]Johansen-Berg and Behrens, 2009;[START_REF] Ghosh | Constrained diffusion kurtosis imaging using ternary quartics & mle[END_REF][START_REF] Ghosh | From second to higher order tensors in diffusion-mri[END_REF].

In particular, the exploration of high order tensors and the development of the theory of tensor invariants as a mathematical framework for computing new bio-markers for dMRI offer exciting and innovative research direction to better characterize white matter abnormalities [START_REF] Ghosh | From second to higher order tensors in diffusion-mri[END_REF][START_REF] Schultz | Higher-order tensors in diffusion imaging[END_REF]. High order tensors will also allow to go beyond simple bio-markers, such as those derived from the classical second order diffusion tensor. Indeed, as of today, DTI based bio-markers such as Mean Diffusivity (MD) and Fractional anisotropy (FA) are those clinically used at this date but are well known to be definitely sensitive to confounding factors such as partial volume and axonal dispersion, severely lack specificity and do not capture the subtle effects that might be early indicators of diseases. New biomakers from high order tensors will provide new insights to better understand the structural connectivity and help characterize important neurodegenerative diseases [START_REF] Papadopoulo | Complete set of invariants of a 4th order tensor: the 12 tasks of hardi from ternary quartics[END_REF][START_REF] Caruyer | On facilitating the use of {HARDI} in population studies by creating rotation-invariant markers[END_REF].

High order dMRI modeling allow also to infer and use micro-structural features, such as axonal diameter and density. Indeed, dMRI which is capable to measure the displacement of water molecules at the scale of tens of microns, is exquisitely sensitive to any change in tissue micro-structure and can provide unique insights and intrinsic parameters to precisely characterize white matter, in particular abnormalities that underpin diseases states. Tissue micro-structure parameters extracted from high order dMRI models could be used as important a-priori information in a unified and efficient micro-structure based tractography algorithm. This clearly pave the way to a new generation of tractography algorithms and opens new perspectives for better reconstruct anatomical brain connectivity and impact on certain neurodegenerative diseases thanks to a better connectivity quantification. [START_REF] Sherbondy | Microtrack: an algorithm for concurrent projectome and microstructure estimation[END_REF][START_REF] Girard | Axtract: Microstructure-driven tractography based on the ensemble average propagator[END_REF].

All these challenging scientific problems, when solved, will lead to direct measurements of important micro-structural based biomakers from high order model providing much greater and new insights to better understand the structural connectivity and help characterize important neurodegenerative disease mechanisms, recovery and development.

Advanced M/EEG for Functional Connectivity

A second important challenge is to complement brain structural connectivity with brain functional connectivity. Indeed, if a better knowledge of the structural wiring of the brain is essential to move ahead, it is clearly not sufficient to understand the functioning of the brain. For instance, if one consider the tiny worm C. elegans, one of the simplest organisms with a nervous system, we are far from fully understanding its brain even though, the wiring diagram of its 302 only neurons has been comprehensively mapped with a connectome shown to be a small-world network. dMRI provides only structural connectivity but, as of today, gives no clue on time-varying functional connectivity. Indeed, existing dMRI methods for mapping structural connectivity have provided numerous valuable insights but still remain inadequate to provide information about directionality. On the other hand, MEG and EEG (jointly denoted as M/EEG) provide (after source reconstruction) time courses of the activity of the cortical areas from which it is possible to compute dynamic connectivity measures. Indeed, M/EEG are two non-invasive neurophysiological techniques which directly record electrical voltages from electrodes placed on the scalp (EEG) and magnetic field from SQUID sensors placed above the head (MEG). Both are fundamentally related through Maxwell's equations to the distribution of dipole moment throughout the brain and head and hence have similar millisecond temporal resolution and are very important functional neuroimaging modalities for studying the dynamics of neural activities MEG and EEG can be measured simultaneously and reveal complementary properties of the electromagnetic fields of the brain [START_REF] Hämäläinen | Magnetoencephalography-theory, instrumentation, and applications to noninvasive studies of the working human brain[END_REF] for studying the dynamics of neural activities. Contrarily to fMRI, which "only" measures an haemodynamic response linked to the metabolic demand, MEG and EEG measure a direct consequence of the electrical activity of the brain.

While the few hundred time signals obtained using M/EEG have a clear clinical interest, the measurements are made outside, or on the surface of the head, and hence only provide partial information on the localization of the sources of activity. With a proper model of the head and the sources of this electromagnetic activity, it is possible to simulate the electrical propagation (forward problem [START_REF] Kybic | A common formalism for the integral formulations of the forward eeg problem[END_REF][START_REF] Gramfort | OpenMEEG: opensource software for quasistatic bioelectromagnetics[END_REF]) and to recover the sources corresponding to measurements using an inverse problem. Solving the inverse problem is the key to identify and localize brain areas of the observed activity [START_REF] Baillet | Electromagnetic brain mapping[END_REF] and developing ground-breaking tools and advanced M/EEG source reconstruction is one of the major objectives to achieve to push far forward the state-of-the-art in M/EEG. 6. Toward sructural and functional brain connectivity mapping dMRI provides a global and systematic view of the long-range structural connectivity within the whole brain. On the opposite, M/EEG measures brain activation over time and provides, after source reconstruction (solving the socalled inverse problem of source reconstruction), time courses of the activity of the cortical areas. Given its very high temporal resolution, the signal of reconstructed sources can be processed to reveal the functional connectivity between the nodes.

Such measures are however limited by source estimation methods that currently do not leverage any information on structural connectivity. Such information could improve functional connectivity by disambiguating multiple solutions corresponding to different hypotheses on how a given brain network operates. While dMRI and M/EEG have been the object of considerable research separately, there have been very few studies on combining the information they provide. Some existing studies deal with the localization of abnormal MEG signals, particularly in the case of epilepsy, and on studying the white matter fibers near the detected abnormal source but to our knowledge there are very few studies merging data coming both from M/EEG and dMRI at the analysis level [START_REF] Philippe | A nested cortex parcellation combining analysis of meg forward problem and diffusion mri tractography[END_REF][START_REF] Belaoucha | Using diffusion mri information in the maximum entropy on mean framework to solve meg/eeg inverse problem[END_REF][START_REF] Philippe | Propagation of epileptic spikes revealed by diffusion-based constrained meg source reconstruction[END_REF].

Combining the structural and functional information provided by dMRI and M/EEG is a difficult problem as the spatial and temporal scales of the two types of measures are extremely different. Nonetheless, combining the measurements obtained by these two techniques can potentially provide a detailed macroscopic view of the functioning brain both in space and time. All these considerations support the need to develop new models and computational tools for dMRI and M/EEG and integrate their information to reconstruct brain connectivities.

To complement the above challenges who's components are mainly concerned with developing sophisticated techniques for acquiring, processing and analyzing complex and multidimensional dMRI & M/EEG data, the different theoretical approaches to be developed above will have to be more closely integrated into a unified computational framework. This could be performed through the construction of a joint dynamical structural-functional network of the brain that will help to identify and characterize brain connectivity and see how far the developed framework can be get close to explaining high-level brain structure and function. Advanced network analysis and machine learning tools will have to be developed to design theoretically well founded methods and metrics for de-veloping, analyzing and validating this joint and dynamical anatomo-functional network, This challenge will lead to the creation of a robust and complete connectivity characterization of patients and controls. The joint structural-functional network developed will define a unique relationship between structure and function thanks to the directionality information and the edgewise quantification of connectivity grounded on fine spatial domain resolution, micro-structure features and high order tensors for dMRI and on exquisite time domain sensitivity and high order models for M/EEG.

Last but not least, the advanced integrated and post-hoc analysis tools for the joint structural-functional connectivity analysis are expected to help elucidate patterns of abnormality in the studied population of patients by providing a measure of abnormal connectivity that can then be correlated with clinical measures of symptom severity, aiding in prognosis.

Conclusion

In this article, I described some of the challenges recently solved and identified some of those still to be solved to achieve new frontiers in brain connectivity mapping through dMRI and M/EEG. Even if the described vision includes objectives for the upcoming years that can be considered as extremely challenging and very ambitious5 . I am confident that the use of structural and functional brain connectomics to identify outcome biomarkers for brain disorders will greatly extend the capabilities of neuroimaging in cognitive and clinical neurosciences and that important progress with significant scientific results and clinical impact will occur in support of these objectives. This will push forward and far beyond the state-of-the-art measurements and modeling approaches that span across the disciplines of computer sciences, applied mathematics, neuroimaging and health science.

http://www.humanconnectomeproject.org/

http://www.nih.gov/science/brain/

https://www.humanbrainproject.eu

http://www.emmanuelcaruyer.com/q-space-sampling.php

The work leading to the objectives listed in this article has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (ERC Advanced Grant agreement No

to start onSept. 1st, 2016.)


V Bukanin
email: vabukanin@nsgsoft.com

A Ivanov

A Zenkov

V Nemkov

Design of induction heating lines using ELTA program

Keywords: induction heating, ELTA program, heating line, computer simulation, forging

come

Design of induction heating lines using ELTA program

V Bukanin, A Ivanov, A Zenkov,

Introduction

In metallurgical and tube industries the induction lines are widely used for heating a variety of long products (rods, billets, tubes, slabs, strips) for heat treating or before hot deformation (rolling, forging, pressing). The parts are moving continuously or in small steps along the rails of being supported by driving rolls. These lines can contain many induction coils (several tens). The coils typically have the same length but their diameters, thermal insulation and turn numbers may be different. Induction coils are often connected in series and parallel group supplied from individual power sources with the same or different frequencies. Depending on techno-economical situation gas furnace may be installed before the induction line for energy savings or electrical resistance furnace may be installed at the end of the line for temperature equalization or for holding and flexibility. Design of such lines is a rather challenging task. It may be divided in two phases. The first phase is determination of the line length, power distribution along the line and optimal frequencies for different parts of the line, which should grand required technological specifications and good techno-economical results. This first stage may be called "the process design". The second phase is a detailed determination of individual coil parameters and their optimization as well as proper selection of the power supply and matching circuitry. It is obvious that computer simulation must be used at both phases of the project development. There are many programs for simulation of induction systems such as Flux 2D and 3D, Ansys Multiphysics, ThermNet, Magnet, Sysweld, Comsol, and others [START_REF] Nemkov | Modeling of induction hardening processes[END_REF]. These programs are mainly based on Finite Element (FE) method that allows the users describe and simulate 2D and 3D systems of complex geometry. They are expensive and require well-trained operator to run them effectively. The coil and part geometry in multi-stage induction lines is simple but the number of stages and the overall space is very big that creates difficulties in using general-purpose FE programs. The ideology of ELTA program is different. It is an engineering program which does not require special knowledge in computer simulation due to simple self-explanatory interface and preinstalled design options and database. Calculations are designed "in series steps", i.e. power and temperature is being calculated for each heating stage (inside the induction coils) and "cooling" stages (spaces between coils, quenching and holding zones). It make possible to simulate systems with very large number of steps, which is typical for induction lines. ELTA is at its best in optimal design of the multistage processes. The second phase (coil design) may be performed by the same program or with more complicate program (2DELTA, FLUX, etc.) if more detailed analysis of the coil parameters is required such as calculation of losses in individual turns, turn spacing variation, electrodynamic forces etc. Induction heating line is an object which may be relatively easy optimized [2 -5]. FDM is effectively used for many tasks of optimization of plane-parallel and cylindrical systems for induction heating. Simple guidelines based on 1D calculations give a good initial approach for the process optimization that must take into account the part geometry, its dimensions, required production rate, available equipment, etc.

Main features of ELTA program

ELTA is based on a combination of one-dimensional numerical (Finite Difference) approach and analytical account for finite lengths of the part and induction coil. Because of that simulation is very fast while providing good accuracy for systems of simple geometry typical for induction lines. ELTA 6.0 program has an option for 2D FD simulation of heating bodies with rectangular cross-section. One-dimensional equations describing electromagnetic field (1) and temperature (2) in cylindrical and flat bodies are:

H j R H R R R   0 ωμμ ρ 1       ) ( , H j x H x   0         ) ( , (1)
w x T x t T C v          ) ( , w R T R R R t T C v         ) λ ( 1 , (2)
were Rradius, xcoordinate, electrical resistivity, H  -magnetic field strength,  -angular frequency,  0magnetic permeability, v C -specific heat, T -temperature, t -time,  -thermal conductivity, w -volumetric power density. Two-dimensional non-linear differential equations for magnetic field H  and temperature T are:

, ωμμ ρ ρ 0 H j y H y x H x              ) ( ) ( (3)
where x, ycoordinates of workpiece cross-section.

. λ λ w y T y x T x t T C v              ) ( ) ( (4)
Final length of the system is taken into account using analytical procedure, called Total Flux method. It is based on composing the magnetic substitution circuit for a system "inductor -workpiece" [START_REF] Nemkov | Theory and calculation of induction Heating Devices[END_REF]. Multiple tests show that it gives good practical results in simulation of 2D systems of simple geometry. In some cases even 3D systems such as heating of slabs in oval or rectangular inductors can be simulated with good accuracy [START_REF] Nemkov | Simulation of induction heating of slabs using ELTA 6[END_REF].

The program can simulate the process of heating using an option "power density on the part surface" without accurate description of the coil design, which is very convenient for the first stage of development (Process Design). In detailed coil design (phase two), the program allows to make simulation of selected coils with different preinstalled power supplying circuits at different regimes (generator or coil power, current or voltage, preinstalled or variable frequency, etc.). Hydraulic calculation of the coil cooling is another useful engineering feature of the program.

Design of induction line for heating square billets

Many induction and combined heating lines have been investigated and designed using ELTA program. All technologies are industrial applications for metallurgical and forging plants. The most interesting variant of real technology is presented bellow. A plant had a forging line with gas furnace that required deep modernization. The problem was to modify the line for billets made of a special steel, sensitive to cracking in temperature range 20 -600 °С. Final temperature must be around 1250 °С, which was higher than the gas furnace could provide for a required production rate of 10 t/hr. It was decided to make a heating line containing induction heater to preheat billets up to 500 -600 °С with controlled temperature gradient and increased rate, gas furnace for further heating up to 800 -900 °С and again induction heating for the final heating. It was very difficult task in comparison with traditional heating of round billets.

Billet dimensions: square 8.3×8.3 cm, length l 2 = 91 -112 cm.

Processing parameters: required heating speed is limited to 1 -1.5 K/sec. This task has many parameters that don't allow us to use optimization methods or intuitively predict design of the line.

Operator-guided design was used with expert's evaluation of the results. The best design was selected after several iterations. After preliminary calculations the most promising variants were selected and more detailed analysis of techno-economical analysis performed. First variant. It is natural to start from traditional type of induction heater, i.e. a continuous heater with induction coils of square configuration (Fig. 1). Billets are moving through the inductors in longitudinal direction on a rolling system. Post-heating induction line may be of the same type as a pre-heating one.

Inductor dimensions: the "window" is 2X 1 ×b 1 =13.0×13.0 cm, length l 1 = 29 cm, coils turn number 12 from copper tubing A×T×d = 2×1×0.2 cm. Thermal insulation -Portland cement concrete 2 cm thickness.

Calculations showed that heating up to mean volume temperature of 563 C (see Fig. 3) with speed 1.45 K/sec requires 45 induction coils and total length of line 19.5 m. Total power and energy consumption are 1185 kW and 100.5 kWhr/t respectively. The values of electrical efficiency are 83 -87 %. Total time of heating in induction line is 374 sec. Main attention was paid to edge effect and temperature distribution in the volume of billet. Using the principle "uniformity in large" we can find temperature differences ΔT = (T aver@edge -T aver@centre ). These values after induction and gas heating are 80 K and 15 K respectively (Fig. 456). In the first variant of longitudinal heating edges of square cross-section are overheated (Fig. 5,6) but ΔT equalizes quickly during transportation and in gas furnace. Second variant. In order to reduce the line length, it was decided to use 3 oval inductors for pre-heating line (Fig. 2). All billets are moving through the inductors semi-continuously in transversal direction by using a pusher. It is more convenient than in the first case. There are 14 billets in each inductor, i.e. total length of the load is 14×8.3 = 116 cm.

Post-heating Gas

A -A

Inductor dimensions: the "window" is 2X 1 ×b 1 =13.0×125 cm, length l 1 = 120 cm, coil turn number 26 from tubing A×T×d = 4×2×0.2 cm. Thermal insulation -Portland cement concrete 2 cm thickness. Distributions of temperature in cross-section of the billet during heating in these lines are shown in Fig. 7 and Fig. 8. Mean volume temperature after 860 sec of induction pre-heating is 556 C. Total power is 1065 kW. Specific energy of induction pre-heating is 99 kWhr/t. The length of line decreases sufficiently from 19 to 4 m. In transverse heating the ends of the billets are overheated more than for the first variant (Fig. 8) and equalization will requires much longer time due to a big distance. In this case the electromagnetic edge effect of the load plays determining role. Simulation of the whole heating line showed that by means of an optimal design of pre-and postfurnace induction lines it is possible to reach required temperature level of 1250 C and uniformity of billets at the end of the whole line using both the first and second variants. All complex induction-gas-induction line, including the second variant of pre-heating induction heater, was successfully designed by V. P. Vologdin Scientific Research Institute of High-Frequency Currents (St. Petersburg) for one of the metallurgical plant. Verification of the calculation results, i.e. integral parameters of inductors and surface temperature, has been performed in laboratory and at the plant. It showed good matching. Distribution of power source and temperature in the cross section of square and rectangular work pieces have been verified also on the base of analytical and numerical methods, as shown for example in [START_REF] Nemkov | Simulation of induction heating of slabs using ELTA 6[END_REF].

Fig. 1 .Fig. 2 .Fig. 3 .Fig. 5 .Fig. 6 .

 12356 Fig. 1. Sketch of longitudinal pre-heating induction line + gas furnace + post-heating induction line

Fig. 6 .

 6 Fig. 6. Temperature variation in the cross-section during the time


B Paya

B De Kepper

00 kW heat recovery test bench to collect induction losses

Keywords: Induction heating, induction losses, energy recovery, Predictive Functional Control, PFC

On industrial induction heating devices, induction coils are water cooled at a high flow rate: this large amount of energy (10 to 60 % of the total energy consumption) is often lost. A heat recovery test bench has been built to collect the coil losses and reuse them when needed. The main parts of the test bench are made up two heat exchangers which can be by-passed according to the control and regulation system. The "recovery heat exchanger" collects the heat to be recovered and transfers it to the reuse, the recovery circuit. The "waste heat exchanger" evacuates the remained energy to ensure a safe use of the induction heating device. These two exchangers are designed for a 100 kW power transfer. A "Predictive Functional Control" (PFC) regulation has been implemented to drive the bench. This system is more efficient compared to a conventional PID controller. Thanks to a precise modeling of the bench parts and the induction device, it is possible to anticipate disturbances such as induction power fluctuations and to adapt the recovering in order to always work in a safe way. It was then possible to recover 66 % of the energy lost in the coil with a water output temperature of 85 to 90 °C.

Introduction: recovering the Joule losses of an induction heating device

On industrial induction heating devices (Fig. 1), induction coils are water cooled at a high flow rate: the water temperature increase is then too small to recover the transported energy. Nothing is done today to recover this large amount of energy that can reach 10 to 60 % of the total energy consumption because of a bad investment payback [START_REF] Ahn | Récupération de chaleur dans les installations des fours à induction à creuset[END_REF].

Fig. 1: Global view of an induction heating device and its cooling system Energy consumption analysis in many industrial sectors show that thermal needs at higher temperature (typically 90 to 140 °C) makes the energy recovery more relevant at this level in order to reduce the energy bill [START_REF]Applications of Industrial Heat Pumps[END_REF]. Preliminary studies on a billet heating configuration proved that a water output temperature of 90°C can be reached in safe conditions and without any efficiency reduction [START_REF] Paya | Investigation on energy recovering in an induction coil for continuous through heating billet[END_REF], [START_REF] Paya | Recovering energy in an induction coil: Impact on the water and coil temperature[END_REF]. In the frame of the French joint research project ISIS with financial support of the French Research National Agency (ANR), EDF has built a heat recovery test bench at industrial scale. This paper presents the main results obtained with this equipment.

Network 50 Hz

Frequency converter Capacitors

Cooler

Water cooling circuit

Inductor and heated piece

Transformer

Heat recovery bench

Network 50 Hz

Frequency converter Capacitors

Cooler

Water cooling circuit

Inductor and heated piece

Transformer

Heat recovery bench Presentation of the heat recovery test bench

The main objective of the heat recovery test bench is to be able to recover as much heat as possible without disturbing hot metal production from the induction heating device. That means that induction losses due to hot metal production have to be evacuated if they cannot be reused or stored. To reach this goal, the main parts of the test bench (Fig. 2) are made up two heat exchangers which can be by-passed according to the control and regulation system. The "recovery heat exchanger" collects the heat to be recovered and transfers it to the reuse, the recovery circuit. The hot water tank can play either the rule of heat storage or heat use by acting on the output valve. The "waste heat exchanger" evacuates the remained energy to ensure a safe use of the induction heating device. These two exchangers are designed for a 100 kW power tranfer. Several sensors are put in the various circuits to control the bench: flow meters, temperature, pressure. The two-way and three-way valves can be controlled either manually or by the regulation system. Also are the pumps of the inductor or recovery circuits. The regulation is realized by a programmable controller Schneider M340. First experiments used the "Proportional Integral Derivate" (PID) controller. But the strong interactions between all the valves made the system difficult to adjust for a smooth behavior. The use of the "Predictive Functional Control" (PFC) strategy appears to be more relevant for a safe control of the bench.

Principles of the "Predictive Functional Control" regulation

"Predictive Functional Control" (PFC) regulation [START_REF] Richalet | Model Predictive Heuristic Control: Application to industrial processes[END_REF] is an advanced method of process control using a model of the pieces of equipment to be controlled to predict the evolution of the controlled variables respect to the control inputs. At each timeslot, the current plant state is sampled. An optimized control strategy is then computed to reach the target in a specific time-horizon. This computation takes into account the behavior of the modeled pieces of equipment to predict the future values of the controlled variables. Finally, the first step the controlled strategy is applied to the system and the resulting plant state is sampled again. The calculations are repeated again with a shifted forward horizon. So, PFC is often called receding horizon control.

Implementation in the heat recovery test bench Modeling the various elements

The heart of the process control is the modeling of the full system; it requires describing the behavior of each component constituting the bench. The two-and three-way valves are characterized using an experimental approach. A set point step is given to the valve and the corresponding response is recorded. Assuming a first order behavior, it is then possible to identify the key parameters: delay, gain and time constant.

The behavior of the heat exchangers is strongly dependant on the temperature and flow of the two fluids exchanging heat. A characterization by experimental approach is not relevant and we preferred a numerical modeling. For each input value (flow rate, input temperature), we calculate first the heat exchange (power transferred, output temperatures) thanks to the thermo-hydraulic equations commonly used for designing this kind of equipment. Considering a first order behavior without delay, we can then evaluate the gain and the time constant for each set of input values. These last parameters are then approximated by polynomial curves.

Comparison with the PID control strategy

Experiments were made during the programming of the PFC regulation to evaluate the interest of this strategy. Fig. 3 shows the difference between the two strategies on the valve controlling the flow through the primary circuit of the recovery heat exchanger. During the first 18 minutes, the valve is controlled by the PID regulation. We observe strong oscillations on the valve and the primary flow. When switching to the PFC regulation, the control softens and reaches a relative stability. The remained oscillations were due to the waste exchanger by-pass which was still PID regulated.

Fig. 3: Comparison between PID and PFC regulation

Experimentations

The objective of the experimentation is to prove the robustness of the bench regards to the events which can occur during a production process. As we use a coil without a heated load, the inverter power data collected directly from the power supply corresponds to the coil Joule losses (blue curves below) but these measures are not precise. For a better accuracy, we preferred to use the calorimetric balance (1) in both inductor and recovery circuit.

𝑃 𝑐𝑖𝑟𝑐𝑢𝑖𝑡 = 𝑚 𝑐𝑖𝑟𝑐𝑢𝑖𝑡 • 𝐶 𝑝 • 𝜃 ℎ𝑜𝑡 -𝜃 𝑐𝑜𝑙𝑑 (1)
Where 𝑚 𝑐𝑖𝑟𝑐𝑢𝑖𝑡 is the flow rate of the circuit, 𝐶 𝑝 , the water specific heat and 𝜃 ℎ𝑜𝑡 -𝜃 𝑐𝑜𝑙𝑑 , the temperature increase at the inductor level or at the secondary of the recovery heat exchanger. The recovery yield is then evaluated as the ratio between these two powers. Two typical scenarios are presented here.

In the first scenario (Fig. 4), the production of the induction heater is constant (100 kW Joule losses) and no use is needed in the workshop. In that case, the losses are stored into the water tank until it is full. The blue curve corresponds to the Joule losses in the coil and the red one, the coil output temperature. The brown curve represents the power collected by the recovery circuit and stored into the tank. During the first 11 minutes, the losses are mainly stored into the inductor circuit and the inductor output temperature increases until its set point. After 11 minutes, the by-pass valve transfers the hot water from the coil to the recovery heat exchanger and the storage begins. After 20 minutes, the storage becomes less efficient and the remained energy is transferred to the waste circuit. After 29 minutes, the storage is full and all the losses are evacuated through the waste heat exchanger. The recovery yield is evaluated to 75 % at the end of the storage phase.

The second scenario (Fig. 5) represents a cycling production of a batch heating with three power steps (100 kW losses, then 80 kW and finally 50 kW corresponding to the blue curve Fig. 5) and a last step without energy consumption simulating the furnace unloading and loading. At the same time, a constant flow rate of hot water is dumped from the water tank. After a period of energy storage (not shown in the graph), we observe that the output temperature of the hot

PID PFC

water coming out of the tank is almost flat (purple curve) at a level of 70 °C. The output inductor temperature (red curve) is under control especially during the high power phases. The safety of the heating device is assured. The mean energy recovery yield along a cycle is estimated to 66 %.

Conclusion: a heat recovery bench ready for industrial use

We have developed a 100 kW heat recovery test bench able to collect the losses coming from an induction heating device. Its recovery yield reaches 66 % during a production cycle and reuse of the collected energy and up to 75 % during the heat storage. The hot water output can reach a level of 90 °C. The PFC regulation implemented to drive the bench makes the regulation robust and easy to adjust. It requires an advanced modeling of the different parts of the unit to make the prediction accurate. The use of two heat exchangers allows a safe work of the induction heater: the losses are either reused or evacuated when the storage is full. This test bench is proposed for insertion in an existing heating facility, for example a 400 to 500 kW steel billet heating device in a forging workshop.

Fig. 2 :

 2 Fig. 2: EDF's heat recovery test bench

Fig. 4 : 5 :

 45 Fig. 4: Energy storage in the water tank Fig. 5: Cycling production and use of recovery hot water

Acknowledgment

This work has been supported by French Research National Agency (ANR) Trough "Éfficacité énergétique et réduction des Émissions de CO 2 dans les Systèmes Industriels" program (project ISIS ANR-09-EESI-004).


A Dietrich
email: dietrich@etp.uni-hannover.de

H Schülbe

B Nacke

Detailed 3-D numerical investigation of continuous longitudinal flux heating for thin strips

Keywords: Induction, longitudinal flux heating, strip heating, numerical simulation

Most of the studies on longitudinal flux heating (LFH) deal with special applications and a strip thickness above 15 mm. The aim of this paper is to present the dependence of different geometric parameters for thinner strips. The influence on the temperature distribution and the efficiency will be systematically analyzed. The length, width and the air gap of the induction coil are subjects of the investigation. For the study a 3D numerical model is used.

Introduction

The industry has a big demand on efficient and flexible heat processing of metal strips. Therefore the LFH becomes more and more interesting in comparison to the conventional heating. It is important to have a deep knowledge of the process, because a homogenous temperature distribution is not always given and various things have to be observed. In past studies different aspects of LFH have been examined. In [START_REF] Rudnev | Handbook of Induction Heating[END_REF] the edge effects and end effects were investigated. The industrial application of LFH was studied in [START_REF] Langejürgen | Comparison of different induction heating approaches for hot strip mills in endless operation by using numerical simulation[END_REF] with the main view on slabs with a thickness of 16 mm and width of 1550 mm. The paper [START_REF] Blinov | Numerical simulation and investigation of induction through-heaters in dynamic operation mode[END_REF] presents the influence of the slab geometry of the temperature distribution and a comparison between LFH and transversal flux heating. All named studies refer to plate thicknesses over 10 mm. The nonlinearity of the ferromagnetic property of steel was discussed in [START_REF] Jestremski | Einfluss stark nichtlinearer Materialeigenschaften auf die induktive Längsfelderwärmung[END_REF] for strips with a thickness under 10 mm.

The numerical model

For the numerical modeling and calculation the commercial software package ANSYS 16 was used, which bases on the finite-element-method (FEM). The calculation algorithm which is used in this study is divided in two sub-steps. In each iteration step a harmonic electromagnetic and a transient thermal calculation is performed. To solve the Maxwell Equation's for the harmonic electromagnetic field the magnetic vector potential is used. The result of this solution is the heat sources distribution in the work piece. In the second sub-step the heat transfer equation is solved for a time interval Δt regarding thermal losses due to radiation and convection. The solution of the electromagnetic calculation is the input value for the heat transfer equation. The output value is the temperature distribution. In the next iteration step the temperature distribution is the input value for the harmonic electromagnetic calculation. The material properties specific electrical conductivity, relative permeability, specific thermal conductivity and the specific heat capacity are adjusted on the new temperature distribution. To realize the moving of the metal strip an algorithm is used, which was developed at the Institute of Electrotechnology. In order to shorten the computation time, the geometry is not shifted, because this would mean to mesh the geometry in every iteration step. The temperature distribution is stored in a matrix and in each new iteration step the temperature is written to the next node. For the calculation of the following results a 3-dimensional model is used. Because of symmetry only the half of the system is modeled to reduced the calculation time. The influence of the slope of the induction coil on the electromagnetic and thermal values is not considered. The induction coil has a rectangle cross section. For an exact calculation it is important to mesh the edges of the strip very fine. The calculations are based on the material properties of C45E steel [START_REF]VDI 3128: Induktive Erwärmung, Physikalische Stoffeigenschaften der Metalle[END_REF]. The reference thickness of the strip is 1 mm, because in many industrial application strips with a thickness in this range are used. The frequency is 38.5 kHz. The width of one winding is 30 mm and the distance between the windings is 8 mm to 10 mm.

Influence of the length of the induction coil on the temperature distribution and the efficiency

The length of the induction coil is an essential factor, because it influences the homogeneity of the temperature distribution. Is there not enough space for the required induction coil length there are other possibilities to take influence on temperature distribution. The length of the induction coil will be varied from 150 mm to 600 mm. The air gap is 40 mm. To ensure the comparability of the results, the power is kept constant in the metal strip by adjusting the current this means the temperature level is constant, too. Fig. 1 shows the temperature distribution from the middle of the strip to the edge of the strip. The four graphs have the same qualitative behavior. With increasing distance from the center of the strip to around 50 mm to the edge the temperature drops down, then the temperature rises to a maximum in 10 mm to the edge and then fall down again. The first drop results in a lower current density in this area. The overheating can be explained as follows: On the one hand the proximity effect works. It means that the currents in the induction coil and in the strip with opposite directions flow as close as possible. Then the inductance of the system is minimal. On the other hand the current in the strip wants to take the shortest way to reduce the electrical resistance. On the basis of these two effects the current is compressed and this results in a high current density. The second drop appears because most of the current does not flow to the edge of the strip. The dependence of the length of the induction coil results in the same qualitative behavior but in different quantitative temperature distribution. With increasing of the length of the induction coil the temperature distribution becomes more homogeneous. Tab. 1 presents the percentage difference between the minimum and maximum temperatures.

The efficiency depends on the length of the induction coil (see Fig. 2) and increases with a longer coil. The current for a 150 mm long induction coil (4 turns) is 1300 A and for a 600 mm long induction coil (16 turns) is 520 A. The electrical resistance increases with the length of the coil. The resistance for 16 turns is approximately 4 times greater than for 4 turns. On the other hand a lower current in the induction coil is required for 16 turns. In an ideal case with a coupling factor of 1.0 the current for the longer coil has to be divided by 2. That means a current of 650 A in a 16 turn coil generates the same power in the induction coil and the work piece as a current of 1300 A in a 4 turn coil. In other words, to cause the same power in the work piece and the coil (600 mm) the magnetic field intensity has to be two times smaller. In the numerical calculation the current is smaller than 650 A. That means the part of the magnetic stray field is smaller and the coupling factor is greater than for a short induction coil. The 3-dimensional power density distribution in the strip is pictured in Fig. 3 and Fig. 4 for two different lengths. The power density distribution is not symmetrical because the moving strip at the end of the induction coil has a higher temperature than at the beginning. As a consequence of this the specific electrical conductivity at the end of the coil is higher which results in different power losses. A comparison of the temperature curves with the power density is in a good agreement, because heat conduction plays in continuous heating processes a minor role, if the velocity of the strip is high enough. But at the edge of the strip is a big difference of the power density distribution in a small area. In this area the temperature distribution is smoother than the power density distribution. Because of a big gradient in this area thermal conduction works to equalize the temperature.

Influence of the distance between the induction coil normal to the strip on the temperature distribution and the efficiency

Another important geometric parameter is the air gap copper to copper normal to the strip surface, because the coupling between induction coil and the strip is influenced. In order to investigate the influence to the temperature profile the air gap is varied among 10 mm to 120 mm. The length of the induction coil (600 mm) and the width (half model 150 mm) of the strip are constant. The power in the strip is kept constant by adjusting the current. The influence of the air gap on the temperature is shown in Fig. 5. An increasing of the air gap results in a more homogeneous temperature distribution and undercooled edges. A smaller air gap results in an undercooled area between 0.025 m to 0.125 m from the strip center and an overheated edge. With a bigger air gap the Proximity-Effect has not so much influence of the current.

The efficiency (see Fig. 6) falls with a larger air gap because the magnetic stray field is bigger. That means a smaller part of the magnetic is coupled with the strip. The advantage of a large air gap is the more homogenous temperature distribution in the middle part of the strip but the edges are undercooled and the efficiency falls. As a solution for the undercooled edges magnetic flux concentrators can be used.

Influence of the width of the induction coil and the strip on the temperature distribution and the efficiency

The length of the induction coil (600 mm) and the air gap (40 mm) are fixed. The end temperature is in every calculation the same. The width of the strip and the coil are varied. A smaller ratio of width to length results in undercooled edges (see Fig. 7). With the increasing of the ratio the curves change their shapes and a temperature maximum arises. The efficiency increases (see Fig. 8) with a bigger ratio (width to length). The significant change of the efficiency has its reason in the power losses of the induction coil. Generally at the corners of the induction coil the power losses are higher as in the straight parts. If the straight part becomes bigger the power losses in the corners are relative smaller to the whole losses in the induction coil.

Conclusion

A 3-dimensional numerical model with a coupled electromagnetic and thermal calculation for LFH was developed. The moving of the strip was realized by shifting the node temperatures. The dependence of changing geometry parameters on the temperature distribution and on the efficiency was showed for different cases. A longer induction coil and a bigger air gap result in a more homogenous temperature distribution. With a longer induction coil the efficiency increases, with a larger airgap the efficiency decreases. Furthermore the ratio of length to width was investigated. The efficiency arises with a wider induction coil. In further investigations the numerical model should be expanded to a full 3D model so that the slope between the winding can be considered. Afterwards the numerical results should be verified by experiments.

1 :Fig. 1 :Fig. 2 :Fig. 3 :Fig. 4 :

 11234 Fig. 1: Temperature distribution normal to the direction of movement in dependence of the length of the induction coil

Fig. 6 :

 6 Fig. 6: Efficiency in dependence of the air gap Fig. 5: Temperature distribution dependence of the air gap

Fig. 7 :

 7 Fig. 7: Temperature distribution in dependence of the width of the induction coil


Emmanuel Abiona

Hongliang Yang
email: hongliang.yang@se.abb.com

Rajneesh Chaudhary

Ravi Kumar Kandasamy

Jan-Erik Eriksson

Development of Plasma Heating and Electromagnetic Stirring in Tundish

Keywords: tundish, induction heating, plasma heating, electromagnetic stirring

Plasma heating is advantageous over induction heating in the aspects of less impact on tundish design and lower maintenance costs.On the other hand, plasma heating has a lower heating efficiency of 60% compared to that of induction heating 90%. The idea of utilizing electromagnetic stirringwas numerically investigated in combination with plasma heating. With plasma heating alone, the flow in the upper part of heating chamber is dampened due to the buoyancy force caused by the large temperature gradient of plasma heating. With EMS, one rotational flow with the velocity of 0.2-0.4 m/s is created in the heating chamber, and the temperature in the heating chamber is homogenized. The temperature response time in the outlet is reduced and heating efficiency of plasma heating can be improved with EMS.

Introduction

Tundish heating has existed for almost 30 years and it shows the capability foraccurately controlling the casting temperaturevariation in the tundish. Historically, there are mainly two differing technologies being utilized for tundish heating, namely induction heating [START_REF] Asea | Tundish heating system brochur[END_REF] [START_REF] Kimura | Innovative Technologies in Continuous Casting Tundish[END_REF] and plasma heating [START_REF] Kittaka | Twin-torch type tundish plasma heater "NS-Plasma II" for continuous caster[END_REF] . Tundish heating is closely related to other aspects in tundish, such as flow control, refractory material, stirring and so on. Tundish heating has not received wide application worldwide except in Japan. A short comparison between plasma heating and tundish heating is listed below -Plasma heating is easy to install in a tundish, especially for the existing tundish where the induction heating is difficult to install.

-

The metallurgical benefits of induction heating and plasma heating are quite comparable; the temperature control can almost reach the same level with the both technologies. -Improvements of steel cleanliness has been reported with both technologies, yet with different principles. Accurate control of temperature and flow pattern in the tundish are crucial for improving steel cleanliness.

-

The heating efficiency with induction heating (90%) is generally higher than with plasma heating (60%~70%). On the other hand, the potential for improving the plasma heating efficiency is promising with appropriate stirring of the molten steel.

In this paper, Plasma heating is chosen to make further developments due to its less impact on tundish design and lower maintenance costs. Numerical simulation was made to investigate the heat transfer in the melt with plasma heating, and how to improve it with electromagnetic stirring. It has been reported of applyinggas stirring or mechanical configuration with dam/weir to improve the heat transfer of plasma heating, electromagnetic stirring is proposed in this paper due to the following advantages: -Better reliability: the stirrer has no contact with tundish melt, and can be operated independently. -Superior repeatability: flow in the tundish can be controlled to have a constant flow pattern, irrespective the melt temperature or the refractory conditions.

Case description

A 30 ton tundish was chosen for the case study, the dimension of which are shown in Fig. 1. The basic parameters are shown in Table 1. The initial inlet temperature is 1565 o C and the inlet temperature is reduced by 0.5 o C/min during one ladle cycle of 50 min. Throughput of the melt steel is 2.5 ton/min, with two outlets. The initial melt temperature in the tundish is set as 1540 o C/min, which is the same temperature as the end of the inlet temperaturefrom the previous ladle cycle. The tundish is a T-type tundish, divided into two parts as inlet chamber and outlet chamber.The weirs (two weirs in the outlet chamber, one weir between theinlet chamber and outlet chamber) form the heating chamber for plasma heating.

The two dams constrain the shortcut flow from the heating chamber to the outlets.The max heating power is set as 800 kW, heating efficiency is treated as constant 60% for all the cases and the heating power is evenly distributed on the surface area of heating chamber. One electromagnetic stirrer is mounted outside the heating chamber. Electromagnetic stirring creates one downward stirring force along the tundish wall.This stirring force agitates a rotational flow in the heating chamber, which in turn homogenizes the temperature and improves the heat transfer in the melt.

The effect of temperature on steel density is considered in the model.This is especially important due to the large temperature gradient caused by plasma heating. The electromagnetic calculation of EMS is made with Vector Field/Opera and the flow simulation is made with ANSYS/Fluent 16. Three cases are simulated as listed in Table 2.

Simulation results

Fig2 shows the flow field in tundish with case 1, 2 and 3 at the end of ladle cycle, i.e. 50 minutes after the start of simulation.

Comparedtocases 1 and 2, it can be seen that flow speed in the upper part of the heat chamber is drastically reduced with plasma heating.The reason can be found in Fig, where it shows that the temperature in the upper part of the heating chamber is very high due to the heating-up with plasma. The highest temperature in the melting surface can be 1900 o C according the numerical simulation. The high temperature leads to the large buoyancy force, which counteracts the convection flow caused by the inlet stream from the ladle. The low flow speed will lead to a slow heat transfer from the plasma to the melt. In practice, it may mean that the heating efficiency of plasma torch can be low, and the erosion of refractory material on top of the heating chamber is increased. With electromagnetic stirring, it can be seen both in Fig. 2 and Fig that a rotational flow is created in the heating chamber.With the speed between 0.2 m/sec and 0.4 m/sec, the vertical temperature gradient in the heating chamber has almost disappeared, but the average temperature outside the heating chamber is higher. The dams stop the shortcutting flow from the heating chamber to the outlet.In practice, this can mean that the heating efficiency of plasma torch can be increased, and that the rated power of plasma heating can be reduced. The surface temperature of the heating chamber can be reduced, which can prolong the lifetime of refractory material.

Another benefit with stirring is that one can reduce the response time of the plasma heating, whichis beneficial for temperature control.

Conclusions

4.1 Plasma heating has advantagesforsimpler installations and lower maintenance costs over induction heating,although the heating efficiency needs to be improved and theconsumption of refractory material in the heating chamber reduced. 4.2 Numerical simulation was carried out to simulate the plasma heating process in the tundish.The results show that, with only plasma heating, the convection flow in the upper part of the heating chamber is dampened due to the buoyancy force caused by the high thermal gradient.This low convection rate in the heating chamber leads to a low heat efficiency ofplasma heating. One electromagnetic stirrer installed outside the tundish wallcreates a rotational flow in the heating chamber, which homogenizes temperature distribution and reduces surface temperature. This will increase the heating efficiency of plasma heating, prolong the lifetime of refractory, and reduce the response time for the melt temperature.

Fig. 1 :

 1 Fig.1: Geometry of a 30 ton tundish with plasma heating

Fig. 2 :

 2 Fig.2:Flow field in tundish with Case 1, 2 and 3

Table 1 :

 1 Basic parameters of numerical simulations

	Tundish capacity	30 ton	Time length of one ladle cycle 50 min
	Throughput from the	2.5 ton/min	Tundish refractory wall	150 mm
	ladle shroud		thickness
	Inner diameter of ladle	150 mm	Tundish stainless plate	12 mm
	shroud		thickness
	Outer diameter of ladle	250 mm	Tundish stainless steel	500 mm X 1200 mm
	shroud		window size
	Melt height	830 mm	Air gap between stirrer and	10 mm
			stainless steel window
	Initial temperature in	1540 o C	Density of steel	𝜌
	tundish			= 7010 + 𝑇 -1470
				• -0.883 𝑘𝑔 𝑚 3
	Starting temperature at	1565 o C	Heat flux at free surface	15 kW/m 2
	ladle shroud
	Temperature decrease at	0.5 o C/min	Heat loss at tundish wall	3.2 kW/m 2
	ladle shroud
	Max thermal input of	800 kW	Power of electromagnetic	60 kW
	plasma torch		stirrer
	Heating efficiency of	60%	Frequency of electromagnetic	4.5 Hz
	plasma torch		stirrer

Table 2 :

 2 Simulation conditions with plasma heating and EMS

		Plasma heating	EMS
	Case 1	No	No
	Case 2	Yes	No
	Case 3	Yes	Yes