Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2016-2018 The VES code team
3 : (see the PEOPLE-VES file at the root of this folder for a list of names)
4 :
5 : See http://www.ves-code.org for more information.
6 :
7 : This file is part of VES code module.
8 :
9 : The VES code module is free software: you can redistribute it and/or modify
10 : it under the terms of the GNU Lesser General Public License as published by
11 : the Free Software Foundation, either version 3 of the License, or
12 : (at your option) any later version.
13 :
14 : The VES code module is distributed in the hope that it will be useful,
15 : but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 : GNU Lesser General Public License for more details.
18 :
19 : You should have received a copy of the GNU Lesser General Public License
20 : along with the VES code module. If not, see <http://www.gnu.org/licenses/>.
21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 :
23 : #include "CoeffsMatrix.h"
24 : #include "CoeffsVector.h"
25 : #include "BasisFunctions.h"
26 :
27 : #include "tools/Tools.h"
28 : #include "core/Value.h"
29 : #include "tools/File.h"
30 : #include "tools/Exception.h"
31 : #include "tools/Random.h"
32 : #include "tools/Communicator.h"
33 :
34 : #include <vector>
35 : #include <cmath>
36 : #include <iostream>
37 : #include <sstream>
38 : #include <cstdio>
39 : #include <cfloat>
40 :
41 :
42 : namespace PLMD {
43 : namespace ves {
44 :
45 0 : CoeffsMatrix::CoeffsMatrix(
46 : const std::string& label,
47 : const std::vector<std::string>& dimension_labels,
48 : const std::vector<unsigned int>& indices_shape,
49 : Communicator& cc,
50 : const bool diagonal,
51 0 : const bool use_iteration_counter):
52 : CoeffsBase(label,dimension_labels,indices_shape,use_iteration_counter),
53 : data(0),
54 : size_(0),
55 : nrows_(0),
56 : ncolumns_(0),
57 : diagonal_(diagonal),
58 : averaging_counter(0),
59 : averaging_exp_decay_(0),
60 0 : mycomm(cc)
61 : {
62 0 : setupMatrix();
63 0 : }
64 :
65 :
66 0 : CoeffsMatrix::CoeffsMatrix(
67 : const std::string& label,
68 : std::vector<Value*>& args,
69 : std::vector<BasisFunctions*>& basisf,
70 : Communicator& cc,
71 : const bool diagonal,
72 0 : const bool use_iteration_counter):
73 : CoeffsBase(label,args,basisf,use_iteration_counter),
74 : data(0),
75 : size_(0),
76 : nrows_(0),
77 : ncolumns_(0),
78 : diagonal_(diagonal),
79 : averaging_counter(0),
80 : averaging_exp_decay_(0),
81 0 : mycomm(cc)
82 : {
83 0 : setupMatrix();
84 0 : }
85 :
86 :
87 0 : CoeffsMatrix::CoeffsMatrix(
88 : const std::string& label,
89 : std::vector<std::vector<Value*> >& argsv,
90 : std::vector<std::vector<BasisFunctions*> >& basisfv,
91 : Communicator& cc,
92 : const bool diagonal,
93 : const bool use_iteration_counter,
94 0 : const std::string& multicoeffs_label):
95 : CoeffsBase(label,argsv,basisfv,use_iteration_counter,multicoeffs_label),
96 : data(0),
97 : size_(0),
98 : nrows_(0),
99 : ncolumns_(0),
100 : diagonal_(diagonal),
101 : averaging_counter(0),
102 : averaging_exp_decay_(0),
103 0 : mycomm(cc)
104 : {
105 0 : setupMatrix();
106 0 : }
107 :
108 :
109 159 : CoeffsMatrix::CoeffsMatrix(
110 : const std::string& label,
111 : CoeffsVector* coeffsVec,
112 : Communicator& cc,
113 159 : const bool diagonal):
114 : CoeffsBase( *(static_cast<CoeffsBase*>(coeffsVec)) ),
115 : data(0),
116 : size_(0),
117 : nrows_(0),
118 : ncolumns_(0),
119 : diagonal_(diagonal),
120 : averaging_counter(0),
121 : averaging_exp_decay_(0),
122 159 : mycomm(cc)
123 : {
124 159 : setLabels(label);
125 159 : setupMatrix();
126 159 : }
127 :
128 :
129 318 : CoeffsMatrix::~CoeffsMatrix() {}
130 :
131 :
132 159 : void CoeffsMatrix::setupMatrix() {
133 159 : nrows_=numberOfCoeffs();
134 159 : ncolumns_=nrows_;
135 159 : if(diagonal_) {
136 159 : size_=nrows_;
137 : }
138 : else {
139 0 : size_=(nrows_*nrows_-nrows_)/2+nrows_;
140 : }
141 159 : clear();
142 159 : }
143 :
144 :
145 408 : size_t CoeffsMatrix::getSize() const {
146 408 : return size_;
147 : }
148 :
149 :
150 1380 : bool CoeffsMatrix::isDiagonal() const {
151 1380 : return diagonal_;
152 : }
153 :
154 :
155 0 : bool CoeffsMatrix::sameShape(CoeffsVector& coeffsvector_in) const {
156 0 : return CoeffsBase::sameShape( *(static_cast<CoeffsBase*>(&coeffsvector_in)) );
157 : }
158 :
159 :
160 0 : bool CoeffsMatrix::sameShape(CoeffsMatrix& coeffsmat_in) const {
161 0 : return CoeffsBase::sameShape( *(static_cast<CoeffsBase*>(&coeffsmat_in)) );
162 : }
163 :
164 :
165 0 : bool CoeffsMatrix::sameShape(CoeffsMatrix& coeffsmat0, CoeffsMatrix& coeffsmat1) {
166 0 : return coeffsmat0.sameShape(coeffsmat1);
167 : }
168 :
169 :
170 0 : bool CoeffsMatrix::sameShape(CoeffsVector& coeffsvec, CoeffsMatrix& coeffsmat) {
171 0 : return coeffsmat.sameShape(coeffsvec);
172 : }
173 :
174 :
175 0 : bool CoeffsMatrix::sameShape(CoeffsMatrix& coeffsmat, CoeffsVector& coeffsvec) {
176 0 : return coeffsmat.sameShape(coeffsvec);
177 : }
178 :
179 :
180 0 : void CoeffsMatrix::sumCommMPI() {
181 0 : mycomm.Sum(data);
182 0 : }
183 :
184 :
185 0 : void CoeffsMatrix::sumCommMPI(Communicator& cc) {
186 0 : cc.Sum(data);
187 0 : }
188 :
189 :
190 0 : void CoeffsMatrix::sumMultiSimCommMPI(Communicator& multi_sim_cc) {
191 0 : if(mycomm.Get_rank()==0) {
192 0 : double nwalkers = static_cast<double>(multi_sim_cc.Get_size());
193 0 : multi_sim_cc.Sum(data);
194 0 : scaleAllValues(1.0/nwalkers);
195 : }
196 0 : mycomm.Bcast(data,0);
197 0 : }
198 :
199 :
200 112554 : size_t CoeffsMatrix::getMatrixIndex(const size_t index1, const size_t index2) const {
201 : size_t matrix_idx;
202 : plumed_dbg_assert(index1<nrows_);
203 : plumed_dbg_assert(index2<ncolumns_);
204 112554 : if(diagonal_) {
205 : // plumed_massert(index1==index2,"CoeffsMatrix: you trying to access a off-diagonal element of a diagonal coeffs matrix");
206 : matrix_idx=index1;
207 : }
208 0 : else if (index1<=index2) {
209 0 : matrix_idx=index2+index1*(nrows_-1)-index1*(index1-1)/2;
210 : }
211 : else {
212 0 : matrix_idx=index1+index2*(nrows_-1)-index2*(index2-1)/2;
213 : }
214 112554 : return matrix_idx;
215 : }
216 :
217 :
218 249 : void CoeffsMatrix::clear() {
219 249 : data.resize(getSize());
220 33735 : for(size_t i=0; i<data.size(); i++) {
221 16743 : data[i]=0.0;
222 : }
223 249 : }
224 :
225 :
226 0 : void CoeffsMatrix::setAllValuesToZero() {
227 0 : for(size_t i=0; i<data.size(); i++) {
228 0 : data[i]=0.0;
229 : }
230 0 : }
231 :
232 :
233 27250 : double CoeffsMatrix::getValue(const size_t index1, const size_t index2) const {
234 54500 : return data[getMatrixIndex(index1,index2)];
235 : }
236 :
237 :
238 0 : double CoeffsMatrix::getValue(const std::vector<unsigned int>& indices1, const std::vector<unsigned int>& indices2) const {
239 0 : return getValue(getIndex(indices1),getIndex(indices2));
240 : }
241 :
242 :
243 0 : void CoeffsMatrix::setValue(const size_t index1, const size_t index2, const double value) {
244 0 : data[getMatrixIndex(index1,index2)]=value;
245 0 : }
246 :
247 :
248 0 : void CoeffsMatrix::setValue(const std::vector<unsigned int>& indices1, const std::vector<unsigned int>& indices2, const double value) {
249 0 : setValue(getIndex(indices1),getIndex(indices2),value);
250 0 : }
251 :
252 :
253 0 : double& CoeffsMatrix::operator()(const size_t index1, const size_t index2) {
254 0 : return data[getMatrixIndex(index1,index2)];
255 : }
256 :
257 :
258 20350 : const double& CoeffsMatrix::operator()(const size_t index1, const size_t index2) const {
259 40700 : return data[getMatrixIndex(index1,index2)];
260 : }
261 :
262 :
263 0 : double& CoeffsMatrix::operator()(const std::vector<unsigned int>& indices1, const std::vector<unsigned int>& indices2) {
264 0 : return data[getMatrixIndex(getIndex(indices1),getIndex(indices2))];
265 : }
266 :
267 :
268 0 : const double& CoeffsMatrix::operator()(const std::vector<unsigned int>& indices1, const std::vector<unsigned int>& indices2) const {
269 0 : return data[getMatrixIndex(getIndex(indices1),getIndex(indices2))];
270 : }
271 :
272 :
273 730 : CoeffsVector operator*(const CoeffsMatrix& coeffs_matrix, const CoeffsVector& coeffs_vector) {
274 730 : CoeffsVector new_coeffs_vector(coeffs_vector);
275 730 : new_coeffs_vector.clear();
276 730 : plumed_massert(coeffs_vector.numberOfCoeffs()==coeffs_matrix.numberOfCoeffs(),"CoeffsMatrix and CoeffsVector are of the wrong size");
277 : size_t numcoeffs = coeffs_vector.numberOfCoeffs();
278 730 : if(coeffs_matrix.isDiagonal()) {
279 41430 : for(size_t i=0; i<numcoeffs; i++) {
280 20350 : new_coeffs_vector(i) = coeffs_matrix(i,i)*coeffs_vector(i);
281 : }
282 : }
283 : else {
284 0 : for(size_t i=0; i<numcoeffs; i++) {
285 0 : for(size_t j=0; j<numcoeffs; j++) {
286 0 : new_coeffs_vector(i) += coeffs_matrix(i,j)*coeffs_vector(j);
287 : }
288 : }
289 : }
290 730 : return new_coeffs_vector;
291 : }
292 :
293 :
294 0 : void CoeffsMatrix::addToValue(const size_t index1, const size_t index2, const double value) {
295 0 : data[getMatrixIndex(index1,index2)]+=value;
296 0 : }
297 :
298 :
299 0 : void CoeffsMatrix::addToValue(const std::vector<unsigned int>& indices1, const std::vector<unsigned int>& indices2, const double value) {
300 0 : addToValue(getIndex(indices1),getIndex(indices2),value);
301 0 : }
302 :
303 :
304 770 : void CoeffsMatrix::scaleAllValues(const double scalef) {
305 61430 : for(size_t i=0; i<data.size(); i++) {
306 30330 : data[i]*=scalef;
307 : }
308 770 : }
309 :
310 :
311 770 : CoeffsMatrix& CoeffsMatrix::operator*=(const double scalef) {
312 770 : scaleAllValues(scalef);
313 770 : return *this;
314 : }
315 :
316 :
317 0 : CoeffsMatrix operator*(const double scalef, const CoeffsMatrix& coeffsmatrix) {
318 0 : return CoeffsMatrix(coeffsmatrix)*=scalef;
319 : }
320 :
321 :
322 0 : CoeffsMatrix operator*(const CoeffsMatrix& coeffsmatrix, const double scalef) {
323 0 : return scalef*coeffsmatrix;
324 : }
325 :
326 :
327 0 : CoeffsMatrix& CoeffsMatrix::operator*=(const CoeffsMatrix& other_coeffsmatrix) {
328 0 : plumed_massert(data.size()==other_coeffsmatrix.getSize(),"Coeffs matrices do not have the same size");
329 0 : for(size_t i=0; i<data.size(); i++) {
330 0 : data[i]*=other_coeffsmatrix.data[i];
331 : }
332 0 : return *this;
333 : }
334 :
335 :
336 0 : CoeffsMatrix CoeffsMatrix::operator*(const CoeffsMatrix& other_coeffsmatrix) const {
337 0 : return CoeffsMatrix(*this)*=other_coeffsmatrix;
338 : }
339 :
340 :
341 0 : void CoeffsMatrix::setValues(const double value) {
342 0 : for(size_t i=0; i<data.size(); i++) {
343 0 : data[i]=value;
344 : }
345 0 : }
346 :
347 :
348 770 : void CoeffsMatrix::setValues(const std::vector<double>& values) {
349 770 : plumed_massert( data.size()==values.size(), "Incorrect size");
350 61430 : for(size_t i=0; i<data.size(); i++) {
351 30330 : data[i]=values[i];
352 : }
353 770 : }
354 :
355 :
356 0 : void CoeffsMatrix::setValues(const CoeffsMatrix& other_coeffsmatrix) {
357 0 : plumed_massert( data.size()==other_coeffsmatrix.getSize(), "Incorrect size");
358 0 : for(size_t i=0; i<data.size(); i++) {
359 0 : data[i]=other_coeffsmatrix.data[i];
360 : }
361 0 : }
362 :
363 :
364 0 : CoeffsMatrix& CoeffsMatrix::operator=(const double value) {
365 0 : setValues(value);
366 0 : return *this;
367 : }
368 :
369 :
370 770 : CoeffsMatrix& CoeffsMatrix::operator=(const std::vector<double>& values) {
371 770 : setValues(values);
372 770 : return *this;
373 : }
374 :
375 :
376 : // CoeffsMatrix& CoeffsMatrix::operator=(const CoeffsMatrix& other_coeffsmatrix) {
377 : // setValues(other_coeffsmatrix);
378 : // return *this;
379 : // }
380 :
381 :
382 0 : CoeffsMatrix CoeffsMatrix::operator+() const {
383 0 : return *this;
384 : }
385 :
386 :
387 0 : CoeffsMatrix CoeffsMatrix::operator-() const {
388 0 : return CoeffsMatrix(*this)*=-1.0;
389 : }
390 :
391 :
392 0 : void CoeffsMatrix::addToValues(const double value) {
393 0 : for(size_t i=0; i<data.size(); i++) {
394 0 : data[i]+=value;
395 : }
396 0 : }
397 :
398 :
399 0 : void CoeffsMatrix::addToValues(const std::vector<double>& values) {
400 0 : plumed_massert( data.size()==values.size(), "Incorrect size");
401 0 : for(size_t i=0; i<data.size(); i++) {
402 0 : data[i]+=values[i];
403 : }
404 0 : }
405 :
406 :
407 0 : void CoeffsMatrix::addToValues(const CoeffsMatrix& other_coeffsmatrix) {
408 0 : plumed_massert( data.size()==other_coeffsmatrix.getSize(), "Incorrect size");
409 0 : for(size_t i=0; i<data.size(); i++) {
410 0 : data[i]+=other_coeffsmatrix.data[i];
411 : }
412 0 : }
413 :
414 :
415 0 : void CoeffsMatrix::subtractFromValues(const double value) {
416 0 : for(size_t i=0; i<data.size(); i++) {
417 0 : data[i]-=value;
418 : }
419 0 : }
420 :
421 :
422 0 : void CoeffsMatrix::subtractFromValues(const std::vector<double>& values) {
423 0 : plumed_massert( data.size()==values.size(), "Incorrect size");
424 0 : for(size_t i=0; i<data.size(); i++) {
425 0 : data[i]-=values[i];
426 : }
427 0 : }
428 :
429 :
430 0 : void CoeffsMatrix::subtractFromValues(const CoeffsMatrix& other_coeffsmatrix) {
431 0 : plumed_massert( data.size()==other_coeffsmatrix.getSize(), "Incorrect size");
432 0 : for(size_t i=0; i<data.size(); i++) {
433 0 : data[i]-=other_coeffsmatrix.data[i];
434 : }
435 0 : }
436 :
437 :
438 0 : CoeffsMatrix& CoeffsMatrix::operator+=(const double value) {
439 0 : addToValues(value);
440 0 : return *this;
441 : }
442 :
443 :
444 0 : CoeffsMatrix operator+(const double value, const CoeffsMatrix& coeffsmatrix) {
445 0 : return coeffsmatrix+value;
446 : }
447 :
448 :
449 0 : CoeffsMatrix operator+(const CoeffsMatrix& coeffsmatrix, const double value) {
450 0 : return CoeffsMatrix(coeffsmatrix)+=value;
451 : }
452 :
453 :
454 0 : CoeffsMatrix& CoeffsMatrix::operator+=(const std::vector<double>& values) {
455 0 : addToValues(values);
456 0 : return *this;
457 : }
458 :
459 :
460 0 : CoeffsMatrix operator+(const std::vector<double>& values, const CoeffsMatrix& coeffsmatrix) {
461 0 : return coeffsmatrix+values;
462 : }
463 :
464 :
465 0 : CoeffsMatrix operator+(const CoeffsMatrix& coeffsmatrix, const std::vector<double>& values) {
466 0 : return CoeffsMatrix(coeffsmatrix)+=values;
467 : }
468 :
469 :
470 0 : CoeffsMatrix& CoeffsMatrix::operator-=(const double value) {
471 0 : subtractFromValues(value);
472 0 : return *this;
473 : }
474 :
475 :
476 0 : CoeffsMatrix operator-(const double value, const CoeffsMatrix& coeffsmatrix) {
477 0 : return -1.0*coeffsmatrix+value;
478 : }
479 :
480 :
481 0 : CoeffsMatrix operator-(const CoeffsMatrix& coeffsmatrix, const double value) {
482 0 : return CoeffsMatrix(coeffsmatrix)-=value;
483 : }
484 :
485 :
486 0 : CoeffsMatrix& CoeffsMatrix::operator-=(const std::vector<double>& values) {
487 0 : subtractFromValues(values);
488 0 : return *this;
489 : }
490 :
491 :
492 0 : CoeffsMatrix operator-(const std::vector<double>& values, const CoeffsMatrix& coeffsmatrix) {
493 0 : return -1.0*coeffsmatrix+values;
494 : }
495 :
496 :
497 0 : CoeffsMatrix operator-(const CoeffsMatrix& coeffsmatrix, const std::vector<double>& values) {
498 0 : return CoeffsMatrix(coeffsmatrix)-=values;
499 : }
500 :
501 :
502 0 : CoeffsMatrix& CoeffsMatrix::operator+=(const CoeffsMatrix& other_coeffsmatrix) {
503 0 : addToValues(other_coeffsmatrix);
504 0 : return *this;
505 : }
506 :
507 :
508 0 : CoeffsMatrix CoeffsMatrix::operator+(const CoeffsMatrix& other_coeffsmatrix) const {
509 0 : return CoeffsMatrix(*this)+=other_coeffsmatrix;
510 : }
511 :
512 :
513 0 : CoeffsMatrix& CoeffsMatrix::operator-=(const CoeffsMatrix& other_coeffsmatrix) {
514 0 : subtractFromValues(other_coeffsmatrix);
515 0 : return *this;
516 : }
517 :
518 :
519 0 : CoeffsMatrix CoeffsMatrix::operator-(const CoeffsMatrix& other_coeffsmatrix) const {
520 0 : return CoeffsMatrix(*this)-=other_coeffsmatrix;
521 : }
522 :
523 :
524 0 : void CoeffsMatrix::averageMatrices(CoeffsMatrix& coeffsmat0, CoeffsMatrix& coeffsmat1) {
525 0 : plumed_massert(sameShape(coeffsmat0,coeffsmat1),"both CoeffsMatrix objects need to have the same shape");
526 0 : for(size_t i=0; i<coeffsmat0.getSize(); i++) {
527 0 : coeffsmat0.data[i] = coeffsmat1.data[i] = 0.5 * (coeffsmat0.data[i]+coeffsmat1.data[i]);
528 : }
529 0 : }
530 :
531 :
532 0 : void CoeffsMatrix::averageMatrices(const std::vector<CoeffsMatrix*>& coeffsmatSet) {
533 0 : double norm_factor = 1.0/static_cast<double>(coeffsmatSet.size());
534 0 : for(unsigned int k=1; k<coeffsmatSet.size(); k++) {
535 0 : plumed_massert(coeffsmatSet[0]->sameShape(*coeffsmatSet[k]),"All CoeffsMatrix objects need to have the same shape");
536 : }
537 0 : for(size_t i=0; i<coeffsmatSet[0]->getSize(); i++) {
538 : double value = 0.0;
539 0 : for(unsigned int k=0; k<coeffsmatSet.size(); k++) {
540 0 : value += coeffsmatSet[k]->data[i];
541 : }
542 0 : value *= norm_factor;
543 0 : for(unsigned int k=0; k<coeffsmatSet.size(); k++) {
544 0 : coeffsmatSet[k]->data[i] = value;
545 : }
546 : }
547 0 : }
548 :
549 :
550 :
551 0 : double CoeffsMatrix::getMinValue() const {
552 : double min_value=DBL_MAX;
553 0 : for(size_t i=0; i<data.size(); i++) {
554 0 : if(data[i]<min_value) {
555 : min_value=data[i];
556 : }
557 : }
558 0 : return min_value;
559 : }
560 :
561 :
562 0 : double CoeffsMatrix::getMaxValue() const {
563 : double max_value=DBL_MIN;
564 0 : for(size_t i=0; i<data.size(); i++) {
565 0 : if(data[i]>max_value) {
566 : max_value=data[i];
567 : }
568 : }
569 0 : return max_value;
570 : }
571 :
572 :
573 0 : void CoeffsMatrix::randomizeValuesGaussian(int randomSeed) {
574 0 : Random rnd;
575 0 : if (randomSeed<0) {randomSeed = -randomSeed;}
576 0 : rnd.setSeed(-randomSeed);
577 0 : for(size_t i=0; i<data.size(); i++) {
578 0 : data[i]=rnd.Gaussian();
579 : }
580 0 : }
581 :
582 :
583 0 : void CoeffsMatrix::resetAveraging() {
584 0 : clear();
585 : resetAveragingCounter();
586 0 : }
587 :
588 :
589 0 : void CoeffsMatrix::addToAverage(const CoeffsMatrix& coeffsmat) {
590 0 : plumed_massert( data.size()==coeffsmat.getSize(), "Incorrect size");
591 : //
592 0 : double aver_decay = 1.0 / ( static_cast<double>(averaging_counter) + 1.0 );
593 0 : if(averaging_exp_decay_>0 && (averaging_counter+1 > averaging_exp_decay_) ) {
594 0 : aver_decay = 1.0 / static_cast<double>(averaging_exp_decay_);
595 : }
596 : //
597 0 : for(size_t i=0; i<data.size(); i++) {
598 0 : data[i]+=(coeffsmat.data[i]-data[i])*aver_decay;
599 : }
600 0 : averaging_counter++;
601 0 : }
602 :
603 :
604 650 : void CoeffsMatrix::writeToFile(OFile& ofile) {
605 650 : writeHeaderToFile(ofile);
606 650 : writeDataToFile(ofile);
607 650 : }
608 :
609 :
610 0 : void CoeffsMatrix::writeToFile(const std::string& filepath, const bool append_file, Action* action_pntr) {
611 0 : OFile file;
612 0 : if(action_pntr!=NULL) {
613 0 : file.link(*action_pntr);
614 : }
615 : else {
616 0 : file.link(mycomm);
617 : }
618 0 : if(append_file) { file.enforceRestart(); }
619 0 : file.open(filepath);
620 0 : writeToFile(file);
621 0 : file.close();
622 0 : }
623 :
624 :
625 650 : void CoeffsMatrix::writeMatrixInfoToFile(OFile& ofile) {
626 650 : std::string field_diagonal = "diagonal_matrix";
627 650 : ofile.addConstantField(field_diagonal).printField(field_diagonal,isDiagonal());
628 650 : }
629 :
630 :
631 650 : void CoeffsMatrix::writeHeaderToFile(OFile& ofile) {
632 650 : ofile.clearFields();
633 650 : if(isIterationCounterActive()) {
634 650 : writeIterationCounterAndTimeToFile(ofile);
635 : }
636 650 : writeCoeffsInfoToFile(ofile);
637 650 : writeMatrixInfoToFile(ofile);
638 650 : }
639 :
640 :
641 650 : void CoeffsMatrix::writeDataToFile(OFile& ofile) {
642 650 : if(diagonal_) {
643 650 : writeDataDiagonalToFile(ofile);
644 : }
645 : else {
646 0 : writeDataFullToFile(ofile);
647 : }
648 650 : }
649 :
650 :
651 650 : void CoeffsMatrix::writeDataDiagonalToFile(OFile& ofile) {
652 : //
653 650 : std::string field_indices_prefix = "idx_";
654 : std::string field_coeffs = getDataLabel();
655 650 : std::string field_index = "index";
656 : //
657 650 : std::string int_fmt = "%8d";
658 650 : std::string str_separate = "#!-------------------";
659 : //
660 650 : char* s1 = new char[20];
661 650 : std::vector<unsigned int> indices(numberOfDimensions());
662 1300 : std::vector<std::string> ilabels(numberOfDimensions());
663 2350 : for(unsigned int k=0; k<numberOfDimensions(); k++) {
664 3400 : ilabels[k]=field_indices_prefix+getDimensionLabel(k);
665 : }
666 : //
667 55150 : for(size_t i=0; i<numberOfCoeffs(); i++) {
668 54500 : indices=getIndices(i);
669 126310 : for(unsigned int k=0; k<numberOfDimensions(); k++) {
670 99060 : sprintf(s1,int_fmt.c_str(),indices[k]);
671 148590 : ofile.printField(ilabels[k],s1);
672 : }
673 81750 : ofile.fmtField(" "+getOutputFmt()).printField(field_coeffs,getValue(i,i));
674 54500 : sprintf(s1,int_fmt.c_str(),i); ofile.printField(field_index,s1);
675 27250 : ofile.printField();
676 : }
677 650 : ofile.fmtField();
678 : // blank line between iterations to allow proper plotting with gnuplot
679 650 : ofile.printf("%s\n",str_separate.c_str());
680 650 : ofile.printf("\n");
681 650 : ofile.printf("\n");
682 650 : delete [] s1;
683 650 : }
684 :
685 :
686 0 : void CoeffsMatrix::writeDataFullToFile(OFile& ofile) {
687 : //
688 0 : std::string field_index_row = "idx_row";
689 0 : std::string field_index_column = "idx_column";
690 : std::string field_coeffs = getDataLabel();
691 : //
692 0 : std::string int_fmt = "%8d";
693 0 : std::string str_separate = "#!-------------------";
694 : //
695 0 : char* s1 = new char[20];
696 : //
697 0 : for(size_t i=0; i<nrows_; i++) {
698 0 : for(size_t j=0; j<ncolumns_; j++) {
699 : sprintf(s1,int_fmt.c_str(),i);
700 0 : ofile.printField(field_index_row,s1);
701 : sprintf(s1,int_fmt.c_str(),j);
702 0 : ofile.printField(field_index_column,s1);
703 0 : ofile.fmtField(" "+getOutputFmt()).printField(field_coeffs,getValue(i,j));
704 0 : ofile.printField();
705 : }
706 : }
707 0 : ofile.fmtField();
708 : // blank line between iterations to allow proper plotting with gnuplot
709 0 : ofile.printf("%s\n",str_separate.c_str());
710 0 : ofile.printf("\n");
711 0 : ofile.printf("\n");
712 0 : delete [] s1;
713 0 : }
714 :
715 :
716 : }
717 4839 : }
|
