8 #ifndef MUELU_QR_INTERFACEEX_DEF_HPP_ 9 #define MUELU_QR_INTERFACEEX_DEF_HPP_ 17 template <
class Scalar,
class LocalOrdinal>
18 void LapackQR(Teuchos::LAPACK<LocalOrdinal,Scalar> &lapack, LocalOrdinal myAggSize,
19 int intFineNSDim, ArrayRCP<Scalar> &localQR, ArrayRCP<Scalar> &tau,
20 ArrayRCP<Scalar> &work, LocalOrdinal &workSize, LocalOrdinal &info)
22 lapack.ORGQR(myAggSize, intFineNSDim, intFineNSDim, localQR.getRawPtr(),
23 myAggSize, tau.getRawPtr(), work.getRawPtr(), workSize, &info );
27 template <
class LocalOrdinal>
28 void LapackQR(Teuchos::LAPACK<LocalOrdinal, std::complex<double> > &lapack,
29 LocalOrdinal myAggSize,
int intFineNSDim, ArrayRCP<std::complex<double> > &localQR,
30 ArrayRCP<std::complex<double> > &tau, ArrayRCP<std::complex<double> > &work,
31 LocalOrdinal &workSize, LocalOrdinal &info)
33 lapack.UNGQR(myAggSize, intFineNSDim, intFineNSDim, localQR.getRawPtr(),
34 myAggSize, tau.getRawPtr(), work.getRawPtr(), workSize, &info );
37 template <
class Scalar,
class LocalOrdinal>
39 tau_ = ArrayRCP<Scalar>(NSDim);
40 work_ = ArrayRCP<Scalar>(NSDim);
43 template <
class Scalar,
class LocalOrdinal>
45 bool bIsZeroColumn =
true;
46 for(LocalOrdinal r = 0; r < myAggSize; ++r) {
47 if(localQR[ myAggSize*nspCol + r ]!=0.0) {
48 bIsZeroColumn =
false;
55 template <
class Scalar,
class LocalOrdinal>
58 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType Magnitude;
63 for (LocalOrdinal k=0; k<myAggSize; ++k) {dtemp += Teuchos::ScalarTraits<Scalar>::magnitude(localQR[k])*Teuchos::ScalarTraits<Scalar>::magnitude(localQR[k]);}
64 dtemp = Teuchos::ScalarTraits<Magnitude>::squareroot(dtemp);
71 for(LocalOrdinal j = 0; j<
NSDim_; ++j) {
88 for(
size_t r = 0; r<myAggSize; ++r) {
106 std::string msg =
"QR_InterfaceEx: dgeqrf (LAPACK QR routine) returned error code " +
Teuchos::toString(
info_);
114 if ( Teuchos::ScalarTraits<Scalar>::magnitude(work_[0]) > workSize_) {
115 workSize_ = Teuchos::as<int>(Teuchos::ScalarTraits<Scalar>::magnitude(work_[0]));
138 template <
class Scalar,
class LocalOrdinal>
141 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType Magnitude;
144 Magnitude dtemp = Teuchos::ScalarTraits<Scalar>::magnitude(localQR[0]);
145 localQR[0] =
tau_[0];
147 for (LocalOrdinal i=0; i<myAggSize; ++i)
156 std::string msg =
"QR_InterfaceEx: dorgqr (LAPACK auxiliary QR routine) returned error code " +
Teuchos::toString(
info_);
162 if ( Teuchos::ScalarTraits<Scalar>::magnitude(
work_[0]) >
workSize_) {
163 workSize_ = Teuchos::as<int>(Teuchos::ScalarTraits<Scalar>::magnitude(
work_[0]));
170 for(
size_t r = 0; r<myAggSize; ++r) {
QR_InterfaceEx(const size_t nullSpaceDimension)
Constructor.
std::string toString(const T &what)
Little helper function to convert non-string types to strings.
LocalOrdinal info_
(out) =0: success; =i, i<0: i-th argument has illegal value
Teuchos::LAPACK< LocalOrdinal, Scalar > lapack_
Teuchos LAPACK wrapper.
LocalOrdinal internalNSDim_
Namespace for MueLu classes and methods.
ArrayRCP< Scalar > tau_
Internal LAPACK variables.
void ExtractQ(LocalOrdinal const &myAggSize, ArrayRCP< Scalar > &localQR)
Calculate the Q factor.
ArrayRCP< Scalar > work_
Temporary work space.
void Compute(LocalOrdinal const &myAggSize, ArrayRCP< Scalar > &localQR)
Compute the QR factorization.
std::vector< LocalOrdinal > nonZeroCols_
void LapackQR(Teuchos::LAPACK< LocalOrdinal, Scalar > &lapack, LocalOrdinal myAggSize, int intFineNSDim, ArrayRCP< Scalar > &localQR, ArrayRCP< Scalar > &tau, ArrayRCP< Scalar > &work, LocalOrdinal &workSize, LocalOrdinal &info)
Non-member templated function to handle extracting Q from QR factorization.
LocalOrdinal NSDim_
Dimension of nullspace.
bool isZeroNspColumn(LocalOrdinal const &myAggSize, ArrayRCP< Scalar > &localQR, LocalOrdinal nspCol)
Exception throws to report errors in the internal logical of the program.
LocalOrdinal workSize_
Length of work vectors. Must be at least dimension of nullspace.
ArrayRCP< Scalar > internalLocalQR_
internal local QR (without zero columns)