Subversion Repositories gelsvn

Rev

Details | Last modification | View Log | RSS feed

Rev Author Line No. Line
464 jab 1 
/*
2 
 *  CSCMatrixBuilder.h
3 
 *  GEL
4 
 *	Added to simplify sequetnial building of sparse matrices in CSC format
5 
 *  At this version this doesn't check if we are outside matrix (no matrix size)
6 
 *  1. Insert the entries column by column (for symmetric matrices only up triangle) by using:
7 
 *				insert_entry(row_number,entry);
8 
 *				next_column()
9 
 *				next_column_nonsort() - can be used instead of next_column where the columns were also filled
10 
 *											sequentially (in their order)
11 
 *		You can also use LoadDense or LoadSymmetricDense if you have dense representation of the matrix in
12 
 *      LinAlg::CMatrix type.
13 
 *  2. Call get_Matrix() to get the representation of the matrix as a ARluSymMatrix<Type>
14 
 *		(you can pass this type directly to the symmetric problem constructor).
15 
 *  Created by Katarzyna Gebal on 18/11/08
16 
 *
17 
 */
18 
#ifndef __MESHEDIT_CSCMATRIXBUILDER_H__
19 
#define __MESHEDIT_CSCMATRIXBUILDER_H__
20 
 
21 
#include <vector>
22 
#include <LinAlg/Matrix.h>
23 
#include "arlsmat.h"
24 
 
25 
template<class T>
26 
class CSCMatrixBuilder
27 
{
28 
	typedef T Type;
29 
	std::vector<Type> vA;
30 
	std::vector<int> vpcol;
31 
	std::vector<int> virow;
32 
 
33 
	ARluSymMatrix<Type> mat;
34 
	int *irow;
35 
	int *pcol;
36 
	Type *A;
37 
	bool mat_created;
38 
 
39 
	int get_nnz()
40 
	{
41 
		return vA.size();
42 
	}
43 
 
44 
	int get_ncol()
45 
	{
46 
		return vpcol.size()-1;
47 
	}
48 
 
49 
	int* get_pcol()
50 
	{
51 
		int *pcol = new int[vpcol.size()];
52 
		std::vector<int>::iterator cit = vpcol.begin();
53 
		for(int i = 0; i < vpcol.size(); i++)
54 
		{
55 
			pcol[i] = *cit;
56 
			++cit;
57 
		}
58 
		return pcol;
59 
	}
60 
 
61 
	int* get_irow()
62 
	{
63 
		int nnz = vA.size();
64 
		std::vector<int>::iterator rit = virow.begin();
65 
		int *irow = new int[nnz];
66 
		for(int i = 0; i < nnz; i++)
67 
		{
68 
			irow[i] = *rit;
69 
			++rit;
70 
		}
71 
		return irow;
72 
	}
73 
 
74 
 
75 
	Type* get_A()
76 
	{
77 
		int nnz = vA.size();
78 
		Type *A = new Type[nnz];
79 
		for(int i = 0; i < nnz; i++)
80 
			A[i] = vA[i];
81 
		return A;
82 
	}
83 
 
84 
	void create_mat()
85 
	{
86 
		if(mat_created)
87 
			return;
88 
		mat_created = true;
89 
		irow = get_irow();
90 
		pcol = get_pcol();
91 
		A = get_A();
92 
		mat = ARluSymMatrix<Type>(get_ncol(), get_nnz(), A, irow, pcol);
93 
 
94 
	}
95 
public:
96 
	ARluSymMatrix<Type>& get_Matrix()
97 
	{
98 
		create_mat();
99 
		return mat;
100 
	}
101 
	CSCMatrixBuilder()
102 
	{
103 
		mat_created = false;
104 
		vpcol.push_back(0);
105 
	}
106 
	~CSCMatrixBuilder()
107 
	{
108 
		if(mat_created)
109 
		{
110 
			delete [] irow;
111 
			delete [] pcol;
112 
			delete [] A;
113 
		}
114 
	}
115 
 
116 
	void sort_entries(int id)
117 
	{
118 
		int start = vpcol[id];
119 
		int stop = vpcol[id+1];
120 
		Type dpom;
121 
		int ipom;
122 
		for(int i = start; i < stop; i++)
123 
			for(int j = i+1; j < stop; j++)
124 
			{
125 
				if(virow[i] > virow[j])
126 
				{
127 
					ipom = virow[i]; virow[i] = virow[j]; virow[j] = ipom;
128 
					dpom = vA[i]; vA[i] = vA[j]; vA[j] = dpom;
129 
				}
130 
			}
131 
	}
132 
 
133 
	void insert_entry(int row, Type entry)
134 
	{
135 
		vA.push_back(entry);
136 
		virow.push_back(row);
137 
	}
138 
 
139 
	void next_column_nonsort()
140 
	{
141 
		vpcol.push_back(vA.size());
142 
	}
143 
 
144 
	void next_column()
145 
	{
146 
		vpcol.push_back(vA.size());
147 
		sort_entries(vpcol.size()-2);
148 
	}
149 
 
150 
	void LoadSymmetricDense(LinAlg::CMatrix Q, Type eps)
151 
	{
152 
		 for(int i = 0; i < Q.Cols(); i++)
153 
		 {
154 
			 for(int j = i; j < Q.Rows(); j++)
155 
			 {
156 
				 if(Q[i][j] > eps || Q[i][j] < -eps)
157 
					 insert_entry(j,Q[i][j]);
158 
			 }
159 
			 next_column();
160 
		 }
161 
	}
162 
 
163 
	void LoadDense(LinAlg::CMatrix Q, Type eps)
164 
	{
165 
		 for(int i = 0; i < Q.Cols(); i++)
166 
		 {
167 
			 for(int j = 0; j < Q.Rows(); j++)
168 
			 {
169 
				 if(Q[i][j] > eps || Q[i][j] < -eps)
170 
					 insert_entry(j,Q[i][j]);
171 
			 }
172 
			 next_column();
173 
		 }
174 
	}
175 
 
176 
};
177 
 
178 
#endif