Line data    Source code

       1             : /* 
       2             :  * File:   IndProbSimplex.cpp
       3             :  * Author: Pantelis Sopasakis
       4             :  * 
       5             :  * Created on January 15, 2016, 2:37 AM
       6             :  * 
       7             :  * ForBES is free software: you can redistribute it and/or modify
       8             :  * it under the terms of the GNU Lesser General Public License as published by
       9             :  * the Free Software Foundation, either version 3 of the License, or
      10             :  * (at your option) any later version.
      11             :  *
      12             :  * ForBES is distributed in the hope that it will be useful,
      13             :  * but WITHOUT ANY WARRANTY; without even the implied warranty of
      14             :  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
      15             :  * GNU Lesser General Public License for more details.
      16             :  *
      17             :  * You should have received a copy of the GNU Lesser General Public License
      18             :  * along with ForBES. If not, see <http://www.gnu.org/licenses/>.
      19             :  */
      20             : 
      21             : #include "IndProbSimplex.h"
      22             : #include "MatrixFactory.h"
      23             : #include <vector>
      24             : #include <algorithm>
      25             : #include <cmath>
      26             : 
      27           3 : IndProbSimplex::IndProbSimplex() {
      28             : 
      29           3 : }
      30             : 
      31           6 : IndProbSimplex::~IndProbSimplex() {
      32           6 : }
      33             : 
      34           3 : int IndProbSimplex::callProx(Matrix& x, double gamma, Matrix& prox) {
      35             : 
      36           3 :     size_t n = x.getNrows();
      37             : 
      38           3 :     std::vector<double> x_hat_vec(x.getData(), x.getData() + n);
      39             : 
      40             :     /* x_hat := rev_sort(x_hat) */
      41           3 :     std::sort(x_hat_vec.rbegin(), x_hat_vec.rend());
      42             : 
      43           3 :     double t = 0.0;
      44           3 :     bool flag = true;
      45           3 :     size_t i = 0;
      46           3 :     double val = 0.0;
      47             : 
      48             :     /* skipping negative data */
      49           3 :     const std::vector<double>::iterator lb = std::lower_bound(x_hat_vec.begin(), x_hat_vec.end(), 0.0, std::greater<double>());
      50             : 
      51           3 :     std::vector<double>::iterator it;
      52           3 :     it += (lb - x_hat_vec.begin());
      53             : 
      54         993 :     for (it = x_hat_vec.begin(); flag && it != x_hat_vec.end(); ++it, ++i) {
      55         990 :         val = std::max(0.0, *it);
      56         990 :         t += val;
      57         990 :         flag = val * (i + 1.0) > t - 1.0;
      58             :     }
      59           3 :     t -= 1.0 + val;
      60             : 
      61           3 :     const double theta = std::max(0.0, t / (i - 1));
      62             : 
      63      500023 :     for (size_t j = 0; j < n; j++) {
      64      500020 :         prox[j] = std::max(0.0, x[j] - theta);
      65             :     }
      66             : 
      67           3 :     return ForBESUtils::STATUS_OK;
      68             : }
      69             : 
      70           1 : int IndProbSimplex::callProx(Matrix& x, double gamma, Matrix& prox, double& f_at_prox) {
      71           1 :     f_at_prox = 0.0;
      72           1 :     return callProx(x, gamma, prox);
      73             : }
      74             : 
      75           0 : int IndProbSimplex::call(Matrix& x, double& f) {
      76             :     /* Implement me! */
      77           0 :     return ForBESUtils::STATUS_UNDEFINED_FUNCTION;
      78             : }
      79             : 
      80           9 : FunctionOntologicalClass IndProbSimplex::category() {
      81           9 :     return FunctionOntologyRegistry::indicator();
      82           3 : }