33#ifndef __H__UG__CPU_ALGEBRA__CORE_SMOOTHERS__
34#define __H__UG__CPU_ALGEBRA__CORE_SMOOTHERS__
105template<
typename Matrix_type,
typename Vector_type>
106void gs_step_LL(
const Matrix_type &A, Vector_type &c,
const Vector_type &d,
const number relaxFactor)
110 typedef typename Matrix_type::value_type matrix_block;
111 typedef typename Matrix_type::const_row_iterator const_row_it;
112 typename Vector_type::value_type
s;
114 const size_t sz = c.size();
115 for (
size_t i = 0; i < sz; ++i)
123 const const_row_it rowEnd = A.end_row(i);
124 const_row_it it = A.begin_row(i);
125 for(; it != rowEnd && it.index() < i; ++it)
130 const matrix_block& A_ii = it.index() == i ? it.value() : matrix_block(0);
149template<
typename Matrix_type,
typename Vector_type>
150void gs_step_UR(
const Matrix_type &A, Vector_type &c,
const Vector_type &d,
const number relaxFactor)
154 typename Vector_type::value_type
s;
156 if(c.size() == 0)
return;
157 size_t i = c.size()-1;
161 typename Matrix_type::const_row_iterator diag = A.get_connection(i, i);
163 typename Matrix_type::const_row_iterator it = diag; ++it;
164 for(; it != A.end_row(i); ++it)
188template<
typename Matrix_type,
typename Vector_type>
189void sgs_step(
const Matrix_type &A, Vector_type &c,
const Vector_type &d,
const number relaxFactor)
197 typename Vector_type::value_type
s;
198 for(
size_t i = 0; i<c.size(); i++)
217template<
typename Matrix_type,
typename Vector_type>
218void diag_step(
const Matrix_type& A, Vector_type& c,
const Vector_type& d,
number damp)
221 UG_ASSERT(c.size() == d.size() && c.size() == A.num_rows(), c <<
", " << d <<
222 " and " << A <<
" need to have same size.");
224 for(
size_t i=0; i < c.size(); i++)
void gs_step_UR(const Matrix_type &A, Vector_type &c, const Vector_type &d, const number relaxFactor)
Performs a backward gauss-seidel-step, that is, solve on the upper right of A. Using gs_step_UR withi...
Definition core_smoothers.h:150
void gs_step_LL(const Matrix_type &A, Vector_type &c, const Vector_type &d, const number relaxFactor)
Gauss-Seidel-Iterations.
Definition core_smoothers.h:106
void diag_step(const Matrix_type &A, Vector_type &c, const Vector_type &d, number damp)
Performs a jacobi-step.
Definition core_smoothers.h:218
void sgs_step(const Matrix_type &A, Vector_type &c, const Vector_type &d, const number relaxFactor)
Performs a symmetric gauss-seidel step. Using sgs_step within a preconditioner-scheme leads to the fa...
Definition core_smoothers.h:189
#define UG_ASSERT(expr, msg)
Definition assert.h:70
double number
Definition types.h:124
bool MatMultAdd(vector_t &dest, const number &alpha1, const vector_t &v1, const number &beta1, const matrix_t &A1, const vector_t &w1)
calculates dest = alpha1*v1 + beta1 * A1 *w1;
Definition operations_mat.h:68
bool MatMult(vector_t &dest, const number &beta1, const matrix_t &A1, const vector_t &w1)
calculates dest = beta1 * A1;
Definition operations_mat.h:59
bool InverseMatMult(number &dest, const double &beta, const TMat &mat, const TVec &vec)
you can implement this function with GetInverse and MatMult