SHOGUN  v3.2.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
RationalApproximationCGMJob.cpp
Go to the documentation of this file.
00001 /*
00002  * This program is free software; you can redistribute it and/or modify
00003  * it under the terms of the GNU General Public License as published by
00004  * the Free Software Foundation; either version 3 of the License, or
00005  * (at your option) any later version.
00006  *
00007  * Written (W) 2013 Soumyajit De
00008  */
00009 
00010 #include <shogun/lib/config.h>
00011 
00012 #ifdef HAVE_EIGEN3
00013 #include <shogun/lib/SGVector.h>
00014 #include <shogun/lib/computation/jobresult/ScalarResult.h>
00015 #include <shogun/mathematics/eigen3.h>
00016 #include <shogun/mathematics/linalg/linop/LinearOperator.h>
00017 #include <shogun/mathematics/linalg/linsolver/CGMShiftedFamilySolver.h>
00018 #include <shogun/mathematics/linalg/ratapprox/logdet/computation/job/RationalApproximationCGMJob.h>
00019 #include <shogun/base/Parameter.h>
00020 
00021 using namespace Eigen;
00022 
00023 namespace shogun
00024 {
00025 
00026 CRationalApproximationCGMJob::CRationalApproximationCGMJob()
00027     : CIndependentJob()
00028 {
00029     init();
00030 }
00031 
00032 CRationalApproximationCGMJob::CRationalApproximationCGMJob(
00033     CStoreScalarAggregator<float64_t>* aggregator,
00034     CCGMShiftedFamilySolver* linear_solver,
00035     CLinearOperator<float64_t>* linear_operator,
00036     SGVector<float64_t> vector,
00037     SGVector<complex128_t> shifts,
00038     SGVector<complex128_t> weights,
00039     float64_t const_multiplier)
00040     : CIndependentJob((CJobResultAggregator*)aggregator)
00041 {
00042     init();
00043 
00044     m_linear_solver=linear_solver;
00045     SG_REF(m_linear_solver);
00046 
00047     m_operator=linear_operator;
00048     SG_REF(m_operator);
00049 
00050     m_vector=vector;
00051 
00052     m_shifts=shifts;
00053     m_weights=weights;
00054     m_const_multiplier=const_multiplier;
00055 }
00056 
00057 void CRationalApproximationCGMJob::init()
00058 {
00059     m_linear_solver=NULL;
00060     m_operator=NULL;
00061     m_const_multiplier=0.0;
00062 
00063     SG_ADD((CSGObject**)&m_linear_solver, "linear_solver",
00064         "Linear solver for complex-shifted system", MS_NOT_AVAILABLE);
00065 
00066     SG_ADD((CSGObject**)&m_operator, "linear_operator",
00067         "Linear operator", MS_NOT_AVAILABLE);
00068 
00069     SG_ADD(&m_vector, "trace_sample",
00070         "Sample vector to apply linear operator on", MS_NOT_AVAILABLE);
00071 
00072     SG_ADD(&m_weights, "complex_shifts",
00073         "Shifts in the linear systems to be solved", MS_NOT_AVAILABLE);
00074 
00075     SG_ADD(&m_weights, "complex_weights",
00076         "Weights to be multiplied to the solution vector", MS_NOT_AVAILABLE);
00077 
00078     SG_ADD(&m_const_multiplier, "constant_multiplier",
00079         "Constant multiplier to be multiplied with the final solution", MS_NOT_AVAILABLE);
00080 }
00081 
00082 CRationalApproximationCGMJob::~CRationalApproximationCGMJob()
00083 {
00084     SG_UNREF(m_linear_solver);
00085     SG_UNREF(m_operator);
00086 }
00087 
00088 void CRationalApproximationCGMJob::compute()
00089 {
00090     SG_DEBUG("Entering\n");
00091 
00092     REQUIRE(m_aggregator, "Job result aggregator is not set!\n");
00093     REQUIRE(m_operator, "Operator is not set!\n");
00094     REQUIRE(m_vector.vector, "Vector is not set!\n");
00095     REQUIRE(m_shifts.vector, "Shifts are not set!\n");
00096     REQUIRE(m_weights.vector, "Weights are not set!\n");
00097     REQUIRE(m_operator->get_dimension()==m_vector.vlen,
00098         "Dimension mismatch! %d vs %d\n", m_operator->get_dimension(), m_vector.vlen);
00099     REQUIRE(m_shifts.vlen==m_weights.vlen,
00100         "Number of shifts and weights are not equal!\n");
00101 
00102     // solve the linear system with the sample vector
00103     SGVector<complex128_t> vec=m_linear_solver->solve_shifted_weighted(
00104         m_operator, m_vector, m_shifts, m_weights);
00105 
00106     // Take negative (see CRationalApproximation for the formula)
00107     Map<VectorXcd> v(vec.vector, vec.vlen);
00108     v=-v;
00109 
00110     // take out the imaginary part of the result before
00111     // applying linear operator
00112     SGVector<float64_t> agg=m_operator->apply(vec.get_imag());
00113 
00114     // perform dot product
00115     Map<VectorXd> map_agg(agg.vector, agg.vlen);
00116     Map<VectorXd> map_vector(m_vector.vector, m_vector.vlen);
00117     float64_t result=map_vector.dot(map_agg);
00118 
00119     result*=m_const_multiplier;
00120 
00121     // form the final result into a scalar result and submit to the aggregator
00122     CScalarResult<float64_t>* final_result=new CScalarResult<float64_t>(result);
00123     SG_REF(final_result);
00124 
00125     m_aggregator->submit_result(final_result);
00126 
00127     SG_UNREF(final_result);
00128 
00129     SG_DEBUG("Leaving\n");
00130 }
00131 
00132 }
00133 #endif // HAVE_EIGEN3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

SHOGUN Machine Learning Toolbox - Documentation