/usr/include/mia-2.0/mia/core/sparse_solver.hh is in libmia-2.0-dev 2.0.13-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 | /* -*- mia-c++ -*-
*
* This file is part of MIA - a toolbox for medical image analysis
* Copyright (c) Leipzig, Madrid 1999-2013 Gert Wollny
*
* MIA is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with MIA; if not, see <http://www.gnu.org/licenses/>.
*
*/
#ifndef mia_core_sparse_solver_hh
#define mia_core_sparse_solver_hh
#include <mia/core/product_base.hh>
#include <mia/core/traits.hh>
NS_MIA_BEGIN
/**
\ingroup misc
\brief solver for sparse systems of equations
This is the templatex base class for solvers for systems of equations
Ax=b if A where the multiplication Ax can be expressed as a convolution
operation with a N dimansional operator.
\tparam F the field this solver works on. The Field type must provide the
types const_iterator and iterator that implement random access iterators.
\remark this class is currently not used at all and it may actually be better to
design it to interface with a library like LAPACK etc.
*/
template <typename F>
class TSparseSolver : public CProductBase{
public:
/// typedef for the field the solver works on
typedef F Field;
/// helper typedef for the plug-in handling
typedef F plugin_data;
/// typeef of the dimension type
typedef typename dim_traits<F>::dimsize_type dimsize_type;
/// helper for the plug-in handling
typedef TSparseSolver<F> plugin_type;
/// plugin type component of plugin path descriptor
static const char * const type_descr;
/// type of the field random access iterator
typedef typename Field::iterator field_iterator;
/// type of the field random access const_iterator
typedef typename Field::const_iterator const_field_iterator;
/// type of the values stored in Field
typedef typename Field::value_type value_type;
/**
\brief A class that implements the multiplication of a cell of the
matrix A with the field x.
This class defines the (sparse) multiplication of the matrix A with a vector x
*/
class A_mult_x: public CProductBase {
public:
/// helper typedef for the plug-in handling
typedef F plugin_data;
/// helper typedef for the plug-in handling
typedef A_mult_x plugin_type;
/// plugin type component of plugin path descriptor
static const char * const type_descr;
/** standard constructor sets the expected fiels size to zero */
A_mult_x():m_size(dimsize_type()){};
/** constructor sets the expected fiels size to given value */
A_mult_x(const dimsize_type& size):m_size(size){};
/** destructor to ensure virtual distruction */
virtual ~A_mult_x() {};
/**
Operator to execute the multiplication
\param ix is the random access iterator to the
central value of the input vector that would be
multiplied by the element on the main diagonal of A
\returns the result for the cell multiplication
\remark For efficiency no tests are done whether the boundary conditions are actually honoured.
this is the job of the solver calling the operator
*/
virtual value_type operator () (const_field_iterator ix) const = 0;
/**
Returns the size of the boundary that one has to take into account
*/
virtual int get_boundary_size() const = 0;
/**
\returns the sizes of the (multidimensional) array that is expected
by the operator
*/
dimsize_type get_size() const {
return m_size;
}
private:
dimsize_type m_size;
};
virtual ~TSparseSolver() {};
/**
Run the solver
\param rhs the right hand side of teh system of linear equations to be solved
\param[in,out] x at input an initial estimate of the solution, on output the solution
\param mult the class implementing the multiplication with matrix A
\returns 0 if all went fine and an algorithm specific error value if something went wrong.
*/
virtual int solve(const Field& rhs, Field& x, const A_mult_x& mult) const = 0;
};
template <typename F>
const char * const TSparseSolver<F>::type_descr = "sparsesolver";
template <typename F>
const char * const TSparseSolver<F>::A_mult_x::type_descr = "amultx";
/**
\brief Class template to implement a matrix - field multiplication.
This class implements the multiplication of a sparse structured matrix by a field representing
the input vector data. This structure is just a stub ans needs to be instanciated for
specific fields.
\tparam F the field that representes the data to be multiplied by the matrix
*/
template <typename F>
struct multiply {
/**
The function to apply the multiplication
\param[out] result location to store the result, ususally it should be allocated and of the same size like x
\param[in] A sparse Matrix
\param[in] x the filed to multiply with the matrix
*/
static void apply(F& MIA_PARAM_UNUSED(result),
const typename TSparseSolver<F>::A_mult_x& MIA_PARAM_UNUSED(A), const F& MIA_PARAM_UNUSED(x)) {
static_assert(sizeof(F)==0, "This operation needs to be specialized");
}
};
/**
Implement the multiplication of a Matrix with an Image,
\param A the matrix
\param x the input image
\returns A * x, (works like a matrix filter)
\remark maybethis should go into the filter section
*/
template <typename F>
F operator * (const typename TSparseSolver<F>::A_mult_x& A, const F& x)
{
F result(x.get_size());
multiply<F>::apply(result, A, x);
return result;
}
NS_MIA_END
#endif
|