OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
rowcol.C
Go to the documentation of this file.
1 /*
2  *
3  * ##### ##### ###### ###### ### ###
4  * ## ## ## ## ## ## ## ### ##
5  * ## ## ## ## #### #### ## # ##
6  * ## ## ## ## ## ## ## ##
7  * ## ## ## ## ## ## ## ##
8  * ##### ##### ## ###### ## ##
9  *
10  *
11  * OOFEM : Object Oriented Finite Element Code
12  *
13  * Copyright (C) 1993 - 2013 Borek Patzak
14  *
15  *
16  *
17  * Czech Technical University, Faculty of Civil Engineering,
18  * Department of Structural Mechanics, 166 29 Prague, Czech Republic
19  *
20  * This library is free software; you can redistribute it and/or
21  * modify it under the terms of the GNU Lesser General Public
22  * License as published by the Free Software Foundation; either
23  * version 2.1 of the License, or (at your option) any later version.
24  *
25  * This program is distributed in the hope that it will be useful,
26  * but WITHOUT ANY WARRANTY; without even the implied warranty of
27  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
28  * Lesser General Public License for more details.
29  *
30  * You should have received a copy of the GNU Lesser General Public
31  * License along with this library; if not, write to the Free Software
32  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
33  */
34 
35 #include "rowcol.h"
36 #include "floatarray.h"
37 #include "intarray.h"
38 #include "mathfem.h"
39 #include "error.h"
40 
41 #include <cstdlib>
42 
43 namespace oofem {
44 RowColumn :: RowColumn(int n, int st)
45 // Constructor. Creates a row-column with number n, starting at index st.
46 {
47  int size;
48 
49  number = n;
50  start = st;
51  diag = 0.;
52  size = number - start;
53  if ( size ) {
54  row = ( double * ) calloc( size, sizeof( double ) );
55  column = ( double * ) calloc( size, sizeof( double ) );
56  if ( !row || !column ) {
57  OOFEM_ERROR("Failed to allocate memory");
58  }
59  } else {
60  row = NULL;
61  column = NULL;
62  }
63 }
64 
65 
67 // Destructor.
68 {
69  free(row);
70  free(column);
71 }
72 
73 
74 #ifdef DEBUG
75 double &RowColumn :: atU(int i)
76 // Returns the i-th coefficient of the column of the receiver. Slow but
77 // safe.
78 {
79  this->checkBounds(i);
80  return column [ i - start ];
81 }
82 #endif
83 
84 
85 #ifdef DEBUG
86 double &RowColumn :: atL(int i)
87 // Returns the i-th coefficient of the row of the receiver. Slow but safe.
88 {
89  this->checkBounds(i);
90  return row [ i - start ];
91 }
92 #endif
93 
94 
96 // Checks that the receiver (k) possesses coefficients (i,k) and (k,i).
97 {
98  if ( i < start || i >= number ) {
99  OOFEM_ERROR("error in bounds of RowColumn %d (start=%d) : no entry %d", number, start, i);
100  }
101 }
102 
103 
105 // If loc points to a coefficient which out of the receiver's range,
106 // enlarges the receiver.
107 {
108  int n, i, j, c, first = start;
109 
110  n = loc.giveSize();
111  if ( !n ) {
112  OOFEM_ERROR("0-size location array");
113  }
114 
115  // get first non-zero coefficient in loc
116  for ( i = 1; i <= n; i++ ) {
117  if ( ( first = loc.at(i) ) ) {
118  break;
119  }
120  }
121 
122  // get min coefficient in loc
123  for ( j = i + 1; j <= n; j++ ) {
124  if ( ( c = loc.at(j) ) ) {
125  first = min(c, first);
126  }
127  }
128 
129  // enlarge the receiver if necessary
130  if ( first < start && first ) {
131  this->growTo(first);
132  }
133 }
134 
135 
136 double RowColumn :: dot(const FloatArray &b, char c, int first, int last)
137 // Returns the dot product of the row (c=='R') or the column (c=='C') of
138 // the receiver with array 'b', from indices 'first' to 'last'.
139 {
140  double answer, *p1;
141  int i;
142 
143  answer = 0.;
144  i = last - first + 1;
145  if ( i > 0 ) {
146  if ( c == 'R' ) {
147  p1 = row + first - start;
148  } else {
149  p1 = column + first - start;
150  }
151 
152  // I do the bad practice, because it speed up debug runs significantly / Mikael
153  const double *p2 = b.givePointer() + first-1 ; // bad practice !
154  //int pos = first;
155  while ( i-- ) {
156  answer += *p1++ * *p2++ ;
157  //answer += * p1++ * b.at(pos++);
158  }
159  }
160 
161  return answer;
162 }
163 
164 
165 void
166 RowColumn :: growTo(int newStart)
167 // Enlarges the receiver (the column upwards, the row leftwards) to
168 // 'start' = newStart.
169 {
170  int newSize, size, i;
171  double *newRow, *newColumn, *p1, *p2;
172 
173  if ( newStart == start ) {
174  return;
175  }
176 
177 # ifdef DEBUG
178  if ( newStart <= 0 || newStart > start ) {
179  OOFEM_ERROR("cannot enlarge RowCol %d (start=%d) to %d", number, start, newStart);
180  }
181 
182 # endif
183 
184  newSize = number - newStart;
185  newRow = ( double * ) calloc( newSize, sizeof( double ) );
186  newColumn = ( double * ) calloc( newSize, sizeof( double ) );
187  if ( !newRow || !newColumn ) {
188  OOFEM_ERROR("Failed to allocate memory");
189  }
190  size = number - start;
191 
192  if ( size ) {
193  p1 = row;
194  p2 = newRow + ( start - newStart );
195  i = size;
196  while ( i-- ) {
197  * p2++ = * p1++;
198  }
199 
200  free(row);
201 
202  p1 = column;
203  p2 = newColumn + ( start - newStart );
204  i = size;
205  while ( i-- ) {
206  * p2++ = * p1++;
207  }
208 
209  free(column);
210  }
211 
212  row = newRow;
213  column = newColumn;
214  start = newStart;
215 }
216 
217 
219 // Prints the receiver on screen.
220 {
221  printf("Row-column %d : start = %d, diag = %.5f\n col : ",
222  number, start, diag);
223  for ( int i = 0; i < number - start; i++ ) {
224  printf(" % .5f", column [ i ]);
225  }
226 
227  printf("\n row : ");
228  for ( int i = 0; i < number - start; i++ ) {
229  printf(" % .5f", row [ i ]);
230  }
231 
232  printf("\n");
233 }
234 
235 
236 void
238 // Sets all coefficients of the receiver to 0.
239 {
240  int size = number - start;
241 
242  diag = 0.;
243  if ( size ) {
244  for ( int i = 0; i < size; i++ ) {
245  column [ i ] = 0.0;
246  row [ i ] = 0.0;
247  }
248  }
249 }
250 
251 RowColumn *
253 {
254  int i;
255  double *newRow, *newColumn, *p1, *p2;
256  int size = number - start;
257 
258  if ( size ) {
259  newRow = ( double * ) malloc( size * sizeof( double ) );
260  newColumn = ( double * ) malloc( size * sizeof( double ) );
261  if ( !newRow || !newColumn ) {
262  OOFEM_ERROR("Failed to allocate memory");
263  }
264 
265  p1 = row;
266  p2 = newRow;
267  i = size;
268  while ( i-- ) {
269  * p2++ = * p1++;
270  }
271 
272  p1 = column;
273  p2 = newColumn;
274  i = size;
275  while ( i-- ) {
276  * p2++ = * p1++;
277  }
278  } else {
279  newRow = newColumn = NULL;
280  }
281 
282  return new RowColumn(this->number, this->start, newRow, newColumn, this->diag);
283 }
284 
285 
286 RowColumn :: RowColumn(int inumber, int istart, double *irow, double *icol, double idiag)
287 {
288  this->number = inumber;
289  this->start = istart;
290  this->row = irow;
291  this->column = icol;
292  this->diag = idiag;
293 }
294 } // end namespace oofem
double * column
Definition: rowcol.h:72
void checkSizeTowards(const IntArray &)
Definition: rowcol.C:104
void printYourself()
Definition: rowcol.C:218
double diag
Definition: rowcol.h:73
This class implements a segment of a unsymmetric matrix stored in segmented form (skyline).
Definition: rowcol.h:65
void checkBounds(int)
Definition: rowcol.C:95
Class implementing an array of integers.
Definition: intarray.h:61
int & at(int i)
Coefficient access function.
Definition: intarray.h:103
void zero()
Definition: rowcol.C:237
#define OOFEM_ERROR(...)
Definition: error.h:61
Class representing vector of real numbers.
Definition: floatarray.h:82
double dot(const FloatArray &, char, int, int)
Definition: rowcol.C:136
double * row
Definition: rowcol.h:71
const double * givePointer() const
Gives the pointer to the raw data, breaking encapsulation.
Definition: floatarray.h:469
double & atU(int i)
Definition: rowcol.h:83
double & atL(int i)
Definition: rowcol.h:84
int min(int i, int j)
Returns smaller value from two given decimals.
Definition: mathfem.h:59
int giveSize() const
Definition: intarray.h:203
void growTo(int)
Definition: rowcol.C:166
the oofem namespace is to define a context or scope in which all oofem names are defined.
RowColumn * GiveCopy()
Definition: rowcol.C:252
RowColumn(int, int)
Definition: rowcol.C:44

This page is part of the OOFEM documentation. Copyright (c) 2011 Borek Patzak
Project e-mail: info@oofem.org
Generated at Tue Jan 2 2018 20:07:31 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011