1 #include "RBGen_ISVDSingle.h"
3 #include "Epetra_Comm.h"
18 const int oldNumCols = numProc_;
19 const int newNumCols = A_->NumVectors();
21 "RBGen::ISVSingle::updateBasis(): number of snapshots to process has exceeded specified maximum humber of columns.");
22 while (numProc_ < oldNumCols+newNumCols) {
38 lup = (lup < lmin_ ? lmin_ : lup);
39 lup = (lup > oldNumCols+newNumCols - numProc_ ? oldNumCols+newNumCols - numProc_ : lup);
40 lup = (lup > lmax_ ? lmax_ : lup);
41 lup = (lup > maxBasisSize_ - curRank_ ? maxBasisSize_ - curRank_ : lup);
52 Aplus = Teuchos::null;
60 for (
int i=0; i<curRank_; i++) {
66 if (comm->
MyPID() == 0 && verbLevel_ >= 1) {
68 <<
"------------- ISVDSingle::updateBasis() -----------" << std::endl
69 <<
"| Current rank: " << curRank_ << std::endl
70 <<
"| Current sigmas: " << std::endl;
71 for (
int i=0; i<curRank_; i++) {
72 std::cout <<
"| " << sigma_[i] << std::endl;
79 void ISVDSingle::makePass() {
82 "RBGen::ISVDSingle::makePass(): Max Num Passes should be 1, but is not.");
84 if (curNumPasses_ > 0)
return;
85 const int numCols = A_->NumVectors();
86 while (numProc_ < numCols) {
102 lup = (lup < lmin_ ? lmin_ : lup);
103 lup = (lup > numCols - numProc_ ? numCols - numProc_ : lup);
104 lup = (lup > lmax_ ? lmax_ : lup);
105 lup = (lup > maxBasisSize_ - curRank_ ? maxBasisSize_ - curRank_ : lup);
114 Unew = Teuchos::null;
115 Aplus = Teuchos::null;
138 int info = ATUlcl.Multiply(
'T',
'N',1.0,*A_,Ulcl,0.0);
140 "RBGen::ISVDMultiCD::makePass(): Error calling Epetra_MultiVector::Multiply for A^T U.");
143 for (
int i=0; i<curRank_; i++) {
147 info = ATUlcl.Multiply(
'N',
'N',-1.0,Vlcl,S,1.0);
149 "RBGen::ISVDMultiCD::computeBasis(): Error calling Epetra_MultiVector::Multiply for V S.");
150 resNorms_.resize(curRank_);
151 ATUlcl.Norm2(&resNorms_[0]);
153 for (
int i=0; i<curRank_; i++) {
154 if (sigma_[i] != 0.0) {
155 resNorms_[i] /= sigma_[i];
163 std::vector<double> errnorms(curRank_);
169 curU(::View,*U_,0,curRank_),
170 curV(::View,*V_,0,curRank_);
174 for (
int i=0; i<curRank_; i++) {
175 curS[i][i] = sigma_[i];
177 info = work.Multiply(
'N',
'N',1.0,curU,curS,0.0);
179 "RBGen::ISVDMultiCD::makePass(): Error calling Epetra_MultiVector::Multiply() for debugging U S.");
180 info = work.Multiply(
'N',
'N',-1.0,*A_,curV,1.0);
182 "RBGen::ISVDMultiCD::makePass(): Error calling Epetra_MultiVector::Multiply() for debugging U S - A V.");
183 work.Norm2(&errnorms[0]);
184 for (
int i=0; i<curRank_; i++) {
185 if (sigma_[i] != 0.0) {
186 errnorms[i] /= sigma_[i];
197 if (comm->
MyPID() == 0 && verbLevel_ >= 1) {
199 <<
"------------- ISVDSingle::makePass() -----------" << std::endl
200 <<
"| Number of passes: " << curNumPasses_ << std::endl
201 <<
"| Current rank: " << curRank_ << std::endl
202 <<
"| Current sigmas: " << std::endl;
203 for (
int i=0; i<curRank_; i++) {
204 std::cout <<
"| " << sigma_[i] << std::endl;
207 std::cout <<
"|DBG US-AV norms: " << std::endl;
208 for (
int i=0; i<curRank_; i++) {
209 std::cout <<
"|DBG " << errnorms[i] << std::endl;
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
ISVDSingle()
Default constructor.
void updateBasis(const Teuchos::RCP< Epetra_MultiVector > &update_ss)
Update the current basis by appending new snapshots.
virtual int MyPID() const =0
void Initialize(const Teuchos::RCP< Teuchos::ParameterList > ¶ms, const Teuchos::RCP< const Epetra_MultiVector > &init, const Teuchos::RCP< RBGen::FileIOHandler< Epetra_Operator > > &fileio=Teuchos::null)
Initialize the method with the given parameter list and snapshot set.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)