Code:
#include<iostream>
using namespace std;
double M[3][3] = { {-2.5, 2, -1},
{1,3,-1},
{1,-1,1}
};
double L;
double H;
int C = 1;
double tol = 1e-5;
int max_passes = 10;
double alpha[3] = {0,0,0};
double b = 0;
int passes = 0;
int num_changed_alphas;
double E[3];
double sum;
int k;
double oldalpha[3];
double eta;
double Max(double N1, double N2)
{
double Maximum = N1;
if( Maximum < N2 ){
Maximum = N2;}
return Maximum;
}
double Min(double N1, double N2)
{
double Minimum = N1;
if( Minimum > N2 ){
Minimum = N2;}
return Minimum;
}
double innerproduct(int index1, int index2){
double out = M[index1][1]*M[index2][1] + M[index1][2]*M[index2][2];
return out;
}
int main(){
while (passes < max_passes){
num_changed_alphas = 0;
for (int i=1; i <= 3; i++){
for (int j=1; j <= 3; j++){
sum += alpha[j]*M[j][3]*innerproduct(j,i);}
E[i] = sum + b - M[i][3];
sum = 0;
if ((M[i][3]*E[i]<-tol && alpha[i]<C) || (M[i][3]*E[i]>tol && alpha[i]>0)){
if (i==3) {k = 1;}
else {k = i++;}
}
for (int m=1; m <= 3; m++){
sum += alpha[m]*M[m][3]*innerproduct(m,k);
}
E[k] = sum + b - M[k][3];
sum = 0;
oldalpha[i] = alpha[i];
oldalpha[k] = alpha[k];
if (M[i][3] != M[k][3]){L = Max(0, alpha[k] - alpha[i]);}
else{ L = Max(0, alpha[i] + alpha[k] - C);}
if (M[i][3] != M[k][3]){H = Min(C,C + alpha[k] - alpha[i]);}
else{ H = Min(C, alpha[i] + alpha[k]);}
if (L==H) {continue;}
eta = 2*innerproduct(i,k) - innerproduct(i,i) - innerproduct(k,k);
if (eta>=0) {continue;}
alpha[k] += M[k][3]*(E[i] - E[k])/eta;
if (alpha[k] > H) {alpha[k] = H;}
else if (alpha[k] < L) {alpha[k] = L;}
if (abs(alpha[k] - oldalpha[k]) < 1e-5) {continue;}
alpha[i] += M[i][3]*M[k][3]*(oldalpha[k] - alpha[k]);
double b1 = b - E[i] - M[i][3]*(alpha[i] - oldalpha[i])*innerproduct(i,i) - M[k][3]*(alpha[k] - oldalpha[k])*innerproduct(i,k);
double b2 = b - E[k] - M[i][3]*(alpha[i] - oldalpha[i])*innerproduct(i,k) - M[k][3]*(alpha[k] - oldalpha[k])*innerproduct(k,k);
if (alpha[i] > 0 && alpha[i] < C) {b = b1;}
else if (alpha[k] > 0 && alpha[k] < C) {b = b2;}
else {b = (b1+b2)/2;}
num_changed_alphas++;
}
}
if (num_changed_alphas==0) {passes++;}
else {passes = 0;}
}
return 0;
}
