Enginee ing Science and Technology, an In e na ional Jou nal 55 (2024) 101721
A ailable online 4 June 2024
2215-0986/© 2024 THE AUTHORS. Published by Else ie BV on behal o Ka abuk Uni e si y. This is an open access a icle unde he CC BY-NC-ND license
(h p://c ea i ecommons.o g/licenses/by-nc-nd/4.0/).
Full Leng h A icle
Full-s a e eedback LQR wi h in eg al gain o con ol o induc ion hea ing
o s eel bille
Sohaibullah Za ghoon
a
, Samuel Emebu
b
,
c
,
e
,
*
, Radek Ma uˇ
sů
b
, Cy il Bela ý
a
, Luk´
aˇ
s Ba alský
a
,
S anisla ˇ
Du iˇ
s
a
, Sabi Husnain
, Cla a Mendoza Ma inez
d
a
Ins i u e o Au oma ion, In o ma iza ion, and Measu emen , Facul y o Mechanical Enginee ing, Slo ak Uni e si y o Technology, N´
am. slobody 17, 812 31 B a isla a,
Slo akia
b
Depa men o Au oma ion and Con ol Enginee ing, Facul y o Applied In o ma ics, Tomas Ba a Uni e si y in Zlín, Nad S ´
anˇ
emi 4511, 760 05 Zlín, Czech Republic
c
Depa men o Chemical Enginee ing, Facul y o Enginee ing, Uni e si y o Benin, PO Box 1154, Benin Ci y, Nige ia
d
Depa men o Ene gy, Lappeen an a-Lah i Uni e si y o Technology LUT, Yliopis onka u 34, FI-53850 Lappeen an a, Finland
e
Depa men o Sepa a ion Science, Lappeen an a-Lah i Uni e si y o Technology LUT, Yliopis onka u 34, FI-53850 Lappeen an a, Finland
Depa men o Mechanical, Mecha onics and Manu ac u ing Enginee ing, Uni e si y o Enginee ing and Technology Laho e, Laho e, Punjab, Pakis an
ARTICLE INFO
Keywo ds:
Induc ion hea ing
Skin e ec s
Ca bon s eel bille
Semi-nume ical me hod
Full-s a e eedback Linea Quad a ic Regula o
(LQR)
ABSTRACT
Induc ion hea ing is widely used in indus ial u naces due o i s apid esponse and ene gy e iciency. Compu e -
aided modelling and simula ion a e necessa y o he design and op imisa ion o hese u naces. This a icle
ocuses on modelling and simula ing induc ion hea ing o a one-dimensional spa ial pa ial di e en ial model o
a ec angula ca bon s eel bille . The simula ion conside ed empe a u e dynamics, including he skin e ec and
hea loss ia adia ion. The pa ial di e en ial model was con e ed in o a s a e-space model o designing a ull-
s a e eedback Linea Quad a ic Regula o (LQR) wi h in eg al ac ion. This con olle was e ec i ely applied o
egula e he bille ’s co e empe a u e om 1000 ◦C o 1200 ◦C, unde di e en inpu , R and s a e, Q, weighing
ma ices o he LQR, as well as in he p esence o dis u bances.
1. In oduc ion
Induc ion hea ing is a widely used hea ing echnique in a ious in-
dus ies due o i s en i onmen al bene i s. T adi ional hea ing me hods,
such as bu ne hea ing u naces, consume mo e ene gy when compa ed
o induc ion hea ing equipmen . As a esul , inco po a ing induc ion
hea ing sys ems in ho - o ging plan s helps educe ene gy consump ion,
as well as ai pollu ion [1]. This me hod e icien ly gene a es hea
di ec ly wi hin he molecula le el o me allic induc i e ma e ials such
as aluminium, b ass, coppe , s eel, and semiconduc ing ma e ials like
silicon ca bide [2], h ough elec omagne ic induc ion. The e o e, in-
duc ion hea ea men has become a globally ecognised echnique o
epai ing, welding bo h mic oscopic componen s and la ge pa s, o -
e ing economic as well as p ac ical ad an ages [3].
In p ac ice, he design o induc ion hea ing u naces o a ious ap-
plica ions o en in ol es a leng hy and expensi e ial-and-e o p ocess.
Measu ing and e alua ing key ac o s such as cu en equency, em-
pe a u e, and hea dis ibu ion can pose challenges. And compu e
modelling is a sui able app oach o such e alua ion, as well as he
design o induc ion hea ing u naces [4]. The modelling o he induc ion
hea ing p ocess is based on he Maxwell equa ions, which cons i u e a
se o Pa ial Di e en ial Equa ions (PDEs) as such a nonlinea model
ha can be di icul o sol e. Howe e , he solu ions can be ob ained o
simpli ied using modelling echniques such as analy ical, nume ical, and
semi-analy ical app oaches [5–8]. De eloping and simula ing models o
induc ion hea ing u naces a e essen ial o op imising and con olling
he p ocess o achie e he desi ed empe a u e o he hea ed bille . This
p ima ily aims o p e en haza ds and ensu e e icien ene gy u ilisa ion
[9]. The induc ion hea ing p ocess can be con olled using ei he
nonlinea o linea con olle s. Linea con olle s a e commonly
employed in indus ial con ol sys ems since many nonlinea p ocesses
can be accu a ely desc ibed by linea models unde speci ic nominal
condi ions h ough app oxima ion me hods like ans e unc ions and
s a e-space models [10]. Addi ionally, linea con olle s ha e a well-
es ablished ounda ion, and iden i ying linea models based on p ocess
da a is ela i ely s aigh o wa d. They p o ide a good app oxima ion in
he icini y o he designed ope a ing condi ions, o e as esponse
imes, and a e cos -e ec i e o implemen [11].
* Co esponding au ho .
E-mail add ess: [email p o ec ed] (S. Emebu).
Con en s lis s a ailable a ScienceDi ec
Enginee ing Science and Technology,
an In e na ional Jou nal
jou nal homepage: www.else ie .com/loca e/jes ch
h ps://doi.o g/10.1016/j.jes ch.2024.101721
Recei ed 20 Decembe 2023; Recei ed in e ised o m 20 Ap il 2024; Accep ed 21 May 2024
Enginee ing Science and Technology, an In e na ional Jou nal 55 (2024) 101721
2
Al hough he impo ance o modelling, simula ion, and con ol has
been discussed, i is wo h no ing ha mos li e a u e s udies ha e p i-
ma ily ocused on he modelling and simula ion o induc ion hea ing
u naces. These s udies ha e o en compa ed he esul s ob ained om
nume ical models, pa icula ly he Fini e Elemen Me hod (FEM), wi h
expe imen al da a, and ha e gene ally ound good ag eemen
[8,12–17]. Howe e , Demido ich e al. [18] no only conside ed FEM
bu also explo ed and compa ed o he nume ical me hods such as ini e-
di e ence me hods (FDM), bounda y elemen me hod (BEM), in ini e
elemen me hod (IEM), as well as hei combina ions, wi hou he use o
expe imen al da a. Fu he mo e, some epo s ha e compa ed analy ical
models wi h expe imen al da a [5,15]. Fo ins ance, Li e al. [5]
de eloped an analy ical model o p edic he empe a u e p o ile in a
plana mo ing induc ion hea ing p ocess. The s udy demons a ed ha
he empe a u e p o ile and cu en in ensi y calcula ed by he analy -
ical model we e simila o he esul s ob ained om he nume ical
model. Impo an ly, he analy ical model exhibi ed highe compu a-
ional e iciency compa ed o he nume ical model. The ad an age o
educed expe imen al and compu a ional complexi y o e ed by he
analy ical model p o ides suppo o s udies such as Pa ida e al.’s[6]
epo on he combina ion o analy ical and nume ical me hods, as well
as A ei ioau ena e al.’s[7] semi-analy ical me hod o modelling he
induc ion hea ing p ocess.
Based on he discussion hus a , he e a e limi ed epo s on he
con ol o he induc ion hea ing p ocess using he de eloped non-linea
PDEs models. Cambe e al. [19] employed he dis ibu ed pa ame e
heo y o he con ol o a con inuous s eel bille induc ion hea e sys em
as lumped-inpu s (i.e., al e na ing cu en ) and dis ibu ed-ou pu s
( empe a u e p o ile o bille ). Goodwin e al. [20] modelled he in-
duc ion hea ing p ocess ia a non-linea PDE, con olled he p ocess
using nonlinea model p edic i e con ol (MPC), as well as compa ed
expe imen al and simula ed esul s. Fu he mo e, Roe ze e al. [21]
also modelled he induc ion hea ing p ocess ia a non-linea PDE,
con olled and compa ed he p ocess wi h di e en ypes o con olle s
(i.e., P opo ional in eg al (PI), and linea -quad a ic Gaussian (LQG)
con olle s). Howe e , as highligh ed by Goodwin e al. [20], he
inco po a ion o models o accoun o skin e ec phenomena, and he
empe a u e dependency o pa ame e s, a e s ill limi ed in li e a u e.
Pa icula ly he con ol o ma e ial wi h he inco po a ion o skin e ec
[22] and he empe a u e dependency o ma e ial p ope ies a he han
cons an alues (such as densi y, emissi i y, elec ical conduc i i y,
he mal conduc i i y, and speci ic hea capaci y) [23], can esul in a
e y complex nonlinea sys em ha is di icul o sol e as well as con-
ol. In p ac ice, his kind o model can be compu a ionally demanding,
slow, and expensi e o con ol using nonlinea con olle s [24]. The e-
o e, o educe cos and inc ease he con olle esponse ime, a ises he
need o explo e linea con olle s such as he LQR, o con ol he in-
duc ion hea ing p ocess using simpli ied app oxima ions o non-linea
PDEs models, wi h bo h conside a ion o skin e ec , and empe a u e
dependency o ma e ial p ope ies. Al hough a emp s ha e been made
o co e his esea ch gas [25,26], howe e esul s we e s ill inadequa e.
Fo ins ance, Asadzadeh e al. [26] applied an analy ical app oxima-
ion, speci ically he in e se model o simula e he hea ing p ocess wi h
empe a u e dependency o he ma e ial p ope ies and a emp ed o
conside skin e ec . Howe e , he in e se model nega ed he spa ial
dis ibu ion o empe a u e along he ma e ial, hence being unable o
e ec i ely po ay he skin e ec on he ma e ial. Fu he mo e, Kapus a
e al. [25] applied a nume ical model o he simula ion o induc ion
hea ing o s eel bille and con ol ia dis ibu ed-inpu and dis ibu ed-
pa ame e -ou pu sys ems. Al hough Kapus a e al. [25] conside ed he
empe a u e dependency o he ma e ial on he model, howe e , he
model may no ha e adequa ely conside ed skin e ec on he spa ial
empe a u e dis ibu ion wi hin he dep h o he ma e ial. The e o e, i
is his esea ch gap ha his wo k a emp s o co e . Hence he no el y
o his wo k en ails he con ol ia linea con olle o acili a e he apid
esponse o he highly nonlinea model o induc ion hea ing o s eel
bille wi h adequa e conside a ion o skin e ec , along wi h empe a u e
dependency o i s ma e ial p ope ies. Speci ically, he LQR is a eliable
linea con olle ha can be applied, because o i s abili y o deli e
op imal as well as apid con ol esponse o he sys em pe he designe ’s
speci ica ions (e.g. ac ua o s cha ac e is ics o limi a ions) [27]. LQR is
howe e mo e sui able o sys ems wi h limi ed s eady-s a e e o on
applica ion o he LQR eedback gain ma ix o he con olled sys em
[28]. The e o e, o sys ems wi h high s eady-s a e e o he LQR is
usually inco po a ed wi h in eg al gain, p opo ional in eg al gain,
p opo ional de i a i e gain, e c [27,29,30]. This addi ional gain o he
LQR gain ma ix se es o elimina e he high s eady-s a e e o , as well
as inc ease he con olle esponse. Howe e , his may esul in an un-
desi able o e shoo on he con olle esponse [31]. The e o e, o limi
he e ec o con olle o e shoo , i is necessa y o op imise he LQR
weighing ma ices [32]. The induc ion hea ing o s eel bille is cha ac-
e ised by high s eady-s a e e o , hence in his s udy, he LQR wi h an
in eg al gain is a po en ial linea con olle ha can be applied. Which o
he bes o he au ho s’knowledge, ha e no been epo ed in li e a u e.
The e o e, his s udy aims o model, simula e, and con ol he in-
duc ion hea ing p ocess o a ec angula s eel bille . This will be ach-
ie ed h ough he ollowing objec i es: Sou ce da a on he p ope ies o
s eel as a unc ion o empe a u e; De elop a model and simula e he
induc ion hea ing p ocess using a semi-nume ical me hod, inco po-
a ing skin e ec phenomena and he p ope ies o s eel as a unc ion o
empe a u e; Linea ise he model using he Jacobian me hod; Apply a
linea con ol scheme (i.e., LQR wi h in eg al ac ion in conjunc ion wi h
op imised weighing ma ices) o egula e he co e empe a u e o he
s eel bille . By accomplishing hese objec i es, he s udy in ends o
p o ide a comp ehensi e unde s anding o he induc ion hea ing p o-
cess and p opose an e ec i e con ol s a egy o main aining he
desi ed empe a u e in he s eel bille . Fu he mo e, he inno a i e
con ibu ion o his s udy would en ail a p edic i e es ima ion and
con ol o he nonlinea spa ial empe a u e a ia ion along he dep h o
s eel bille because o skin e ec , as well as he non-uni o m ene gy o
cu en dis ibu ion, especially ega ding powe losses associa ed wi h
his phenomenon.
2. Me hodologies
2.1. P ocess desc ip ion and model de elopmen
The Fig. 1a illus a es an induc ion hea ing u nace, in which he coil
ca ies high- equency al e na ing cu en (AC) and gene a es ime-
a ying elec omagne ic ields. These a ying ields, in u n, induce
eddy cu en s in he bille , leading o hea gene a ion h ough he Joule
e ec .
The go e ning model o desc ibe he induc ion hea ing can be
deduced om Ampe e’s law, Equa ions (1) a subse o Maxwell’s law
[33] in ela ion o Equa ions (2), (3), and (4) o magne ic lux densi y, B
(N.A
-1
m
−1
), elec ic lux densi y, D (C.m
−2
), and conduc ion cu en
densi y, J(A.m
−2
) espec i ely.
∇ × H=J+
∂
D
∂
(1)
D=
ε
E (2)
B=
μ
H (3)
J=
σ
E (4)
In sol ing o elec ic ield in ensi y, E (N.C
-1
), and magne ic ield in-
ensi y, H (A.m
−1
), i is mo e con enien o in oduce he magne ic
ec o po en ial, A(N.A
-1
), and elec ic scala po en ial, φ(V), as gi en
by Equa ion (5) –(6). The pa ame e s
ε
,
μ
,and
σ
a e espec i ely
emissi i y, magne ic pe meabili y (H.m
−1
), and elec ical conduc i i y
(S.m
−1
) o he bille .
S. Za ghoon e al.
Enginee ing Science and Technology, an In e na ional Jou nal 55 (2024) 101721
3
B= ∇ × A(5)
E= ∇φ+
∂
A
∂
(6)
Subs i u ing Equa ions (2), (3), and (4) in o (1) yields Equa ion (7d).
No ing ha
∂ε
/
∂
=0, since
ε
is a cons an , and (
σ
E+
ε∂
E/
∂
) ≈
σ
E.
1
μ
∇B=
σ
E+
∂
(
ε
E)
∂
(7a)
1
μ
∇B=
σ
E+
ε∂
E
∂
+E
∂ε
∂
(7b)
1
μ
∇B=
σ
E+
ε∂
E
∂
(7c)
1
μ
∇B=
σ
E (7d)
Fu he mo e, subs i u ing Equa ions (5) and (6) in o (7d), yields Equa-
ion (8b) upon ea anging. No ing ha J = −
σ
∇φ.
1
μ
∇2A=
σ
(∇φ+
∂
A
∂
)=
σ
∇φ+
σ∂
A
∂
(8a)
1
μ
∇2A= − J+
σ∂
A
∂
(8b)
Assuming hea ing o he s eel bille is implemen ed by al e na ing cu -
en , ypically a sinusoidal ha monic, so ha he magne ic ec o po-
en ial, A, and cu en densi y, J, can be desc ibed as A=Aoej
ω
and
J=Jsej
ω
espec i ely. No ing ha Aois he maximum magne ic ec o
po en ial and Jsis he maximum cu en densi y on bille su aces [8].
The e o e, Equa ion (8b) can be exp essed in e ms o he sinusoidal
ha monic, Equa ion (9c). In his equa ion, j implies an imagina y s a e,
ω
=2
π
, is he angula equency and, is he ime o obse a ion. In
sol ing his Equa ion (9c), he Di ichle (Ao=0) and Neumann
bounda y condi ions a e espec i ely used in he co e and ou e su aces
o he bille .
1
μ
∇2Aoej
ω
= − Jsej
ω
+
σ∂
(Aoej
ω
)
∂
(9a)
1
μ
∇2Aoej
ω
= − Jsej
ω
+
σ
Aoj
ω
ej
ω
(9b)
1
μ
∇2Ao= − Js+j
ωσ
Ao(9c)
The Equa ion (9c) can be modi ied o inco po a e he skin e ec , i.e., he
endency o he sou ce cu en o be g ea es a he su ace. This implies
ha a highe cu en is induced a he su ace o he bille bu dec eases
owa ds i s co e as desc ibed by Equa ion (10) [4,7] and Fig. 2.
In Equa ion (10),dis he dis ance o he co e o he bille o i s
su ace, δis he skin dep h, Equa ion (11) [4,7], i.e., he dep h below he
su ace o he bille a which he cu en densi y has allen by he in e se
exponen o he maximum sou ce cu en a he bille ou e su ace, Jso.
In Equa ion (11), (Hz) is he equency o he sou ce cu en , and
μ
is
he ela i e magne ic pe meabili y o he bille o acuum.
Js=Jsoe−d/δ(10)
δ=
1/(
πμσ
)
√≅503 
1/(
σμ
)
√(11)
The gene ic magne ic ec o po en ial on he bille su ace, Ao, can be
deduced om Equa ion (9c) –(11), and subsequen ly he sou ce hea
Fig. 1. Illus a ion o induc ion hea ing o a ec angula s eel bille .
Fig. 2. Illus a ion o he skin e ec o a ec angula bille .
S. Za ghoon e al.
Enginee ing Science and Technology, an In e na ional Jou nal 55 (2024) 101721
4
lux, q, can be calcula ed as gi en by Equa ion (12) in e ms o he
conduc ion cu en densi y, J, which can also be exp essed in e ms o
Aoas gi en by Equa ion (13) [8]. In his equa ion, Rsigni ies he eal
s a e, and (AoAo)*indica es a p oduc o complex conjuga e.
J=
σ
E= − j
ωσ
Ao(12)
q=[R(J)]2
σ
=1
2
ω
2
σ
(AoAo)*(13)
The hea di usion equa ion o induc ion hea ing on he bille can be
exp essed by Equa ion (14). In he equa ion,
ρ
(Kg.m
−3
) is he densi y,
Cp(J. Kg
−1
. K
−1
) is he speci ic hea capaci y, and k (W.m
−1
. k
−1
) is he
he mal conduc i i y o he bille .
ρ
Cp
∂
T
∂
=k∇2T+q (14)
The Equa ions (9c) and (14) a e sol ed using he o dina y di e en ial
equa ion sol e in MATLAB (speci ically, ode15s) in collabo a ion wi h
Equa ions (10), (11), and (13). In sol ing hese equa ions he spa ial
empe a u e dis ibu ion o pa ial di e en ial, ∇, was app oxima ed
using he Cen al scheme o he Fini e Di e ence Me hod (FDM). This
esul s in N
h
numbe o equa ions ha ep esen s each space poin , i =
1, 2, 3…, Nxas shown in Fig. 1b. Fu he mo e, in Equa ion (14) he
Neumann bounda y condi ion is applied o he co e, and o he ou e
su ace o he bille , adia ion hea loss, i.e., Equa ion (15), is applied.
The Equa ion (15) ela es he en i onmen o su ounding ai empe -
a u e, TE, o he ou e empe a u e o he bille , T, ia he S e an-
Bol zmann cons an ,
σ
*, he mal conduc i i y, and emissi i y.
k∇T=
σ
*
ε
(T4−T4
E)(15)
These bounda y condi ions applied a e based on he ope a ing condi ion
o indus ial induc ion hea ing u naces as well as popula ly epo ed in
li e a u e, o which bo h con ec ion and adia i e hea losses a e
conside ed [34,35]. Howe e , conside ing ha he adia i e hea loss is
≫con ec i e hea loss[36], only adia i e hea loss was conside ed in
his wo k.
3. P ocess con ol heo y
To e ec i ely con ol he sys em, i is necessa y o ini ially in es i-
ga e he sys em p ope ies such as he sys em’s dynamics, s abili y, and
pe o mance (i.e., he di e ence in s eady-s a e e o ). P io obse a-
ions o induc ion hea ing sys ems om li e a u e [7,15] indica e ha
he sys em is gene ally s able bu may exhibi signi ican s eady-s a e
e o . To mi iga e his s eady-s a e e o , a ull-s a e eedback Linea
Quad a ic Regula o (LQR) can be employed. In gene al, he LQR con ol
aims o con ol linea sys em s a e a iables using op imal con ol cos
unc ion in ini e ho izon o in ini e ho izon. I aims o deduce gain
ma ices ha minimises a cos unc ion o ensu e adequa e sys em s a-
bili y and con ol pe o mance. The e o e, he choice o weighing
ma ices (i.e., s a e weighing and inpu weighing ma ices) is c ucial o
adequa ely design he LQR [37].
The ull-s a e eedback LQR me hodology is based on he concep o a
cos unc ion, deno ed as, j, whe ein he objec i e is o de e mine an
op imal inpu o ac ua o e o , u, ha changes he s a e ec o , x,
h ough minimisa ion o Equa ion (16) by adjus ing he s a e weighing
ma ix, Q. The inpu weighing ma ix, R, is u ilised o penalise he
ac ua o e o and achie e he desi ed con ol pe o mance.
j=∫(xTQx +uTRu)d (16)
Q and R a e diagonal ma ices, which mus be posi i e-de ini e o ensu e
ha when mul iplied by he s a e, x(i.e., empe a u e), and inpu , u(i.e.,
sou ce hea lux), he column ec o s yield xTQx ≥0 and uTRu >0. The
dimension Q and R ma ices a e de ined by he numbe o s a es, s, and
inpu s, w, as gi en by Equa ion (17).
Q=⎡
⎢
⎢
⎣
Q1,10⋯0
0 Q2,2⋯0
0 0 ⋱ ⋮
0 0 ⋯Qs,s
⎤
⎥
⎥
⎦
&R=⎡
⎢
⎢
⎣
R1,10⋯0
0 R2,2⋯0
0 0 ⋱ ⋮
0 0 ⋯Rw,w
⎤
⎥
⎥
⎦(17)
Typically, he LQR using he in ini e ho izon cos unc ion, Equa ion
(16), can be s uc u ed as gi en by Equa ion (18). This ollows he
adi ional ze o- e minal quad a ic cos unc ion, applied in conjunc ion
wi h dec e ised s a e-space model, Equa ion (22) –(23) as cons ain s.
min
uj=1
2∑
N
k=1
(x(k )*−x(k ))TQ(x(k )*−x(k ))+u(k )TRu(k )(18)
Subjec o:
x(k +1) = Ax(k ) + Bu(k )
y(k ) = Cx(k )
In Equa ion (18), he ime scala o he p ocess model is di ided in o
disc e e poin s, k ∈ [1 N ], and x(k)*is he e e ence alue o ou pu s.
The solu ion o Equa ion (18) is equi alen o sol ing he Ricca i
Fig. 3. Illus a ion o a ull-s a e eedback sys em wi h in eg al ac ion.
Table 1
Simula ion pa ame e s u ilised o he induc i e hea ing u nace.
Pa ame e , symbol Value, uni
Du a ion o simula ion, 5000 s
Thickness o bille ,ϑ0.12 m
Numbe o disc e ised spa ial poin s,Nx5
Numbe o disc e ised ime poin s,N 100
Di e en ial dis ance be ween disc e ised spa ial poin s, Δxϑ/2Nx
S e an-Bol zmann cons an ,
σ
*5.6704 ×10
8
Cu ie empe a u e o s eel,Tc750 ℃
Pe meabili y ac o ,n 0.5
Magne ic pe meabili y o acuum,
μ
o4
π
×10
-7
H.m
−1
Rela i e pe meabili y o ca bon s eel a oom empe a u e,
μ
o 10H.m
−1
F equency o al e na ing cu en , [17] 100 Hz
Maximum cu en densi y on bille su ace,Jso 50 ×10
5
A.m
−2
Su ounding ai empe a u e,TE35 ℃
Ini ial empe a u e o bille ,Tin 1000 ℃
Ini ial maximum magne ic ec o po en ial,Ao0 N.A
-1
S. Za ghoon e al.
Enginee ing Science and Technology, an In e na ional Jou nal 55 (2024) 101721
5
equa ion [38]. In sol ing he Ricca i equa ion, o deduce he inpu , u,
Equa ion (21), he so-called “s a e ansi ion ma ix”, P, mus be
deduced ia Equa ion (19), o es ima e he con olle gain, K, Equa ion
(20).
P=ATP+PA−PBR−1BTP+Q (19)
K=R−1(BTP)(20)
u= − R−1(BTPx)(21)
The ull-s a e eedback LQR wi h in eg al ac ion wi hin he eed o wa d
pa hway connec ing he e o compa a o and he sys em is depic ed in
Fig. 3. In summa y o he diag am shown in Fig. 3, is he desi ed se -
poin , u is he p edic ed con olle inpu , A, B,and C a e ma ices ep-
esen ing he p ocess model, Kxis he LQR eedback gain ma ix and Kiis
Fig. 4. P ope ies o ca bon s eel as a unc ion o empe a u e.
S. Za ghoon e al.

Enginee ing Science and Technology, an In e na ional Jou nal 55 (2024) 101721
6
he in eg al gain. The addi ion o in eg al ac ion o he con ol sys em
allows he con olle o ack he se poin (ins ead o elying on he
classical null s abilisa ion) and, as such ensu e a be e pe u ba ion
obus ness [39]. The aim is o s abilise he con olled sys em as well as
in eg a e he s eady-s a e e o . Howe e , he in eg al ac ion can lead o
con olle o e shoo s, he e o e o esol e his issue[40], adequa e
measu es mus be aken o op imise he LQR h ough i s weigh ing
ma ices, i.e. Q and R. Adequa e design o uning o hese wo ma ices
a e he mos impo an pa ame e s o he LQR which can be edious and
challenging[41]. To apply he p ocess model o he LQR, i s s a e and
ou pu dynamics a e usually exp essed as s a e-space models, Equa ion
(22) and (23). The s a e-space model de ined by ma ices, A, B,and C
can be deduced om he linea isa ion o he dynamic models, Equa ion
(9c) and (14). No e ha C is a ma ix o ze os and ones ha de ines he
ou pu s om he s a es.
˙
x=Ax +Bu (22)
y=Cx (23)
˙
I= −y= −Cx (24)
The esul ing e o om he sys em, ˙
I, can be deduced om he di e -
ence be ween he desi ed se poin , , and ou pu , y, o he sys em.
Equa ions (22) and (24), can be used o design a new s a e-space model
wi h in eg al ac ion as gi en by he s a e-space ma ices, Equa ion (25),
Fig. 5. Illus a ion o empe a u e dynamics o he disc e ised bille su aces.
Table 2
E alua ion me ic o he compa ison o de eloped model
and da a o ca bon s eel p ope ies.
Ca bon s eel p ope ies R
2
Densi y,
ρ
0.9985
Speci ic hea capaci y,Cp0.9997
The mal conduc i i y, k 0.9998
Emissi i y,
ε
1.0000
Elec ical conduc i i y,
σ
0.9950
Fig. 6. Illus a ion o he esul ing magne ic ec o po en ial om he cu en densi y on bille su aces.
S. Za ghoon e al.
Enginee ing Science and Technology, an In e na ional Jou nal 55 (2024) 101721
7
which can be u he simpli ied in o Equa ion (26).
[˙
x
˙
I]=[A 0
−C 0][x
I]+[B
0]u+[0
1] (25)
˙
x=AIxI+BIu+B (26)
The con olle inpu , u, o he ull-s a e eedback LQR can be exp essed
as gi en in Equa ion (27). Subs i u ing Equa ion (27) in o (26) yields
Equa ion (28). The gain ma ices o he con olle a e gi en as K =
[Kx⋮Ki]. Typically Kxis a single- ow ma ix whose column size is he
same as he numbe o s a e a iables, x, i.e. he i e empe a u es o
each disc e ised poin , Kx=[Kx1, Kx2, Kx3, Kx4, Kx5] and Kia single
alue. The e o e, K =[Kx1, Kx2, Kx3, Kx4, Kx5, Ki].
u= − KxI(27)
˙
x= (AI−KBI)xI+B (28)
In he a o emen ioned equa ions, he gain ma ix, K, can be de e mined
h ough he u ilisa ion o he “lq ” unc ion in MATLAB R2021b by
sol ing he Ricca i equa ion, Equa ion (19), and he con ol sys em is
implemen ed wi h SIMULINK R2021b in he same design amewo k
gi en by Fig. 3.
4. Resul s and discussion
4.1. Simula ion o induc ion hea ing
To simula e he induc ion hea ing p ocess, Equa ions (9c) and (14)
a e conside ed as one-dimensional spa ial di e en ial, ∇, along he x-
axis, and disc e ised in o i e segmen s (i.e., Nx=5), Fig. 1b. The cen al
di e ence app oxima ions a e used o de i e i e equa ions o ep esen
Equa ions (9c) and (14) each. These equa ions a e hen sol ed simul-
aneously wi h Equa ions (10), (11), and (13). The simula ion esul s,
ob ained using he pa ame e s speci ied in Table 1, and Equa ions (29) –
(34), a e p esen ed in Figs. 4 –5.
The Equa ion (29) –(34) shows he empe a u e dependence (i.e., T
in kel ins, K) o he densi y,
ρ
[42,43], speci ic hea capaci y, Cp
[42,44–48], he mal conduc i i y, k [42–48], emissi i y,
ε
[49], elec-
ical conduc i i y,
σ
[50], and ela i e magne ic pe meabili y,
μ
[34,51,52], o ca bon s eel. The Equa ions (29) –(33) we e cu e- i ed
as a unc ion o empe a u e, using da a om he espec i e highligh ed
li e a u e. Howe e , Equa ion (34) is as epo ed in li e a u e. The
esul ing compa ison o he de eloped models and da a is shown in
Fig. 4, wi h hei espec i e R-squa ed (R
2
) e alua ion me ics gi en in
Table 2, and he esul showed adequa e i (R
2
>0.95).
ε
=⎧
⎨
⎩
T≤663→0.28
653 ≤T≤793→
T≥793→0.69
0.6999exp(−(T−848.7278
203.6169 )2)(32)
μ
=
μ
μ
o
=⎧
⎪
⎨
⎪
⎩
T≤1023→1+ (
μ
o −1)(1−(T−273
TC(◦C))n)
T≥1023→1
(34)
In Fig. 5a, upon he applica ion o sou ce cu en o he bille su ace, Js,
he ou e mos su ace (i.e., Nx5) ini ially expe iences a dec ease in
empe a u e due o adia ion hea loss, Equa ion (15) o he su ounding
ai ilm. This hea loss occu s because o he sha p empe a u e con as
ρ
(kg.m3)=⎧
⎪
⎪
⎪
⎪
⎪
⎨
⎪
⎪
⎪
⎪
⎪
⎩
273 ≤T≤1673K→7849exp(−(T−146.6
4858 )2)+367.2exp(−(T−1914
468 )2)
1673 ≤T≤1803K→12420exp( − 0.0003T)
1803 ≤T≤1973K→7120exp(1.635E−8T)
(29)
Cp(J.kg−1K−1)=⎧
⎪
⎪
⎪
⎪
⎪
⎨
⎪
⎪
⎪
⎪
⎪
⎩
273 ≤T≤993K→4780exp(−(T−996.9
13.39 )2)+710.2exp(−(T−2046
2627 )2)
993 ≤T≤1973K→4412exp(−(T−990.6
42.22 )2)+603.7exp(−(T−21320
356500 )2)(30)
k(W.m−1K−1)=⎧
⎪
⎪
⎪
⎪
⎪
⎨
⎪
⎪
⎪
⎪
⎪
⎩
273 ≤T≤1673K→12.13exp(−(T−349.3
384.4)2)+3.292E+15exp(−(T−641800
113100 )2)
1673 ≤T≤1773K→23.32exp( − 0.0002T) + 1.706E−13exp(0.0197T)
1773 ≤T≤1973K→259.7exp(7.473E−7T)
(31)
σ
(S.m−1) = ⎧
⎪
⎪
⎪
⎪
⎪
⎨
⎪
⎪
⎪
⎪
⎪
⎩
273 ≤T≤1773K→5364015.153exp(−(T−4977.4071
4746.0594 )2)+239071.5586exp(−(T−781.0628
317.3960 )2)
1773 ≤T≤1803K→1.4860exp(0.0074T)
1803 ≤T≤1973K→2.2269T2−8614.6218T+9202966.681
(33)
S. Za ghoon e al.
Enginee ing Science and Technology, an In e na ional Jou nal 55 (2024) 101721
8
be ween he bille su ace (1000 ℃) and he su ounding ai empe a-
u e (35 ℃). This hea loss akes place e en hough he skin e ec he-
o em p edic s ha Jsshould be maximum a Nx5as such should
expe ience an ins an aneous inc ease in empe a u e (i.e., in he absence
o hea loss). The e o e, due o his coun e ac ing e ec , a e his ini ial
hea loss, he su ace begins o s eadily inc ease in empe a u e. In
con as o he obse a ion a he ou e mos su ace, he inne su aces
(i.e., Nx4 o Nx1) expe ience an ins an aneous inc ease in empe a u e.
Du ing he ini ial phase o hea ing (≤205 s), he middle su ace, Nx3,
expe iences he highes hea ing a e, ollowed by Nx2and Nx1. This is
because he skin e ec has a g ea e in luence on he Nx3su ace
compa ed o adia ion hea loss. On he o he hand, he Nx4su ace
expe iences he lowes hea ing a e due o i s p oximi y o he ou e
su ace, esul ing in mo e adia ion hea loss o he su ounding su ace,
Nx5. Howe e , as he hea ing p ocess con inues (>205 s), he e is mo e
hea accumula ion owa d he co e, Nx1, o he bille due o negligible
adia ion hea loss o he su oundings. This con inues e en hough,
acco ding o he skin e ec heo em, Jsdec eases owa ds he co e, Nx1.
Consequen ly, as he hea ing p og esses, when all su aces each hei
s eady-s a e empe a u es a app oxima ely 4400 s, he co e empe a-
u e inally eaches he highes empe a u e. The p eceding su ace
empe a u es ollow he o de : Nx1(≅1259 ◦C) >Nx2(≅1258 ◦C) >
Nx3(≅1246 ◦C) >Nx4(≅1217 ◦C) >Nx5(≅1164 ◦C). Conside ing ha
he spa ial disc e isa ion o s eel bille was linea , i.e., Δx=ϑ/2Nx. I
would ha e been expec ed ha he spa ial dis ibu ion o empe a u e
along he dep h o he s eel bille should be co espondingly linea . The
esul s, (Tx1−Tx2≅1◦C, Tx2−Tx3≅12 ◦C, Tx3−Tx4≅29 ◦C,
Tx4−Tx5≅53 ◦C) and Fig. 5a, howe e shows o he wise. This highly
nonlinea spa ial empe a u e dis ibu ion can be a ibu ed o he skin
e ec , as well as he nonlinea bounda y condi ions (i.e., he Di ichle
bounda y condi ion ia adia i e hea loss on Nx5, and he Neumann
bounda y on Nx1).
Fu he mo e, Fig. 5b shows a su ace plo illus a ing he dynamic
changes in empe a u e ac oss he disc e ised su aces as ime p o-
g esses. I isually ep esen s how he empe a u e e ol es h oughou
he bille . In addi ion, Fig. 6 showcases he esul ing magne ic ec o
po en ial, A, on he bille su aces. As an icipa ed, A, dec eases owa ds
he co e o he bille , ollowing he in luence o he skin e ec on Js. This
demons a es he dis ibu ion and beha io o he magne ic ield wi hin
he induc ion hea ing p ocess.
The esul s o he simula ion o his model hus a con o m o hose
epo ed in li e a u e. Luozzo e al.[34] o ca bon s eel ube, Hansson
and Fisk [35] o s ainless s eel ube, Jelicic [53] and Baldan [54] o
solid s eel bille , all epo ed simila esul s ha e i y he esul o his
wo k. Conside ing ha bo h adia ion and con ec ion hea losses on he
su ace we e conside ed, he inne su ace also expe ienced highe
empe a u es han he ou e su ace, as ound in his s udy. Howe e , in
his wo k, he empe a u e di e ence was highe compa ed o hose
epo ed in li e a u e. The di e ence can be a ibu ed o he ob ious
di e ence in ini ial condi ion (i.e. 1000 ◦C in his wo k as compa ed o
abou 20 o 35 ◦C o o he wo k), emissi i y, as well as he assump ion
o con ec i e hea ans e conside ed. Fu he mo e Jankowski e al.
[55] al hough conside ed a cylind ical s eel bille , he end o he esul
also showed simila i y wi h ha o his wo k, e en hough he spa ial
Fig. 8. LQR empe a u e con ol o bille co e su ace based on changes in inpu weighing ma ix.
Fig. 7. LQR empe a u e con ol o bille co e su ace based on changes in s a e weighing ma ix.
S. Za ghoon e al.
Enginee ing Science and Technology, an In e na ional Jou nal 55 (2024) 101721
9
empe a u e dis ibu ion wasn’ conside ed in hei wo k. Janu ien˙
e and
Maˇ
zeika [56] howe e epo ed highe empe a u es on he bille su -
ace han a he co e o a solid cylind ical ca bon s eel bille , con a y o
he esul o his wo k. This disc epancy is a esul o negligible hea loss
assumed a he su ace, as such due o he p e ailing skin e ec , he
esul is expec ed o such an assump ion.
5. P ocess con ol o induc ion hea ing
To con ol he induc ion hea ing p ocess, he s a e-space model
desc ibed by Equa ion (22) –(28) is u ilised. To deduce he A and B
ma ices equi ed o his pu pose, he dynamic models, Equa ions (9c)
and (14) in conjunc ion wi h Equa ions (10), (11), and (13) a e line-
a ised using Taylo ’s se ies o Jacobian linea isa ion heo em, as
ex ensi ely desc ibed in li e a u e [11]. The alue o Jso gi en in Table 1
was used as he s eady-s a e alue o he maximum cu en densi y, Jso,
o he linea isa ion p ocess, along wi h he p elimina y s eady-s a e
empe a u e alues o su aces Nx1 o Nx5i.e., T1=1259 ◦C,T2=
1258 ◦C, T3=1246 ◦C, T4=1217 ◦C andT5=1164 ◦C. These empe a u e
alues a e used o he p elimina y calcula ion o ca bon s eel p ope ies
as highligh ed in Equa ion (29) –(34). The esul ing ma ices ob ained
om he linea isa ion p ocess a e gi en by Equa ion (35), wi h he
co esponding ou pu ma ix, C = [ 10000], and po en ially he
dis u bance o di ec ansmission ma ix, D =0. The ypical applica ion
o he ma ices A, B, and C o he s a e-space model, Equa ion (22) –
(23), is gi en in Equa ion (30b).
To e i y he ea lie assump ion o sys em s abili y, he open-loop
s abili y o he induc ion hea ing sys em was e alua ed based on he
eigen alues, λ, ob ained om he de e minan o he ma ix A and i s
iden i y ma ix, I, i.e, |A−λI|= 0. The sys em s abili y is con i med i all
eigen alues ha e nega i e eal pa s, i.e., λ<0. The esul ing eigen-
alues (−0.0017, −0.0308, −0.0954, −0.1850, −0.1599) a i m he
sys em is s able. Wi h he ma ices A and B de e mined, and s abili y
con i med, he LQR con ol sys em can be designed. The esul ing AI,
and BIma ices a e gi en by Equa ion (38), wi h he p e iously speci ied
ma ix, B =[ 0 1 ]T. The LQR, i espec i e o he speci ic alues
assigned o he Q and R pa ame e s, se es as design elemen s o
penalise bo h he s a e a iables and he con ol signals [57,58]. In his
s udy, he Q and R ma ices we e es ed unde wo case s udies. In case-
s udy-1, R was kep cons an while Q was a ied, and in case-s udy-2, Q
was kep cons an while R was assigned di e en alues. The objec i e o
hese case s udies was o e alua e he sys em’s ou pu when a emp ing
o each a se poin o 1200 ℃a he co e o he bille .
In case-s udy-1, he ollowing alues o Q and R ma ices we e
espec i ely applied, i.e., Q1=diag (1,1,1,1,1,100), Q2=diag
(1,1,1,1,1,150), Q3=diag (1,1,1,1,1,200), Q4=diag (1,1,1,1,1,300),
and Q5=diag (1,1,1,1,1,350), Q6=diag (1,1,1,1,1,600) as well as R =
2, i.e., he R ma ix is kep cons an . Fig. 7a illus a es he esponse o he
sys em upon he implemen a ion o a ull-s a e eedback LQR, wi h he
esul ing gain ma ices gi en by Equa ion (39). The obse ed beha iou
indica es ha as he alue o Q inc eases, he sys em app oaches he
se poin wi h minimal o e shoo and an imp essi e se ling ime bo h
be o e and a e he occu ence o dis u bance. Gene ally, in a con ol
sys em, a as e con olle esponse o se ling ime o each he se poin
is desi able, while o e shoo is conside ed undesi able. Among he
es ed alues o Q, Q3demons a es he bes pe o mance as i esul s in
minimal o e shoo and an adequa e se ling ime. Inc easing he alue
o Q gene ally leads o as e se ling ime bu a he expense o highe
o e shoo ing due o he ini ial o e eac ion o he con olle [59], which
may ha e esul ed in o ampli ied by he in eg al ac ion added o he
LQR. This can be obse ed in he esponse o Q6(i.e., he highes Q
alue), which achie es he as es se ling ime bu wi h he highes
o e shoo , while Q1(i.e., he lowes Q alue) exhibi s he opposi e
beha iou . O e all, he con olle p o es o be e ec i e e en in he
p esence o dis u bance a 4000 s, as he sys em quickly e u ns o i s
se poin . The esul ing maximum cu en densi y, Jso, o achie ing
empe a u e con ol o he bille co e su ace is depic ed in Fig. 7b.
A=⎡
⎢
⎢
⎢
⎢
⎣
−0.0908 0.0908 0.0000 0.0000 0.0000
0.0454 −0.0908 0.0454 0.0000 0.0000
0.0000 0.0455 −0.0909 0.0455 0.0000
0.0000 0.0000 0.0456 −0.0911 0.0456
0.0000 0.0000 0.0000 0.0915 −0.1091
⎤
⎥
⎥
⎥
⎥
⎦
&B
=⎡
⎢
⎢
⎢
⎢
⎣
5.1363E −8
1.7614E −7
3.2884E −7
4.6366E −7
5.3544E −7
⎤
⎥
⎥
⎥
⎥
⎦
(35)
˙
x=⎡
⎢
⎢
⎢
⎢
⎣
dT1/d
dT2/d
dT3/d
dT4/d
dT5/d
⎤
⎥
⎥
⎥
⎥
⎦
=⎡
⎢
⎢
⎢
⎢
⎣
−0.0908 0.0908 0.0000 0.0000 0.0000
0.0454 −0.0908 0.0454 0.0000 0.0000
0.0000 0.0455 −0.0909 0.0455 0.0000
0.0000 0.0000 0.0456 −0.0911 0.0456
0.0000 0.0000 0.0000 0.0915 −0.1091
⎤
⎥
⎥
⎥
⎥
⎦
⏟⏞⏞⏟
A
⎡
⎢
⎢
⎢
⎢
⎣
T1
T2
T3
T4
T5
⎤
⎥
⎥
⎥
⎥
⎦
⏟⏞⏞⏟
x
+⎡
⎢
⎢
⎢
⎢
⎣
5.1363E −8
1.7614E −7
3.2884E −7
4.6366E −7
5.3544E −7
⎤
⎥
⎥
⎥
⎥
⎦
⏟⏞⏞⏟
B
q
⏟⏞⏞⏟
u
(36)
Fig. 9. LQR empe a u e con ol o bille co e su ace based on speci ic s a e and inpu weighing ma ix.
S. Za ghoon e al.