solver_ddpqiluc2.hpp 4.87 KB
Newer Older
Kirill Terekhov's avatar
Fixes  
Kirill Terekhov committed
1

Kirill Terekhov's avatar
Kirill Terekhov committed
2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19
#ifndef __SOLVER_DDPQILUC2__
#define __SOLVER_DDPQILUC2__
#include <iomanip>

#include "inmost_solver.h"
#include "solver_prototypes.hpp"
class ILUC_preconditioner : public Method
{
	typedef struct Interval_t
	{
		INMOST_DATA_ENUM_TYPE first, last;
		INMOST_DATA_ENUM_TYPE Size() { return last - first; }
	} Interval;
	typedef std::pair<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> coord;
	typedef std::pair<INMOST_DATA_REAL_TYPE, coord > wgt_coord;
	typedef std::vector< wgt_coord > wgt_coords;
	typedef struct row_col_t
	{
Kirill Terekhov's avatar
Kirill Terekhov committed
20
		Sparse::Row row, col;
Kirill Terekhov's avatar
Kirill Terekhov committed
21 22 23 24
		INMOST_DATA_REAL_TYPE diag;
	} row_col;
	typedef dynarray<INMOST_DATA_ENUM_TYPE,256> levels_t;
	//result of multilevel preconditioner
Kirill Terekhov's avatar
Kirill Terekhov committed
25
	//        |LDU  F |
Kirill Terekhov's avatar
Kirill Terekhov committed
26 27 28
	//  A  =  |       |
	//        | E   C |
	// LU decomposition is stored in skyline format with reversed direction
Kirill Terekhov's avatar
Kirill Terekhov committed
29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46
	//
	// For the next step of multilevel factorization we have to compute
	// S = C - (E U) D (L F)
	//
	// U is stored by rows
	// L is stored by columns
	// F is stored by rows
	// E is stored by rows
	// C is stored by rows
	//
	// During LF calculation F is traversed transposed (by columns) and
	// each column is solved with L. 
	// LF is obtained by columns.
	//
	// During EU calculating E is traversed by rows and
	// each row is solved with U. Then re
	//
	// For the faster multiplication of 
Kirill Terekhov's avatar
Kirill Terekhov committed
47
	std::vector<Sparse::Row::entry> LU_Entries, B_Entries;
Kirill Terekhov's avatar
Kirill Terekhov committed
48 49
	interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_REAL_TYPE> LU_Diag;
	interval<INMOST_DATA_ENUM_TYPE, Interval> U_Address, L_Address, B_Address;
Kirill Terekhov's avatar
Kirill Terekhov committed
50
	std::vector<Sparse::Row::entry> E_Entries, F_Entries;
Kirill Terekhov's avatar
Kirill Terekhov committed
51 52 53 54 55 56 57
	std::vector<interval<INMOST_DATA_ENUM_TYPE, Interval> *> E_Address;
	interval<INMOST_DATA_ENUM_TYPE, Interval> F_Address;
	interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_REAL_TYPE> temp; // temporal place for solve phase
	levels_t level_size; //remember size of each level
	std::vector<Interval> level_interval;
	//reordering information
	interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE > ddP,ddQ;
Kirill Terekhov's avatar
Kirill Terekhov committed
58 59
	INMOST_DATA_ENUM_TYPE reorder_nnz, ddpq_tau_adapt, estimator;
	INMOST_DATA_REAL_TYPE ddpq_tau, iluc2_tau;
Kirill Terekhov's avatar
Kirill Terekhov committed
60 61
	INMOST_DATA_REAL_TYPE tau, eps;
	INMOST_DATA_ENUM_TYPE sciters;
Kirill Terekhov's avatar
Kirill Terekhov committed
62
	Sparse::Matrix * Alink;
Kirill Terekhov's avatar
Kirill Terekhov committed
63 64
	Solver::OrderInfo * info;
	bool init;
Kirill Terekhov's avatar
Kirill Terekhov committed
65
	void DumpMatrix(interval<INMOST_DATA_ENUM_TYPE, Interval> & Address, 
Kirill Terekhov's avatar
Kirill Terekhov committed
66
									std::vector<Sparse::Row::entry> & Entries,
Kirill Terekhov's avatar
Kirill Terekhov committed
67 68
									INMOST_DATA_ENUM_TYPE wmbeg, INMOST_DATA_ENUM_TYPE wmend,
									std::string file_name);
Kirill Terekhov's avatar
Kirill Terekhov committed
69
	void SwapEntries(interval<INMOST_DATA_ENUM_TYPE, Interval> & Address, 
Kirill Terekhov's avatar
Kirill Terekhov committed
70
					 std::vector<Sparse::Row::entry> & Entries, 
Kirill Terekhov's avatar
Kirill Terekhov committed
71
					 INMOST_DATA_ENUM_TYPE rbeg, INMOST_DATA_ENUM_TYPE rend, 
Kirill Terekhov's avatar
Kirill Terekhov committed
72 73 74
					 INMOST_DATA_ENUM_TYPE k, INMOST_DATA_ENUM_TYPE j);
	void SwapLine(interval<INMOST_DATA_ENUM_TYPE, Interval> & Line, INMOST_DATA_ENUM_TYPE i, INMOST_DATA_ENUM_TYPE j);
	void SwapE(INMOST_DATA_ENUM_TYPE i, INMOST_DATA_ENUM_TYPE j);
Kirill Terekhov's avatar
Kirill Terekhov committed
75 76 77 78 79 80 81
	
