if(new_target < old_target){    KKTerror = 0;    //      target_count=0;    if(parameters->descend < old_target - new_target){      target_count=0;    }    else{      convError=1;    };    if(parameters->verbosity>=5){      cout<<"descend = "<<old_target-new_target<<"  ("<<old_target<<" --> "<<new_target<<")"<<endl;    };  }  else{    // nothing we can do    KKTerror=0;    convError=1;    if(parameters->verbosity>=5){      cout<<"WARNING: no descend ("<<old_target<<" -> "<<new_target<<"), stopping"<<endl;    };  };  if(convError){    target_count++;    if(new_target>=old_target){      for(i=0;i<working_set_size;i++){	primal[i] = qp.A[i]*all_alphas[working_set[i]];      };      };                                                                            if(parameters->verbosity>=5){	      cout<<"WARNING: Convergence error, setting sigfig = "<<sigfig_max<<endl;    };  };  if(target_count>50){    // non-recoverable numerical error    feasible_epsilon=1;    convergence_epsilon*=2;    if(parameters->verbosity>=1)      cout<<"WARNING: reducing KKT precision to "<<convergence_epsilon<<endl;    target_count=0;  };  if(parameters->verbosity>=5){	    cout<<"Resulting values:"<<endl;    for(i=0;i<working_set_size;i++){      cout<<i<<": "<<primal[i]<<endl;    };  };  time_optimize += get_time() - time_start;};void svm_nu_regression_c::print_special_statistics(){  // calculate tube size epsilon  SVMFLOAT b = examples->get_b();  SVMFLOAT epsilon_pos = 0;  SVMFLOAT epsilon_neg = 0;  SVMINT pos_count = 0;  SVMINT neg_count = 0;  SVMINT i;  for(i=0;i<examples_total;i++){    if((all_alphas[i] > is_zero) && (all_alphas[i]-Cpos<-is_zero)){      epsilon_neg += all_ys[i]-sum[i]-b;      neg_count++;    }    else if((all_alphas[i] <- is_zero) && (all_alphas[i]+Cneg>+is_zero)){      epsilon_pos += -all_ys[i]+sum[i]+b;      pos_count++;    };  };  if((parameters->Lpos == parameters->Lneg) ||     (pos_count == 0) ||     (neg_count == 0)){    // symmetrical    epsilon_pos += epsilon_neg;    pos_count += neg_count;    if(pos_count>0){      epsilon_pos /= (SVMINT)pos_count;      cout<<"epsilon = "<<epsilon_pos<<endl;    }    else{      cout<<"ERROR: could not calculate epsilon."<<endl;      cout<<pos_count<<"\t"<<neg_count<<endl; // @@@@@@@    };  }  else{    // asymmetrical    epsilon_pos /= (SVMINT)pos_count;    cout<<"epsilon+ = "<<epsilon_pos<<endl;    epsilon_neg /= (SVMINT)neg_count;    cout<<"epsilon- = "<<epsilon_pos<<endl;  };};/** * * svm_nu_pattern_c * **/SVMFLOAT svm_nu_pattern_c::nabla(const SVMINT i){  if(all_ys[i] > 0){    return( sum[i]);  }  else{    return(-sum[i]);  };};void svm_nu_pattern_c::init(kernel_c* new_kernel, parameters_c* new_parameters){  new_parameters->realC = 1;  svm_nu_regression_c::init(new_kernel,new_parameters);};void svm_nu_pattern_c::init_optimizer(){  // Cs are dived by examples_total in init_optimizer  svm_nu_regression_c::init_optimizer();  SVMINT i;  for(i=0;i<working_set_size;i++){    qp.l[i] = 0;  };};void svm_nu_pattern_c::update_working_set(){  svm_c::update_working_set();  SVMINT i;  for(i=0;i<working_set_size;i++){    if(qp.A[i]>0){      qp.c[i] += all_ys[working_set[i]];    }    else{      qp.c[i] -= all_ys[working_set[i]];    };  };};void svm_nu_pattern_c::init_working_set(){  // calculate nu-sum   if(examples->initialised_alpha()){    project_to_constraint();  };  sum_alpha_nu=0;  SVMFLOAT the_nu_sum = 0;  SVMFLOAT the_sum=0;  SVMINT pos_count=0;  SVMINT neg_count=0;  SVMINT ni;  for(ni=0;ni<examples_total;ni++){    the_sum += all_alphas[ni];    the_nu_sum += abs(all_alphas[ni]);    if(is_alpha_neg(ni)> 0){      neg_count++;    }    else{      pos_count++;    };  };  if((abs(the_sum) > is_zero) || (abs(the_nu_sum-nu) > is_zero)){    // set initial feasible point    // neg alpha: -nu/2n    // pos alpha:  nu/2p    if((nu*(SVMFLOAT)examples_total>2*(SVMFLOAT)pos_count) ||       (nu*(SVMFLOAT)examples_total>2*(SVMFLOAT)neg_count)){      nu = 2*((SVMFLOAT)pos_count)/((SVMFLOAT)examples_total);      if(nu > 2*((SVMFLOAT)neg_count)/((SVMFLOAT)examples_total)){	nu = 2*((SVMFLOAT)neg_count)/((SVMFLOAT)examples_total);      };      nu -= is_zero; // just to make sure      cout<<"ERROR: nu too large, setting nu = "<<nu<<endl;    };    for(ni=0;ni<examples_total;ni++){      if(is_alpha_neg(ni)> 0){	examples->put_alpha(ni,nu/(2*(SVMFLOAT)neg_count));      }      else{	examples->put_alpha(ni,-nu/(2*(SVMFLOAT)pos_count));      };    };    examples->set_initialised_alpha();  };  svm_c::init_working_set();};void svm_nu_pattern_c::print_special_statistics(){  // calculate margin rho  SVMFLOAT b = examples->get_b();  SVMFLOAT rho_pos = 0;  SVMFLOAT rho_neg = 0;  SVMINT pos_count = 0;  SVMINT neg_count = 0;  SVMINT i;  for(i=0;i<examples_total;i++){    if((all_alphas[i] > is_zero) && (all_alphas[i]-Cpos<-is_zero)){      rho_neg += sum[i]+b;      neg_count++;    }    else if((all_alphas[i] <- is_zero) && (all_alphas[i]+Cneg>+is_zero)){      rho_pos += -sum[i]-b;      pos_count++;    };  };  if((parameters->Lpos == parameters->Lneg) ||     (pos_count == 0) ||     (neg_count == 0)){    // symmetrical    rho_pos += rho_neg;    pos_count += neg_count;    if(pos_count>0){      rho_pos /= (SVMINT)pos_count;      cout<<"margin = "<<rho_pos<<endl;    }    else{      cout<<"ERROR: could not calculate margin."<<endl;    };  }  else{    // asymmetrical    rho_pos /= (SVMINT)pos_count;    cout<<"margin+ = "<<rho_pos<<endl;    rho_neg /= (SVMINT)neg_count;    cout<<"margin- = "<<rho_pos<<endl;  };};/** * * svm_distribution_c * **/int svm_distribution_c::is_alpha_neg(const SVMINT i){  // variable i is alpha*  return 1;};SVMFLOAT svm_distribution_c::nabla(const SVMINT i){  return( sum[i]);};SVMFLOAT svm_distribution_c::lambda(const SVMINT i){  // size lagrangian multiplier of the active constraint  SVMFLOAT alpha;  SVMFLOAT result = 0;  alpha=all_alphas[i];  if(alpha>is_zero){    // alpha*    if(alpha-Cneg >= - is_zero){      // upper bound active      result = -lambda_eq-sum[i];    }    else{      result = -abs(sum[i]+lambda_eq);    };  }  else{    // lower bound active    result = sum[i] + lambda_eq;  };  return result;};int svm_distribution_c::feasible(const SVMINT i){  // is direction i feasible to minimize the target function  // (includes which_alpha==0)  if(at_bound[i] >= shrink_const){ return 0; };  SVMFLOAT alpha;  SVMFLOAT result;  alpha=all_alphas[i];  if(alpha-Cneg >= - is_zero){    // alpha* at upper bound    result = -lambda_eq - sum[i];    if(result>=-feasible_epsilon){      return 0;     };  }  else if(alpha<=is_zero){    // lower bound active    result = sum[i]+lambda_eq;    if(result>=-feasible_epsilon){      return 0;     };  }  else{    // not at bound    result= abs(sum[i]+lambda_eq);    if(result<=feasible_epsilon){      return 0;     };  };  return 1;};int svm_distribution_c::feasible(const SVMINT i, SVMFLOAT* the_nabla, SVMFLOAT* the_lambda, int* atbound){  // is direction i feasible to minimize the target function  // (includes which_alpha==0)  int is_feasible=1;  if(at_bound[i] >= shrink_const){ is_feasible = 0; };  SVMFLOAT alpha;  alpha=all_alphas[i];  *the_nabla = sum[i];  if(alpha >= Cneg){ //alpha-Cneg >= - is_zero){    // alpha* at upper bound    *atbound = 1;    *the_lambda = -lambda_eq - *the_nabla; //sum[i] + 1;    if(*the_lambda >= 0){      at_bound[i]++;      if(at_bound[i] == shrink_const) to_shrink++;    }    else{      at_bound[i] = 0;    };  }  else if(alpha <= 0){    // lower bound active    *atbound = -1;    *the_lambda = lambda_eq + *the_nabla; //sum[i] + 1;    if(*the_lambda >= 0){      at_bound[i]++;      if(at_bound[i] == shrink_const) to_shrink++;    }    else{      at_bound[i] = 0;    };  }  else{    // not at bound    *atbound = 0;    *the_lambda = -abs(*the_nabla+lambda_eq);    at_bound[i] = 0;  };  if(*the_lambda >= feasible_epsilon){    is_feasible = 0;   };  return is_feasible;};void svm_distribution_c::init(kernel_c* new_kernel, parameters_c* new_parameters){  new_parameters->realC = 1;  nu = new_parameters->nu;  convergence_epsilon = 1e-4;  svm_pattern_c::init(new_kernel,new_parameters);  //  is_pattern = 1;};void svm_distribution_c::init_optimizer(){  // Cs are dived by examples_total in init_optimizer  svm_pattern_c::init_optimizer();};void svm_distribution_c::project_to_constraint(){  SVMINT total = 0;  SVMFLOAT alpha_sum=sum_alpha-nu;  SVMFLOAT alpha=0;  SVMINT i;  for(i=0;i<examples_total;i++){    alpha = all_alphas[i];    alpha_sum += alpha;    if((alpha>is_zero) && (alpha-Cneg < -is_zero)){      total++;    };  };  if(total>0){    // equality constraint violated    alpha_sum /= (SVMFLOAT)total;    for(i=0;i<examples_total;i++){      if((alpha>is_zero) && (alpha-Cneg < -is_zero)){	all_alphas[i] -= alpha_sum;      };    };  };};int svm_distribution_c::convergence(){  long time_start = get_time();  SVMFLOAT the_lambda_eq = 0;  SVMINT total = 0;  SVMFLOAT alpha_sum=0;  SVMFLOAT alpha=0;  SVMINT i;  int result=1;  // actual convergence-test  total = 0; alpha_sum=0;  //  cout<<Cneg<<"\t"<<nu<<"\t"<<all_alphas[0]<<endl;  for(i=0;i<examples_total;i++){    alpha = all_alphas[i];    alpha_sum += alpha;    if((alpha>is_zero) && (alpha-Cneg < -is_zero)){      // alpha^* = - nabla      the_lambda_eq += -sum[i];      total++;    };  };  if(parameters->verbosity>= 4){    cout<<"lambda_eq = "<<(the_lambda_eq/total)<<endl;  };  if(total>0){    lambda_eq = the_lambda_eq / total;  }  else{    // keep WS lambda_eq    lambda_eq = lambda_WS;    if(parameters->verbosity>= 4){      cout<<"*** no SVs in convergence(), lambda_eq = "<<lambda_eq<<"."<<endl;    };  };  if(target_count>2){    if(target_count>20){      // desperate!      lambda_eq = ((40-target_count)*lambda_eq + (target_count-20)*lambda_WS)/20;      if(parameters->verbosity>=5){	cout<<"Re-Re-calculated lambda from WS: "<<lambda_eq<<endl;      };      if(target_count>40){	// really desperate, kick one example out!	i = working_set[target_count%working_set_size];	lambda_eq = -sum[i];	if(parameters->verbosity>=5){	  cout<<"set lambda_eq to nabla("<<i<<"): "<<lambda_eq<<endl;	};      };    }    else{      lambda_eq = lambda_WS;      if(parameters->verbosity>=5){	cout<<"Re-calculated lambda_eq from WS: "<<lambda_eq<<endl;      };    };  };  // check linear constraint  if(abs(alpha_sum+sum_alpha-nu) > convergence_epsilon){    // equality constraint violated    if(parameters->verbosity>= 4){      cout<<"No convergence: equality constraint violated: |"<<(alpha_sum+sum_alpha)<<"| >> 0"<<endl;    };    project_to_constraint();    result = 0;    };  i=0;  while((i<examples_total) && (result != 0)){    if(lambda(i)>=-convergence_epsilon){      i++;    }    else{      result = 0;    };  };  time_convergence += get_time() - time_start;  return result;};void svm_distribution_c::init_working_set(){  // calculate sum  SVMINT i,j;  if(nu>1){    cout<<"ERROR: nu too large, setting nu to 1"<<endl;    nu = 1-is_zero;  };  SVMFLOAT the_sum=0;  for(i=0; i<examples_total;i++){    the_sum += all_alphas[i];  };  if(abs(the_sum-nu) > is_zero){    for(i=0; i<examples_total;i++){      examples->put_alpha(i,nu/((SVMFLOAT)examples_total));    };    examples->set_initialised_alpha();  };  if(parameters->verbosity >= 3){    cout<<"Initialising variables, this may take some time."<<endl;  };  for(i=0; i<examples_total;i++){    all_ys[i] = 1;    sum[i] = 0;    at_bound[i] = 0;    for(j=0; j<examples_total;j++){      sum[i] += all_alphas[j]*kernel->calculate_K(i,j);    };  };  calculate_working_set();  update_working_set();};void svm_distribution_c::print_special_statistics(){  // calculate margin rho  SVMFLOAT rho = 0;  SVMINT count = 0;  SVMFLOAT norm_x;  SVMFLOAT max_norm_x=-infinity;  //  SVMFLOAT xi_i;  //  SVMINT estim_loo=examples_total;  //  SVMINT estim_loo2=examples_total;  SVMINT svs=0;  SVMINT i;  for(i=0;i<examples_total;i++){    if((all_alphas[i] > is_zero) && (all_alphas[i]-Cpos<-is_zero)){      rho += sum[i];      count++;    };    if(all_alphas[i] != 0){      svs++;      norm_x = kernel->calculate_K(i,i);      if(norm_x>max_norm_x){	max_norm_x = norm_x;      };    };  };  if(count == 0){    cout<<"ERROR: could not calculate margin."<<endl;  }  else{    // put -rho as b (same decision function)    rho /= (SVMINT)count;    examples->put_b(-rho);    cout<<"margin = "<<rho<<endl;  };  cout<<"examples in distribution support : "<<count<<" ("<<((SVMINT)(10000.0*(SVMFLOAT)count/((SVMFLOAT)examples_total)))/100.0<<"%)."<<endl;};