44 #include <Epetra_Operator.h> 45 #include <Epetra_CrsMatrix.h> 46 #include <Epetra_LinearProblem.h> 47 #include <Amesos_BaseSolver.h> 48 #include <Epetra_MultiVector.h> 50 #include "shylu_local_schur_operator.h" 54 ShyLU_Local_Schur_Operator::ShyLU_Local_Schur_Operator(
59 Epetra_LinearProblem *LP, Amesos_BaseSolver *solver,
60 Ifpack_Preconditioner *ifSolver, Epetra_CrsMatrix *C,
61 Epetra_Map *localDRowMap,
int nvectors)
71 localDRowMap_ = localDRowMap;
72 orig_lhs_ = ssym_->OrigLP->GetLHS();
73 orig_rhs_ = ssym_->OrigLP->GetRHS();
74 ResetTempVectors(nvectors);
82 matvec_time_ = RCP<Time>(
new Time(
"First 2 Matvec time in prober: "));
83 localize_time_ = RCP<Time>(
new Time(
"Localize time in prober: "));
84 trisolve_time_ = RCP<Time>(
new Time(
"Triangular solve time in prober: "));
85 dist_time_ = RCP<Time>(
new Time(
"Distribute time in prober: "));
86 matvec2_time_ = RCP<Time>(
new Time(
"Second 1 matvec time in prober: "));
87 apply_time_ = RCP<Time>(
new Time(
"Apply time in prober: "));
88 update_time_ = RCP<Time>(
new Time(
"Update time in prober: "));
92 int ShyLU_Local_Schur_Operator::SetUseTranspose(
bool useTranspose)
97 int ShyLU_Local_Schur_Operator::Apply(
const Epetra_MultiVector &X,
98 Epetra_MultiVector &Y)
const 101 apply_time_->start();
111 matvec_time_->start();
114 err = G_->Multiply(
false, X, *temp2);
116 err = C_->Multiply(
false, X, *temp);
120 matvec_time_->stop();
125 int nrows = C_->RowMap().NumMyElements();
126 cout <<
"DEBUG MODE" << endl;
127 ASSERT((nrows == localDRowMap_->NumGlobalElements()));
129 int gids[nrows], gids1[nrows];
130 C_->RowMap().MyGlobalElements(gids);
131 localDRowMap_->MyGlobalElements(gids1);
133 for (
int i = 0; i < nrows; i++)
135 ASSERT((gids[i] == gids1[i]));
140 localize_time_->start();
144 localize_time_->stop();
145 trisolve_time_->start();
149 if (config_->amesosForDiagonal)
152 ifSolver_->ApplyInverse(*temp, *localX);
155 trisolve_time_->stop();
161 matvec2_time_->start();
164 R_->Multiply(
false, *localX, Y);
167 matvec2_time_->stop();
168 update_time_->start();
171 err = Y.Update(1.0, *temp2, -1.0);
176 update_time_->stop();
184 void ShyLU_Local_Schur_Operator::PrintTimingInfo()
187 cout << matvec_time_->name()<< matvec_time_->totalElapsedTime() << endl;
188 cout << localize_time_->name()<< localize_time_->totalElapsedTime() << endl;
189 cout << trisolve_time_->name()<< trisolve_time_->totalElapsedTime() << endl;
190 cout << dist_time_->name() <<dist_time_->totalElapsedTime() << endl;
191 cout << matvec2_time_->name() <<matvec2_time_->totalElapsedTime() << endl;
192 cout << update_time_->name() <<update_time_->totalElapsedTime() << endl;
193 cout << apply_time_->name() <<apply_time_->totalElapsedTime() << endl;
197 void ShyLU_Local_Schur_Operator::ResetTempVectors(
int nvectors)
200 nvectors_ = nvectors;
203 (
new Epetra_MultiVector(View, *(ssym_->Drhs), 0, nvectors));
205 (
new Epetra_MultiVector(View, *(ssym_->Dlhs), 0, nvectors));
207 (
new Epetra_MultiVector(View, *(ssym_->Gvec), 0, nvectors));
208 ssym_->OrigLP->SetLHS(localX.
getRawPtr());
210 ssym_->ReIdx_LP->fwd();
214 int ShyLU_Local_Schur_Operator::ApplyInverse(
const Epetra_MultiVector &X,
215 Epetra_MultiVector &Y)
const 220 double ShyLU_Local_Schur_Operator::NormInf()
const 225 const char* ShyLU_Local_Schur_Operator::Label()
const 227 return "Shylu Local Schur Operator (Wide Separator)";
230 bool ShyLU_Local_Schur_Operator::UseTranspose()
const 235 bool ShyLU_Local_Schur_Operator::HasNormInf()
const 240 const Epetra_Comm& ShyLU_Local_Schur_Operator::Comm()
const 245 const Epetra_Map& ShyLU_Local_Schur_Operator::OperatorDomainMap()
const 248 return G_->DomainMap();
251 const Epetra_Map& ShyLU_Local_Schur_Operator::OperatorRangeMap()
const 254 return G_->RangeMap();