SIMULATION OF SURFACE HEATING FOR ARBITRARY SHAPE ’ S MOVING BODIES / SOURCES BY USING R-FUNCTIONS

The purpose of this article is to propose an efficient algorithm for determining the place of an action of a heat source with a given motion law for a body of an arbitrary shape using methods of analytical geometry. The solution to this problem is an important part of a modeling of a laser, plasma, ion beam treatment. In addition, it can also be used for mass transfer problems, such as simulation of coating, sputtering, painting etc. The problem is solved by the method of R-functions to define the shape of the test body and the heat source and the analytical determination zone shadowing. As an example, we consider the problem of using the method of ion cleaning parameters optimization considering temperature limitations. Application of the R-functions can significantly reduce the amount of computation with usage of the ray tracing algorithm. The numerical realization of the proposed method requires an accurate creation of a numerical mesh. The best results in terms of accuracy of determination the scope of the source can be expected when applying adaptive tunable meshes. In case of integration of the R-functions into the CAD system, the use of the proposed method would be simple enough. The proposed method allows to determine the range of the source by the expression, which is constructed only once for the body and the source of arbitrary geometric shapes moving in any law. This distinguishes the proposed approach against all known algorithms for ray tracing. The proposed method can also be used for time-dependent multisource with arbitrary shapes, which move in different directions.


Introduction
Moving body/source heating analysis has an application in several manufacturing processes [1,2].In most well-known studies, in this formulation, the problem is considered for a circular heat source, which moves in a straight line [3,4].However, there are papers that outline the temperature distribution in a half-space, because of the complexity of the shape of the moving heat source [5], and because of the fact that a heat source moves according to a more complex laws [6].
The problem of heating bodies of finite dimensions by moving heat sources is considered in fewer studies.The objects of the study are likely to be the bodies of simple shapes (cylinder or parallelepipeds).For example, a number of studies were carried out by [7], to investigate the temperature distribution in a rotating cylinder heated with the laser heat source.This problem may be associated with the calculation of laser hardening regimes, laser-assisted machining, etc.
At the same time, it is necessary to calculate the temperature of an arbitrary shaped body caused by heating of a moving heat source.The law of motion and the shape of the heat source can be quite complex.Such problem is typical for cases, when body heat happens because of an energy flux action (radiation, heating, laser, ion or electron beam processing).
In this case, the definition of the heated area's borders is a complicated problem.In this paper, we propose an analytical method for the solution of this problem by using the R-functions.The geometric shapes of the body, the source and the law and their relative motion can be arbitrary.After definition of the action zone of the heat source , the further temperature calculation was determined numerically by a finite element method.

