45 #include "Thyra_VectorStdOps.hpp" 46 #include "Thyra_ProductVectorSpaceBase.hpp" 47 #include "Thyra_get_Epetra_Operator.hpp" 48 #include "Thyra_TpetraLinearOp.hpp" 50 #include "Tpetra_CrsMatrix.hpp" 55 template <
typename ScalarT,
typename LocalOrdinalT,
typename GlobalOrdinalT,
typename NodeT>
63 bool x_matches=
false, f_matches=
false, dxdt_matches=
false;
66 RCP<const VectorSpaceBase<ScalarT> > range = get_A()->range();
67 RCP<const VectorSpaceBase<ScalarT> > domain = get_A()->domain();
70 x_matches = range->isCompatible(*get_x()->space());
75 dxdt_matches = range->isCompatible(*get_dxdt()->space());
80 f_matches = range->isCompatible(*get_f()->space());
84 else if(get_x()!=null && get_dxdt()!=null) {
86 x_matches = get_x()->space()->isCompatible(*get_dxdt()->space());
87 dxdt_matches = x_matches;
90 f_matches = x_matches = dxdt_matches =
true;
93 return x_matches && dxdt_matches && f_matches;
96 template <
typename ScalarT,
typename LocalOrdinalT,
typename GlobalOrdinalT,
typename NodeT>
101 using Thyra::PhysicallyBlockedLinearOpBase;
102 using Thyra::ProductVectorSpaceBase;
104 using Teuchos::rcp_dynamic_cast;
106 if(get_x()!=Teuchos::null) Thyra::assign<ScalarT>(x.ptr(),0.0);
107 if(get_dxdt()!=Teuchos::null) Thyra::assign<ScalarT>(get_dxdt().ptr(),0.0);
108 if(get_f()!=Teuchos::null) Thyra::assign<ScalarT>(get_f().ptr(),0.0);
109 if(get_A()!=Teuchos::null) {
110 RCP<PhysicallyBlockedLinearOpBase<ScalarT> > Amat
111 = rcp_dynamic_cast<PhysicallyBlockedLinearOpBase<ScalarT> >(get_A(),
true);
112 RCP<const ProductVectorSpaceBase<ScalarT> > range = Amat->productRange();
113 RCP<const ProductVectorSpaceBase<ScalarT> > domain = Amat->productDomain();
116 for(
int i=0;i<range->numBlocks();i++) {
117 for(
int j=0;j<domain->numBlocks();j++) {
118 RCP<LinearOpBase<ScalarT> > block = Amat->getNonconstBlock(i,j);
119 if(block!=Teuchos::null) {
120 RCP<Tpetra::Operator<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> > t_block =
121 rcp_dynamic_cast<Thyra::TpetraLinearOp<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> >(block,
true)->getTpetraOperator();
123 RCP<const MapType> map_i = t_block->getRangeMap();
124 RCP<const MapType> map_j = t_block->getDomainMap();
126 RCP<Tpetra::CrsMatrix<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> > mat =
127 rcp_dynamic_cast<Tpetra::CrsMatrix<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> >(t_block,
true);
130 mat->setAllToScalar(0.0);
131 mat->fillComplete(map_j,map_i);
138 template <
typename ScalarT,
typename LocalOrdinalT,
typename GlobalOrdinalT,
typename NodeT>
143 using Thyra::PhysicallyBlockedLinearOpBase;
144 using Thyra::ProductVectorSpaceBase;
146 using Teuchos::rcp_dynamic_cast;
148 if(get_A()!=Teuchos::null) {
149 RCP<PhysicallyBlockedLinearOpBase<ScalarT> > Amat
150 = rcp_dynamic_cast<PhysicallyBlockedLinearOpBase<ScalarT> >(get_A(),
true);
151 RCP<const ProductVectorSpaceBase<ScalarT> > range = Amat->productRange();
152 RCP<const ProductVectorSpaceBase<ScalarT> > domain = Amat->productDomain();
155 for(
int i=0;i<range->numBlocks();i++) {
156 for(
int j=0;j<domain->numBlocks();j++) {
157 RCP<LinearOpBase<ScalarT> > block = Amat->getNonconstBlock(i,j);
158 if(block!=Teuchos::null) {
159 RCP<Tpetra::Operator<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> > t_block =
160 rcp_dynamic_cast<Thyra::TpetraLinearOp<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> >(block,
true)->getTpetraOperator();
163 RCP<const MapType> map_i = t_block->getRangeMap();
164 RCP<const MapType> map_j = t_block->getDomainMap();
166 RCP<Tpetra::CrsMatrix<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> > mat =
167 rcp_dynamic_cast<Tpetra::CrsMatrix<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> >(t_block,
true);
170 mat->setAllToScalar(
value);
171 mat->fillComplete(map_j,map_i);
178 template <
typename ScalarT,
typename LocalOrdinalT,
typename GlobalOrdinalT,
typename NodeT>
182 set_x(Teuchos::null);
183 set_dxdt(Teuchos::null);
184 set_f(Teuchos::null);
185 set_A(Teuchos::null);
188 template <
typename ScalarT,
typename LocalOrdinalT,
typename GlobalOrdinalT,
typename NodeT>
193 using Thyra::PhysicallyBlockedLinearOpBase;
194 using Thyra::ProductVectorSpaceBase;
196 using Teuchos::rcp_dynamic_cast;
198 if(get_A()!=Teuchos::null) {
199 RCP<PhysicallyBlockedLinearOpBase<ScalarT> > Amat
200 = rcp_dynamic_cast<PhysicallyBlockedLinearOpBase<ScalarT> >(get_A(),
true);
201 RCP<const ProductVectorSpaceBase<ScalarT> > range = Amat->productRange();
202 RCP<const ProductVectorSpaceBase<ScalarT> > domain = Amat->productDomain();
205 for(
int i=0;i<range->numBlocks();i++) {
206 for(
int j=0;j<domain->numBlocks();j++) {
207 RCP<LinearOpBase<ScalarT> > block = Amat->getNonconstBlock(i,j);
208 if(block!=Teuchos::null) {
209 RCP<Tpetra::Operator<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> > t_block =
210 rcp_dynamic_cast<Thyra::TpetraLinearOp<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> >(block,
true)->getTpetraOperator();
212 RCP<Tpetra::CrsMatrix<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> > mat =
213 rcp_dynamic_cast<Tpetra::CrsMatrix<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> >(t_block,
true);
222 template <
typename ScalarT,
typename LocalOrdinalT,
typename GlobalOrdinalT,
typename NodeT>
227 using Thyra::PhysicallyBlockedLinearOpBase;
228 using Thyra::ProductVectorSpaceBase;
230 using Teuchos::rcp_dynamic_cast;
232 if(get_A()!=Teuchos::null) {
233 RCP<PhysicallyBlockedLinearOpBase<ScalarT> > Amat
234 = rcp_dynamic_cast<PhysicallyBlockedLinearOpBase<ScalarT> >(get_A(),
true);
235 RCP<const ProductVectorSpaceBase<ScalarT> > range = Amat->productRange();
236 RCP<const ProductVectorSpaceBase<ScalarT> > domain = Amat->productDomain();
239 for(
int i=0;i<range->numBlocks();i++) {
240 for(
int j=0;j<domain->numBlocks();j++) {
241 RCP<LinearOpBase<ScalarT> > block = Amat->getNonconstBlock(i,j);
242 if(block!=Teuchos::null) {
243 RCP<Tpetra::Operator<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> > t_block =
244 rcp_dynamic_cast<Thyra::TpetraLinearOp<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> >(block,
true)->getTpetraOperator();
246 RCP<const MapType> map_i = t_block->getRangeMap();
247 RCP<const MapType> map_j = t_block->getDomainMap();
249 RCP<Tpetra::CrsMatrix<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> > mat =
250 rcp_dynamic_cast<Tpetra::CrsMatrix<ScalarT,LocalOrdinalT,GlobalOrdinalT,NodeT> >(t_block,
true);
252 mat->fillComplete(map_j,map_i);
bool checkCompatibility() const
Make sure row and column spaces match up.
void initializeMatrix(ScalarT value)
Put a particular scalar in the matrix.
virtual void initialize()