Spiral bevel gears transmit torque between intersecting shafts with a contact ratio and motion-error characteristic that few other gear forms can match, which is why they dominate final drives in aircraft engine and landing-gear transmissions, automotive rear axles and differentials, wind-turbine yaw drives, and heavy earthmoving equipment. In each of these applications the gear is hardened before finishing, and grinding is the operation that establishes the final flank geometry and surface integrity. The same operation, however, imposes a mechanical-thermal loading cycle on a thin surface layer: contact pressure from the abrasive grains produces plastic compression, while the heat flux generated in the contact zone produces thermal expansion followed by contraction during cooling. The self-equilibrated stress field that remains after the wheel disengages is therefore a direct product of the grinding parameters, and it controls fatigue initiation, wear resistance and dimensional stability in service. A compressive subsurface peak retards crack initiation; an excessive tensile zone beneath it promotes crack growth and, in the presence of a stress concentrator, premature failure. Quantifying how the parameters of the grinding cycle map onto this stress profile is a prerequisite for process planning rather than a post-hoc quality check.
Tooth Generation Kinematics: From the Imaginary Crown Gear to the Ground Flank
The geometry of a spiral bevel gear is defined by a generating process rather than by a direct profile. Both cutting and grinding rely on the concept of an imaginary flat-top (crown) gear: a virtual gear concentric with the machine cradle, whose tooth surfaces are swept by the cutting edges of a tool head mounted eccentrically on the cradle and swinging with it. The virtual gear meshes without backlash with the workpiece blank, and the envelope of the successive tool positions removes the tooth space. Grinding reproduces this kinematics with a cup wheel whose active surface plays the role of the cutter blade; the difference between milling and grinding lies in the tool, not in the generating principle.
Two kinematic families are used. In generating grinding the blank rotates about its own axis while the wheel head rotates about the cradle axis, so that the two surfaces remain in line contact along an instantaneous contact line; the superposition of these contact lines over the cycle envelopes the flank. In form grinding the wheel axis keeps a fixed spatial relationship with the gear axis, and the wheel profile is a direct replica of the required tooth profile; the spiral is produced by feeding the wheel along the tooth space. Generating wheels are further classified by the face cone angle of the generating gear: a plane generating wheel has a 90° face cone angle with the cutter axis parallel to the cradle axis, whereas a conical generating wheel has a non-90° face cone angle and a non-parallel axis arrangement.
Abrasive–Workpiece Interaction: Rubbing, Plowing and Cutting
Each active grain acts as a micro-cutting tool with a large negative rake angle. The interaction with the flank proceeds through three regimes. In the rubbing stage the grain contacts the surface elastically and slides; friction dissipates most of the energy as heat, raising the local temperature and, in the limit, causing thermal damage. In the plowing stage the grain penetrates and displaces material laterally, forming a groove with raised ridges on both sides; the subsurface undergoes severe plastic strain and work hardening, and a residual stress field begins to establish. In the cutting stage the undeformed chip thickness exceeds a critical value, shear localization occurs ahead of the grain and a chip flows over the rake face. Because the grains protrude over a distribution of heights, only a small fraction reaches the cutting regime at any instant; the remainder rub and plow, which is the origin of the high specific energy of grinding and of the strongly coupled mechanical-thermal character of the resulting stress state.
Wheel Configurations Used in Production
Three wheel arrangements are in industrial use for spiral bevel gears, and their differences are essentially differences in contact geometry. The straight cup wheel mirrors the blade profile of the milling cutter and contacts the tooth space over the full face width; it is applicable to gears with pitch cone angles up to about 50°, but the full-width contact raises the normal force and restricts coolant access, so the fillet region is prone to grinding burn. For this reason it is rarely selected for form grinding of the wheel member of a pair. The flaring cup wheel uses a 30° cone angle arranged so that the outer surface curvature radius is smaller than that of the concave slot and the inner surface radius is larger than that of the convex slot. The curvature mismatch enforces single-point contact, which opens a coolant channel into the contact zone, lowers the instantaneous thermal load and suppresses burn; multi-axis interpolation sweeps the wheel along the slot to reproduce the tooth space.
The Wguri principle dispenses with the special cup profile. The wheel spindle is mounted on an eccentric mechanism with an offset of roughly 0.1 mm, and the orbital motion of the spindle relative to the tooth space produces an intermittent, clearance-type contact that alternately addresses the concave and convex flanks as the spindle feeds axially. Coolant and chips are evacuated through the additional clearance, and a standard cylindrical wheel replaces the expensive flaring cup. The trade-off is machine-side: a precision eccentric spindle assembly must be added, increasing structural complexity and cost.
Force and Heat Inputs to the Analysis
For a predictive model the grinding action must be reduced to equivalent mechanical and thermal loads. The normal and tangential force components per unit width are each decomposed into a sliding, a plowing and a chip-formation contribution:
$$F’_t = F’_{ts} + F’_{tp} + F’_{tch}, \qquad F’_n = F’_{ns} + F’_{np} + F’_{nch}$$
The sliding term is frictional in origin, with the tangential component proportional to the normal one through the grain–workpiece friction coefficient, and the normal component scaling with the fourth root of the ratio of feed speed to the product of wheel speed and deformed wheel diameter, modulated by a sliding coefficient and the effective contact area. The plowing terms are related linearly through a plowing constant, and the chip-formation terms follow an energy balance in which the specific chip-formation energy, the feed speed and the actual depth of cut are divided by the wheel speed:
$$F’_{np} = k_{p1} F’_{tp}, \qquad F’_{tch} = \frac{u_{ch}\, v_w\, a}{v_s}, \qquad F’_{nch} = k_{ch} F’_{tch}$$
Heat is treated as a band source moving over a semi-infinite body at the wheel speed. The heat flux entering the workpiece is
$$q = R_w \frac{F_t \,(v_s \pm v_w)}{b\, l_k}$$
where the sign is negative for down grinding and positive for up grinding, b is the contact width and lk the contact arc length. The energy partition ratio Rw is obtained from the ratio of the thermophysical properties (conductivity, density, specific heat) of the abrasive and the workpiece, together with the grain geometry factor, the effective number of active grains per unit area and the real contact area; a coefficient of 1.06 appears in the denominator of the standard partition expression.
Finite-Element Realization
The analysis targets the 44-tooth wheel member of a spiral bevel gear pair in AISI 4340. Tooth flank coordinates were computed numerically and exported as a point cloud, from which continuous concave, convex and root surface patches were lofted and stitched into a single tooth slot; the slot was then trimmed from the blank, rotated into position and arrayed about the axis to give the full gear. The main geometry is listed in Table 1.
| Parameter | Value | Parameter | Value |
|---|---|---|---|
| Number of teeth | 44 | Face cone angle | 73.39° |
| Module | 3 mm | Root cone angle | 67° |
| Pressure angle | 20° | Outer cone distance | 69.73 mm |
| Mean spiral angle | 35° | Mounting distance | 27 mm |
| Hand | Left | Addendum | 1.52 mm |
| Pitch diameter | 132 mm | Dedendum | 4.15 mm |
| Shaft angle | 90° | Face width | 20.2 mm |
Because a grain is orders of magnitude smaller than the gear, the contact problem was reduced to a representative patch of the flank. A single grain was idealized as a rigid cone 0.07 mm high with a 90° included angle and a 0.01 mm edge radius. The workpiece was treated as an elastic–plastic body governed by the Johnson–Cook constitutive relation, which couples strain hardening, strain-rate sensitivity and thermal softening:
$$\sigma = \left(A + B\varepsilon^{n}\right)\left(1 + C\ln\frac{\dot{\varepsilon}}{\dot{\varepsilon}_0}\right)\left[1 – \left(\frac{T – T_0}{T_m – T_0}\right)^{m}\right]$$
Damage was evaluated with the Johnson–Cook fracture criterion, in which the equivalent plastic strain increment is accumulated against a failure strain that depends on stress triaxiality, strain rate and homologous temperature:
$$\varepsilon_f = \left(D_1 + D_2 e^{D_3 p/q}\right)\left(1 + D_4 \ln\frac{\dot{\varepsilon}}{\dot{\varepsilon}_0}\right)\left(1 + D_5\frac{T – T_0}{T_m – T_0}\right), \qquad D = \sum \frac{\Delta\varepsilon}{\varepsilon_f}$$
Material data are collected in Tables 2 and 3.
| Property | AISI 4340 workpiece | Abrasive grain |
|---|---|---|
| Density | 7830 kg/m³ | 3480 kg/m³ |
| Elastic modulus | 200 GPa | 706 GPa |
| Poisson ratio | 0.29 | 0.15 |
| Yield strength | 792 MPa | — |
| Tensile strength | 980 MPa | — |
| Thermal conductivity | 38 W/(m·°C) | 79.54 W/(m·°C) |
| A (MPa) | B (MPa) | n | C | m | D₁ | D₂ | D₃ | D₄ | D₅ |
|---|---|---|---|---|---|---|---|---|---|
| 792 | 510 | 0.26 | 0.014 | 1.03 | 0.05 | 3.44 | −2.12 | 0.002 | 0.61 |
Hexahedral coupled temperature–displacement elements were used for the workpiece and tetrahedral elements for the rigid grain. Contact was described by a penalty formulation with a friction coefficient of 0.3, with frictional heating and heat conduction active at the interface; the base of the workpiece was fixed. Each case was solved in two stages: the grinding pass by an explicit solver, followed by a cooling stage in the implicit solver, in which the tool, its material and the contact definition were removed and convective boundary conditions were applied. The stress component along the grinding direction (S11) was extracted layer by layer, and 98 equally spaced nodes were averaged at each depth to suppress edge artefacts.
A single grain cannot represent the interaction of neighbouring grains, so a multi-grain patch was also constructed. The wheel surface was partitioned into a regular lattice of cells; a virtual-lattice algorithm placed one grain per cell and allowed a random offset of the grain centre within each cell according to a probability density function, giving macroscopic uniformity and microscopic randomness simultaneously. Grain protrusion heights were distributed between 0.05 and 0.07 mm. The constitutive, contact and extraction settings were kept identical to the single-grain case so that the two sets of results are directly comparable.
Predicted Stress Profile and Parameter Sensitivity
The computed profiles share a consistent shape. The as-ground surface carries compressive stress; the maximum compression develops not at the surface but in the subsurface layer, where the mechanical and thermal contributions overlap most strongly; with increasing depth the stress crosses zero and becomes tensile, then decays toward zero as the elastically recovered zone is left behind. The affected layer extends to roughly 0.3 mm, beyond which the stress is negligible.
Single-factor simulations were run over four levels of each parameter (Table 4). Increasing wheel speed from 10 to 25 m/s raises the residual stress. The mechanism is not a simple force effect: a higher wheel speed increases the number of grains engaged per unit time, so the undeformed chip thickness per grain and the overall grinding force fall, but the heat generation rate rises, the surface and subsurface temperatures increase, and the resulting thermal stress contribution dominates. Increasing depth of cut from 0.01 to 0.04 mm produces the largest change, because the contact area between grain and flank grows approximately linearly, raising both the normal force and the plastic strain penetration depth, and the volume removed per unit time rises sharply. Increasing feed speed from 3 to 4.5 m/min raises the stress only modestly: the material removal rate and the instantaneous temperature gradient increase, but the contact time per unit length shortens, so heat accumulation and sustained plastic deformation are partly offset.
| Level | Wheel speed vs (m/s) | Depth of cut ap (mm) | Feed speed vw (m/min) |
|---|---|---|---|
| 1 | 10 | 0.01 | 3.0 |
| 2 | 15 | 0.02 | 3.5 |
| 3 | 20 | 0.03 | 4.0 |
| 4 | 25 | 0.04 | 4.5 |
Ranking the three factors by their influence on the stress state gives depth of cut first, wheel speed second and feed speed third — the same order that emerges from the statistical treatment described below.
Single-Grain versus Multi-Grain Fields
Under the same nominal conditions (25 m/s, 0.04 mm, 4.5 m/min) the multi-grain patch reproduces the depth profile of the single-grain model but with a lower peak compressive stress. The reason is superposition: a grain compresses the material directly beneath it, while adjacent grains impose tensile fields on the same region, so part of the compressive field is cancelled. The similarity of the profile shapes nevertheless establishes a usable mapping between single-grain and multi-grain behaviour, which is what makes the much cheaper single-grain model a legitimate surrogate for trend studies.
Parameter Optimization: Orthogonal Array, Grey Relational Analysis, Particle Swarm and Entropy Weighting
Two responses are of interest and they conflict: the peak compressive stress should be as large as possible, and the peak tensile stress as small as possible. A three-factor, four-level orthogonal array of 16 simulation runs was used as the design (Table 5), with compressive magnitude treated as a larger-the-better index and tensile magnitude as a smaller-the-better index.
| Run | vs (m/s) | ap (mm) | vw (m/min) | Max. compressive S₁ (MPa) | Max. tensile S₂ (MPa) |
|---|---|---|---|---|---|
| 1 | 10 | 0.01 | 3.0 | −235.9 | 12.4 |
| 2 | 10 | 0.02 | 3.5 | −266.5 | 12.7 |
| 3 | 10 | 0.03 | 4.0 | −295.4 | 14.2 |
| 4 | 10 | 0.04 | 4.5 | −363.6 | 26.7 |
| 5 | 15 | 0.01 | 3.5 | −249.1 | 12.2 |
| 6 | 15 | 0.02 | 3.0 | −285.2 | 11.9 |
| 7 | 15 | 0.03 | 4.5 | −332.9 | 22.7 |
| 8 | 15 | 0.04 | 4.0 | −356.9 | 20.9 |
| 9 | 20 | 0.01 | 4.0 | −287.4 | 11.6 |
| 10 | 20 | 0.02 | 4.5 | −321.7 | 18.3 |
| 11 | 20 | 0.03 | 3.0 | −335.2 | 17.1 |
| 12 | 20 | 0.04 | 3.5 | −357.6 | 20.6 |
| 13 | 25 | 0.01 | 4.5 | −279.3 | 9.1 |
| 14 | 25 | 0.02 | 4.0 | −347.6 | 22.2 |
| 15 | 25 | 0.03 | 3.5 | −347.9 | 23.8 |
| 16 | 25 | 0.04 | 3.0 | −342.3 | 20.5 |
Grey relational analysis converts the two responses into a single grade. Each response sequence is first normalized — the compressive sequence by a larger-the-better transformation and the tensile sequence by a smaller-the-better transformation — and the grey relational coefficient is then computed from the absolute deviation between the normalized comparison sequence and the reference sequence, with a resolution coefficient of 0.5; the grade is the mean coefficient over the two responses. Averaging the grade at each factor level (Table 6) identifies the best combination as 25 m/s, 0.04 mm and 4.5 m/min, and the range of the level means ranks depth of cut (0.0776) well above wheel speed (0.0472) and feed speed (0.0317).
| Factor | Level 1 | Level 2 | Level 3 | Level 4 | Range |
|---|---|---|---|---|---|
| Wheel speed | 0.5776 | 0.5890 | 0.6117 | 0.6248 | 0.0472 |
| Depth of cut | 0.6033 | 0.5764 | 0.5728 | 0.6504 | 0.0776 |
| Feed speed | 0.5842 | 0.5920 | 0.6109 | 0.6159 | 0.0317 |
A second-order response surface fitted to the grade data, after removal of terms with p > 0.05, is statistically sound: F = 35.73 with p < 0.0001, R² = 95.97% and adjusted R² = 93.29%. Equivalent surfaces fitted to the two stress responses give R² = 97.11% (adjusted 94.58%) for the compressive peak and R² = 97.18% (adjusted 93.95%) for the tensile peak, and the residuals of both lie close to the fitted line with no systematic structure.
The same two responses were then optimized directly with a multi-objective particle swarm algorithm, using a population of 100, 500 iterations, cognitive and social coefficients of 2, and a linearly decreasing inertia weight from 0.9 to 0.4 to shift the balance from exploration to exploitation. Five non-dominated solutions were obtained (Table 7), all located near the upper bound of depth of cut and at moderate wheel speeds of 20.7–22.6 m/s and feeds of 3.37–3.45 m/min.
| Candidate | vs (m/s) | ap (mm) | vw (m/min) | S₁ (MPa) | S₂ (MPa) | Composite score | Rank |
|---|---|---|---|---|---|---|---|
| Grey relational optimum | 25.00 | 0.040 | 4.50 | −411.0 | 38.6 | 0.899 | 1 |
| PSO solution 1 | 22.55 | 0.0395 | 3.43 | −357.3 | 22.3 | 0.143 | 2 |
| PSO solution 2 | 21.41 | 0.0400 | 3.45 | −355.9 | 21.1 | 0.128 | 4 |
| PSO solution 3 | 21.79 | 0.0400 | 3.43 | −356.0 | 21.4 | 0.129 | 3 |
| PSO solution 4 | 20.71 | 0.0400 | 3.44 | −353.8 | 20.1 | 0.101 | 6 |
| PSO solution 5 | 22.09 | 0.0400 | 3.37 | −354.2 | 21.0 | 0.102 | 5 |
The six candidates were finally ranked by the entropy weight method, which assigns each response a weight proportional to the dispersion of its normalized values and thereby avoids subjective weighting. The grey relational optimum scores 0.899 against 0.101–0.143 for the five swarm solutions. That margin should be read with care. The compressive peak varies across the candidates over a range of about 57 MPa, whereas the tensile peak varies over about 18 MPa, so after normalization the compressive response dominates the composite score; the winning scheme is the one that maximizes compression, and it does so at the cost of a tensile peak of 38.6 MPa, higher than any of the sixteen orthogonal runs. The competition between a deep compressive layer and a suppressed tensile tail is therefore not resolved by the weighting but inherited from it, and a process planner whose failure mode is tensile-driven cracking would legitimately prefer one of the swarm solutions.
Two further observations follow from Table 7. The selected optimum lies on the boundary of the search space in both wheel speed and depth of cut, which indicates that both responses are monotonic in the ranges explored; the true optimum, subject to burn and wheel-life constraints, lies outside the simulated box and the box should be widened in any extension of the work. In addition, the corner (25 m/s, 0.04 mm, 4.5 m/min) is not represented by any of the sixteen array runs, so the predicted values of −411.0 MPa and 38.6 MPa are extrapolations of the fitted surface rather than interpolations, which is precisely why an independent confirmation was necessary.
Experimental Verification
Grinding trials were carried out on an eight-axis, five-linkage CNC spiral bevel gear grinder equipped with an eccentric grinding spindle, using the optimized parameters: wheel speed 25 m/s, depth of cut 0.04 mm, feed speed 4.5 m/min. Each finished gear was sectioned by wire electrical discharge machining into single-tooth specimens, with the cut located well away from the flank and cooling applied to limit thermally induced stress relaxation. Residual stress on the flank was measured by X-ray diffraction with an Xstress 3000G3 analyzer, which was chosen over mechanical hole-drilling and slotting methods because those relieve stress by destroying the surface being characterized. Three measurement locations were marked on the flank and are designated as experimental groups 1–3. The stress profile with depth was obtained by successive electrolytic polishing and re-measurement.
The measured profiles follow the predicted shape: compression at the surface, a subsurface maximum, a transition to tension and a decay to zero within about 0.3 mm. The measured magnitudes are consistently larger than the simulated ones. Two explanations are consistent with the data. First, the real grain density on the wheel is likely higher than the density assumed in the multi-grain patch, so overlapping stress fields from adjacent grains reinforce rather than cancel; second, the simulation assumes an ideal environment, whereas the actual process involves ambient temperature fluctuation, machine vibration and non-uniform coolant distribution, all of which perturb the thermal cycle.
| Source | S₁ measured (MPa) | S₁ predicted (MPa) | Error | S₂ measured (MPa) | S₂ predicted (MPa) | Error |
|---|---|---|---|---|---|---|
| Finite-element simulation | −396.5 | −411.0 | 3.67% | 36.9 | 38.6 | 4.52% |
| Experiment, group 1 | −450.2 | −411.0 | 8.77% | 42.9 | 38.6 | 10.02% |
| Experiment, group 2 | −463.8 | −411.0 | 11.37% | 43.4 | 38.6 | 11.16% |
| Experiment, group 3 | −454.5 | −411.0 | 9.55% | 42.4 | 38.6 | 8.94% |
Against the regression prediction, the simulation shows 3.67% and 4.52% relative error for the compressive and tensile peaks, and the three experimental locations show 8.77–11.37% and 8.94–11.16% respectively. Comparing experiment directly with simulation gives 11.9–15.01% across the three locations. All discrepancies remain below 16%, which supports the use of the model for parameter ranking and for interpolative prediction within the calibrated window, but the consistent sign of the bias — measurement always exceeding prediction — indicates that the multi-grain patch systematically under-represents the compressive field and would need a density calibration before being used as an absolute predictor.
Implications for Process Planning
Three conclusions carry over to production planning. Depth of cut is the dominant lever on the residual stress field, with wheel speed a secondary and feed speed a weak third; a change in depth of cut of 0.03 mm shifts the compressive peak by more than 100 MPa, which is larger than the spread produced by any combination of the other two factors in the tested range. The stress profile is self-compensating in shape but not in magnitude: the same parameters that deepen the subsurface compression also raise the tensile tail, so specifying a target compressive peak implicitly specifies a tensile penalty and both must be quoted. And the single-grain model, though it overestimates the compressive peak relative to the multi-grain case, reproduces the profile shape and the parameter ranking and is therefore adequate for design screening; absolute values require the multi-grain patch plus the kind of experimental calibration reported here.
The model has recognized limits. The grain was treated as a rigid cone of fixed geometry, so the influence of rake angle and grain wear is absent; the wheel topography was idealized by a virtual lattice rather than measured; and no microstructural phase transformation, such as the retained-austenite conversion that accompanies severe tempering in hardened bearing steels, was represented. Extending the parameter space beyond the boundary optima found here, coupling the model to measured wheel topography, and replacing the deterministic single-grain assumption with a distribution of grain shapes are the natural next steps if the approach is to move from trend prediction to absolute specification.