binary_library.cpp
Go to the documentation of this file.
00001 // This file is part of Eigen, a lightweight C++ template library
00002 // for linear algebra.
00003 //
00004 // Copyright (C) 2009 Benoit Jacob <jacob.benoit.1@gmail.com>
00005 //
00006 // Eigen is free software; you can redistribute it and/or
00007 // modify it under the terms of the GNU Lesser General Public
00008 // License as published by the Free Software Foundation; either
00009 // version 3 of the License, or (at your option) any later version.
00010 //
00011 // Alternatively, you can redistribute it and/or
00012 // modify it under the terms of the GNU General Public License as
00013 // published by the Free Software Foundation; either version 2 of
00014 // the License, or (at your option) any later version.
00015 //
00016 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
00017 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
00018 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
00019 // GNU General Public License for more details.
00020 //
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License and a copy of the GNU General Public License along with
00023 // Eigen. If not, see <http://www.gnu.org/licenses/>.
00024 
00025 // This C++ file compiles to binary code that can be linked to by your C program,
00026 // thanks to the extern "C" syntax used in the declarations in binary_library.h.
00027 
00028 #include "binary_library.h"
00029 
00030 #include <Eigen/Core>
00031 
00032 using namespace Eigen;
00033 
00034 /************************* pointer conversion methods **********************************************/
00035 
00037 
00038 inline MatrixXd& c_to_eigen(C_MatrixXd* ptr)
00039 {
00040   return *reinterpret_cast<MatrixXd*>(ptr);
00041 }
00042 
00043 inline const MatrixXd& c_to_eigen(const C_MatrixXd* ptr)
00044 {
00045   return *reinterpret_cast<const MatrixXd*>(ptr);
00046 }
00047 
00048 inline C_MatrixXd* eigen_to_c(MatrixXd& ref)
00049 {
00050   return reinterpret_cast<C_MatrixXd*>(&ref);
00051 }
00052 
00053 inline const C_MatrixXd* eigen_to_c(const MatrixXd& ref)
00054 {
00055   return reinterpret_cast<const C_MatrixXd*>(&ref);
00056 }
00057 
00059 
00060 inline Map<MatrixXd>& c_to_eigen(C_Map_MatrixXd* ptr)
00061 {
00062   return *reinterpret_cast<Map<MatrixXd>*>(ptr);
00063 }
00064 
00065 inline const Map<MatrixXd>& c_to_eigen(const C_Map_MatrixXd* ptr)
00066 {
00067   return *reinterpret_cast<const Map<MatrixXd>*>(ptr);
00068 }
00069 
00070 inline C_Map_MatrixXd* eigen_to_c(Map<MatrixXd>& ref)
00071 {
00072   return reinterpret_cast<C_Map_MatrixXd*>(&ref);
00073 }
00074 
00075 inline const C_Map_MatrixXd* eigen_to_c(const Map<MatrixXd>& ref)
00076 {
00077   return reinterpret_cast<const C_Map_MatrixXd*>(&ref);
00078 }
00079 
00080 
00081 /************************* implementation of classes **********************************************/
00082 
00083 
00085 
00086 
00087 C_MatrixXd* MatrixXd_new(int rows, int cols)
00088 {
00089   return eigen_to_c(*new MatrixXd(rows,cols));
00090 }
00091 
00092 void MatrixXd_delete(C_MatrixXd *m)
00093 {
00094   delete &c_to_eigen(m);
00095 }
00096 
00097 double* MatrixXd_data(C_MatrixXd *m)
00098 {
00099   return c_to_eigen(m).data();
00100 }
00101 
00102 void MatrixXd_set_zero(C_MatrixXd *m)
00103 {
00104   c_to_eigen(m).setZero();
00105 }
00106 
00107 void MatrixXd_resize(C_MatrixXd *m, int rows, int cols)
00108 {
00109   c_to_eigen(m).resize(rows,cols);
00110 }
00111 
00112 void MatrixXd_copy(C_MatrixXd *dst, const C_MatrixXd *src)
00113 {
00114   c_to_eigen(dst) = c_to_eigen(src);
00115 }
00116 
00117 void MatrixXd_copy_map(C_MatrixXd *dst, const C_Map_MatrixXd *src)
00118 {
00119   c_to_eigen(dst) = c_to_eigen(src);
00120 }
00121 
00122 void MatrixXd_set_coeff(C_MatrixXd *m, int i, int j, double coeff)
00123 {
00124   c_to_eigen(m)(i,j) = coeff;
00125 }
00126 
00127 double MatrixXd_get_coeff(const C_MatrixXd *m, int i, int j)
00128 {
00129   return c_to_eigen(m)(i,j);
00130 }
00131 
00132 void MatrixXd_print(const C_MatrixXd *m)
00133 {
00134   std::cout << c_to_eigen(m) << std::endl;
00135 }
00136 
00137 void MatrixXd_multiply(const C_MatrixXd *m1, const C_MatrixXd *m2, C_MatrixXd *result)
00138 {
00139   c_to_eigen(result) = c_to_eigen(m1) * c_to_eigen(m2);
00140 }
00141 
00142 void MatrixXd_add(const C_MatrixXd *m1, const C_MatrixXd *m2, C_MatrixXd *result)
00143 {
00144   c_to_eigen(result) = c_to_eigen(m1) + c_to_eigen(m2);
00145 }
00146 
00147 
00148 
00150 
00151 
00152 C_Map_MatrixXd* Map_MatrixXd_new(double *array, int rows, int cols)
00153 {
00154   return eigen_to_c(*new Map<MatrixXd>(array,rows,cols));
00155 }
00156 
00157 void Map_MatrixXd_delete(C_Map_MatrixXd *m)
00158 {
00159   delete &c_to_eigen(m);
00160 }
00161 
00162 void Map_MatrixXd_set_zero(C_Map_MatrixXd *m)
00163 {
00164   c_to_eigen(m).setZero();
00165 }
00166 
00167 void Map_MatrixXd_copy(C_Map_MatrixXd *dst, const C_Map_MatrixXd *src)
00168 {
00169   c_to_eigen(dst) = c_to_eigen(src);
00170 }
00171 
00172 void Map_MatrixXd_copy_matrix(C_Map_MatrixXd *dst, const C_MatrixXd *src)
00173 {
00174   c_to_eigen(dst) = c_to_eigen(src);
00175 }
00176 
00177 void Map_MatrixXd_set_coeff(C_Map_MatrixXd *m, int i, int j, double coeff)
00178 {
00179   c_to_eigen(m)(i,j) = coeff;
00180 }
00181 
00182 double Map_MatrixXd_get_coeff(const C_Map_MatrixXd *m, int i, int j)
00183 {
00184   return c_to_eigen(m)(i,j);
00185 }
00186 
00187 void Map_MatrixXd_print(const C_Map_MatrixXd *m)
00188 {
00189   std::cout << c_to_eigen(m) << std::endl;
00190 }
00191 
00192 void Map_MatrixXd_multiply(const C_Map_MatrixXd *m1, const C_Map_MatrixXd *m2, C_Map_MatrixXd *result)
00193 {
00194   c_to_eigen(result) = c_to_eigen(m1) * c_to_eigen(m2);
00195 }
00196 
00197 void Map_MatrixXd_add(const C_Map_MatrixXd *m1, const C_Map_MatrixXd *m2, C_Map_MatrixXd *result)
00198 {
00199   c_to_eigen(result) = c_to_eigen(m1) + c_to_eigen(m2);
00200 }


libicr
Author(s): Robert Krug
autogenerated on Mon Jan 6 2014 11:32:30