GRASS Programmer's Manual
6.4.2(2012)
|
#include <math.h>
#include <unistd.h>
#include <stdio.h>
#include <string.h>
#include "grass/N_pde.h"
#include "solvers_local_proto.h"
Go to the source code of this file.
Functions | |
int | N_solver_pcg (N_les *L, int maxit, double err, int prec) |
The iterative preconditioned conjugate gradients solver for symmetric positive definite matrices. | |
int | N_solver_cg (N_les *L, int maxit, double err) |
The iterative conjugate gradients solver for symmetric positive definite matrices. | |
int | N_solver_bicgstab (N_les *L, int maxit, double err) |
The iterative biconjugate gradients solver with stabilization for unsymmetric non-definite matrices. | |
void | N_matrix_vector_product (N_les *L, double *x, double *result) |
Calculates the matrix - vector product of matrix L->A and vector x. | |
void | N_sparse_matrix_vector_product (N_les *L, double *x, double *result) |
Calculates the matrix - vector product of sparse matrix L->Asp and vector x. | |
int | check_symmetry (N_les *L) |
Check if the matrix in les is symmetric. | |
N_les * | N_create_diag_precond_matrix (N_les *L, int prec) |
Compute a diagonal preconditioning matrix for krylov space solver. |
int check_symmetry | ( | N_les * | L | ) |
Check if the matrix in les is symmetric.
L | N_les* – the linear equation system |
Definition at line 823 of file N_solvers_krylov.c.
References N_les::A, N_les::Asp, N_spvector::cols, G_debug(), G_warning(), N_spvector::index, N_SPARSE_LES, N_les::quad, N_les::rows, SYMM_TOLERANCE, N_les::type, and N_spvector::values.
Referenced by N_solver_cg(), N_solver_cholesky(), and N_solver_pcg().
Compute a diagonal preconditioning matrix for krylov space solver.
L | N_les* prec int – the preconditioner which should be choosen N_DIAGONAL_PRECONDITION, N_ROWSUM_PRECONDITION |
Definition at line 909 of file N_solvers_krylov.c.
References N_les::A, N_les::Asp, N_spvector::cols, N_les::cols, gui_modules.psmap_dialogs::cols, N_spvector::index, N_add_spvector_to_les(), N_alloc_les_A(), N_alloc_spvector(), N_DIAGONAL_PRECONDITION, N_NORMAL_LES, N_ROWSCALE_ABSSUMNORM_PRECONDITION, N_ROWSCALE_EUKLIDNORM_PRECONDITION, N_SPARSE_LES, N_les::rows, N_les::type, and N_spvector::values.
Referenced by N_solver_pcg().
void N_matrix_vector_product | ( | N_les * | L, |
double * | x, | ||
double * | result | ||
) |
Calculates the matrix - vector product of matrix L->A and vector x.
The result is written to vector named result. This function only works with regular quadratic matrices.
L | N_les * |
x | double * |
result | double * |
Definition at line 594 of file N_solvers_krylov.c.
References N_les::A, N_les::cols, and N_les::rows.
Referenced by N_les_integrate_dirichlet_2d(), N_les_integrate_dirichlet_3d(), N_solver_bicgstab(), N_solver_cg(), and N_solver_pcg().
int N_solver_bicgstab | ( | N_les * | L, |
int | maxit, | ||
double | err | ||
) |
The iterative biconjugate gradients solver with stabilization for unsymmetric non-definite matrices.
The result is written to the vector L->x of the les. This iterative solver works with sparse matrices and regular quadratic matrices.
The parameter maxit specifies the maximum number of iterations. If the maximum is reached, the solver will abort the calculation and writes the current result into the vector L->x. The parameter err defines the error break criteria for the solver.
L | N_les * – the linear equatuin system |
maxit | int – the maximum number of iterations |
err | double – defines the error break criteria |
Definition at line 407 of file N_solvers_krylov.c.
References b, N_les::b, G_free(), G_message(), G_warning(), N_matrix_vector_product(), N_SPARSE_LES, N_sparse_matrix_vector_product(), N_les::quad, r, N_les::rows, gui_modules.psmap_dialogs::s, N_les::type, vectmem(), and N_les::x.
int N_solver_cg | ( | N_les * | L, |
int | maxit, | ||
double | err | ||
) |
The iterative conjugate gradients solver for symmetric positive definite matrices.
This iterative solver works with symmetric positive definite sparse matrices and symmetric positive definite regular quadratic matrices.
The parameter maxit specifies the maximum number of iterations. If the maximum is reached, the solver will abort the calculation and writes the current result into the vector L->x. The parameter err defines the error break criteria for the solver.
L | N_les * – the linear equatuin system |
maxit | int – the maximum number of iterations |
err | double – defines the error break criteria |
Definition at line 242 of file N_solvers_krylov.c.
References b, N_les::b, check_symmetry(), G_free(), G_message(), G_warning(), N_matrix_vector_product(), N_SPARSE_LES, N_sparse_matrix_vector_product(), N_les::quad, r, N_les::rows, gui_modules.psmap_dialogs::s, N_les::type, vectmem(), and N_les::x.
int N_solver_pcg | ( | N_les * | L, |
int | maxit, | ||
double | err, | ||
int | prec | ||
) |
The iterative preconditioned conjugate gradients solver for symmetric positive definite matrices.
This iterative solver works with symmetric positive definite sparse matrices and symmetric positive definite regular quadratic matrices.
The parameter maxit specifies the maximum number of iterations. If the maximum is reached, the solver will abort the calculation and writes the current result into the vector L->x. The parameter err defines the error break criteria for the solver.
L | N_les * – the linear equatuin system |
maxit | int – the maximum number of iterations |
err | double – defines the error break criteria |
prec | int – the preconditioner which shoudl be used N_DIAGONAL_PRECONDITION, N_ROWSUM_PRECONDITION |
Definition at line 64 of file N_solvers_krylov.c.
References b, N_les::b, check_symmetry(), G_free(), G_message(), G_warning(), N_create_diag_precond_matrix(), N_matrix_vector_product(), N_SPARSE_LES, N_sparse_matrix_vector_product(), N_les::quad, r, N_les::rows, gui_modules.psmap_dialogs::s, N_les::type, vectmem(), and N_les::x.
void N_sparse_matrix_vector_product | ( | N_les * | L, |
double * | x, | ||
double * | result | ||
) |
Calculates the matrix - vector product of sparse matrix L->Asp and vector x.
The result is written to vector named result. This function only works with sparse matrices matrices.
L | N_les * |
x | double * |
result | double * |
Definition at line 622 of file N_solvers_krylov.c.
References N_les::Asp, N_spvector::cols, N_spvector::index, N_les::rows, and N_spvector::values.
Referenced by N_les_integrate_dirichlet_2d(), N_les_integrate_dirichlet_3d(), N_solver_bicgstab(), N_solver_cg(), and N_solver_pcg().