Mathematical formulation
The problem with determining the action zone of the moving heat source is very similar to the problem with determining the shadow zones location that is typical to computer graphics.. Beginning with the studies [8,9] such problems are solved by different algorithms of ray tracing.It required creating a lot of rays from points of the body surface in the direction of the radiation source.
These tasks require large amount of computations.At the same time, for complex shape bodies, the precise definition of the shaded area is associated with considerable difficulties in case of the moving sources, which change the intensity of the radiation.
We assume that the geometry of the source and the body is known.It can be independent of time, or change according to the known law.The source intensity is able to change arbitrary on time and on its cross section.The law of the source movement is also known.
Solution of the problem is reduced to the construction of a switch-function, which would be equal to 1 on the heated surface and to 0 on the rest of the surface of the part.
For this purpose, it is enough to create a function that has the following properties: ω(x, y, z, t) > 0 inside the heated surface part, ω(x, y, z, t) = 0 on the border of the heated surface part, ω(x, y, z, t) < 0 ouside the heeated surface part. ( Such relation can be described using the Rfunctions [10].The R-function method was developed as an improvement of Ritz methods for solving boundary-value problems.It is known for its utilization of solving the heat conduction problems in [11].However, in this paper it will be used only as a tool of analytical geometry.
The R-functions are functions of continuous real arguments.Their sign is determined by the sign of the argument.If, instead the sign of the "+" and "−", we use the values 0 and 1, R-functions are equivalent to some Boolean logic functions.Boolean function equivalent to a particular R-function, called the companion function.
Almost in every CAD system, a geometric shape of a complex object is created by using Boolean operations with geometric primitives.These primitives can be defined by algebraic equations or inequalities.For example, inequality ω 1 = (R 2 − x 2 − y 2 ≥ 0) in the plane X0Y defines a circle of radius R centered at the origin.
It is proved by [12] that if the geometry of the complex object is created by using Boolean operations (∨, ∧, ¬) with geometric primitives which are described by the inequalities ω i ≥ 0, the replacement of companion Boolean function with R-functions, al- lows to obtain the inequality ω = f (ω 1 , . . ., ω n ) ≥ 0, with the properties (1).The simplest complete system of these R-functions is the following one: Let it be required to create an expression ω(x, y) ≥ 0 for the area which is shown in Figure 1.As geometric primitives, the following items are selected: ω 1 = 1 − |x/a| − |y/b| ≥ 0 -the part of the plane inside a diamond with vertices (±a, 0), (0, ±b); ω 2 = (x 2 + y 2 − r 2 ≥ 0) -the part of the plane outside of the circle with radius r centered at the origin.
The shape of the region is determined by ω = ω 1 ∧ ω 2 , or after replacing the Boolean operations with R-functions: Further, for simplicity, the symbols ∧ R , ∨ R will be used instead of the expanded form of R-function.Method of analytical determination of the coverage of the source considers the example of two-dimensional problem.We propose dependences that describe the geometric shape of the heat source ω source ≥ 0 and body ω body = R(ω i ) ≥ 0, where R(ω i )-system of R-functions; ω i , i = 1, . . ., N -geometric primitives (Figure 2).

We propose the expression
which is equal to zero in the surface regions Γ i , primitives fall into the domain ω source and less than zero at all other points.We also use the expression where H (Heaviside function) is equal to one for all points on the surface of the body, lying on the Γ i , and zero for the remaining points.It means that ray tracing are needed only for points, where Ψ i is equal to one.
Without loss of generality, we assume that the action of the source is directed along the axis y source .We use lines f i = x source − x point for points on the surface lying on Γ i .Expression H(−f 2 i ) × Ψ i × |y source | determines the distance between the source and the point Γ i , which lies on f i .Point of the line f i , which lies closest to the source, will have coordinate y i,min source = min body H(−f 2 i ) × Ψ i × |y source | .Seeking function equal to 1 in areas which are affected by the source and 0 for the rest of the Γ i can be written as: With the known law of motion of the body, connecting the coordinates (x body , y body ) and (x source , y source ), functions H(g body ) and Φ clearly defines the scope of the source for the body of arbitrary shape given by the expression ω body .
The method of creating Equation (2) may seems cumbersome.However, this expression is created only once.This distinguishes the proposed approach against all known algorithms for ray tracing.The proposed method can also be used for time-dependent multisource with arbitrary shapes which move in different directions.
At solving such problems using the finite element method, quality of domain finite element mesh has a serious influence on the accuracy of the source locale.If the finite element mesh of the object is not performed accurately enough, the mesh generator determines the coordinates of the surface nodes with an error.In this case, the elements heated by source on the surface of the body are defined incorrectly.
You can see such influence in the following test problem.Determination of the moving heat source site of action for complex shape body was solved as a test problem.A rectangular shape source with a cut performs a reciprocating motion along the axis x body , and rotates around the axis y source (Figure 3).
For this problem, the coordinate system of the source and body are connected by expressions: x source = (x body + l − 2l sin φ1 t) cos φ2 t +z body sin φ2 t, y source = y body , z source = −(x body + l − 2l sin φ1 t) sin φ2 t +z body cos φ2 t.
Figure 4 shows the results of numerical calculations for two cases: when using automatic generation of unstructured mesh without additional conditions and with the mesh obtained by the forced association of mesh nodes to the body surface.Area of the source is defined incorrectly on curved surfaces and planes parallel to it for case one.

Numerical simulation result and discussion
Simulation of gas turbine engine blade heating during ion sputtering has been done using the proposed method.The process is performed on an ion-plasma device and is a preparatory stage before coating.Sputtering is a slow process, so to improve performance, a set of parts was treated at the same time.So batch of blades is set on a turntable which rotates with an angular velocity φ1 (Figure 5).For more uniform sputtering and coating, the blade also rotates around its own axis with an angular velocity φ2 .Energy of the ion beam, sputtering time, the rotation speed of the table and of the blades must be assigned so as to provide a predetermined quality of treatment and to prevent the blade material from overheating above the phase transition temperature.For the heat resistant alloys, which are most commonly used in a turbine blade manufacturing, this temperature is equal to 1270 K.The intersection of the ion beam and the blade surface causes the blade heating.It was assumed that the ion flux occupies a cylindrical shape area with radius r, located perpendicularly to the plane Z0X, displaced along the axis z from the plane of the table on a distance h.Movement of ions occurs along a positive direction of the axis y (Figure 5).
The blade heating occurs near the ion source.Therefore, the shape of the ion beam was set by the expression ω source = (r 2 − x 2 − (z − h) 2 ) ∧ R −y.Coordinate system of the body and the source are related by the expressions: x source = (x body + R cos φ1 t) cos(− φ2 )t +(y body + R sin φ1 t) sin(− φ2 )t, Blade surface geometry was defined by point cloud.
The heat flux is determined on the basis of the energy balance at the surface between input heat Q in and heat loss Q out .
For plasma processes Q in described by the expression [13]: Here J rad,in is the heat radiation towards the surface; J ch is the power transferred by charge carriers (electrons and ions); J n is the contribution of neutral species of the background gas and the neutral particles; J ads and J react,in are energies released by an absorption or a condensation and the reaction energy of exothermic processes including molecular surface recombination; J ext,in is an input heat flux by external sources that influences the thermal balance of the substrate.
Clearing by ion sputtering is carried out at low pressure with using a high-purity neutral gases, without additional heating and cooling.In this case, the expression for determining Q in can be simplified [13]: where ε is the spectral emittance of the radiation source at a temperature T rad ; ε s represents the spectral absorbance of the substrate surface at a temperature T s ; σ denotes the Stefan-Boltzmann constant; j i is the ion flux density of the surface; E ion is the ionization potential of the incident ion; k c is the energy transfer coefficient; M i , M s are masses ratio of the colliding particles (ion and surface atom); θ is the angle of incidence; e i is the ion charge; V bias is the sum of the plasma potential and the substrate potential.
The heat loss Q out of the substrate during plasma processing consist of the following terms [13]: Here J rad,out is the energy radiated from the substrate at a temperature T s ; J particle is the energy transport from the substrate due to sputtering of surface atoms and the secondary electron emission; J des is the energy sink due to the desorption of particles into the gas phase and the diffusion into the solid bulk; J react,out is the reaction energy of exothermic processes, including molecular surface recombination; J ext,out is the heat loss caused by an external cooling.
Subject to the terms of the ion cleaning, the expression for determining Q out can be written as: where ε env is the emissivity of the environment; T env is the environmental temperature (reactor walls, etc.); k sput is the sputtering coefficient dependent on the ion energy and angle of incidence; E coh is the cohesive energy of sputtered material's atoms.The temperature field in the clearing blade was calculated during the simulation.For this, the transient heat conduction equation was solved.The heat flux through surface of the blade was specified by the expression: where Φ is a switch-function, which is determined by the expression (2).Dependence of the heat flux value on the angle between the direction of ions flux and the blade surface is a feature of the process.This feature required calculating the directional cos of the normal to the surface.Using R-functions makes it simple to solve this problem.If the equation ω i ≥ 0 for geometric primitives is known, and for a complex domain the equation ω(ω i ) ≥ 0 is constructed by using R-functions, than for the function ω * i ≥ 0, where , the expressions ∂ω(ω * i )/∂x, ∂ω(ω * i )/∂y, ∂ω(ω * i )/∂z will determine the direction cosines of the interior normal to the corresponding coordinate axes [12].
Parameters of sputtering (such as value of the ions energy, rotation speed of the turntable and the blade) have been determined with respect to the temperature limit, so the maximum surface temperature should not exceed the temperature of the material phase transition.Additionally, the condition , in which the minimum value of the sputtered material layer should not be less than a predetermined value throughout surface of the blade, was applied.
The graphs of the maximum temperature change in a checkpoint of the blade at a various ion beam energy is shown in Figure 6.Based on the simulation results, recommendations on the choice parameters of the sputtering of the blades have been made.Criterion of choice was the minimum time of the sputtering while respecting the temperature limitations.

Conclusions
The use of R-functions can significantly reduce the amount of tracing computations.The abovementioned effect is due to the analytical determination of points falling within the scope of the source.The point clouds produced by the 3D scanners can be used for this purpose.
The numerical realization of the proposed method requires an accurate build numerical mesh.The best results in terms of accuracy of determination the scope of the source can be expected, when applying adaptive tunable meshes.
The proposed approach can be applied to cases with an arbitrary number of heat sources considering that the law of body motion is known.It can also be used in mass transfer problems, such as simulation coating, painting or clearing bodies with complex geometric shapes.In case of integration of the R-functions in CAD system, the use of the proposed method would be simple enough.

Figure 3 .
Figure 3. Design scheme for the test problem.

Figure 4 .
Figure 4.The results of numerical simulations using unstructured mesh without conditions (left) and mesh obtained by the forced association of mesh nodes from the body surface (right).

Figure 5 .
Figure 5.The equipment for sputtering and coating -turntable with blades and design scheme for the task of the blade heating under ion sputtering.

Figure 6 .
Figure 6.Maximum blade surface temperature versus the sputtering time.