smoothed_aggregation.h

Go to the documentation of this file.

/*
 *  Copyright 2008-2014 NVIDIA Corporation
 *
 *  Licensed under the Apache License, Version 2.0 (the "License");
 *  you may not use this file except in compliance with the License.
 *  You may obtain a copy of the License at
 *
 *      http://www.apache.org/licenses/LICENSE-2.0
 *
 *  Unless required by applicable law or agreed to in writing, software
 *  distributed under the License is distributed on an "AS IS" BASIS,
 *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *  See the License for the specific language governing permissions and
 *  limitations under the License.
 */

#pragma once

#include <cusp/detail/config.h>
#include <cusp/detail/multilevel.h>

#include <cusp/array1d.h>
#include <cusp/complex.h>
#include <cusp/precond/aggregation/detail/sa_view_traits.h>

#include <thrust/execution_policy.h>
#include <thrust/detail/use_default.h>

#include <vector> // TODO replace with host_vector

namespace cusp
{
namespace precond
{
namespace aggregation
{

/* \cond */
template<typename MatrixType>
struct sa_level
{
    public:

    typedef typename MatrixType::index_type IndexType;
    typedef typename MatrixType::value_type ValueType;
    typedef typename MatrixType::memory_space MemorySpace;
    typedef typename cusp::norm_type<ValueType>::type NormType;

    MatrixType A_;                                                // matrix
    MatrixType T;                                                 // matrix
    cusp::array1d<IndexType,MemorySpace> aggregates;      // aggregates
    cusp::array1d<ValueType,MemorySpace> B;               // near-nullspace candidates

    size_t   num_iters;
    NormType rho_DinvA;

    sa_level(void) : num_iters(1), rho_DinvA(0) {}

    template<typename SALevelType>
    sa_level(const SALevelType& L)
      : A_(L.A_),
        aggregates(L.aggregates),
        B(L.B),
        num_iters(L.num_iters),
        rho_DinvA(L.rho_DinvA)
    {}
};
/* \endcond */

template <typename IndexType,
          typename ValueType,
          typename MemorySpace,
            typename SmootherType = thrust::use_default,
            typename SolverType   = thrust::use_default,
          typename Format       = thrust::use_default>
class smoothed_aggregation :
    public cusp::multilevel<IndexType,ValueType,MemorySpace,Format,SmootherType,SolverType>::container
{
  private:

    typedef typename detail::select_sa_matrix_type<IndexType,ValueType,MemorySpace>::type   SetupMatrixType;
    typedef typename cusp::multilevel<IndexType,ValueType,MemorySpace,Format,SmootherType,SolverType>::container ML;

  public:

    /* \cond */
    std::vector< sa_level<SetupMatrixType> > sa_levels;
    /* \endcond */

    smoothed_aggregation(void) : ML() {};

    template <typename MatrixType>
    smoothed_aggregation(const MatrixType& A);

    template <typename MatrixType,typename ArrayType>
    smoothed_aggregation(const MatrixType& A, const ArrayType& B,
                         typename thrust::detail::enable_if_convertible<typename ArrayType::format,cusp::array1d_format>::type* = 0);

    template <typename MemorySpace2,typename SmootherType2,typename SolverType2,typename Format2>
    smoothed_aggregation(const smoothed_aggregation<IndexType,ValueType,MemorySpace2,SmootherType2,SolverType2,Format2>& M);

    template <typename MatrixType>
    void initialize(const MatrixType& A);

    template <typename MatrixType, typename ArrayType>
    void initialize(const MatrixType& A, const ArrayType& B,
                    typename thrust::detail::enable_if_convertible<typename MatrixType::format,cusp::known_format>::type* = 0);

    /* \cond */
    template <typename DerivedPolicy, typename MatrixType>
    void initialize(const thrust::detail::execution_policy_base<DerivedPolicy> &exec,
                    const MatrixType& A);

    template <typename DerivedPolicy, typename MatrixType, typename ArrayType>
    void initialize(const thrust::detail::execution_policy_base<DerivedPolicy> &exec,
                    const MatrixType& A, const ArrayType& B);
    /* \endcond */

protected:

    /* \cond */
    template <typename DerivedPolicy, typename MatrixType>
    void extend_hierarchy(const thrust::detail::execution_policy_base<DerivedPolicy> &exec,
                          const MatrixType& A);
    /* \endcond */
};
} // end namespace aggregation
} // end namespace precond
} // end namespace cusp

#include <cusp/precond/aggregation/detail/smoothed_aggregation.inl>