How to use the casadi.jacobian function in casadi
To help you get started, we’ve selected a few casadi examples, based on popular ways it is used in public projects.
Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
pymoca / pymoca / src / pymoca / backends / casadi / model.py View on Github 
elif value.is_one():
value = one
value = value if value.numel() != 1 else ca.repmat(value, *variable.symbol.size())
attribute_lists[attribute_list_index].append(value)
expr = ca.horzcat(*[ca.veccat(*attribute_list) for attribute_list in attribute_lists])
if len(self.parameters) > 0 and isinstance(expr, ca.MX):
f = ca.Function('f', [in_var], [expr])
# NOTE: This is not a complete list of operations that can be
# handled in an affine expression. That said, it should
# capture the most common ways variable attributes are
# expressed as a function of parameters.
allowed_ops = {ca.OP_INPUT, ca.OP_OUTPUT, ca.OP_CONST,
ca.OP_SUB, ca.OP_ADD, ca.OP_SUB, ca.OP_MUL, ca.OP_DIV, ca.OP_NEG}
f_ops = {f.instruction_id(k) for k in range(f.n_instructions())}
contains_unallowed_ops = not f_ops.issubset(allowed_ops)
zero_hessian = ca.jacobian(ca.jacobian(expr, in_var), in_var).is_zero()
if contains_unallowed_ops or not zero_hessian:
is_affine = False
out.append(expr)
if len(self.parameters) > 0 and is_affine:
# Rebuild variable metadata as a single affine expression, if all
# subexpressions are affine.
in_var_ = ca.MX.sym('in_var', in_var.shape)
out_ = []
for o in out:
Af = ca.Function('Af', [in_var], [ca.jacobian(o, in_var)])
bf = ca.Function('bf', [in_var], [o])
A = Af(0)
A = ca.sparsify(A)
b = bf(0)self.__jac_u = ca.Function('jac_x', [x, u, p],
[ca.jacobian(ode_casadi(x,u,p), u)])
# Jacobian of discrete RK4 system
self.__discrete_rk4_jac_x = ca.Function('jac_x', [x, u, p],
[ca.jacobian(self.rk4(x,u,p), x)])
self.__discrete_rk4_jac_u = ca.Function('jac_x', [x, u, p],
[ca.jacobian(self.rk4(x,u,p), u)])
# Jacobian of exact discretization
self.__discrete_jac_x = ca.Function('jac_x', [x, u, p],
[ca.jacobian(self.Integrator(x0=x,
p=ca.vertcat(u,p))['xf'], x)])
self.__discrete_jac_u = ca.Function('jac_u', [x, u, p],
[ca.jacobian(self.Integrator(x0=x,
p=ca.vertcat(u,p))['xf'], u)])elif expr.is_symbolic():
if expr.name() in states:
return der_states[expr.name()].symbol
elif expr.name() in alg_states:
# This algebraic state must now become a differentiated state.
states[expr.name()] = alg_states.pop(expr.name())
der_sym = ca.MX.sym('der({})'.format(expr.name()))
der_states[expr.name()] = Variable(der_sym, float)
return der_sym
else:
return 0.0
else:
# Differentiate using CasADi and chain rule
orig_deps = ca.symvar(expr)
deps = ca.vertcat(*orig_deps)
J = ca.Function('J', [deps], [ca.jacobian(expr, deps)])
J_sparsity = J.sparsity_out(0)
der_deps = [get_derivative(dep) if J_sparsity.has_nz(0, j) else ca.DM.zeros(dep.size()) for j, dep in enumerate(orig_deps)]
return ca.mtimes(J(deps), ca.vertcat(*der_deps))def T2WJ(T,p):
"""
w_101 = T2WJ(T_10,p).diff(p,t)
"""
R = T2R(T)
RT = R.T
temp = []
for i,k in [(2,1),(0,2),(1,0)]:
#temp.append(c.mul(c.jacobian(R[:,k],p).T,R[:,i]).T)
temp.append(c.mul(RT[i,:],c.jacobian(R[:,k],p)))
return c.vertcat(temp)if len(self.delay_arguments) > 0:
self.delay_arguments = self._substitute_delay_arguments(self.delay_arguments, variables, values)
if options['reduce_affine_expression']:
logger.info("Collapsing model into an affine expression")
for equation_list in ['equations', 'initial_equations']:
equations = getattr(self, equation_list)
if len(equations) > 0:
states = ca.veccat(*self._symbols(itertools.chain(self.states, self.der_states, self.alg_states, self.inputs)))
constants = ca.veccat(*self._symbols(self.constants))
parameters = ca.veccat(*self._symbols(self.parameters))
equations = ca.veccat(*equations)
Af = ca.Function('Af', [states, constants, parameters], [ca.jacobian(equations, states)])
bf = ca.Function('bf', [states, constants, parameters], [equations])
# Work around CasADi issue #172
if len(self.constants) == 0 or not ca.depends_on(equations, constants):
constants = 0
else:
logger.warning('Not all constants have been eliminated. As a result, the affine DAE expression will use a symbolic matrix, as opposed to a numerical sparse matrix.')
if len(self.parameters) == 0 or not ca.depends_on(equations, parameters):
parameters = 0
else:
logger.warning('Not all parameters have been eliminated. As a result, the affine DAE expression will use a symbolic matrix, as opposed to a numerical sparse matrix.')
A = Af(0, constants, parameters)
b = bf(0, constants, parameters)
# Replace veccat'ed states with brand new state vectors so as to avoid the value copy operations induced by veccat.self.Integrator = ca.integrator('Integrator', 'cvodes', dae, options)
#TODO: Fix discrete DAE model
if alg is None:
""" Create discrete RK4 model """
ode_casadi = ca.Function("ode", [x, u, p], [ode(x,u,z,p)])
k1 = ode_casadi(x, u, p)
k2 = ode_casadi(x + dt/2*k1, u, p)
k3 = ode_casadi(x + dt/2*k2, u, p)
k4 = ode_casadi(x + dt*k3,u, p)
xrk4 = x + dt/6*(k1 + 2*k2 + 2*k3 + k4)
self.rk4 = ca.Function("ode_rk4", [x, u, p], [xrk4])
# Jacobian of continuous system
self.__jac_x = ca.Function('jac_x', [x, u, p],
[ca.jacobian(ode_casadi(x,u,p), x)])
self.__jac_u = ca.Function('jac_x', [x, u, p],
[ca.jacobian(ode_casadi(x,u,p), u)])
# Jacobian of discrete RK4 system
self.__discrete_rk4_jac_x = ca.Function('jac_x', [x, u, p],
[ca.jacobian(self.rk4(x,u,p), x)])
self.__discrete_rk4_jac_u = ca.Function('jac_x', [x, u, p],
[ca.jacobian(self.rk4(x,u,p), u)])
# Jacobian of exact discretization
self.__discrete_jac_x = ca.Function('jac_x', [x, u, p],
[ca.jacobian(self.Integrator(x0=x,
p=ca.vertcat(u,p))['xf'], x)])
self.__discrete_jac_u = ca.Function('jac_u', [x, u, p],
[ca.jacobian(self.Integrator(x0=x,der_dep = _new_mx("der({})".format(dep.name()), dep.size())
der_dep._modelica_shape = \
self.nodes[self.current_class][dep.name()]._modelica_shape
self.derivative[dep.name()] = der_dep
self.nodes[self.current_class][der_dep.name()] = der_dep
return der_dep[slice_info['start']:slice_info['stop']:slice_info['step']]
else:
return self.derivative[dep.name()][slice_info['start']:slice_info['stop']:slice_info['step']]
# Case 4: s is an expression that requires differentiation, e.g. MX(x2 * x2)
# Need to do this sort of expansion: der(x1 * x2) = der(x1) * x2 + x1 * der(x2)
else:
# Differentiate expression using CasADi
orig_deps = ca.symvar(s)
deps = ca.vertcat(*orig_deps)
J = ca.Function('J', [deps], [ca.jacobian(s, deps)])
J_sparsity = J.sparsity_out(0)
der_deps = [self.get_derivative(dep) if J_sparsity.has_nz(0, j) else ca.DM.zeros(dep.size()) for j, dep in enumerate(orig_deps)]
return ca.mtimes(J(deps), ca.vertcat(*der_deps))V0 = ca.vertcat([
pl.ones(3), \
pl.zeros(N), \
ydata[0,:].T
])
sol = nlpsolver(x0 = V0)
p_est_single_shooting = sol["x"][:3]
tstart_Sigma_p = time()
J_s = ca.jacobian(r, V)
F_s = ca.mul([J_s.T, Sigma, J_s])
beta = (ca.mul([r.T, Sigma, r]) / (r.size() - V.size()))
Sigma_p_s = beta * ca.solve(F_s, ca.MX.eye(F_s.shape[0]), "csparse")
beta_fcn = ca.MXFunction("beta_fcn", [V], [beta])
print beta_fcn([sol["x"]])[0]
Sigma_p_s_fcn = ca.MXFunction("Sigma_p_s_fcn", \
[V] , [Sigma_p_s])
Cov_p = Sigma_p_s_fcn([sol["x"]])[0][:3, :3]
tend_Sigma_p = time()casadi
CasADi -- framework for algorithmic differentiation and numeric optimization
