IAES
Inter
national
J
our
nal
of
Robotics
and
A
utomation
(IJRA)
V
ol.
15,
No.
2,
June
2026,
pp.
307
∼
318
ISSN:
2722-2586,
DOI:
10.11591/ijra.v15i2.pp307-318
❒
307
Design
and
implementation
of
NMPC
f
or
a
tw
o-DOF
r
obotic
arm
using
CasADi
Lahcen
Boulbalah
1
,
F
aiza
Dib
1
,
Nabil
Benaya
1
,
Khaddouj
Ben
Meziane
2
1
Department
of
Ph
ysics,
F
aculty
of
Sciences
and
T
echniques,
Abdelmalek
Essaadi
Uni
v
ersity
,
T
etouan,
Morocco
2
Laboratory
of
Inno
v
ation
in
Management
and
Engineering
(LIMIE),
Higher
Institute
of
Engineering
and
Business
(ISGA),
Fez,
Morocco
Article
Inf
o
Article
history:
Recei
v
ed
No
v
14,
2025
Re
vised
Feb
1,
2026
Accepted
Apr
19,
2026
K
eyw
ords:
CasADi
Dynamic
model
Model
predicti
v
e
control
Nonlinear
control
T
w
o-link
robot
arm
ABSTRA
CT
Achie
ving
ac
curate
joint-space
tracking
in
multi-link
robotic
arms
is
compli-
cated
by
strong
conguration-dependent
nonlinearities
and
mandatory
actuator
limits
that
classical
controllers
are
structurally
unable
to
enforce.
This
paper
presents
a
nonlinear
model
predicti
v
e
control
(NMPC)
scheme
for
a
tw
o-de
gree-
of-freedom
(2-DOF)
serial
robotic
arm,
implemented
within
the
CasADi
sym-
bolic
computing
en
vironment
to
le
v
erage
automatic
dif
ferentiation
and
sparse
interior
-point
solving.
The
complete
set
of
Lagrangian
equations
of
motion-
inertia,
Coriolis,
and
gra
vity
terms-is
incorporated
directly
into
the
optimizer’
s
prediction
model
through
fourth-order
Runge-K
utta
(RK4)
inte
gration,
eliminat-
ing
the
need
for
linearization.
T
orque,
v
elocity
,
and
angle
bounds
are
imposed
as
nati
v
e
hard
inequality
constraints
at
e
v
ery
step
of
the
nite-horizon
optimization.
Systematic
simulations
pit
the
proposed
NMPC
ag
ainst
a
Zie
gler
-Nichols-tuned
decentralized
PID
at
tw
o
distinct
sampling
periods.
The
NMPC
achie
v
ed
a
95%
reduction
in
peak
tracking
error
relati
v
e
to
PID
(0.0058
rad
vs.
0.1347
rad
for
Joint
1),
with
mean
error
decreases
of
64.65%
and
57.58%
for
Joi
nts
1
and
2
respecti
v
ely
,
at
an
a
v
erage
solv
er
time
of
0.053
s-comfortably
within
the
0.1
s
control
c
ycle.
The
ndings
demonstrate
that
online
NMPC
with
unabridged
nonlinear
dynamics
is
computationally
practical
for
real-time
joint
control
on
standard
computing
hardw
are.
This
is
an
open
access
article
under
the
CC
BY
-SA
license
.
Corresponding
A
uthor:
Lahcen
Boulbalah
Department
of
Ph
ysics,
F
aculty
of
Sciences
and
T
echniques,
Abdelmalek
Essaadi
Uni
v
ersity
T
etouan,
Morocco
Email:
boulbalah.lahcen@gmail.com
1.
INTR
ODUCTION
High-delity
joint-space
tracking
in
robotic
manipulators
hinges
on
the
ability
to
handle
coupled,
conguration-dependent
nonlinearities
while
simultaneously
respecting
actuator
saturation
limits.
Proportional-
inte
gral-deri
v
ati
v
e
(PID)
controllers
are
widely
adopted
in
practice
due
to
their
simplicity
,
b
ut
the
per
-joint,
decoupled
structure
means
the
y
cannot
compensate
Coriolis
and
inertia
cross-coupling,
nor
can
the
y
enforce
torque
or
v
elocity
bounds
in
a
principled
manner
[1],
[2].
Model
predicti
v
e
control
reframes
the
actuation
problem
as
a
constrained
nite-horizon
optimization
e
x
ecuted
at
e
v
ery
sampling
step.
Linear
model
predic-
ti
v
e
control
(MPC)
v
ariants
are
computationally
lightweight
b
ut
introduce
linearization
mismatch
for
plants
with
strong
nonlinear
beha
vior
.
Nonlinear
model
predicti
v
e
control
(NMPC)
elim
inates
this
mismatch
by
em-
bedding
the
complete
no
nl
inear
plant
model
inside
the
optimizer;
recent
progress
in
sparse
NLP
s
olv
ers
and
symbolic
computat
ion
tools
has
brought
online
NMPC
within
reach
of
both
embedded
processors
and
re-
J
ournal
homepage:
http://ijr
a.iaescor
e
.com
Evaluation Warning : The document was created with Spire.PDF for Python.
308
❒
ISSN:
2722-2586
search
platforms
[3],
[4].
A
wide
range
of
alternati
v
e
controllers
has
been
studied
for
robotic
manipulation.
Neural-netw
ork
path
planners
for
mobile
robots
[5],
fuzzy-logic
re
gulators
[6],
[7],
and
adapti
v
e
controllers
[8]
each
impro
v
e
on
basic
PID
performance,
yet
none
enforces
ph
ysical
constraints
as
an
inte
gral
design
feature.
MPC
applied
to
robotic
arms
[9]-[11]
substantiates
the
benets
of
prediction-based
actuation,
though
man
y
re-
ported
implementations
rely
on
linearized
or
reduced-order
plant
models.
Recent
contrib
utions
ha
v
e
pushed
the
boundary
of
intelligent
h
ybrid
control:
a
three-link
manipulator
controller
combining
adapti
v
e
fuzzy
logi
c
with
sliding-mode
disturbance
rejection
w
as
de
v
eloped
in
[12];
a
chattering-free
fuzzy
super
-twisting
sliding-mode
scheme
for
2-DOF
arms
under
parametric
uncertainty
w
as
reported
in
[13];
and
a
neural-netw
ork-enhanced
nonlinear
PID
with
backstepping
for
wheeled
mobile
robots
w
as
presented
in
[14].
Notwithstanding
these
adv
ances,
hard
constraint
satisf
action
is
a
structural
property
of
MPC
that
is
absent
from
sliding-mode
and
neural-netw
ork
paradigms
by
construction.
Laguerre-function
MPC
[15]-[17]
and
P
article
sw
arm
optimization
(PSO)-optimized
predicti
v
e
controllers
[18]
demonstrate
ef
cienc
y
g
ains
in
po
wer
electronics
and
autonomous
v
ehicle
domains;
ho
we
v
er
,
full-order
NMPC
from
the
complete
Lagrangian
model
remains
the
gold
standard
for
arm
s
with
pronounced
conguration-dependent
nonlinearities
[19],
[20].
The
specic
contrib
utions
of
this
w
ork
are:
−
A
fully
symbolic
NMPC
architecture
for
a
2-DOF
serial
arm
implemented
in
CasADi,
using
RK4
inte
gra-
tion
of
the
complete
Lagrangian
model
and
e
x
ecuting
at
10
Hz
online
−
Rigorous
comparison
ag
ainst
a
Zie
gler
-Nichols-tuned
PID
baseline,
yielding
95%
peak-error
reduction
and
mean
impro
v
ements
abo
v
e
57%
for
both
joints
−
Structural
guarantee
of
torque
and
v
elocity
constraint
satisf
action
at
e
v
ery
time
step-a
property
not
achie
v-
able
with
the
PID
baseline
−
Demonstration
that
CasADi/IPOPT
-based
NMPC
operates
within
the
0.1
s
control
period
on
of
f-the-shelf
laptop
hardw
are
Section
2
of
this
paper
deri
v
es
the
full
Lagrangian
dynamic
model
of
the
tw
o-link
arm.
Section
3
formulates
the
NMPC
problem
and
describes
the
PID
benchmark.
Section
4
presents
simulation
results
and
a
comparati
v
e
analysis.
Section
5
summarizes
the
w
ork
and
identies
directions
for
future
research.
2.
PR
OBLEM
FORMULA
TION
AND
METHODS
2.1.
Dynamic
equations
The
plant
under
st
udy
is
a
planar
tw
o-link
robotic
arm.
F
or
each
link
i
∈
{
1
,
2
}
,
θ
i
denotes
the
joint
angle,
L
i
the
link
length,
and
M
i
the
link
mass;
g
is
the
gra
vitat
ional
acceleration.
The
ph
ysical
conguration
is
illustrated
in
Figure
1.
Figure
1.
T
w
o-link
robot
arm
IAES
Int
J
Rob
&
Autom,
V
ol.
15,
No.
2,
June
2026:
307–318
Evaluation Warning : The document was created with Spire.PDF for Python.
IAES
Int
J
Rob
&
Autom
ISSN:
2722-2586
❒
309
The
equations
of
motion
are
deri
v
ed
using
the
Lagrangian
ener
gy
method
[21].
This
approach
con-
structs
the
scalar
Lagrangian
L
=
T
−
V
from
the
total
kinetic
and
potential
ener
gies,
then
applies
the
Euler
-
Lagrange
operator
to
obtain
the
go
v
erning
dif
ferential
equations
for
each
generalized
coordinate.
2.1.1.
Lagrangian
f
ormulation
W
ith
the
joint
angles
θ
1
and
θ
2
as
generalized
coordinates,
the
kinetic
ener
gy
T
and
potential
ener
gy
V
are
e
v
aluated
and
combined
to
yield
the
Lagrangian
L
.
The
total
kinetic
ener
gy
T
is
gi
v
en
by:
T
=
1
2
M
1
˙
x
2
1
+
1
2
M
2
(
˙
x
2
1
+
˙
x
2
2
+
2
˙
x
1
˙
x
2
cos(
θ
2
))
,
(1)
where
˙
x
1
and
˙
x
2
are
the
v
elocities
of
the
centers
of
mass
of
the
links,
and
M
1
and
M
2
are
the
masses
of
the
rst
and
second
links,
respecti
v
ely
.
Gra
vitational
potential
ener
gy
stored
in
both
links
summed
o
v
er
their
respecti
v
e
center
-of-mass
height
s
gi
v
es:
V
=
(
M
1
g
L
1
cos(
θ
1
))
+
(
M
2
g
(
L
1
cos(
θ
1
)
+
L
2
cos(
θ
1
+
θ
2
)))
,
(2)
where
g
denotes
gra
vitational
acceleration.
The
Lagrangian
is
therefore:
L
=
T
−
V
.
(3)
Substituting
L
into
the
Euler
-Lagrange
operator
for
each
generalized
coordinate
yields
the
go
v
erning
torque
equations:
d
dt
∂
L
∂
˙
θ
i
−
∂
L
∂
θ
i
=
τ
i
,
i
=
1
,
2
,
(4)
with
τ
1
and
τ
2
being
the
joint
torques
that
serv
e
as
control
inputs.
2.1.2.
Inertia
matrix
Dif
ferentiating
L
twice
with
respect
to
the
generalized
v
elocities
produces
the
conguration-dependent
inertia
matrix
M
(
θ
2
)
:
M
(
θ
2
)
=
D
1
D
2
D
2
D
4
,
(5)
where:
D
1
=
(
M
1
+
M
2
)
L
2
1
+
M
2
L
2
2
+
2
M
2
L
1
L
2
cos(
θ
2
)
,
(6)
D
2
=
M
2
L
2
2
+
M
2
L
1
L
2
cos(
θ
2
)
,
(7)
and
D
4
=
M
2
L
2
2
.
(8)
2.1.3.
Coriolis
matrix
Inter
-joint
v
elocity
interactions
generate
Coriolis
forces
that
are
collected
in
the
matrix
C
(
˙
θ
1
,
˙
θ
2
,
θ
2
)
:
C
(
˙
θ
1
,
˙
θ
2
,
θ
2
)
=
C
1
C
2
,
(9)
where:
C
1
=
−
M
2
L
1
L
2
2
˙
θ
1
˙
θ
2
+
˙
θ
2
2
sin(
θ
2
)
,
(10)
and
C
2
=
−
M
2
L
1
L
2
˙
θ
1
˙
θ
2
sin(
θ
2
)
.
(11)
Design
and
implementation
of
NMPC
for
a
two-DOF
r
obotic
arm
using
CasADi
(Lahcen
Boulbalah)
Evaluation Warning : The document was created with Spire.PDF for Python.
310
❒
ISSN:
2722-2586
2.1.4.
Gra
vity
v
ector
The
reaction
torques
at
each
joint
due
to
gra
vity
are
assembled
into
the
v
ector:
G
(
θ
1
,
θ
2
)
=
G
1
G
2
,
(12)
where:
G
1
=
−
(
M
1
+
M
2
)
g
L
1
sin(
θ
1
)
−
M
2
g
L
2
sin(
θ
1
+
θ
2
)
,
(13)
and
G
2
=
−
M
2
g
L
2
sin(
θ
1
+
θ
2
)
.
(14)
2.1.5.
State
equations
Joint
accelerations
¨
θ
are
reco
v
ered
by
in
v
erting
the
inertia
matrix:
¨
θ
=
M
−
1
(
τ
−
C
−
G
)
,
(15)
where
τ
=
[
τ
1
,
τ
2
]
⊤
collects
the
joint
torques
acting
as
control
inputs,
and
C
,
G
are
the
v
elocity-coupling
and
gra
vitational
terms
deri
v
ed
abo
v
e.
2.1.6.
F
orward
Kinematics
The
Cartesian
coordinates
(
x,
y
)
of
the
end-ef
fector
are
obtained
from
the
joint
angles
through
the
standard
planar
chain
geometry:
x
=
L
1
cos(
θ
1
)
+
L
2
cos(
θ
1
+
θ
2
)
,
(16)
and
y
=
L
1
sin(
θ
1
)
+
L
2
sin(
θ
1
+
θ
2
)
.
(17)
These
geometric
e
xpressions
permit
joint-space
tracking
metrics
to
be
directly
interpreted
as
end-ef
fector
po-
sitioning
accurac
y
in
the
w
orkspace.
3.
CONTR
OL
STRA
TEGY
At
each
control
step,
the
proposed
NMPC
formulat
es
and
solv
es
a
constrained
nonlinear
program
(NLP)
that
selects
joint
torques
so
as
to
jointly
minimize
trajectory
tracking
de
viation
and
actuator
ener
gy
consumption
o
v
er
a
nite
prediction
horizon.
Future
states
within
that
horizon
are
propag
ated
using
a
fourth-
order
Runge-K
utta
(RK4)
inte
gration
scheme
[22],
which
preserv
es
the
arm’
s
conguration-dependent
inertial
and
Coriolis
characteristic
s
without
approximation.
Bounds
on
joint
position,
v
elocity
,
and
torque
are
imposed
directly
inside
the
NLP
as
hard
inequality
constraints,
so
the
optimizer
can
only
produce
solutions
that
are
ph
ysically
realizable.
3.1.
Cost
function
The
NLP
objecti
v
e
accumulates,
o
v
er
the
N
-step
prediction
windo
w
,
a
quadratic
penalty
on
state
de
viation
from
the
reference
trajectory
and
a
quadratic
penalty
on
the
applied
torques,
together
with
a
terminal
stage
cost:
J
=
N
−
1
X
k
=0
(
x
k
−
x
ref
)
⊤
Q
(
x
k
−
x
ref
)
+
u
⊤
k
Ru
k
+
x
⊤
N
Q
N
x
N
(18)
Here
x
k
∈
R
4
collects
the
joint
angles
and
angular
v
elocities
at
prediction
step
k
,
x
ref
is
the
reference
trajectory
to
be
track
ed,
and
u
k
∈
R
2
holds
the
tw
o
joint
torques.
The
positi
v
e
semi-denite
weight
Q
balances
angle
and
v
elocity
tracking
accurac
y
,
R
penalizes
actuator
ef
fort
to
pre
v
ent
e
xcessi
v
e
torque
usage,
and
the
terminal
weight
Q
N
pro
vides
a
L
yapuno
v-lik
e
endpoint
penalty
that
promotes
closed-loop
stability
.
IAES
Int
J
Rob
&
Autom,
V
ol.
15,
No.
2,
June
2026:
307–318
Evaluation Warning : The document was created with Spire.PDF for Python.
IAES
Int
J
Rob
&
Autom
ISSN:
2722-2586
❒
311
3.2.
P
arameter
selection
The
prediction
horizon
N
and
the
weight
matrices
were
selected
fol
lo
wing
the
systematic
tuning
frame
w
ork
of
[23],
[24].
A
grid
search
o
v
er
N
∈
{
3
,
5
,
8
,
10
}
sho
wed
that
N
=
5
pro
vides
the
optimal
cost-computation
tradeof
f:
hori
zons
shorter
than
5
produced
near
-sighted
torque
commands
with
persistently
ele
v
ated
s
teady-state
error
,
whereas
e
xtending
to
N
>
5
reduced
RMSE
by
under
2%
at
a
disproportionate
com-
putational
cost
[4].
Initial
diagonal
entries
of
Q
were
set
according
to
Bryson’
s
in
v
erse-square
rule
[24],
using
the
reciprocal
of
the
squared
admissible
de
viation
for
each
state.
This
w
as
follo
wed
by
iterati
v
e
simulation-
based
renement,
yielding
Q
=
di
ag
(55
,
55
,
0
.
1
,
0
.
1)
and
R
=
diag
(0
.
001
,
0
.
001)
.
The
ele
v
ated
Q
11
/
Q
22
entries
prioritize
angular
accurac
y;
the
smaller
Q
33
/
Q
44
v
alues
accept
moderate
v
elocity
de
viations;
and
the
small
R
entries
allo
w
adequat
e
torque
authority
for
tracking
demanding
trajectories
without
unnecessary
ener
gy
consumption
[23].
3.3.
System
dynamics
Dening
the
state
v
ector
as
x
=
[
θ
1
,
θ
2
,
˙
θ
1
,
˙
θ
2
]
⊤
and
the
input
v
ector
as
u
=
[
τ
1
,
τ
2
]
⊤
,
the
arm’
s
continuous-time
dynamics
tak
e
the
compact
form:
˙
x
=
f
(
x
,
u
)
,
(19)
where
the
right-hand
side
f
encapsulates
the
complete
Lagrangian
model,
co
v
ering
conguration-v
arying
iner
-
tia,
Coriolis
v
elocity
coupling,
and
gra
vity
loading
as
deri
v
ed
in
section
2.
Each
predicted
state
in
the
horizon
is
adv
anced
by
one
sampling
step
through
the
classical
fourth-
o
r
der
Runge-K
utta
(RK4)
inte
grator
[22]:
x
k
+1
=
x
k
+
∆
t
6
(
k
1
+
2
k
2
+
2
k
3
+
k
4
)
,
(20)
where
k
1
=
f
(
x
k
,
u
k
)
,
k
2
=
f
(
x
k
+
∆
t
2
k
1
,
u
k
)
,
k
3
=
f
(
x
k
+
∆
t
2
k
2
,
u
k
)
,
k
4
=
f
(
x
k
+∆
t
k
3
,
u
k
)
.
3.4.
Constraints
At
each
of
the
N
shooting
nodes
k
∈
{
0
,
.
.
.
,
N
−
1
}
the
NLP
enforces
three
constraint
groups:
bilat-
eral
bounds
on
the
state
v
ector
x
min
≤
x
k
≤
x
max
,
bilateral
bounds
on
the
control
input
u
min
≤
u
k
≤
u
max
,
and
continuity
equations
x
k
+1
=
f
(
x
k
,
u
k
)
that
couple
successi
v
e
nodes.
Bound
v
alues
are
listed
in
T
able
1.
T
able
1.
Simulation
parameters
and
conguration
P
arameter
Symbol
V
alue
Link
masses
M
1
,
M
2
1
.
0
kg
Link
lengths
L
1
,
L
2
1
.
0
m
Gra
vity
g
9
.
81
m
/
s
2
Simulation
time
–
50
s
Sampling
time
T
0
.
1
s
Prediction
horizon
N
5
State
weight
Q
d
i
a
g
(55
,
55
,
0
.
1
,
0
.
1)
Control
weight
R
diag
(0
.
001
,
0
.
001)
PID
g
ains
(
K
p
,
K
i
,
K
d
)
–
500
,
10
0
,
20
3.5.
Implementation
details
The
complete
NLP
is
assembled
using
CasADi’
s
symbolic
API
within
MA
TLAB
[3],
[25],
which
constructs
e
xact
rst-
and
second-order
deri
v
ati
v
es
of
both
the
objecti
v
e
and
the
constraints
via
automatic
dif-
ferentiation.
The
resulting
sparse
der
i
v
ati
v
e
structures
are
passed
directly
to
IPOPT
[26],
an
interior
-point
solv
er
that
e
xploits
this
sparsity
for
ef
cient
matrix
f
actorization.
State
discretization
uses
the
direct
multiple-
shooting
transcription
[27],
[28]:
each
node’
s
state
x
k
is
introduced
as
an
independent
optimization
v
ariable
and
coupled
to
its
successor
through
the
RK4
equality
constraint,
which
yields
a
well-conditioned
problem
compared
with
single-shooting
alternati
v
es.
At
each
control
c
ycle,
only
the
leading
torque
v
ector
of
the
com-
puted
optim
al
sequence
U
opt
is
dispatched
to
the
plant;
the
whole
procedure
then
repeats
from
the
updated
state
measurement.
Design
and
implementation
of
NMPC
for
a
two-DOF
r
obotic
arm
using
CasADi
(Lahcen
Boulbalah)
Evaluation Warning : The document was created with Spire.PDF for Python.
312
❒
ISSN:
2722-2586
3.5.1.
PID
contr
oller
The
baseline
controller
is
a
per
-joint
decentralized
PID
[1],
[2],
in
which
each
joint
independently
computes
its
torque
as
τ
i
(
t
)
=
K
pi
e
i
(
t
)
+
K
ii
R
t
0
e
i
dξ
+
K
di
˙
e
i
(
t
)
,
with
e
i
=
θ
i,
ref
−
θ
i
the
signed
angular
tracking
error
for
joint
i
.
3.5.2.
PID
tuning
Initial
g
ains
were
determined
by
the
Zie
gler
-Nichols
ultimate-g
ain
procedure
[29],
[30]
and
further
rened
through
simulation
trials
for
each
of
the
tw
o
sampling
periods
tested
(dt
=
0.01
s
and
dt
=
0.1
s);
the
nal
v
alues
are
compiled
in
T
able
1.
The
requirement
for
independent
re-tuning
a
t
each
sampling
rate,
com-
bined
with
the
controller’
s
inherent
inability
to
account
for
Coriolis
cross-coupling,
underlines
a
fundamental
limitation
of
the
decentralized
PID
architecture.
4.
SIMULA
TION
SETUP
AND
RESUL
TS
4.1.
Simulation
parameters
The
NMPC
frame
w
ork
w
as
v
alidated
on
a
tw
o-link
arm
with
parameters
and
conguration
sum-
marised
in
T
able
1.
Constraints
bound
joint
angles
to
±
π
rad,
v
elocities
to
±
10
rad/s,
and
torques
to
±
20
Nm.
4.2.
Contr
ol
system
ar
chitectur
e
Both
controllers
operate
within
the
closed-loop
architecture
depicted
in
Figure
2,
where
the
plant
output
(measured
joint
angles
and
v
elocities)
feeds
back
to
the
controller
at
e
v
ery
sampling
instant.
Figure
2.
Closed-loop
NMPC
4.3.
Simulation
r
esults
Full-range
sinusoidal
reference
signals
spanning
the
entire
joint
w
orkspace
were
issued
to
both
con-
trollers,
creating
a
stringent
e
v
aluation
of
transient
agility
and
long-run
ste
ady-state
performance
under
realistic
operating
conditions.
4.3.1.
P
erf
ormance
analysis
The
quantitati
v
e
adv
antage
of
the
NMPC
visible
in
T
able
2
originates
from
three
structural
properties
that
distinguish
it
from
an
y
reacti
v
e
controller
.
−
Anticipatory
torque
generation:
The
NMPC
computes
an
optimal
torque
sequence
by
minimizing
the
cost
J
o
v
er
a
rolling
v
e-step
horizon,
so
the
act
u
a
tor
is
adjusted
in
adv
ance
of
predicted
trajectory
curv
ature
rather
than
after
it
produces
an
error
.
This
pre-empti
v
e
capability
is
architecturally
absent
from
e
v
ery
reacti
v
e
control
scheme.
−
Coupled
dynamics
compensation:
Because
the
prediction
model
incorporates
the
full
Coriolis
matrix
C
(
˙
θ
1
,
˙
θ
2
,
θ
2
)
,
the
cross-joint
v
elocity
interactions
are
accounted
for
when
computing
optimal
torques.
The
PID
handles
each
joint
independently;
Coriolis
ef
fects
then
appear
as
unmeasured
disturbances
that
gro
w
with
joint
speed
and
manifest
as
persistent
tracking
of
fset.
−
Built-in
constraint
satisf
action:
The
inequality
constraints
|
τ
i
|
≤
20
Nm
,
|
θ
i
|
≤
π
rad
,
and
|
˙
θ
i
|
≤
10
rad
/
s
are
embedded
inside
the
NLP
,
so
the
optimizer
is
mathematically
prohibited
from
returning
limit-violating
torques.
The
PID
relies
on
an
e
xternal
saturation
element
which,
when
acti
v
e,
interacts
with
the
inte
grator
state
and
can
erode
stability
mar
gins.
IAES
Int
J
Rob
&
Autom,
V
ol.
15,
No.
2,
June
2026:
307–318
Evaluation Warning : The document was created with Spire.PDF for Python.
IAES
Int
J
Rob
&
Autom
ISSN:
2722-2586
❒
313
The
statistical
tracking
metrics
and
computation
times
across
both
sampling
rates
are
compiled
in
T
able
2.
At
dt
=
0.01
s,
the
NMPC
recorded
mean
errors
of
0.0035
rad
and
0.0014
rad
for
Joints
1
and
2,
ag
ainst
0.0099
rad
and
0.0033
rad
for
the
PID,
corresponding
to
reductions
of
64.65%
and
57.58%
respecti
v
ely
.
Peak
de
viations
dropped
from
0.1347
rad
to
0.0058
rad
for
Joint
1
(
−
95.69%)
and
from
0.0510
rad
to
0.0019
rad
f
o
r
Joint
2
(
−
96.27%).
At
the
coarser
dt
=
0.1
s
rate,
analogous
impro
v
ements
of
56.00%
and
63.89%
were
obtained
with
no
adjustment
to
the
NMPC
weight
matrices.
The
mean
solv
e
time
remained
between
0.053
s
and
0.056
s
across
both
rates,
safely
inside
the
0.1
s
b
udget.
T
able
2.
NMPC
and
PID
performance
comparison
dt
P
arameters
J
oint
1
J
oint
2
T
iming
(s)
(lr)3-5
(lr)6-8
(lr)9-10
Mean
Max
RMS
Mean
Max
RMS
Mean
Max
0
.
010
NMPC
0.0035
0.0058
0.0038
0.0014
0.0019
0.0015
0.053
0.261
PID
0.0099
0.1347
0.0170
0.0033
0.0510
0.0050
3
.
61
×
10
−
6
0.001595
Impr
o
v
ement
(%)
64.65
95.69
77.65
57.58
96.27
70.00
-
-
(lr)2-11
NMPC
Q
=
500
R
=
0
.
001
Q
=
55
R
=
0
.
0001
PID
Kp
=
400
Ki
=
100
Kd
=
20
0
.
100
NMPC
0.0055
0.0076
0.0055
0.0013
0.0022
0.0015
0.056
0.246
PID
0.0125
0.1374
0.0177
0.0036
0.0493
0.0056
1
.
08
×
10
−
5
0.001315
Impr
o
v
ement
(%)
56.00
94.47
68.93
63.89
95.54
73.21
-
-
(lr)2-11
F
or
NMPC
Q
=
500
R
=
0
.
001
Q
=
55
R
=
0
.
0001
PID
Kp
=
30
Ki
=
15
Kd
=
1.2
4.3.2.
Real-time
perf
ormance
and
parameter
sensiti
vity
All
tests
ran
on
an
Intel
Core
i7-8550U
laptop
under
MA
TLAB
R2023a.
A
v
erage
NMPC
solv
e
ti
mes
were
0.053
s
at
dt
=
0.01
s
and
0.056
s
at
dt
=
0.1
s,
both
k
eeping
comfortably
inside
the
0.1
s
control
b
udget.
The
single
peak
of
0.261
s
occurred
only
at
the
rst
iteration
due
to
a
cold-start
initialization;
once
a
w
arm-start
guess
w
as
a
v
ailable,
e
v
ery
subsequent
solv
e
sat
well
belo
w
the
b
udgeted
0.1
s.
The
horizon
sensiti
vity
study
o
v
er
N
∈
{
3
,
5
,
8
,
10
}
v
eried
that
N
=
5
strik
es
the
best
balance
bet
ween
accurac
y
and
computation
load:
smaller
windo
ws
led
to
myopic
torque
strate
gies
with
noticeable
steady-state
of
fset,
while
N
=
10
trimmed
RMSE
by
no
more
than
2%
at
a
1.8
×
increase
in
a
v
erage
solv
e
time.
The
weight
matrices
selected
during
tuning
remained
insensiti
v
e
to
changes
in
trajectory
shape
because
the
receding-horizon
correction
mechanism
automatically
compensates
for
plant-model
discrepancies
at
each
control
step,
eliminating
t
he
need
for
g
ain
scheduling.
The
w
orkspace
trajectory
traced
by
the
end-ef
fector
under
both
controllers
is
sho
wn
in
Figure
3.
The
NMPC
path
closely
follo
ws
the
reference
locus,
whereas
the
PID
trajectory
e
xhibits
notable
de
viations,
particularly
at
trajectory
turning
points.
Figure
3.
Comparison
of
the
end
ef
fector
trajectory
of
NMPC
vs
PID
Figure
4
presents
the
angle
time
histories
for
each
joint
under
both
controllers.
The
NMPC
t
racks
the
sinusoidal
reference
with
near
-zero
residual
error
throughout
the
full
50
s
run;
proacti
v
e
torque
planning
k
eeps
the
controller
ahead
of
reference
curv
ature
changes.
The
PID
e
v
entually
con
v
er
ges
b
ut
sho
ws
noticeable
phase
lag
during
f
ast
transitions
and
persistent
oscillation
on
Joint
2,
which
is
a
characteristic
symptom
of
uncompensated
Coriolis
coupling
from
the
other
joint.
Design
and
implementation
of
NMPC
for
a
two-DOF
r
obotic
arm
using
CasADi
(Lahcen
Boulbalah)
Evaluation Warning : The document was created with Spire.PDF for Python.
314
❒
ISSN:
2722-2586
Figure
4.
Comparison
of
joint
trajectories:
PID
vs
NMPC
The
error
time
series
in
Figure
5
mak
e
the
performance
dif
ference
quantitati
v
ely
clear
.
After
a
brief
initial
transient,
the
NMPC
error
settles
to
a
narro
w
,
near
-stationary
band,
whereas
the
PID
error
sho
ws
recur
-
ring
oscillatory
spik
es
whose
ampl
itude
and
frequenc
y
track
v
ariations
in
trajectory
curv
ature.
The
NMPC’
s
consistently
lo
w
error
is
direct
e
vidence
that
anticipatory
torque
planning
pre
v
ents
the
error
b
ursts
that
are
una
v
oidable
with
purely
reacti
v
e,
error
-dri
v
en
control.
Figure
5.
T
racking
error
comparison:
PID
vs
NMPC
IAES
Int
J
Rob
&
Autom,
V
ol.
15,
No.
2,
June
2026:
307–318
Evaluation Warning : The document was created with Spire.PDF for Python.
IAES
Int
J
Rob
&
Autom
ISSN:
2722-2586
❒
315
Figure
6
display
the
torque
w
a
v
eforms
deli
v
ered
to
Joint
1.
The
NMPC
torque
stays
smooth
and
remains
strictly
within
±
20
Nm
at
all
times
because
the
constraint
is
b
uilt
into
the
NLP;
the
optimizer
spreads
correcti
v
e
ef
fort
across
the
whole
horizon
rather
than
issuing
single
lar
ge
impulses.
The
PID
signal
is
lar
ger
at
trajectory
onset,
e
xhibits
abrupt
step
changes
at
curv
ature
inections,
and
requires
an
e
xternal
saturation
block
to
clip
limit
violations-a
reacti
v
e
strate
gy
that
is
fundamentally
incapable
of
pre
v
enting
cons
traint
breaches
before
the
y
occur
.
Figure
6.
Control
torque
for
joint
1:
NMPC
vs
PID
Figure
7
e
xtend
the
torque
comparison
to
Joint
2.
The
NMPC
optimizes
τ
1
and
τ
2
simultaneously
as
a
coupled
tw
o-dimensional
control
problem,
so
the
Joint
2
torque
is
naturally
coordinated
with
Joint
1
to
share
Coriolis
loads
between
the
tw
o
joints.
The
output
is
smooth
and
ener
gy-ef
cient.
W
ithout
access
to
Joint
1’
s
state,
the
PID
Joint
2
loop
produces
more
irre
gular
torque;
the
oscillation
frequenc
y
in
its
output
tracks
Joint
1
angular
speed
directly
,
conrming
that
unmodeled
Coriolis
forces
are
dri
ving
the
disturbance.
Figure
7.
Control
torque
for
joint
2:
NMPC
vs
PID
Design
and
implementation
of
NMPC
for
a
two-DOF
r
obotic
arm
using
CasADi
(Lahcen
Boulbalah)
Evaluation Warning : The document was created with Spire.PDF for Python.
316
❒
ISSN:
2722-2586
5.
CONCLUSION
This
w
ork
de
v
eloped
and
v
alidated
an
NMPC
frame
w
ork
for
joint-space
trajectory
tracking
of
a
2-DOF
serial
robotic
arm,
b
uilt
within
the
CasADi
symbolic
en
vironment
and
solv
ed
at
each
control
step
by
IPOPT’
s
interior
-point
algorithm.
The
prediction
model
embeds
the
complete
Lagrangian
dynamics-
conguration-dependent
inertia,
Coriolis
coupling,
and
gra
vitational
loading-discretized
via
RK4
inte
gration,
and
actuator
and
state
limits
are
imposed
as
hard
NLP
inequalities
rather
than
post-hoc
saturations.
Compared
with
a
Zie
gler
-Nichols-tuned
decentralized
PID
across
tw
o
sampling
rates,
the
NMPC
trimmed
peak
angular
error
by
95.69%
for
Joint
1
and
96.27%
for
Joint
2,
and
cut
mean
error
by
64.65%
and
57.58%
respecti
v
ely
.
A
v
erage
solv
e
times
of
0.053–0.056
s
conrmed
that
the
optimizer
consistently
ts
within
the
0.1
s
control
windo
w
on
a
standard
laptop,
establishing
real-time
suitability
without
needing
specialized
embedded
hardw
are.
Three
structural
features
of
the
NMPC
formulation-look-ahead
torque
planning,
inte
grated
coupling
compensation,
and
mathematically
guaranteed
constraint
adherence-collecti
v
ely
account
for
the
measured
performance
adv
antage
o
v
er
the
decentralized
PID.
Critically
,
none
of
these
attrib
utes
requires
problem-specic
engineering
ef
fort
be
yond
pro
viding
the
system
description
to
CasADi.
Planned
e
xtensions
include:
i)
Hardw
are-
in-the-loop
v
erication
on
a
ph
ysical
tw
o-link
a
rm
testbed
to
assess
disturbance
rejection
and
model-mismatch
rob
ustness
under
real
operating
conditions.
ii)
Scaling
to
higher
-DOF
manipulators
with
machine-learning
surrog
ate
dynamics
to
reduce
online
computational
demand.
iii)
Reinforcement-learning-based
automated
tuning
of
the
weighting
m
atrices
Q
and
R
to
remo
v
e
the
manual
calibration
step.
and
i
v)
Exploration
of
adapti
v
e
and
rob
ust
NMPC
formulations
capa
b
l
e
of
handling
time-v
arying
parameters
such
as
payload
mass
changes.
A
CKNO
WLEDGMENTS
The
authors
declare
that
no
specic
ackno
wledgments
are
required
for
this
w
ork.
FUNDING
INFORMA
TION
The
authors
state
no
funding
in
v
olv
ed.
A
UTHOR
CONTRIB
UTIONS
ST
A
TEMENT
This
journal
uses
the
Contrib
utor
Roles
T
axonomy
(CRediT)
to
recognize
indi
vidual
author
contrib
u-
tions,
reduce
authorship
disputes,
and
f
acilitate
collaboration.
Name
of
A
uthor
C
M
So
V
a
F
o
I
R
D
O
E
V
i
Su
P
Fu
Lahcen
Boulbalah
✓
✓
✓
✓
✓
✓
✓
✓
✓
F
aiza
Dib
✓
✓
✓
Nabil
Benaya
✓
✓
✓
Khaddouj
Ben
Meziane
✓
✓
✓
C
:
C
onceptualization
I
:
I
n
v
estig
ation
V
i
:
V
i
sualization
M
:
M
ethodology
R
:
R
esources
Su
:
Su
pervision
So
:
So
ftw
are
D
:
D
ata
Curation
P
:
P
roject
Administrati
on
V
a
:
V
a
lidation
O
:
Writing
-
O
riginal
Draft
Fu
:
Fu
nding
Acquisition
F
o
:
F
o
rmal
Analysis
E
:
Writing
-
Re
vie
w
&
E
diting
CONFLICT
OF
INTEREST
ST
A
TEMENT
Authors
state
no
conict
of
interest.
REFERENCES
[1]
K.
J.
Astr
om
and
R.
M.
Murray
,
Feedback
systems:
an
introduction
for
scientists
and
engineers,
2nd
ed.
Princeton
Uni
v
ersity
Press,
2021.
[2]
K.
J
.
Astr
om
and
T
.
H
agglund,
Adv
anced
PID
control.
ISA-The
Instrumentation,
Systems,
and
Automation
Society
,
2006.
IAES
Int
J
Rob
&
Autom,
V
ol.
15,
No.
2,
June
2026:
307–318
Evaluation Warning : The document was created with Spire.PDF for Python.