[Trilinos-Users] [EXTERNAL] Re: using teko

[Hoang Giang Bui]

Sounds good but it doesn't compile again. The configure output contains:
...
-- Processing enabled TPL: Boost
--   TPL_Boost_INCLUDE_DIRS='/opt/boost/include'
...
-- Setting TeuchosCore_ENABLE_Boost=ON since TPL_ENABLE_Boost=ON
...
Is that Boost enabled as you mention?

BR
Giang Bui


On Fri, Jan 4, 2013 at 6:52 PM, Bartlett, Roscoe A. <bartlettra at ornl.gov>wrote:

> **
> The boost::shared_ptr is handled correctly when passed to the RCP. Just
> enable support for boost.
>
> -Ross
>
> Sent from my android phone. Excuse the terse message.
>
>
>
> -----Original Message-----
> *From: *Cyr, Eric C [eccyr at sandia.gov]
> *Sent: *Friday, January 04, 2013 12:22 PM Eastern Standard Time
> *To: *Hoang Giang Bui; trilinos-users at software.sandia.gov
> *Subject: *Re: [Trilinos-Users] [EXTERNAL] Re: using teko
>
> Giang Bui
>
> I think that the rcp's associated with the teko linear operator are going
> out of scope. Thus when the 'a' (boost::shared_ptr) is accessed the
> underlying
> epetra operator has been deleted by the Teuchos::RCP. This is easily
> resolved by ensuring the RCP does not delete the operator when it goes
> out of scope and gets deconstructed. You have a few options in this regard.
>
> I like this one (though others don't seem to use it as much as me)
>
> Teko::LinearOp block_a  = Thyra::epetraLinearOp(Teuchos::rcpFromRef(*a));
> Teko::LinearOp block_b1 = Thyra::epetraLinearOp(Teuchos::rcpFromRef(*b1));
> Teko::LinearOp block_b2 = Thyra::epetraLinearOp(Teuchos::rcpFromRef(*b2));
> Teko::LinearOp block_c  = Thyra::epetraLinearOp(Teuchos::rcpFromRef(*c));
>
> This builds an RCP that understands that it _does_not_ own the underlying
> pointer
> and thus should not try to delete it when the reference count goes to
> zero. It
> does however imply that the user is responsible for make sure that the
> reference
> stays in scope (is not deallocated) for as long as the RCP is required. An
> option with
> a different statement but with the same meaning is to do
>
> Teko::LinearOp block_a  = Thyra::epetraLinearOp(Teuchos::rcp(&(*a),false));
> Teko::LinearOp block_b1 =
> Thyra::epetraLinearOp(Teuchos::rcp(&(*b1),false));
> Teko::LinearOp block_b2 =
> Thyra::epetraLinearOp(Teuchos::rcp(&(*b2),false));
> Teko::LinearOp block_c  = Thyra::epetraLinearOp(Teuchos::rcp(&(*c),false));
>
> This says the rcp does not own the memory. A final option, which I have
> never
> seen used but think it might be fun to try is
>
> Teko::LinearOp block_a  = Thyra::epetraLinearOp(Teuchos::rcp(a));
> Teko::LinearOp block_b1 = Thyra::epetraLinearOp(Teuchos::rcp(b1));
> Teko::LinearOp block_b2 = Thyra::epetraLinearOp(Teuchos::rcp(b2));
> Teko::LinearOp block_c  = Thyra::epetraLinearOp(Teuchos::rcp(c));
>
> This uses the boost::shared_ptr directly, not sure if the pointer is shared
> between the RCP and boost::shared_ptr in some intelligent way or not.
>
> Anyway, I see that you are looking at Multi Phase, I would be interested
> to see what you are doing, if you wanted to discuss it or had a
> preprint/publication to send on. Again, any questions with teko please
> doesn't hesitate to contact me.
>
> Best,
> Eric Cyr
>
>  From: Hoang Giang Bui <hgbk2008 at gmail.com>
> Date: Fri, 4 Jan 2013 12:01:51 +0100
> To: Eric Cyr <eccyr at sandia.gov>, "trilinos-users at software.sandia.gov" <
> trilinos-users at software.sandia.gov>
> Subject: [EXTERNAL] Re: [Trilinos-Users] using teko
>
>
> Hi
>
> This seems working. The preconditioner by Teko could compute the Schur
> complement and the Aztec solver works. However, there is one issue which I
> did not fully understand:
>
> I created Teko linear operator for each of my blocks by:
>
> Teko::LinearOp block_a  = Thyra::epetraLinearOp(Teuchos::rcp(&(*a)));
> Teko::LinearOp block_b1 = Thyra::epetraLinearOp(Teuchos::rcp(&(*b1)));
> Teko::LinearOp block_b2 = Thyra::epetraLinearOp(Teuchos::rcp(&(*b2)));
> Teko::LinearOp block_c  = Thyra::epetraLinearOp(Teuchos::rcp(&(*c)));
>
> Where a, b1, b2, c are of boost::shared_ptr< Epetra_FECrsMatrix >
>
> creating the block:
>
> Teuchos::RCP<Teko::Epetra::EpetraOperatorWrapper> P = Teuchos::rcp(new
> Teko::Epetra::EpetraOperatorWrapper(Teko::block2x2(block_a,block_b1,block_b2,block_c)));
>
> creating the factory:
> Teuchos::RCP<Teko::BlockPreconditionerFactory> precFact = Teuchos::rcp(new
> MultiphaseBlock2x2PreconditionerFactory(inverseA, inverseS));
>
> build the preconditioner:
>
> Teko::Epetra::EpetraBlockPreconditioner prec(precFact);
> prec.buildPreconditioner(P);
>
>
> pass to aztec and solve:
> ....
>            aztec_solver.SetPrecOperator(&prec);
>            aztec_solver.Iterate(mmax_iter, mtol);
> ....
>
> However after solving, I tried to probe the value of matrix 'a' by
> std::cout << (*a)[0][0] leads to a segmentation fault error. This prevents
> me to propagate to next analysis step. I'm worry that the declaration
> Teko::LinearOp block_a  = Thyra::epetraLinearOp(Teuchos::rcp(&(*a)) has
> some implications. Can you advise me if there is anything wrong with the
> scheme above ?
>
> BR
> Giang Bui
>
>
>
>
>
>
>
>
> On Thu, Jan 3, 2013 at 11:52 PM, Cyr, Eric C <eccyr at sandia.gov> wrote:
>
>>  Giang Bui,
>>
>> Thanks for the question. I'll try to explain what is going on.
>>
>>    const Teko::LinearOp invA = Teko::buildInverse(*invAfact_, A);
>>    const Teko::LinearOp invS = Teko::buildInverse(*invSfact_, S);
>>
>>  These operations create some approximation inverse of 'A' and 'S'. None
>> of the standard inverse factories create a concrete Epetra_CrsMatrix (they
>> are possibly Epetra_Operator objects). Now when you call
>>
>>    const Teko::LinearOp tmpP = Teko::explicitMultiply(invS,
>> Teko::explicitMultiply(B2, invA));
>>
>>  What this says is that you want to form the product of invS*B2*invA
>> explicitly. That means
>> it will try to do first multiply B2*invA and store it as a CRS matrix.
>> Generally this operation is not possible
>> since B2 is probably a sparse operator and invA is likely an
>> approximation of a dense
>> inverse (using Aztec for instance). Teko will not allow this operation
>> because it would be
>> very, very expensive. What Teko does allow is "implicit" multiply
>> operations. That is you
>> can compute the action of invS*B2*invA by chaining together their
>> application. To do that
>> the pertinent line is
>>
>>     const Teko::LinearOp tmpP = Teko::multiply(invS, Teko::multiply(B2,
>> invA));
>>
>>  So what is explicitMultiply and explicitAdd for then? Well sometimes
>> you know you have
>> a concrete sparse matrix and you want to explicitly form an operator so
>> you can apply
>> multigrid or ILU to it as an approximate inverse. For example, in the
>> SIMPLE preconditioner
>> for Navier-Stokes you need to form a Schur complement that looks like
>>
>>    S=C-B*H*Bt
>>
>> Where C, B, and Bt are sparse CRS matrices and H is a diagonal matrix.
>> Thus the Teko code
>> for combining those is (not the most efficient version mind you)
>>
>>    const Teko::LinearOp S =
>> Teko::explicitAdd(C,Teko::scale(-1.0,Teko::explicitMulitiply(B,H,Bt)));
>>
>> This gives a concrete sparse CRS matrix that you can then apply ML,
>> IFPACK or Amesos to. The
>> efficient construction of these approximations is what makes block
>> preconditioning challenging,
>> but this is an issue that has to be resolved and understood at the
>> mathematical level.
>>
>> Please let me know if some I said is not clear. I appreciate your
>> interest in Teko and would be eager
>> to hear about your project and its goals, so please feel free to contact
>> me.
>>
>> Eric
>>
>> From: "Heroux, Michael A" <maherou at sandia.gov>
>> Date: Thu, 3 Jan 2013 15:30:26 -0700
>> To: Eric Cyr <eccyr at sandia.gov>, Hoang Giang Bui <hgbk2008 at gmail.com>
>> Subject: Re: [Trilinos-Users] using teko
>>
>>   Giang Bui,
>>
>> I have copied Eric Cyr on this message.  He is the principal author of
>> Teko.
>>
>> Eric, Do you have any advice?
>>
>> Thanks.
>>
>> Mike
>>
>> From: Hoang Giang Bui <hgbk2008 at gmail.com>
>> Date: Tuesday, December 25, 2012 7:42 AM
>> To: "trilinos-users at software.sandia.gov" <
>> trilinos-users at software.sandia.gov>
>> Subject: [EXTERNAL] [Trilinos-Users] using teko
>>
>>     Hi
>>
>> When using Teko to build the preconditioner operator for Aztec solver I
>> faced this error:
>>
>> I tried to create the LinearOp of invS * B2 * invA. The invS is the
>> approximation of S^{-1} and invA is approximation of A^{-1}. I created invS
>> and invA by this code
>>
>> const Teko::LinearOp invA = Teko::buildInverse(*invAfact_, A);
>> const Teko::LinearOp invS = Teko::buildInverse(*invSfact_, S);
>>
>> Then
>>
>> const Teko::LinearOp tmpP = Teko::explicitMultiply(invS,
>> Teko::explicitMultiply(B2, invA));
>>
>> At this point, the program crashed and thow the error:
>>
>> Throw number = 1
>>
>> Throw test that evaluated to true: true
>>
>> dyn_cast<Thyra::EpetraLinearOp>(Thyra::LinearOpBase<double>) : Error, the
>> object with the concrete type 'Thyra::DefaultInverseLinearOp<double>'
>> (passed in through the interface type 'Thyra::LinearOpBase<double>')  does
>> not support the interface 'Thyra::EpetraLinearOp' and the dynamic cast
>> failed!
>>
>> As I assumed, the type of invA and invS is
>> Thyra::DefaultInverseLinearOp<double> and not inherited from
>> Thyra::EpetraLinearOp, hence the dynamic_cast is failed. I also think that
>> the InverseLinearOp does not try to create Crs_Matrix so it is not
>> compatible with EpetraLinearOp.
>>
>> Please help me on this case. It's important for me to create a block
>> preconditioner for my problem.
>>
>> BR
>> Giang Bui
>>
>>
>
>
> --
> With Best Regards !
> Giang Bui
> To learn and to excel
>
>


-- 
With Best Regards !
Giang Bui
To learn and to excel
-------------- next part --------------
An HTML attachment was scrubbed...
URL: https://software.sandia.gov/pipermail/trilinos-users/attachments/20130105/ceabdcb5/attachment.html