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 conguration-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, conguration-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 benets 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 conguration-dependent nonlinearities [19], [20]. The specic 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 identies 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 conguration 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 conguration-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 conguration-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-denite 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 renement, 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 Dening 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 conguration-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 conguration 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 + 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 rened 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 conguration 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 eried 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 inections, 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 , conrming 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- conguration-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 conrmed 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-specic engineering ef fort be yond pro viding the system description to CasADi. Planned e xtensions include: i) Hardw are- in-the-loop v erication 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 specic 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 conict 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.