How to use casadi - 10 common examples
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.
adbuerger / casiopeia / test / test_casadi_interface.py View on Github 
def setUp(self):
self.f = ca.MX.sym("f")
self.g = ca.MX.sym("g")# - m: representing the ball of the mass in kg
# - L: the length of the pendulum bar in meters
# - g: the gravity constant in m/s^2
# - psi: the actuation angle of the manuver in radians, which stays
# constant for this problem
m = 1.0
L = 3.0
g = 9.81
# psi = pl.pi / 2.0
psi = pl.pi / (180.0 * 2)
# System
x = ca.MX.sym("x", 2)
p = ca.MX.sym("p", 1)
u = ca.MX.sym("u", 1)
# f = ca.vertcat([x[1], p[0]/(m*(L**2))*(u-x[0]) - g/L * pl.sin(x[0])])
f = ca.vertcat(x[1], p[0]/(m*(L**2))*(u-x[0]) - g/L * x[0])
phi = x
system = cp.system.System(x = x, u = u, p = p, f = f, phi = phi)
data = pl.loadtxt('data_pendulum.txt')
time_points = data[:500, 0]
numeas = data[:500, 1]
wmeas = data[:500, 2]
N = time_points.size
ydata = pl.array([numeas,wmeas])
udata = [psi] * (N-1)#
# - m: representing the ball of the mass in kg
# - L: the length of the pendulum bar in meters
# - g: the gravity constant in m/s^2
# - psi: the actuation angle of the manuver in radians, which stays
# constant for this problem
m = 1.0
L = 3.0
g = 9.81
# psi = pl.pi / 2.0
psi = pl.pi / (180.0 * 2)
# System
x = ca.MX.sym("x", 2)
p = ca.MX.sym("p", 1)
u = ca.MX.sym("u", 1)
# f = ca.vertcat([x[1], p[0]/(m*(L**2))*(u-x[0]) - g/L * pl.sin(x[0])])
f = ca.vertcat(x[1], p[0]/(m*(L**2))*(u-x[0]) - g/L * x[0])
phi = x
system = cp.system.System(x = x, u = u, p = p, f = f, phi = phi)
data = pl.loadtxt('data_pendulum.txt')
time_points = data[:500, 0]
numeas = data[:500, 1]
wmeas = data[:500, 2]
N = time_points.size
ydata = pl.array([numeas,wmeas])def setUp(self):
self.x = ca.MX.sym("x", 4)
self.p = ca.MX.sym("p", 6)
self.u = ca.MX.sym("u", 2)
self.f = ca.vertcat( \
self.x[3] * np.cos(self.x[2] + self.p[0] * self.u[0]),
self.x[3] * np.sin(self.x[2] + self.p[0] * self.u[0]),
self.x[3] * self.u[0] * self.p[1],
self.p[2] * self.u[1] \
- self.p[3] * self.u[1] * self.x[3] \
- self.p[4] * self.x[3]**2 \
- self.p[5] \
- (self.x[3] * self.u[0])**2 * self.p[1] * self.p[0])
self.phi = self.x
data = np.array(np.loadtxt("test/data_2d_vehicle_pe.dat", \report.write("\n\n repetitions = " + str(repetitions))
report.write("\n sigma = " + str(sigma))
report.write("\n\n p_true = " + str(ca.DM(ptrue)))
report.write("\n\n p_mean = " + str(ca.DM(p_mean)))
report.write("\n phat_last_exp = " + \
str(ca.DM(pe_test.estimated_parameters)))
report.write("\n\n p_sd = " + str(ca.DM(p_std)))
report.write("\n sd_from_covmat = " \
+ str(ca.diag(ca.sqrt(pe_test.covariance_matrix))))
report.write("\n beta = " + str(pe_test.beta))
report.write("\n\n delta_abs_sd = " + str(ca.fabs(ca.DM(p_std) - \
ca.diag(ca.sqrt(pe_test.covariance_matrix)))))
report.write("\n delta_rel_sd = " + str(ca.fabs(ca.DM(p_std) - \
ca.diag(ca.sqrt(pe_test.covariance_matrix))) / ca.DM(p_std)) \
+ "\n")
report.close()
try:
os.system("rst2pdf " + fname)
except:
print("Generating PDF report failed, is rst2pdf installed correctly?")(in1,v1,x[[1,0],0],sparsify(DM([[0,1],[1,0]]))),
(in1,v1,w,sparsify(DM([[1,0],[0,2]]))),
(in1,v1,w2,blockcat([[1,MX(1,1)],[x[1],x[0]]])),
(in1,v1,ww,2*c.diag(x)),
(in1,v1,wwf,vertcat(*[x[[1,0]].T,x[[1,0]].T])),
(in1,v1,yy[:,0],DM.eye(2)),
(in1,v1,yy2[:,0],2*c.diag(x)),
(in1,v1,yyy[:,0],sparsify(DM([[0,1],[1,0]]))),
(in1,v1,mtimes(y,x),y),
(in1,v1,mtimes(x.T,y.T),y),
(in1,v1,mac(y,x,DM.zeros(Sparsity.triplet(2,1,[1],[0]))),y[Sparsity.triplet(2,2,[1,1],[0,1])]),
(in1,v1,mac(x.T,y.T,DM.zeros(Sparsity.triplet(2,1,[1],[0]).T)),y[Sparsity.triplet(2,2,[1,1],[0,1])]),
(in1,v1,mtimes(y[Sparsity.triplet(2,2,[0,1,1],[0,0,1])],x),y[Sparsity.triplet(2,2,[0,1,1],[0,0,1])]),
(in1,v1,mtimes(x.T,y[Sparsity.triplet(2,2,[0,1,1],[0,0,1])].T),y[Sparsity.triplet(2,2,[0,1,1],[0,0,1])]),
(in1,v1,mtimes(y,x**2),y*2*vertcat(*[x.T,x.T])),
(in1,v1,sin(x),c.diag(cos(x))),
(in1,v1,sin(x**2),c.diag(cos(x**2)*2*x)),
(in1,v1,x*y[:,0],c.diag(y[:,0])),
(in1,v1,x*y.nz[[0,1]],c.diag(y.nz[[0,1]])),
(in1,v1,x*y.nz[[1,0]],c.diag(y.nz[[1,0]])),
(in1,v1,x*y[[0,1],0],c.diag(y[[0,1],0])),
(in1,v1,x*y[[1,0],0],c.diag(y[[1,0],0])),
(in1,v1,c.dot(x,x),(2*x).T),
(in1,v1,c.dot(x**2,x),(3*x**2).T),
#(in1,v1,c.det(horzcat(*[x,DM([1,2])])),DM([-1,2])), not implemented
(in1,v1,f1.call(in1)[1],y),
(in1,v1,f1.call([x**2,y])[1],y*2*vertcat(*[x.T,x.T])),
(in1,v1,f2.call(in1)[0],DM.zeros(0,2)),
(in1,v1,f2(x**2,y),DM.zeros(0,2)),
(in1,v1,f3.call(in1)[0],DM.zeros(0,2)),
(in1,v1,f3.call([x**2,y])[0],DM.zeros(0,2)),
(in1,v1,f4.call(in1)[0],DM.zeros(0,2)),pe_test.print_estimation_results()
# Generate report
print("\np_mean = " + str(ca.DM(p_mean)))
print("phat_last_exp = " + str(ca.DM(pe_test.estimated_parameters)))
print("\np_sd = " + str(ca.DM(p_std)))
print("sd_from_covmat = " + str(ca.diag(ca.sqrt(pe_test.covariance_matrix))))
print("beta = " + str(pe_test.beta))
print("\ndelta_abs_sd = " + str(ca.fabs(ca.DM(p_std) - \
ca.diag(ca.sqrt(pe_test.covariance_matrix)))))
print("delta_rel_sd = " + str(ca.fabs(ca.DM(p_std) - \
ca.diag(ca.sqrt(pe_test.covariance_matrix))) / ca.DM(p_std)))
fname = os.path.basename(__file__)[:-3] + ".rst"
report = open(fname, "w")
report.write( \
'''Concept test: covariance matrix computation
===========================================
Simulate system. Then: add gaussian noise N~(0, sigma^2), estimate,
store estimated parameter, repeat.
.. code-block:: python
y_randn = sim_true.simulation_results + sigma * \
(np.random.randn(*sim_true.estimated_parameters.shape))def test_general_convex_sparse(self):
self.message("Convex sparse QP with solvers: " + str([conic for conic,options,aux_options in conics]))
H = c.diag([2,1,0.2,0.7,1.3])
H[1,2]=0.1
H[2,1]=0.1
G = DM([-2,-6,1,0,0])
A = DM([[1, 0,0.1,0.7,-1],[0.1, 2,-0.3,4,0.1]])
A = sparsify(A)
LBA = DM([-inf])
UBA = DM([2, 2])
LBX = DM([0]*5)
UBX = DM([inf]*5)
for conic, qp_options, aux_options in conics:def test_diag_sparse(self):
self.message("diag sparse")
for n in [[0,1,0,0,2,3,4,5,6,0],[1,2,3,0],[0,1,2,3]]:
d = DM(n)
D = DM(n)
d = sparsify(d)
m = c.diag(d)
M = sparsify(c.diag(D))
self.checkarray(m.sparsity().colind(),M.sparsity().colind())
self.checkarray(m.sparsity().row(),M.sparsity().row())self.assertEqual(len(A.parameters), len(B.parameters))
if not isinstance(A, CachedModel) and not isinstance(B, CachedModel):
self.assertEqual(len(A.equations), len(B.equations))
self.assertEqual(len(A.initial_equations), len(B.initial_equations))
for f_name in ['dae_residual', 'initial_residual', 'variable_metadata']:
this = getattr(A, f_name + '_function')
that = getattr(B, f_name + '_function')
np.random.seed(0)
args_in = []
for i in range(this.n_in()):
sp = this.sparsity_in(0)
r = ca.DM(sp, np.random.random(sp.nnz()))
args_in.append(r)
this_out = this.call(args_in)
that_out = that.call(args_in)
# N.B. Here we require that the order of the equations in the two models is identical.
for i, (a, b) in enumerate(zip(this_out, that_out)):
for j in range(a.size1()):
for k in range(a.size2()):
if a[j, k].is_regular() or b[j, k].is_regular():
test = float(ca.norm_2(ca.vec(a[j, k] - b[j, k]))) <= tol
if not test:
print(j)
print(k)
print(a[j,k])
print(b[j,k])casadi
CasADi -- framework for algorithmic differentiation and numeric optimization
