00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042 #include "Epetra_config.h"
00043 #include "Epetra_MultiVector.h"
00044 #include "LOCA_Epetra_LeftPreconditionedOp.H"
00045
00046 LOCA::Epetra::LeftPreconditionedOp::LeftPreconditionedOp(
00047 const Teuchos::RCP<Epetra_Operator>& jacOperator,
00048 const Teuchos::RCP<Epetra_Operator>& precOperator) :
00049 label("LOCA::Epetra::LeftPreconditionedOp"),
00050 J(jacOperator),
00051 M(precOperator),
00052 useTranspose(false)
00053 {
00054 }
00055
00056 LOCA::Epetra::LeftPreconditionedOp::~LeftPreconditionedOp()
00057 {
00058 }
00059
00060 int
00061 LOCA::Epetra::LeftPreconditionedOp::SetUseTranspose(bool UseTranspose)
00062 {
00063 useTranspose = UseTranspose;
00064 int res_1 = J->SetUseTranspose(UseTranspose);
00065 int res_2 = M->SetUseTranspose(UseTranspose);
00066
00067 return res_1 + res_2;
00068 }
00069
00070 int
00071 LOCA::Epetra::LeftPreconditionedOp::Apply(const Epetra_MultiVector& Input,
00072 Epetra_MultiVector& Result) const
00073 {
00074
00075 Epetra_MultiVector tmp(Input);
00076 int res_1, res_2;
00077
00078 if (!useTranspose) {
00079
00080
00081 res_1 = J->Apply(Input, tmp);
00082
00083
00084 res_2 = M->ApplyInverse(tmp, Result);
00085
00086 }
00087 else {
00088
00089
00090 res_1 = M->ApplyInverse(Input, tmp);
00091
00092
00093 res_2 = J->Apply(tmp, Result);
00094
00095 }
00096
00097 return res_1 + res_2;
00098 }
00099
00100 int
00101 LOCA::Epetra::LeftPreconditionedOp::ApplyInverse(
00102 const Epetra_MultiVector& Input,
00103 Epetra_MultiVector& Result) const
00104 {
00105
00106 Epetra_MultiVector tmp(Input);
00107 int res_1, res_2;
00108
00109 if (!useTranspose) {
00110
00111
00112 res_1 = M->Apply(Input, tmp);
00113
00114
00115 res_2 = J->ApplyInverse(tmp, Result);
00116
00117 }
00118 else {
00119
00120
00121 res_1 = J->ApplyInverse(Input, tmp);
00122
00123
00124 res_2 = M->Apply(tmp, Result);
00125
00126 }
00127
00128 return res_1 + res_2;
00129 }
00130
00131 double
00132 LOCA::Epetra::LeftPreconditionedOp::NormInf() const
00133 {
00134 return J->NormInf() + 1.0/M->NormInf();
00135 }
00136
00137
00138 const char*
00139 LOCA::Epetra::LeftPreconditionedOp::Label () const
00140 {
00141 return const_cast<char*>(label.c_str());
00142 }
00143
00144 bool
00145 LOCA::Epetra::LeftPreconditionedOp::UseTranspose() const
00146 {
00147 return useTranspose;
00148 }
00149
00150 bool
00151 LOCA::Epetra::LeftPreconditionedOp::HasNormInf() const
00152 {
00153 return J->HasNormInf() && M->HasNormInf();
00154 }
00155
00156 const Epetra_Comm &
00157 LOCA::Epetra::LeftPreconditionedOp::Comm() const
00158 {
00159 return J->Comm();
00160 }
00161 const Epetra_Map&
00162 LOCA::Epetra::LeftPreconditionedOp::OperatorDomainMap() const
00163 {
00164 return J->OperatorDomainMap();
00165 }
00166
00167 const Epetra_Map&
00168 LOCA::Epetra::LeftPreconditionedOp::OperatorRangeMap() const
00169 {
00170 return M->OperatorDomainMap();
00171 }