NeoPZ
pztransfer.cpp
Go to the documentation of this file.
1 
6 #include "pztransfer.h"
7 #include "pzfmatrix.h"
8 #include <stdlib.h>
9 
10 using namespace std;
11 
12 template<class TVar>
14 fNTVarVar(0),fRowBlock(),fColBlock(),fColPosition(0),fNumberofColumnBlocks(0),fColumnBlockPosition(0),
15 fColumnBlockNumber(0),fColumnBlockLastUsed(0),fDoubleValues(0),fDoubleValLastUsed(0) {
16 }
17 
18 // the sparse matrix blocks are defined by row, col
19 template<class TVar>
20 TPZTransfer<TVar>::TPZTransfer(TPZBlock<TVar> &row, TPZBlock<TVar> &col,int nvar, int nrowblocks, int ncolblocks) :
22  TPZMatrix<TVar>(), fNTVarVar(nvar), fRowBlock(), fColBlock(),
26 {
27  SetBlocks(row,col,nvar,nrowblocks,ncolblocks);
28 }
29 
30 template<class TVar>
31 void TPZTransfer<TVar>::Print(const char *name,ostream &out,const MatrixOutputFormat form) const {
32  if(form == EFormatted) {
33  if(name) {
34  out << name << endl;
35  } else {
36  out << "TPZTransfer<TVar>::Print : ";
37  }
38  out << "rows : " << this->Rows() << " cols : " << this->Cols() << endl;
39  } else {
40  out << this->Rows() << ' ' << this->Cols() << endl;
41  }
42  int irb;
43  for(irb=0 ; irb < fRowBlock.MaxBlockSize(); irb++) {
44  if(form == EFormatted) {
45  out << "block row number : " << irb << endl;
46  out << "number of column blocks : " << fNumberofColumnBlocks[irb] << endl;
47  out << "position of first column block : " << fColPosition[irb] << endl;
48  }
49  int colpos = fColPosition[irb];
50  int numcolbl = fNumberofColumnBlocks[irb];
51  int icbcounter;
52  for(icbcounter=0; icbcounter < numcolbl; icbcounter++) {
53  int icb = fColumnBlockNumber[colpos+icbcounter];
54  TVar *locval = &fDoubleValues[fColumnBlockPosition[colpos+icbcounter]];
55  if(form == EFormatted) {
56  out << "column block counter : " << icbcounter <<
57  " column block number " << icb << endl;
58  }
59  int rownumber = fRowBlock.Position(irb);
60  int colnumber = fColBlock.Position(icb);
61  int rowsize = fRowBlock.Size(irb);
62  int colsize = fColBlock.Size(icb);
63  if(form == EFormatted) {
64  out << "row position number : " << rownumber <<
65  " column position number " << colnumber << endl;
66  out << "block sizes : row : " << rowsize << " col " << colsize << endl;
67  }
68  TPZFMatrix<TVar> loc(rowsize,colsize,locval,rowsize*colsize);
69  if(form == EFormatted) {
70  loc.Print(0,out,form);
71  } else {
72  int ir,ic;
73  for(ir=0; ir<rowsize; ir++) {
74  for(ic=0; ic<colsize; ic++) {
75  out << (rownumber+ir) << ' ' << (colnumber+ic) << ' ' << loc(ir,ic) << endl;
76  }
77  }
78  }
79  }
80  }
81 }
82 
83 //void TPZTransfer<TVar>::SetBlocks(TPZBlock<TVar> &row,TPZBlock<TVar> &col,int nvar, int nrowblocks, int ncolblocks){
84 template<class TVar>
85 void TPZTransfer<TVar>::SetBlocks(TPZBlock<TVar> &row,TPZBlock<TVar> &col,int nvar, int nrowblocks, int ncolblocks){
86  // this operation will reset the matrix to zero
87  // with no rows defined
88  fRowBlock = row;
89  fColBlock = col;
90  fRowBlock.SetNBlocks(nrowblocks);
91  fColBlock.SetNBlocks(ncolblocks);
92  fNTVarVar = nvar;
93  if(nvar!=1) {
94  int i;
95  for(i=0;i<nrowblocks;i++) fRowBlock.Set(i,row.Size(i)/nvar);
96  for(i=0;i<ncolblocks;i++) fColBlock.Set(i,col.Size(i)/nvar);
97  fRowBlock.Resequence();
98  fColBlock.Resequence();
99  }
100  this->fRow = fRowBlock.Dim()*nvar;
101  this->fCol = fColBlock.Dim()*nvar;
103  fColPosition.Fill(-1,0);
110  fDoubleValLastUsed = 0;
111 }
112 
113 template<class TVar>
115  // returns 1 if the row is defined (i.e. has column entries)
116  if(fColPosition[row] == -1) return 0;
117  return 1;
118 }
119 
120 template<class TVar>
122  if(HasRowDefinition(row)) {
123  cout << "TPZTransfer:SetBlocks called for an already defined row = " <<
124  row << endl;
125  DebugStop();
126  }
128  fNumberofColumnBlocks[row] = colnumbers.NElements();
129  // will specify the sparsity pattern of row
130  ExpandColumnVectorEntries(colnumbers.NElements());
131  int ic= fColumnBlockLastUsed, lastic = ic+colnumbers.NElements();
132  for(;ic < lastic; ic++) {
133  fColumnBlockNumber[ic] = colnumbers[ic-fColumnBlockLastUsed];
134  }
135  fColumnBlockLastUsed = lastic;
136 }
137 
138 template<class TVar>
141  int nextsize = fColumnBlockNumber.NAlloc()+1000;
142  fColumnBlockNumber.Expand(nextsize);
143  fColumnBlockPosition.Expand(nextsize);
144  }
149 }
150 
151 
152 template<class TVar>
154  // sets the row,col block equal to matrix mat
155  // if row col was not specified by AddBlockNumbers, an error
156  // will be issued and exit
157  // find row,col
158  int colpos = fColPosition[row];
159  int numcolblocks = fNumberofColumnBlocks[row];
160  if(colpos == -1 || numcolblocks == -1) {
161  cout << "TPZTransfer<TVar>::SetBlockMatrix called for ilegal parameters : "
162  << " row = " << row << " col = " << col << " colpos = " << colpos <<
163  " numcolblocks = " << numcolblocks << endl;
164  DebugStop();
165  }
166  int ic = colpos, lastic = colpos+numcolblocks;
167  for(;ic<lastic;ic++) {
168  if(fColumnBlockNumber[ic] == col) break;
169  }
170  if(ic == lastic) {
171  cout << "TPZTransfer<TVar>::SetBlockMatrix column not found for row = " << row <<
172  " col = " << col << endl;
173  DebugStop();
174  }
175  int nblrows = fRowBlock.Size(row);
176  int nblcols = fColBlock.Size(col);
177  if(nblrows != mat.Rows() || nblcols != mat.Cols()) {
178  cout << "TPZTransfer<TVar>::SetBlockMatrix matrix has incompatible dimensions : "
179  " nblrows = " << nblrows << " nblcols = " << nblcols << " mat.rows = " <<
180  mat.Rows() << " mat.cols " << mat.Cols() << endl;
181  DebugStop();
182  }
183  ExpandDoubleValueEntries(nblrows*nblcols);
184  TPZFMatrix<TVar> bl(nblrows,nblcols,&fDoubleValues[fDoubleValLastUsed],nblrows*nblcols);
185  bl = mat;
187  fDoubleValLastUsed += nblrows*nblcols;
188 }
189 
190 template<class TVar>
194  // Print("After the double value expansion");
195  }
196  int nextsize = fDoubleValues.NElements()+num;
197  fDoubleValues.Resize(nextsize);
198 }
199 
200 
201 template<class TVar>
203  TVar alpha, TVar beta, int opt) const{
204  // multiplies the transfer matrix and puts the result in z
205  if ((!opt && this->Cols() != x.Rows()) || (opt && this->Rows() != x.Rows()))
206  this->Error( "TPZTransfer<TVar>::MultAdd <matrices with incompatible dimensions>" );
207  if(x.Cols() != y.Cols() || x.Cols() != z.Cols()) {
208  this->Error ("TPZTransfer<TVar>::MultiplyAdd incompatible dimensions\n");
209  }
210  int rows = fRowBlock.MaxBlockSize();
211  int xcols = x.Cols();
212  int ic, c, r;
213  this->PrepareZ(y,z,beta,opt);
214  if (fNTVarVar == 1) {
215  MultAddScalar(x, y, z, alpha, beta, opt);
216  }
217  else
218  {
219  int nrc = x.Rows();
220  int ncc = x.Cols();
221  int thisr = this->Rows()/fNTVarVar;
222  int thisc = this->Cols()/fNTVarVar;
223  for(int iv=0; iv<fNTVarVar; iv++) {
224  TPZFMatrix<TVar> tempcoarse(thisc,ncc), tempfine(thisr,ncc);
225  for (int i=0; i<thisr; i++) {
226  for (int c=0; c<ncc; c++) {
227  tempcoarse(i,c) = x.GetVal(iv+i*fNTVarVar,c);
228  }
229  }
230  MultAddScalar(tempcoarse, tempfine, tempfine, alpha, 0., opt);
231  for (int i=0; i<thisr; i++) {
232  for (int c=0; c<ncc; c++) {
233  z(iv+i*fNTVarVar,c) = tempfine.GetVal(i,c);
234  }
235  }
236  }
237  }
238 }
239 
240 template<class TVar>
242  TVar alpha, TVar beta, int opt) const{
243  // multiplies the transfer matrix and puts the result in z
244  if ((!opt && this->Cols() != x.Rows()) || (opt && this->Rows() != x.Rows()))
245  this->Error( "TPZTransfer<TVar>::MultAdd <matrices with incompatible dimensions>" );
246  if(x.Cols() != y.Cols() || x.Cols() != z.Cols()) {
247  this->Error ("TPZTransfer<TVar>::MultiplyAdd incompatible dimensions\n");
248  }
249  int rows = fRowBlock.MaxBlockSize();
250  int xcols = x.Cols();
251  int ic, c, r;
252  this->PrepareZ(y,z,beta,opt);
253  for (ic = 0; ic < xcols; ic++) {
254  if(!opt) {
255  for ( r=0; r < rows; r++) {
256  int rowblockpos = fRowBlock.Position(r);
257  int rowblocksize = fRowBlock.Size(r);
258  if(!rowblocksize) continue;
259  int colpos = fColPosition[r];
260  if(colpos == -1) continue;
261  int numcolbl = fNumberofColumnBlocks[r];
262  if(numcolbl == 0) continue;
263  TPZFMatrix<TVar> zloc(rowblocksize,1,&z(rowblockpos,ic),rowblocksize);
264  for( c=0; c < numcolbl; c++) {
265  int col = fColumnBlockNumber[colpos+c];
266  int colblockpos = fColBlock.Position(col);
267  int colblocksize = fColBlock.Size(col);
268  TPZFMatrix<TVar> xloc(colblocksize,1,&x.g(colblockpos,ic), colblocksize);
269  TPZFMatrix<TVar> aloc(rowblocksize,colblocksize,
270  &fDoubleValues[fColumnBlockPosition[colpos+c]],rowblocksize*colblocksize);
271  aloc.MultAdd(xloc,zloc,zloc,alpha,1.,opt);
272  }
273  }
274  } else {
275  for ( r=0; r < rows; r++) {
276  int rowblockpos = fRowBlock.Position(r);
277  int rowblocksize = fRowBlock.Size(r);
278  if(!rowblocksize) continue;
279  TPZFMatrix<TVar> xloc(rowblocksize,1,&x.g(rowblockpos,ic),rowblocksize);
280  int colpos = fColPosition[r];
281  if(colpos == -1) continue;
282  int numcolbl = fNumberofColumnBlocks[r];
283  if(numcolbl == 0) continue;
284  for( c=0; c < numcolbl; c++) {
285  int col = fColumnBlockNumber[colpos+c];
286  int colblockpos = fColBlock.Position(col);
287  int colblocksize = fColBlock.Size(col);
288  TPZFMatrix<TVar> zloc(colblocksize,1,&z(colblockpos,ic),colblocksize);
289  TPZFMatrix<TVar> aloc(rowblocksize,colblocksize,
290  &fDoubleValues[fColumnBlockPosition[colpos+c]],rowblocksize*colblocksize);
291  aloc.MultAdd(xloc,zloc,zloc,alpha,1.,opt);
292  col = rowblockpos-1;
293  rowblockpos = col+1;
294  }
295  }
296  }
297  }
298 }
299 
304 template<class TVar>
306  int iv;
307  int nrf = finesol.Rows();
308  int ncf = finesol.Cols();
309  int nrc = coarsesol.Rows();
310  int ncc = coarsesol.Cols();
311  int thisr = this->Rows();
312  int thisc = this->Cols();
313  if(nrf != this->Rows()*fNTVarVar || ncf != ncc) {
314  nrf = this->Rows()*fNTVarVar;
315  ncf = ncc;
316  finesol.Redim(nrf,ncf);
317  }
318  if (fNTVarVar == 1) {
319  MultiplyScalar(coarsesol,finesol,0);
320  }
321  else
322  {
323  for(iv=0; iv<fNTVarVar; iv++) {
324  TPZFMatrix<TVar> tempcoarse(thisc,ncc), tempfine(thisr,ncf);
325  for (int i=0; i<thisr; i++) {
326  for (int c=0; c<ncf; c++) {
327  tempcoarse(i,c) = coarsesol.GetVal(iv+i*fNTVarVar,c);
328 
329  }
330  }
331  MultiplyScalar(tempcoarse, tempfine, 0);
332  for (int i=0; i<thisr; i++) {
333  for (int c=0; c<ncf; c++) {
334  finesol(iv+i*fNTVarVar,c) = tempfine.GetVal(i,c);
335 
336  }
337  }
338  }
339  }
340 }
341 
346 template<class TVar>
348  int iv;
349  int nrf = fine.Rows();
350  int ncf = fine.Cols();
351  int nrc = coarse.Rows();
352  int ncc = coarse.Cols();
353  int thisr = this->Rows();
354  int thisc = this->Cols();
355  if(ncc != ncf || nrc != this->Cols()*fNTVarVar) {
356  ncc = ncf;
357  nrc = this->Cols()*fNTVarVar;
358  coarse.Redim(nrc,ncf);
359  }
360  if (fNTVarVar == 1) {
361  MultiplyScalar(fine,coarse,1);
362  }
363  else
364  {
365  for(iv=0; iv<fNTVarVar; iv++) {
366  TPZFMatrix<TVar> tempcoarse(thisr,ncc), tempfine(thisc,ncf);
367  for (int i=0; i<thisc; i++) {
368  for (int c=0; c<ncf; c++) {
369  tempfine(i,c) = fine.GetVal(iv+i*fNTVarVar,c);
370 
371  }
372  }
373  MultiplyScalar(tempfine, tempcoarse, 1);
374  for (int i=0; i<thisr; i++) {
375  for (int c=0; c<ncf; c++) {
376  coarse(iv+i*fNTVarVar,c) = tempcoarse(i,c);
377 
378  }
379  }
380  }
381  }
382 
383 }
384 
385 template<class TVar>
387  if ((opt==0 && this->Cols() != A.Rows()) || (opt ==1 && this->Rows() != A.Rows()))
388  {
389  this->Error( "TPZTransfer<TVar>::Multiply incompatible dimensions" );
390  }
391  if(!opt && (B.Rows()*fNTVarVar != this->Rows() || B.Cols() != A.Cols())) {
392  B.Redim(this->Rows()/fNTVarVar,A.Cols()/fNTVarVar);
393  }
394  else if (opt && (B.Rows()*fNTVarVar != this->Cols() || B.Cols() != A.Cols())) {
395  B.Redim(this->Cols()/fNTVarVar,A.Cols()/fNTVarVar);
396  }
397  MultAddScalar( A, B, B, 1.0, 0.0, opt);
398 }
399 
400 template class TPZTransfer<float>;
401 template class TPZTransfer<double>;
402 template class TPZTransfer<long double>;
403 
404 template class TPZTransfer<std::complex<double> >;
405 template class TPZTransfer<std::complex<float> >;
407 
TPZManVector< int > fColumnBlockPosition
Vector indicating the starting point of each column block.
Definition: pztransfer.h:125
int ClassId() const override
Define the class id associated with the class.
Definition: pztransfer.h:138
void SetBlockMatrix(int row, int col, TPZFMatrix< TVar > &mat)
Sets the row,col block equal to matrix mat if row col was not specified by AddBlockNumbers, an error will be issued and exit.
Definition: pztransfer.cpp:153
TPZManVector< TVar > fDoubleValues
Storage space for the matrix blocks.
Definition: pztransfer.h:131
static int Error(const char *msg, const char *msg2=0)
Returns error messages.
Definition: pzmatrix.cpp:1402
TPZBlock< TVar > fColBlock
Block sizes of the columns.
Definition: pztransfer.h:119
void AddBlockNumbers(int row, TPZVec< int > &colnumbers)
Will specify the sparsity pattern of row.
Definition: pztransfer.cpp:121
MatrixOutputFormat
Defines output format.
Definition: pzmatrix.h:55
void SetBlocks(TPZBlock< TVar > &row, TPZBlock< TVar > &col, int nvar, int nrowblocks, int ncolblocks)
This operation will reset the matrix to zero with no rows defined.
Definition: pztransfer.cpp:85
int64_t fRow
Number of rows in matrix.
Definition: pzmatrix.h:779
int fNTVarVar
Number of variables associated with each shape function.
Definition: pztransfer.h:113
void PrepareZ(const TPZFMatrix< TVar > &y, TPZFMatrix< TVar > &z, const TVar beta, const int opt) const
Is an auxiliar method used by MultiplyAdd.
Definition: pzmatrix.cpp:135
int fDoubleValLastUsed
Indicates the next free position of fDoubleValues.
Definition: pztransfer.h:133
virtual void Print(const char *name=NULL, std::ostream &out=std::cout, const MatrixOutputFormat form=EFormatted) const override
It prints the matrix data in a MatrixFormat Rows X Cols.
Definition: pztransfer.cpp:31
int NAlloc() const
Returns number of elements allocated for this object.
Definition: pzmanvector.h:92
void Expand(const int64_t newsize)
Expands the allocated storage to fit the newsize parameter.
Definition: pzmanvector.h:367
virtual void Resize(const int64_t newsize, const T &object)
Resizes the vector object.
Definition: pzmanvector.h:426
Contains the TPZTransfer class which implements a rectangular sparse block matrix.
void ExpandColumnVectorEntries(int numcol)
Increases the storage allocated int fColPosition to include numcol more values.
Definition: pztransfer.cpp:139
virtual void Resize(const int64_t newsize, const T &object)
Resizes the vector object reallocating the necessary storage, copying the existing objects to the new...
Definition: pzvec.h:373
TPZVec< int > fNumberofColumnBlocks
Vector indicating the number of column blocks associated with each row.
Definition: pztransfer.h:123
void MultiplyScalar(const TPZFMatrix< TVar > &A, TPZFMatrix< TVar > &B, int opt) const
Definition: pztransfer.cpp:386
int HasRowDefinition(int row)
Returns 1 if the row is defined (i.e. has column entries)
Definition: pztransfer.cpp:114
Contains TPZMatrixclass which implements full matrix (using column major representation).
#define DebugStop()
Returns a message to user put a breakpoint in.
Definition: pzerror.h:20
TPZTransfer()
Default constructor.
Definition: pztransfer.cpp:13
void TransferResidual(const TPZFMatrix< TVar > &fine, TPZFMatrix< TVar > &coarse)
Will transfer the residual, taking into acount there may be more than one TVar variable.
Definition: pztransfer.cpp:347
int64_t Rows() const
Returns number of rows.
Definition: pzmatrix.h:803
TPZVec< int > fColPosition
Vector indicating the starting column block for each row.
Definition: pztransfer.h:121
Full matrix class. Matrix.
Definition: pzfmatrix.h:32
int Redim(const int64_t newRows, const int64_t newCols) override
Redimension a matrix and ZERO your elements.
Definition: pzfmatrix.h:616
void TransferSolution(const TPZFMatrix< TVar > &coarsesol, TPZFMatrix< TVar > &finesol)
Will transfer the solution, taking into acount there may be more than one TVar variable.
Definition: pztransfer.cpp:305
int64_t fCol
Number of cols in matrix.
Definition: pzmatrix.h:781
int fColumnBlockLastUsed
Indicates the next free position.
Definition: pztransfer.h:129
Implements block matrices. Matrix utility.
TPZBlock< TVar > fRowBlock
Block sizes of the rows.
Definition: pztransfer.h:116
int Size(const int block_diagonal) const
Returns block dimension.
Definition: pzblock.h:171
void Fill(const T &copy, const int64_t from=0, const int64_t numelem=-1)
Will fill the elements of the vector with a copy object.
Definition: pzvec.h:460
virtual void MultAdd(const TPZFMatrix< TVar > &x, const TPZFMatrix< TVar > &y, TPZFMatrix< TVar > &z, const TVar alpha=1., const TVar beta=0., const int opt=0) const override
It computes z = beta * y + alpha * opt(this)*x but z and x can not overlap in memory.
Definition: pzfmatrix.cpp:757
int64_t Cols() const
Returns number of cols.
Definition: pzmatrix.h:809
void MultAdd(const TPZFMatrix< TVar > &x, const TPZFMatrix< TVar > &y, TPZFMatrix< TVar > &z, const TVar alpha, const TVar beta, const int opt=0) const override
Multiplies the transfer matrix and puts the result in z.
Definition: pztransfer.cpp:202
virtual void Print(std::ostream &out) const
Definition: pzmatrix.h:253
int64_t NElements() const
Returns the number of elements of the vector.
Definition: pzvec.h:190
int MaxBlockSize() const
Returns the max number of blocks on diagonal.
Definition: pzblock.h:163
void ExpandDoubleValueEntries(int numval)
Increases the storage space available in the fDoubleValues vector to include numval entries...
Definition: pztransfer.cpp:191
Implements rectangular matrix which extends a solution vector of the coarse mesh to a solution vector...
Definition: pzcmesh.h:33
TPZManVector< int > fColumnBlockNumber
Vector indicating the number of the column corresponding to the block.
Definition: pztransfer.h:127
const TVar & GetVal(const int64_t row, const int64_t col) const override
Get values without bounds checking This method is faster than "Get" if DEBUG is defined.
Definition: pzfmatrix.h:566
TVar & g(const int64_t row, const int64_t col) const
Definition: pzfmatrix.h:594
void MultAddScalar(const TPZFMatrix< TVar > &x, const TPZFMatrix< TVar > &y, TPZFMatrix< TVar > &z, const TVar alpha, const TVar beta, const int opt=0) const
Multiplies the transfer matrix and puts the result in z.
Definition: pztransfer.cpp:241
Root matrix class (abstract). Matrix.
Definition: pzmatrix.h:60