Line data    Source code

       1             : /* 
       2             :  * File:   Norm2.cpp
       3             :  * Author: Pantelis Sopasakis
       4             :  * 
       5             :  * Created on October 30, 2015, 6:54 PM
       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 "Norm2.h"
      22             : #include <cmath>
      23             : 
      24             : double vecNorm2(Matrix& x);
      25             : 
      26           8 : Norm2::Norm2() : Norm() {
      27           8 :     m_mu = 1.0;
      28           8 : }
      29             : 
      30           3 : Norm2::Norm2(double mu) : Norm(), m_mu(mu) {
      31           3 : }
      32             : 
      33          21 : Norm2::~Norm2() {
      34          21 : }
      35             : 
      36         220 : inline double vecNorm2(Matrix& x) {
      37         220 :     double norm_x = 0.0;
      38         562 :     for (size_t i = 0; i < x.getNrows(); i++) {
      39         342 :         norm_x += std::pow(x[i], 2);
      40             :     }
      41         220 :     norm_x = std::sqrt(norm_x);
      42         220 :     return norm_x;
      43             : }
      44             : 
      45         176 : int Norm2::call(Matrix& x, double& f) {
      46             :     //LCOV_EXCL_START
      47             :     if (!x.isColumnVector()) {        
      48             :         throw std::invalid_argument("x must be a column-vector");
      49             :     }
      50             :     //LCOV_EXCL_STOP
      51         176 :     f = m_mu * vecNorm2(x);
      52         176 :     return ForBESUtils::STATUS_OK;
      53             : }
      54             : 
      55          41 : int Norm2::dualNorm(Matrix& x, double& norm) {
      56             :     //LCOV_EXCL_START
      57             :     if (!x.isColumnVector()) {
      58             :         throw std::invalid_argument("x must be a column-vector");
      59             :     }   
      60             :     //LCOV_EXCL_STOP
      61          41 :     norm = vecNorm2(x) / m_mu;
      62          41 :     return ForBESUtils::STATUS_OK;
      63             : }
      64             : 
      65           1 : int Norm2::callProx(Matrix& x, double gamma, Matrix& prox) {
      66             :     //LCOV_EXCL_START
      67             :     if (!x.isColumnVector()) {
      68             :         throw std::invalid_argument("x must be a column-vector");
      69             :     }
      70             :     //LCOV_EXCL_STOP
      71           1 :     double norm_x = vecNorm2(x);
      72           1 :     double gm = gamma*m_mu;
      73           1 :     if (norm_x > gm) {
      74           1 :         double s = 1 - gm / norm_x;
      75           1 :         prox = s * x;
      76             :     }
      77           1 :     return ForBESUtils::STATUS_OK;
      78             : }
      79             : 
      80           2 : int Norm2::callProx(Matrix& x, double gamma, Matrix& prox, double& f_at_prox) {
      81             :     //LCOV_EXCL_START
      82             :     if (!x.isColumnVector()) {
      83             :         throw std::invalid_argument("x must be a column-vector");
      84             :     }
      85             :     //LCOV_EXCL_STOP
      86           2 :     double norm_x = vecNorm2(x);    
      87           2 :     double gm = gamma*m_mu;
      88           2 :     if (norm_x > gm) {
      89           1 :         double s = 1 - gm / norm_x;
      90          11 :         for (size_t i = 0; i < x.getNrows(); i++) {
      91          10 :             prox.set(i, 0, s * x[i]); // prox = s * x;
      92             :         }
      93           1 :         f_at_prox = m_mu * s * norm_x;
      94             :     } else {
      95           1 :         f_at_prox = 0.0;
      96             :     }
      97           2 :     return ForBESUtils::STATUS_OK;
      98             : }
      99             : 
     100           3 : FunctionOntologicalClass Norm2::category() {
     101           3 :     return FunctionOntologyRegistry::norm();
     102           9 : }
     103             : 
     104             : 
     105             : 
     106             : 
     107             : 
     108             :