You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
I'm writing here with the hope that someone could address my problem with the Jacobian computation using CasADi (maybe for you it is very easy to find a proper answer).
The point is that i'm implementing a neural state-space model (MLP with 1 layer) in an MPC loop and i want to compute the jacobian of the output of the model. I implemented the code in pyTorch and it works fine. Now, i'm trying to shift to CasADi, because the next step should be to compute everything in continous time. But, i'm still stucked in the computation of the Jacobian in descrete time. The point is that the values obtained are quite far from thoose computed with pyTorch (verified manually): I'm doing something wrong with the CasADi implementation? It is normal that there are these discrepancies between pytorch and CasADi?
MATLAB code:
clear all
addpath('yourCasadiFolderPath/casadi-3.6.5-windows64-matlab2018b')
import casadi.*
%Init model params
rng(1,"twister")
A=rand(100,5)
C=rand(3,100)
b=rand(100,1)
d=rand(3,1)
scal= [1.306500000000000,0.212830000000000,0.066820000000000,0.349070000000000,13];
scalOut=[0.039047000000000,0.042727000000000,0.013325000000000];
N=5
n_theta=903;%n. of parameters of the model
u_opt=rand(N,2);
ptheta=rand(n_theta,n_theta);
%Define symbolic variables
x = MX.sym('x',6);
% Define parameters
Ac = MX.sym('Ac', size(A,1), size(A,2));
Cc = MX.sym('Cc', size(C,1), size(C,2));
bc = MX.sym('bc', size(b,1));
dc = MX.sym('dc', size(d,1));
SCALIN=MX.sym('scalin', 5,1)
SCAL= MX.sym('scal', 6,1)
SCALOUT= MX.sym('scalout', 3,1)
DT= MX.sym('Dt', 1)
% Define control input sequence
U = MX.sym('U', [2,N]);
% Initialize state and Jacobian storage
x_curr = x;
X_sim = MX.zeros(6, N); % To store simulated states
J_Ac_list = []; % To store Jacobians for Ac
J_Cc_list = []; % To store Jacobians for Cc
J_bc_list = []; % To store Jacobians for bc
J_dc_list = []; % To store Jacobians for dc
% Loop over each time step
for k = 1:N
% Define the current input
U_k = [U(1,k) U(2,k)];
%EULER
input_scaled = [x_curr(1)/SCALIN(1); x_curr(2)/SCALIN(2); x_curr(3)/SCALIN(3); U_k(1)/SCALIN(4); U_k(2)/SCALIN(5)];
xdot = [(Cc(1,:) * tanh(Ac * input_scaled + bc) + dc(1)) .* SCALOUT(1);
(Cc(2,:) * tanh(Ac * input_scaled + bc) + dc(2)) .* SCALOUT(2);
(Cc(3,:) * tanh(Ac * input_scaled + bc) + dc(3)) .* SCALOUT(3);
x(3);
x(1) * cos(x(4)) - x(2) * sin(x(4));
x(1) * sin(x(4)) + x(2) * cos(x(4))];
x_next = (x_curr + DT *(xdot))./SCAL
% Store the state
X_sim(:, k) = x_next;
% Compute Jacobian of next state with respect to parameters
J_Ac = jacobian(x_next, Ac);
J_Cc = jacobian(x_next, Cc);
J_bc = jacobian(x_next, bc);
J_dc = jacobian(x_next, dc);
% Store the Jacobians
J_Ac_list = [J_Ac_list full(J_Ac)];
J_Cc_list = [J_Cc_list full(J_Cc)];
J_bc_list = [J_bc_list full(J_bc)];
J_dc_list = [J_dc_list full(J_dc)];
% Update current state
x_curr = x_next.*SCAL;
end
% Create functions for simulation and Jacobian computation
X_sim_func = Function('X_sim_func', {x, U, Ac, Cc, bc, dc,SCAL,SCALOUT,DT,SCALIN}, {X_sim});
J_Ac_func = Function('J_Ac_func', {x, U, Ac, Cc, bc, dc,SCAL,SCALOUT,DT,SCALIN}, {J_Ac_list});
J_Cc_func = Function('J_Cc_func', {x, U, Ac, Cc, bc, dc,SCAL,SCALOUT,DT,SCALIN}, {J_Cc_list});
J_bc_func = Function('J_bc_func', {x, U, Ac, Cc, bc, dc,SCAL,SCALOUT,DT,SCALIN}, {J_bc_list});
J_dc_func = Function('J_dc_func', {x, U, Ac, Cc, bc, dc,SCAL,SCALOUT,DT,SCALIN}, {J_dc_list});
% Init initial state and inputs
% x0_val = [vu_opt(1); vv_opt(1); vr_opt(1); rho_opt(1); X_opt(1); Y_opt(1)];
x0_val = [0.7 0 0 0 0 0];
u_val = u_opt';
A_val = A;
C_val = C;
b_val = b;
d_val = d;
SCAL_val=[scal(1:3)'; 1; 1; 1];
SCALOUT_val=[scalOut']
DT_val=2;
SCALIN_val=scal';
% Compute the state trajectory
X_traj = X_sim_func(x0_val, u_val, A_val, C_val, b_val, d_val,SCAL_val,SCALOUT_val,DT_val,SCALIN_val);
% Compute the Jacobians at each time step
J_Ac_traj = J_Ac_func(x0_val, u_val, A_val, C_val, b_val, d_val,SCAL_val,SCALOUT_val,DT_val,SCALIN_val);
J_Cc_traj = J_Cc_func(x0_val, u_val, A_val, C_val, b_val, d_val,SCAL_val,SCALOUT_val,DT_val,SCALIN_val);
J_bc_traj = J_bc_func(x0_val, u_val, A_val, C_val, b_val, d_val,SCAL_val,SCALOUT_val,DT_val,SCALIN_val);
J_dc_traj = J_dc_func(x0_val, u_val, A_val, C_val, b_val, d_val,SCAL_val,SCALOUT_val,DT_val,SCALIN_val);
JVu = zeros(N, n_theta);
JVv = zeros(N, n_theta);
JVr = zeros(N, n_theta);
JX = zeros(N, n_theta);
JY = zeros(N, n_theta);
for i = 1:N
J_Ac_val = full(J_Ac_traj);
J_Cc_val = full(J_Cc_traj);
J_bc_val = full(J_bc_traj);
J_dc_val = full(J_dc_traj);
% Combine the Jacobians for each parameter
J_combined = [J_Ac_val(:,(i-1)*numel(A)+1:(i-1)*numel(A)+numel(A)), J_bc_val(:,(i-1)*numel(b)+1:(i-1)*numel(b)+numel(b)), J_Cc_val(:,(i-1)*numel(C)+1:(i-1)*numel(C)+numel(C)), J_dc_val(:,(i-1)*numel(d)+1:(i-1)*numel(d)+numel(d))];
% J_combined = [JAc((i-1)*6+1:(i-1)*6+6,:), JCc((i-1)*6+1:(i-1)*6+6,:),Jbc((i-1)*6+1:(i-1)*6+6,:),Jdc((i-1)*6+1:(i-1)*6+6,:)];
% Assign the Jacobians for each state variable
JVu(i,:) = J_combined(1,:);
JVv(i,:) = J_combined(2,:);
JVr(i,:) = J_combined(3,:);
JX(i,:) = J_combined(5,:);
JY(i,:) = J_combined(6,:);
end
reacted with thumbs up emoji reacted with thumbs down emoji reacted with laugh emoji reacted with hooray emoji reacted with confused emoji reacted with heart emoji reacted with rocket emoji reacted with eyes emoji
-
Dear all,
I'm writing here with the hope that someone could address my problem with the Jacobian computation using CasADi (maybe for you it is very easy to find a proper answer).
The point is that i'm implementing a neural state-space model (MLP with 1 layer) in an MPC loop and i want to compute the jacobian of the output of the model. I implemented the code in pyTorch and it works fine. Now, i'm trying to shift to CasADi, because the next step should be to compute everything in continous time. But, i'm still stucked in the computation of the Jacobian in descrete time. The point is that the values obtained are quite far from thoose computed with pyTorch (verified manually): I'm doing something wrong with the CasADi implementation? It is normal that there are these discrepancies between pytorch and CasADi?
MATLAB code:
Thanks.
J
Beta Was this translation helpful? Give feedback.
All reactions