1function [dx rec_out events] = equations_of_motion(t, x, t_all, x_all, rec, N, dt, simcons) 
 2%#codegen 
 3q = x(1:9); 
 4dq = zeros(9,1); 
 5dq = x(10:18); 
 6s = x(19:52); 
 7
 8[Tj Tm Tp dTj dTm dTp ddTm] = Tmatrices(q, dq, simcons.p_in); 
 9
10F = zeros(9, 1); 
11ds = zeros(34, 1); 
12rec_out = zeros(100, 1); 
13
14events = zeros(1, 1); 
15[Fc, ~, rec_out([1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36])] = contactforce_control(t, N, dt, q, dq, t_all, x_all, simcons.control.c1, s, rec, Tj, Tm, Tp, dTj, dTm, dTp, simcons.M); 
16F = F + Fc; 
17[Fc, ~, ~] = gravity_control(t, N, dt, q, dq, t_all, x_all, simcons.control.c2, s, rec, Tj, Tm, Tp, dTj, dTm, dTp, simcons.M); 
18F = F + Fc; 
19[Fc, ~, rec_out([37 38 39 40 41 42])] = ligaments_control(t, N, dt, q, dq, t_all, x_all, simcons.control.c3, s, rec, Tj, Tm, Tp, dTj, dTm, dTp, simcons.M); 
20F = F + Fc; 
21[Fc, ds([1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30]), rec_out([43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92])] = reflex_control(t, N, dt, q, dq, t_all, x_all, simcons.control.c4, s, rec, Tj, Tm, Tp, dTj, dTm, dTp, simcons.M); 
22F = F + Fc; 
23[Fc, ds([31 32 33 34]), rec_out([93 94])] = phase_detect_control(t, N, dt, q, dq, t_all, x_all, simcons.control.c5, s, rec, Tj, Tm, Tp, dTj, dTm, dTp, simcons.M); 
24F = F + Fc; 
25[Fc, ~, rec_out([95 96 97 98 99 100])] = ff_gaitpattern_control(t, N, dt, q, dq, t_all, x_all, simcons.control.c6, s, rec, Tj, Tm, Tp, dTj, dTm, dTp, simcons.M); 
26F = F + Fc; 
27
28ddq = (dTm' * simcons.M * dTm) \ (F - dTm' * simcons.M * (ddTm * dq)); 
29dx = [dq; ddq; ds]; 
30
31events(1) = walker_falling_event(t, N, dt, q, dq, t_all, x_all, s, rec, Tj, Tm, Tp, dTj, dTm, dTp, simcons.M);