#include "MatrixSolver.h"
#include <math.h>
#include <malloc.h>
#include <string>
#include <ctime>
#include "../VariousFunctions/variousFunctions.h"
#include "../Exceptions/error.h"
#include "../Lapack/Include/cpplapack.h"
#include "../GMRES/itlin.h"
#include <iostream>
Go to the source code of this file.
Defines | |
#define | EE 1.e-19 |
#define | I(l, k) (l)*n+(k) |
Functions | |
void | matvec (int m_size, double *x, double *b) |
Vector to matrix multiplication for GMRES. | |
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. | |
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-diffusion matrix (tridiagonal). | |
void | AddBoundary (DiagMatrix &Matrix_A, string type, int in, int id1, double dh) |
Supportive sub-function to add boundary conditions to model matrix. | |
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. | |
void | jacobi (int m_size, double *z, double *w) |
Jacobi (or diagonal) preconditioner:. | |
void | SOR (int m_size, double *z, double *w) |
SOR preconditioner:. | |
void | for_norm_A (int n, double *x, double *b) |
Simple preconditioner. | |
void | gmres_wrapout (CalculationMatrix &A, Matrix1D< double > &B, Matrix1D< double > &X, int maxiter, int i_max, double max_err) |
GMRES inversion. | |
void | Lapack (DiagMatrix &A, Matrix1D< double > &B, Matrix1D< double > &X) |
Lapack inversion. | |
void | over_relaxation_diag (DiagMatrix &A, Matrix1D< double > &B, Matrix1D< double > &X, int max_steps, double EE) |
Over relaxation iteration method. | |
void | GetDerivativeVector (string derivativeType, int &dL, int &dPc, int &dAlpha) |
Get change in indexes according to the derivatives direction. | |
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. | |
double | max2 (double a, double b) |
Function returns maximum between a and b. | |
void | mult_vect2 (double *af, double *vf, double *wf, int nf) |
Multiplocation between vector and matrix. | |
void | gauss_solve (double *a, double *b, int n) |
Gauss inversion. | |
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[]. |
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 |
#define I | ( | l, | |||
k | ) | (l)*n+(k) |
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 267 of file MatrixSolver.cpp.
Referenced by MakeModelMatrix_3D().
void for_norm_A | ( | int | n, | |
double * | x, | |||
double * | b | |||
) |
void gauss_solve | ( | double * | a, | |
double * | b, | |||
int | n | |||
) |
Gauss inversion.
Definition at line 1030 of file MatrixSolver.cpp.
Referenced by PSD::Diffusion_pc_alpha().
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 933 of file MatrixSolver.cpp.
Referenced by SecondDerivativeApproximation_3D().
void gmres | ( | int | n, | |
double * | y, | |||
MATVEC * | matvec, | |||
PRECON * | preconr, | |||
PRECON * | preconl, | |||
double * | b, | |||
struct ITLIN_OPT * | opt, | |||
struct ITLIN_INFO * | info | |||
) |
Definition at line 300 of file gmres.c.
References CheckEachIter, CheckOnRestart, ITLIN_DATA::codeid, DATA, ITLIN_OPT::datafile, ITLIN_OPT::datalevel, DATALEVEL, ERROR, ITLIN_OPT::errorfile, ITLIN_OPT::errorlevel, ERRORLEVEL, False, ITLIN_DATA::Final, FITER, FMISC, FRES, ITLIN_DATA::GMRES, gmres_norm2(), gmres_qrfact(), gmres_qrsolv(), i, ITLIN_OPT::i_max, ITLIN_DATA::Initial, ITLIN_DATA::Intermediate, ITLIN_INFO::iter, ITLIN_OPT::iterfile, itlin_dataout(), itlin_noprecon(), itlin_parcheck_and_print(), matvec(), ITLIN_OPT::maxiter, MAXITER_DEFAULT, ITLIN_OPT::miscfile, ITLIN_DATA::mode, MONITOR, ITLIN_OPT::monitorfile, ITLIN_OPT::monitorlevel, MONITORLEVEL, ITLIN_INFO::nomatvec, ITLIN_INFO::noprecl, ITLIN_INFO::noprecr, ITLIN_DATA::normdx, ITLIN_INFO::precision, preconr(), RCODE, ITLIN_DATA::res, ITLIN_OPT::resfile, ITLIN_DATA::residuum, ITLIN_DATA::Solution, ITLIN_DATA::t, ITLIN_DATA::tau, ITLIN_OPT::termcheck, ITLIN_OPT::tol, True, zibnum_fwalloc(), and zibnum_pfwalloc().
Referenced by gmres_wrapout().
double gmres_norm2 | ( | int | n, | |
double * | v | |||
) |
void gmres_wrapout | ( | CalculationMatrix & | A, | |
Matrix1D< double > & | B, | |||
Matrix1D< double > & | X, | |||
int | maxiter, | |||
int | i_max, | |||
double | max_err | |||
) |
GMRES inversion.
Potentially with bugs.
A * X = B - equation
Definition at line 656 of file MatrixSolver.cpp.
References CheckOnRestart, ITLIN_OPT::datafile, ITLIN_OPT::datalevel, Output::echo(), ITLIN_OPT::errorfile, ITLIN_OPT::errorlevel, gmres(), ITLIN_OPT::i_max, ITLIN_INFO::iter, ITLIN_OPT::iterfile, matvec(), VF::max(), ITLIN_OPT::maxiter, ITLIN_OPT::miscfile, ITLIN_OPT::monitorfile, ITLIN_OPT::monitorlevel, ITLIN_INFO::nomatvec, None, ITLIN_INFO::noprecl, ITLIN_INFO::noprecr, ITLIN_INFO::precision, ITLIN_INFO::rcode, ITLIN_OPT::resfile, ITLIN_OPT::termcheck, ITLIN_OPT::tol, and Verbose.
Referenced by PSD::Diffusion_pc_alpha().
void jacobi | ( | int | m_size, | |
double * | z, | |||
double * | w | |||
) |
void Lapack | ( | DiagMatrix & | A, | |
Matrix1D< double > & | B, | |||
Matrix1D< double > & | X | |||
) |
Lapack inversion.
A * X = B - equation
Definition at line 791 of file MatrixSolver.cpp.
References Output::echo(), i, and VF::max().
Referenced by PSD::Diffusion_pc_alpha().
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 56 of file MatrixSolver.cpp.
References B(), VF::G(), and i.
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 287 of file MatrixSolver.cpp.
References AddBoundary(), VF::alc(), VF::bounce_time(), CalculationMatrix::change_ind, CalculationMatrix::index1d(), VF::Kfunc(), and SecondDerivativeApproximation_3D().
Referenced by PSD::Diffusion_alpha(), PSD::Diffusion_L(), PSD::Diffusion_pc(), and PSD::Diffusion_pc_alpha().
void matvec | ( | int | n, | |
double * | x, | |||
double * | b | |||
) |
Vector to matrix multiplication for GMRES.
Definition at line 559 of file MatrixSolver.cpp.
Referenced by gmres(), and gmres_wrapout().
double max2 | ( | double | a, | |
double | b | |||
) |
void mult_vect2 | ( | double * | af, | |
double * | vf, | |||
double * | wf, | |||
int | nf | |||
) |
void over_relaxation_diag | ( | DiagMatrix & | A, | |
Matrix1D< double > & | B, | |||
Matrix1D< double > & | X, | |||
int | max_steps, | |||
double | EE | |||
) |
Over relaxation iteration method.
Definition at line 884 of file MatrixSolver.cpp.
References Output::echo(), i, and VF::max().
Referenced by PSD::Diffusion_pc_alpha().
void preconr | ( | int | n, | |
double * | x, | |||
double * | b | |||
) |
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 960 of file MatrixSolver.cpp.
References GetDerivativeVector(), and CalculationMatrix::index1d().
Referenced by MakeModelMatrix_3D().
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-diffusion matrix (tridiagonal).
Solver for 1d general equation.
Definition at line 226 of file MatrixSolver.cpp.
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 621 of file MatrixSolver.cpp.
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.
Definition at line 1080 of file MatrixSolver.cpp.
Referenced by PSD::Diffusion_alpha(), PSD::Diffusion_L(), PSD::Diffusion_pc(), SolveMatrix(), and steady_state().