00001
00002 #include <iostream>
00003 #include <Eigen/Dense>
00004
00005
00006 #include <ros/common.h>
00007 #if ROS_VERSION_GE(ROS_VERSION_MAJOR, ROS_VERSION_MINOR, ROS_VERSION_PATCH, 1, 11, 10)
00008 namespace Eigen
00009 {
00010 namespace internal
00011 {
00012 float abs(float a) { return std::abs(a); }
00013 float sqrt(float v) { return ::sqrt(v); }
00014 }
00015 }
00016 #endif
00017
00018
00019
00020 #include "eiquadprog.hpp"
00021
00022 using namespace Eigen;
00023
00024 #define print(var) \
00025 std::cout<<#var"= "<<std::endl<<var<<std::endl
00026
00027 int flag = -1 ;
00028
00029 extern "C" {
00030 int get_constraints_check_flag(){
00031 return flag ;
00032 }
00033 }
00034
00035 extern "C" {
00036 int check_constraints(double* CE, double* ce0, double* CI, double* ci0,
00037 double* x, int x_len, int ce_len, int ci_len, double eqthre, double* ce_err,
00038 double* ci_err) {
00039 int ret = ce_len + ci_len ;
00040 for (int i = 0; i < ce_len; i++) {
00041 for (int j = 0; j < x_len; j++) {
00042 ce_err[i] += CE[i * x_len + j] * x[j];
00043 }
00044 ce_err[i] += ce0[i];
00045 if ( ce_err[i] < eqthre && ce_err[i] > -eqthre ) ret-- ;
00046 }
00047 for (int i = 0; i < ci_len; i++) {
00048 for (int j = 0; j < x_len; j++) {
00049 ci_err[i] += CI[i * x_len + j] * x[j];
00050 }
00051 ci_err[i] += ci0[i];
00052 if ( ci_err[i] > -eqthre ) ret-- ;
00053 }
00054 return ret;
00055 }
00056 }
00057
00058
00059 extern "C" {
00060 double* solve_eiquadprog(double* G, double* g0, double* CE, double* ce0, double* CI,
00061 double* ci0, double* x,
00062 int x_len, int ce_len, int ci_len,
00063 double eqthre,
00064 int debug,
00065 double* ret_buf, double* ce_err, double* ci_err) {
00066 MatrixXd G_buf(x_len, x_len);
00067 Eigen::VectorXd g0_buf(x_len);
00068 MatrixXd CE_buf(x_len, ce_len);
00069 Eigen::VectorXd ce0_buf(ce_len);
00070 MatrixXd CI_buf(x_len, ci_len);
00071 Eigen::VectorXd ci0_buf(ci_len);
00072 Eigen::VectorXd x_buf(x_len);
00073
00074 for ( int i=0 ; i<x_len ; i++ ){
00075 for ( int j=0 ; j<x_len ; j++ ){
00076 G_buf(i,j) = G[i+ x_len*j] ;
00077
00078
00079
00080 if ( i == 0 ){
00081 g0_buf(j) = g0[j] ;
00082 x_buf(j) = x[j] ;
00083 }
00084 }
00085 for ( int j=0 ; j<ce_len ; j++ ){
00086 CE_buf(i,j) = CE[i+ x_len*j] ;
00087 if ( i == 0 ){
00088 ce0_buf(j) = ce0[j] ;
00089 }
00090 }
00091 for ( int j=0 ; j<ci_len ; j++ ){
00092 CI_buf(i,j) = CI[i+ x_len*j] ;
00093 if ( i == 0 ){
00094 ci0_buf(j) = ci0[j] ;
00095 }
00096 }
00097 }
00098
00099
00100 if (debug>1) {
00101 print(G_buf);
00102 print(g0_buf);
00103 print(CE_buf);
00104 print(ce0_buf);
00105 print(CI_buf);
00106 print(ci0_buf);
00107 }
00108
00109 ret_buf[0] = solve_quadprog(G_buf, g0_buf, CE_buf, ce0_buf, CI_buf,
00110 ci0_buf, x_buf) ;
00111 for (int i = 0; i < x_buf.size(); i++)
00112 x[i] = x_buf(i) ;
00113
00114 flag = check_constraints(CE, ce0, CI, ci0, x, x_len, ce_len, ci_len, eqthre, ce_err, ci_err);
00115
00116 if (debug>0) {
00117 std::cout << "[eus-eiquadprog]" << std::endl;
00118 std::cout << " :minimized-object "
00119 << ret_buf[0] << std::endl;
00120 std::cout << " :optimal-state [";
00121 for (int i = 0; i < x_buf.size(); i++)
00122 std::cout << x[i] << ' ';
00123 std::cout << "]" << std::endl;
00124
00125 std::cout << " :eq-constraint || " ;
00126 for ( int i=0 ; i <ce_len ; i++ ) std::cout << ce_err[i] << " " ;
00127 std::cout << "|| < " << eqthre ;
00128 std::cout << std::endl ;
00129 std::cout << " :iq-constraint [" ;
00130 for ( int i=0 ; i <ci_len ; i++ ) std::cout << ci_err[i] << " " ;
00131 std::cout << "] > " << -eqthre ;
00132 std::cout << std::endl ;
00133 std::cout << " :constraint-check " ;
00134 std::cout << flag ;
00135 std::cout << std::endl ;
00136 }
00137
00138 return x;
00139 }
00140 }