	void ReorderEF(INMOST_DATA_ENUM_TYPE mobeg, 
					INMOST_DATA_ENUM_TYPE cbeg,
					INMOST_DATA_ENUM_TYPE cend,
					INMOST_DATA_ENUM_TYPE wend, 
					interval<INMOST_DATA_ENUM_TYPE, bool> & donePQ, 
					interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> & localP,
Kirill Terekhov's avatar
Kirill Terekhov committed
82
					interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> & localQ);
Kirill Terekhov's avatar
Kirill Terekhov committed
83 84 85 86 87
	void inversePQ(INMOST_DATA_ENUM_TYPE wbeg,
				   INMOST_DATA_ENUM_TYPE wend, 
				   interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> & localP,
				   interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> & localQ,
				   interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> & invP,
Kirill Terekhov's avatar
Kirill Terekhov committed
88
				   interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> & invQ);
Kirill Terekhov's avatar
Kirill Terekhov committed
89 90 91 92 93
	void applyPQ(INMOST_DATA_ENUM_TYPE wbeg,
				 INMOST_DATA_ENUM_TYPE wend,
				 interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> & localP,
				 interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> & localQ,
				 interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> & invP,
Kirill Terekhov's avatar
Kirill Terekhov committed
94
				 interval<INMOST_DATA_ENUM_TYPE, INMOST_DATA_ENUM_TYPE> & invQ);
Kirill Terekhov's avatar
Kirill Terekhov committed
95
public:
Kirill Terekhov's avatar
Kirill Terekhov committed
96 97 98 99 100 101 102 103 104 105
	INMOST_DATA_ENUM_TYPE & EnumParameter(std::string name);
	INMOST_DATA_REAL_TYPE & RealParameter(std::string name);
	void Copy(const Method * other);
	ILUC_preconditioner(const ILUC_preconditioner & other);
	ILUC_preconditioner & operator =(ILUC_preconditioner const & other);
	ILUC_preconditioner(Solver::OrderInfo & info);
	bool isInitialized();
	bool isFinalized();
	bool Initialize();
	bool Finalize();
Kirill Terekhov's avatar
Kirill Terekhov committed
106 107 108 109 110 111 112 113 114
	void Multiply(int level, Sparse::Vector & input, Sparse::Vector & output);
	void ApplyLU(int level, Sparse::Vector & input,Sparse::Vector & output);
	void ApplyB(int level, double alpha, Sparse::Vector & x, double beta, Sparse::Vector & y);
	int Descend(int level, Sparse::Vector & input, Sparse::Vector & output);
	int Ascend(int level, Sparse::Vector & input, Sparse::Vector & output);
	bool Solve(Sparse::Vector & input, Sparse::Vector & output);
	bool ReplaceMAT(Sparse::Matrix & A);
	bool ReplaceSOL(Sparse::Vector & x);
	bool ReplaceRHS(Sparse::Vector & b);
Kirill Terekhov's avatar
Kirill Terekhov committed
115 116
	Method * Duplicate();
	~ILUC_preconditioner();
Kirill Terekhov's avatar
Kirill Terekhov committed
117 118 119
};

#endif //__SOLVER_DDPQILUC2__