00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #include "BonQuadCut.hpp"
00011
00012 namespace Bonmin {
00013
00014 QuadCut::QuadCut():
00015 OsiRowCut(),
00016 c_(0),
00017 Q_(),
00018 type_(Upper)
00019 {}
00020
00021 QuadCut::QuadCut(const QuadCut & other):
00022 OsiRowCut(other),
00023 c_(other.c_),
00024 Q_(other.Q_),
00025 type_(other.type_){
00026 }
00027
00028 QuadCut &
00029 QuadCut::operator=(const QuadCut & rhs){
00030 if(this != &rhs){
00031 OsiRowCut::operator=(rhs);
00032 c_ = rhs.c_;
00033 Q_ = rhs.Q_;
00034 type_ = rhs.type_;
00035 }
00036 return *this;
00037 }
00038
00039 OsiRowCut *
00040 QuadCut::clone() const {
00041 return new QuadCut(*this);
00042 }
00043
00044 QuadCut::~QuadCut(){
00045 }
00046
00048 double
00049 QuadCut::violated(const double * solution) const{
00050 double rhs = c_;
00051 rhs += row().dotProduct(solution);
00052 const int * indice = Q_.getIndices();
00053 const double * val = Q_.getElements();
00054 const int * start = Q_.getVectorStarts();
00055 const int * length = Q_.getVectorLengths();
00056 int n = Q_.getMajorDim();
00057
00058 for(int i = 0 ; i < n ; i++){
00059 int s=start[i];
00060 int l=length[i];
00061 for(int k = s ; k < s+l ; k++){
00062 if(i == indice[k]) rhs += solution[i] * solution[indice[k]] * val[k];
00063 else rhs += 2*solution[i] * solution[indice[k]] * val[k];
00064 }
00065 }
00066
00067 #if 0
00068 for(int i = 0 ; i < n ; i++){
00069 for(int k = 0 ; k < length[i] ; k++){
00070 if(i == indice[k]) rhs += solution[i] * solution[indice[k]] * val[k];
00071 else rhs += 2*solution[i] * solution[indice[k]] * val[k];
00072 }
00073 start += length[i];
00074 }
00075 #endif
00076 return std::max(lb() - rhs, rhs - ub());
00077 }
00078
00079
00080
00081 void
00082 QuadCut::print() const{
00083 std::cout<<"Quadratic cut has lower bound "<<lb()<<" and upper bound "<<ub()
00084 <<std::endl;
00085
00086 std::cout<<"Linear part has "<<row().getNumElements()<<" non zeroes:"
00087 <<std::endl;
00088
00089 const int& nElem = row().getNumElements();
00090 const int * indices = row().getIndices();
00091 const double * elements = row().getElements();
00092
00093 for(int i = 0 ; i < nElem ; i++){
00094 if(i > 0 && elements[i] > 0.)
00095 std::cout<<"+ ";
00096 std::cout<< elements[i] <<" x["<<indices[i]<<"]\t";
00097 if(i > 0 && i % 5 == 0) std::cout<<std::endl;
00098 }
00099 std::cout<<std::endl;
00100 if(Q_.getNumElements()){
00101 std::cout<<"Quadratic part is given by the matrix:"<<std::endl;
00102 Q_.dumpMatrix();
00103 }
00104 }
00105 Cuts::Cuts():
00106 OsiCuts(),
00107 quadCuts_(0){
00108 }
00109
00110 Cuts::Cuts(const Cuts & other):
00111 OsiCuts(other),
00112 quadCuts_(other.quadCuts_.size()){
00113 for(unsigned int i = 0 ; i < quadCuts_.size() ; i++){
00114 quadCuts_[i] = new QuadCut(*other.quadCuts_[i]);
00115 }
00116 }
00117
00118 Cuts &
00119 Cuts::operator=(const Cuts & rhs){
00120 if(this != &rhs){
00121 OsiCuts::operator=(rhs);
00122 for(unsigned int i = 0 ; i < quadCuts_.size() ; i++)
00123 {
00124 delete quadCuts_[i];
00125 }
00126 quadCuts_.resize(rhs.quadCuts_.size());
00127 for(unsigned int i = 0 ; i < quadCuts_.size() ; i++){
00128 quadCuts_[i] = new QuadCut(*rhs.quadCuts_[i]);
00129 }
00130 }
00131 return *this;
00132 }
00133
00134 Cuts::~Cuts(){
00135 for(unsigned int i = 0 ; i < quadCuts_.size() ; i++)
00136 {
00137 delete quadCuts_[i];
00138 }
00139 }
00140
00141 void
00142 Cuts::printCuts() const {
00143 OsiCuts::printCuts();
00144 std::cout<<quadCuts_.size()<<" quadratic cuts."<<std::endl;
00145 for(unsigned int i = 0 ; i < quadCuts_.size() ; i++){
00146 quadCuts_[i]->print();
00147 }
00148 }
00149
00150 }