VERB_code_2.2
2
|
#include "MatrixSolver.h"
#include <math.h>
#include <malloc.h>
#include <string>
#include <ctime>
#include "../Parameters/Parameters.h"
#include "../VariousFunctions/variousFunctions.h"
#include "../Exceptions/error.h"
#include "../Lapack/Include/cpplapack.h"
#include "../GMRES/itlin.h"
#include <iostream>
#include "../GMRES2/GMRES2.h"
Go to the source code of this file.
Classes | |
class | Preconditioner_jacobi |
class | Preconditioner_SOR |
Macros | |
#define | EE 1.e-19 |
#define | I(l, k) (l)*n+(k) |
Functions | |
void | matvec (int n, double *x, double *b) |
void | preconr (int n, double *x, double *b) |
double | gmres_norm2 (int n, double *v) |
void | gmres (int n, double *y, MATVEC *matvec, PRECON *preconr, PRECON *preconl, double *b, struct ITLIN_OPT *opt, struct ITLIN_INFO *info) |
bool | MakeMatrix (double *A, double *B1, double *C, double *R, double *x, double tau, double n, double f_lower, double f_upper, int nx, double dt, double *Dxx, double taulc, double alc, int g_flag, string lower_border_condition_type, string upper_border_condition_type, string approximationMethod) |
Make matrix for 1d diffusion. More... | |
bool | SolveMatrix (double *f, double *A, double *B1, double *C, double *R, double f_lower, double f_upper, int nx, double dt) |
Solver for 1d general equation. More... | |
void | AddBoundary (DiagMatrix &Matrix_A, string type, int in, int id1, double dh) |
bool | MakeModelMatrix_3D (CalculationMatrix &matr_A, CalculationMatrix &matr_B, CalculationMatrix &matr_C, Matrix3D< double > &L, Matrix3D< double > &pc, Matrix3D< double > &alpha, int L_size, int pc_size, int alpha_size, Matrix2D< double > &L_lowerBoundaryCondition, Matrix2D< double > &L_upperBoundaryCondition, Matrix2D< double > &pc_lowerBoundaryCondition, Matrix2D< double > &pc_upperBoundaryCondition, Matrix2D< double > &alpha_lowerBoundaryCondition, Matrix2D< double > &alpha_upperBoundaryCondition, string L_lowerBoundaryCondition_calculationType, string L_upperBoundaryCondition_calculationType, string pc_lowerBoundaryCondition_calculationType, string pc_upperBoundaryCondition_calculationType, string alpha_lowerBoundaryCondition_calculationType, string alpha_upperBoundaryCondition_calculationType, Matrix3D< double > &DLL, Matrix3D< double > &Dpcpc, Matrix3D< double > &DpcpcLpp, Matrix3D< double > &Daa, Matrix3D< double > &DaaLpp, Matrix3D< double > &Dpca, Matrix3D< double > &DpcaLpp, Matrix3D< double > &Jacobian, double dt, double Lpp, double tau, double tauLpp) |
void | jacobi (int m_size, double *z, double *w) |
void | SOR (int m_size, double *z, double *w) |
void | for_norm_A (int n, double *x, double *b) |
void | gmres_wrapout (CalculationMatrix &A, Matrix1D< double > &B, Matrix1D< double > &X, Parameters_structure::PSD::GMRES_parameters_structure GMRES_parameters) |
void | gmres2_wrapout (CalculationMatrix &A, Matrix1D< double > &B, Matrix1D< double > &X, Parameters_structure::PSD::GMRES_parameters_structure GMRES_parameters) |
void | Lapack (DiagMatrix &A, Matrix1D< double > &B, Matrix1D< double > &X) |
void | over_relaxation_diag (DiagMatrix &A, Matrix1D< double > &B, Matrix1D< double > &X, int max_steps, double EE) |
void | GetDerivativeVector (string derivativeType, int &dL, int &dPc, int &dAlpha) |
void | SecondDerivativeApproximation_3D (CalculationMatrix &matr_A, int il, int im, int ia, string FirstDerivative, string SecondDerivative, Matrix3D< double > &L, Matrix3D< double > &pc, Matrix3D< double > &alpha, Matrix3D< double > &Dxx, Matrix3D< double > &Jacobian, double multiplicator) |
double | max2 (double a, double b) |
Function returns maximum between a and b. More... | |
void | mult_vect2 (double *af, double *vf, double *wf, int nf) |
Multiplocation between vector and matrix. More... | |
void | gauss_solve (double *a, double *b, int n) |
Gauss inversion. More... | |
bool | tridag (double a[], double b[], double c[], double r[], double u[], long n) |
Solve the AU=R system of equations, where A - tridiagonal matrix nxn with diagonals a[], b[], c[]. More... | |
bool | MakeModelMatrix_3D_KC (CalculationMatrix &matr_A, CalculationMatrix &matr_B, CalculationMatrix &matr_C, Matrix3D< double > &L, Matrix3D< double > &pc, Matrix3D< double > &alpha, int L_size, int pc_size, int alpha_size, Matrix2D< double > &L_lowerBoundaryCondition, Matrix2D< double > &L_upperBoundaryCondition, Matrix2D< double > &pc_lowerBoundaryCondition, Matrix2D< double > &pc_upperBoundaryCondition, Matrix2D< double > &alpha_lowerBoundaryCondition, Matrix2D< double > &alpha_upperBoundaryCondition, string L_lowerBoundaryCondition_calculationType, string L_upperBoundaryCondition_calculationType, string pc_lowerBoundaryCondition_calculationType, string pc_upperBoundaryCondition_calculationType, string alpha_lowerBoundaryCondition_calculationType, string alpha_upperBoundaryCondition_calculationType, Matrix3D< double > &DLL, Matrix3D< double > &Dpcpc, Matrix3D< double > &DpcpcLpp, Matrix3D< double > &Daa, Matrix3D< double > &DaaLpp, Matrix3D< double > &Dpca, Matrix3D< double > &DpcaLpp, Matrix3D< double > &Jacobian, double dt, double Lpp, double tau, double tauLpp) |
void | SecondDerivativeApproximation_3D_KC_Diagonal (CalculationMatrix &matr_A, int il, int im, int ia, string FirstDerivative, string SecondDerivative, Matrix3D< double > &L, Matrix3D< double > &pc, Matrix3D< double > &alpha, Matrix3D< double > &Dxx, Matrix3D< double > &Jacobian, double multiplicator) |
void | SecondDerivativeApproximation_3D_KC_Mixed (CalculationMatrix &matr_A, int il, int im, int ia, string FirstDerivative, string SecondDerivative, Matrix3D< double > &L, Matrix3D< double > &pc, Matrix3D< double > &alpha, Matrix3D< double > &Dxx, Matrix3D< double > &Jacobian, double multiplicator) |
Variables | |
static CalculationMatrix | A_copy |
Making model matrixes, solving model matrixes.
Matrix form of linear equations: A*X[t+1] = B*X[t], where A - model matrix, B - RHS, X[t] - known values of function (PSD), X[t+1] - unknown values of function.
In that file there are procedures for making model matrix for 1d-diffusion, 2d-diffusion, some ideas of 3d-diffusion and mixed terms. Solver for tridiagonal matrix, solver by gauss method and iteration method (upper relaxation).
This file is under development, and has a lot of commented code, unfinished etc code.
There is a checked version for 1d diffusion (or split method of 2d, 3d diffusions) - 1d_universal_solver.
Definition in file MatrixSolver.cpp.
#define EE 1.e-19 |
Definition at line 1186 of file MatrixSolver.cpp.
#define I | ( | l, | |
k | |||
) | (l)*n+(k) |
Definition at line 1187 of file MatrixSolver.cpp.
void matvec | ( | int | m_size, |
double * | x, | ||
double * | b | ||
) |
Vector to matrix multiplication for GMRES
Definition at line 578 of file MatrixSolver.cpp.
References A_copy.
void preconr | ( | int | n, |
double * | x, | ||
double * | b | ||
) |
double gmres_norm2 | ( | int | n, |
double * | v | ||
) |
void gmres | ( | int | n, |
double * | y, | ||
MATVEC * | matvec, | ||
PRECON * | preconr, | ||
PRECON * | preconl, | ||
double * | b, | ||
struct ITLIN_OPT * | opt, | ||
struct ITLIN_INFO * | info | ||
) |
bool MakeMatrix | ( | double * | A, |
double * | B1, | ||
double * | C, | ||
double * | R, | ||
double * | x, | ||
double | tau, | ||
double | n, | ||
double | f_lower, | ||
double | f_upper, | ||
int | nx, | ||
double | dt, | ||
double * | Dxx, | ||
double | taulc, | ||
double | alc, | ||
int | g_flag, | ||
string | lower_border_condition_type, | ||
string | upper_border_condition_type, | ||
string | approximationMethod | ||
) |
Make matrix for 1d diffusion.
Make matrix for 1D diffusion.
Definition at line 71 of file MatrixSolver.cpp.
bool SolveMatrix | ( | double * | f, |
double * | A, | ||
double * | B1, | ||
double * | C, | ||
double * | R, | ||
double | f_lower, | ||
double | f_upper, | ||
int | nx, | ||
double | dt | ||
) |
Solver for 1d general equation.
Solver for 1d-diffusion matrix (tridiagonal).
Definition at line 241 of file MatrixSolver.cpp.
References tridag().
void AddBoundary | ( | DiagMatrix & | Matrix_A, |
string | type, | ||
int | in, | ||
int | id1, | ||
double | dh | ||
) |
Supportive sub-function to add boundary conditions to model matrix
Definition at line 282 of file MatrixSolver.cpp.
bool MakeModelMatrix_3D | ( | CalculationMatrix & | matr_A, |
CalculationMatrix & | matr_B, | ||
CalculationMatrix & | matr_C, | ||
Matrix3D< double > & | L, | ||
Matrix3D< double > & | pc, | ||
Matrix3D< double > & | alpha, | ||
int | L_size, | ||
int | pc_size, | ||
int | alpha_size, | ||
Matrix2D< double > & | L_lowerBoundaryCondition, | ||
Matrix2D< double > & | L_upperBoundaryCondition, | ||
Matrix2D< double > & | pc_lowerBoundaryCondition, | ||
Matrix2D< double > & | pc_upperBoundaryCondition, | ||
Matrix2D< double > & | alpha_lowerBoundaryCondition, | ||
Matrix2D< double > & | alpha_upperBoundaryCondition, | ||
string | L_lowerBoundaryCondition_calculationType, | ||
string | L_upperBoundaryCondition_calculationType, | ||
string | pc_lowerBoundaryCondition_calculationType, | ||
string | pc_upperBoundaryCondition_calculationType, | ||
string | alpha_lowerBoundaryCondition_calculationType, | ||
string | alpha_upperBoundaryCondition_calculationType, | ||
Matrix3D< double > & | DLL, | ||
Matrix3D< double > & | Dpcpc, | ||
Matrix3D< double > & | DpcpcLpp, | ||
Matrix3D< double > & | Daa, | ||
Matrix3D< double > & | DaaLpp, | ||
Matrix3D< double > & | Dpca, | ||
Matrix3D< double > & | DpcaLpp, | ||
Matrix3D< double > & | Jacobian, | ||
double | dt, | ||
double | Lpp, | ||
double | tau, | ||
double | tauLpp | ||
) |
Create matrix form of the Fokker-Planck equation: matr_A * PSD(t+1) = matr_B * PSD(t) + matr_C
Calculation matr_A, matr_B, and matr_C, i.e. numerical approximation of the derivatives Takes boundary conditions and diffusion coefficients as an input. L, pc, alpha should be orthogonal for 3D!!!
Definition at line 302 of file MatrixSolver.cpp.
References AddBoundary(), VF::alc(), VF::bounce_time_new(), CalculationMatrix::change_ind, CalculationMatrix::index1d(), VC::pi, and SecondDerivativeApproximation_3D().
void jacobi | ( | int | m_size, |
double * | z, | ||
double * | w | ||
) |
Jacobi (or diagonal) preconditioner:
Definition at line 620 of file MatrixSolver.cpp.
References A_copy.
void SOR | ( | int | m_size, |
double * | z, | ||
double * | w | ||
) |
SOR preconditioner:
P = (D/rel + L) w = P^-1 * z; P * w = z (D/rel + L) * w = z
Definition at line 640 of file MatrixSolver.cpp.
References A_copy.
void for_norm_A | ( | int | n, |
double * | x, | ||
double * | b | ||
) |
Simple preconditioner
Definition at line 662 of file MatrixSolver.cpp.
void gmres_wrapout | ( | CalculationMatrix & | A, |
Matrix1D< double > & | B, | ||
Matrix1D< double > & | X, | ||
Parameters_structure::PSD::GMRES_parameters_structure | GMRES_parameters | ||
) |
GMRES inversion. Potentially with bugs.
A * X = B - equation
Definition at line 676 of file MatrixSolver.cpp.
References A_copy, B(), Output::echo(), for_norm_A(), gmres(), jacobi(), matvec(), VF::max(), Parameters_structure::PSD::GMRES_parameters_structure::preconditioner_type, Parameters_structure::PSD::GMRES_parameters_structure::SOL_i_max, Parameters_structure::PSD::GMRES_parameters_structure::SOL_max_iter_err, Parameters_structure::PSD::GMRES_parameters_structure::SOL_maxiter, SOR(), and Parameters_structure::PSD::GMRES_parameters_structure::use_normalization.
void gmres2_wrapout | ( | CalculationMatrix & | A, |
Matrix1D< double > & | B, | ||
Matrix1D< double > & | X, | ||
Parameters_structure::PSD::GMRES_parameters_structure | GMRES_parameters | ||
) |
Definition at line 871 of file MatrixSolver.cpp.
References B(), Output::echo(), Parameters_structure::PSD::GMRES_parameters_structure::preconditioner_type, Matrix1D< T >::size_x, Parameters_structure::PSD::GMRES_parameters_structure::SOL_i_max, Parameters_structure::PSD::GMRES_parameters_structure::SOL_max_iter_err, and Parameters_structure::PSD::GMRES_parameters_structure::SOL_maxiter.
void Lapack | ( | DiagMatrix & | A, |
Matrix1D< double > & | B, | ||
Matrix1D< double > & | X | ||
) |
Lapack inversion.
A * X = B - equation
Definition at line 953 of file MatrixSolver.cpp.
References Output::echo(), and VF::max().
void over_relaxation_diag | ( | DiagMatrix & | A, |
Matrix1D< double > & | B, | ||
Matrix1D< double > & | X, | ||
int | max_steps, | ||
double | EE | ||
) |
Over relaxation iteration method.
Definition at line 1046 of file MatrixSolver.cpp.
References Output::echo(), and VF::max().
void GetDerivativeVector | ( | string | derivativeType, |
int & | dL, | ||
int & | dPc, | ||
int & | dAlpha | ||
) |
Get change in indexes according to the derivatives direction.
Used in approximation. If it gets derivativeType == DT_ALPHA_left, for example, that means we have alpha-left-derivative, which is ( f(alpha) - f(alpha-1) ) / ( delta alpha ). So it returns dAlpha = -1, as a derivative direction vector.
Definition at line 1095 of file MatrixSolver.cpp.
void SecondDerivativeApproximation_3D | ( | CalculationMatrix & | matr_A, |
int | il, | ||
int | im, | ||
int | ia, | ||
string | FirstDerivative, | ||
string | SecondDerivative, | ||
Matrix3D< double > & | L, | ||
Matrix3D< double > & | pc, | ||
Matrix3D< double > & | alpha, | ||
Matrix3D< double > & | Dxx, | ||
Matrix3D< double > & | Jacobian, | ||
double | multiplicator | ||
) |
Second derivative approximation, returns coefficients to be putted into the model matrix.
Samarskiy, page 261
Returns coefficients to be put into model matrix for an approximation of a second derivative.
Definition at line 1122 of file MatrixSolver.cpp.
References GetDerivativeVector(), and CalculationMatrix::index1d().
double max2 | ( | double | a, |
double | b | ||
) |
Function returns maximum between a and b.
Definition at line 1190 of file MatrixSolver.cpp.
void mult_vect2 | ( | double * | af, |
double * | vf, | ||
double * | wf, | ||
int | nf | ||
) |
Multiplocation between vector and matrix.
Definition at line 1194 of file MatrixSolver.cpp.
void gauss_solve | ( | double * | a, |
double * | b, | ||
int | n | ||
) |
bool tridag | ( | double | a[], |
double | b[], | ||
double | c[], | ||
double | r[], | ||
double | u[], | ||
long | n | ||
) |
Solve the AU=R system of equations, where A - tridiagonal matrix nxn with diagonals a[], b[], c[].
Solver for a system of equations with 3-diagonal matrix.
Au = r Where A = diag(a, b, c),
a[] | - array, diagonal '-1' of the matrix |
b[] | - array, diagonal '0' of the matrix |
c[] | - array, diagonal '+1' of the matrix |
r[] | - array, r-vector |
u[] | - array, result |
n | - size of the matrix |
Definition at line 1256 of file MatrixSolver.cpp.
bool MakeModelMatrix_3D_KC | ( | CalculationMatrix & | matr_A, |
CalculationMatrix & | matr_B, | ||
CalculationMatrix & | matr_C, | ||
Matrix3D< double > & | L, | ||
Matrix3D< double > & | pc, | ||
Matrix3D< double > & | alpha, | ||
int | L_size, | ||
int | pc_size, | ||
int | alpha_size, | ||
Matrix2D< double > & | L_lowerBoundaryCondition, | ||
Matrix2D< double > & | L_upperBoundaryCondition, | ||
Matrix2D< double > & | pc_lowerBoundaryCondition, | ||
Matrix2D< double > & | pc_upperBoundaryCondition, | ||
Matrix2D< double > & | alpha_lowerBoundaryCondition, | ||
Matrix2D< double > & | alpha_upperBoundaryCondition, | ||
string | L_lowerBoundaryCondition_calculationType, | ||
string | L_upperBoundaryCondition_calculationType, | ||
string | pc_lowerBoundaryCondition_calculationType, | ||
string | pc_upperBoundaryCondition_calculationType, | ||
string | alpha_lowerBoundaryCondition_calculationType, | ||
string | alpha_upperBoundaryCondition_calculationType, | ||
Matrix3D< double > & | DLL, | ||
Matrix3D< double > & | Dpcpc, | ||
Matrix3D< double > & | DpcpcLpp, | ||
Matrix3D< double > & | Daa, | ||
Matrix3D< double > & | DaaLpp, | ||
Matrix3D< double > & | Dpca, | ||
Matrix3D< double > & | DpcaLpp, | ||
Matrix3D< double > & | Jacobian, | ||
double | dt, | ||
double | Lpp, | ||
double | tau, | ||
double | tauLpp | ||
) |
kckim test Create matrix form of the Fokker-Planck equation: matr_A * PSD(t+1) = matr_B * PSD(t) + matr_C
Calculation matr_A, matr_B, and matr_C, i.e. numerical approximation of the derivatives Takes boundary conditions and diffusion coefficients as an input. L, pc, alpha should be orthogonal for 3D!!!
Definition at line 1299 of file MatrixSolver.cpp.
References AddBoundary(), CalculationMatrix::change_ind, CalculationMatrix::index1d(), SecondDerivativeApproximation_3D_KC_Diagonal(), and SecondDerivativeApproximation_3D_KC_Mixed().
void SecondDerivativeApproximation_3D_KC_Diagonal | ( | CalculationMatrix & | matr_A, |
int | il, | ||
int | im, | ||
int | ia, | ||
string | FirstDerivative, | ||
string | SecondDerivative, | ||
Matrix3D< double > & | L, | ||
Matrix3D< double > & | pc, | ||
Matrix3D< double > & | alpha, | ||
Matrix3D< double > & | Dxx, | ||
Matrix3D< double > & | Jacobian, | ||
double | multiplicator | ||
) |
kckim test Second derivative approximation, returns coefficients to be putted into the model matrix.
Samarskiy, page 261
Returns coefficients to be put into model matrix for an approximation of a second derivative.
Definition at line 1657 of file MatrixSolver.cpp.
References GetDerivativeVector(), and CalculationMatrix::index1d().
void SecondDerivativeApproximation_3D_KC_Mixed | ( | CalculationMatrix & | matr_A, |
int | il, | ||
int | im, | ||
int | ia, | ||
string | FirstDerivative, | ||
string | SecondDerivative, | ||
Matrix3D< double > & | L, | ||
Matrix3D< double > & | pc, | ||
Matrix3D< double > & | alpha, | ||
Matrix3D< double > & | Dxx, | ||
Matrix3D< double > & | Jacobian, | ||
double | multiplicator | ||
) |
Definition at line 1720 of file MatrixSolver.cpp.
References GetDerivativeVector(), and CalculationMatrix::index1d().
|
static |
Definition at line 573 of file MatrixSolver.cpp.