当前位置: 首页>>代码示例>>C++>>正文


C++ LAPACK::GEQRF方法代码示例

本文整理汇总了C++中teuchos::LAPACK::GEQRF方法的典型用法代码示例。如果您正苦于以下问题:C++ LAPACK::GEQRF方法的具体用法?C++ LAPACK::GEQRF怎么用?C++ LAPACK::GEQRF使用的例子?那么恭喜您, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在teuchos::LAPACK的用法示例。


在下文中一共展示了LAPACK::GEQRF方法的1个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的C++代码示例。

示例1: K

Stokhos::MonoProjPCEBasis<ordinal_type, value_type>::
MonoProjPCEBasis(
   ordinal_type p,
   const Stokhos::OrthogPolyApprox<ordinal_type, value_type>& pce,
   const Stokhos::Quadrature<ordinal_type, value_type>& quad,
   const Stokhos::Sparse3Tensor<ordinal_type, value_type>& Cijk,
   bool limit_integration_order_) :
  RecurrenceBasis<ordinal_type, value_type>("Monomial Projection", p, true),
  limit_integration_order(limit_integration_order_),
  pce_sz(pce.basis()->size()),
  pce_norms(pce.basis()->norm_squared()),
  a(pce_sz), 
  b(pce_sz),
  basis_vecs(pce_sz, p+1),
  new_pce(p+1)
{
  // If the original basis is normalized, we can use the standard QR
  // factorization.  For simplicity, we renormalize the PCE coefficients
  // for a normalized basis
  Stokhos::OrthogPolyApprox<ordinal_type, value_type> normalized_pce(pce);
  for (ordinal_type i=0; i<pce_sz; i++) {
    pce_norms[i] = std::sqrt(pce_norms[i]);
    normalized_pce[i] *= pce_norms[i];
  }

  // Evaluate PCE at quad points
  ordinal_type nqp = quad.size();
  Teuchos::Array<value_type> pce_vals(nqp);
  const Teuchos::Array<value_type>& weights = quad.getQuadWeights();
  const Teuchos::Array< Teuchos::Array<value_type> >& quad_points =
    quad.getQuadPoints();
  const Teuchos::Array< Teuchos::Array<value_type> >& basis_values =
    quad.getBasisAtQuadPoints();
  for (ordinal_type i=0; i<nqp; i++) {
    pce_vals[i] = normalized_pce.evaluate(quad_points[i], basis_values[i]);
  }

  // Form Kylov matrix up to order pce_sz
  matrix_type K(pce_sz, pce_sz);

  // Compute matrix
  matrix_type A(pce_sz, pce_sz);
  typedef Stokhos::Sparse3Tensor<ordinal_type, value_type> Cijk_type;
  for (typename Cijk_type::k_iterator k_it = Cijk.k_begin();
       k_it != Cijk.k_end(); ++k_it) {
    ordinal_type k = index(k_it);
    for (typename Cijk_type::kj_iterator j_it = Cijk.j_begin(k_it); 
	 j_it != Cijk.j_end(k_it); ++j_it) {
      ordinal_type j = index(j_it);
      value_type val = 0;
      for (typename Cijk_type::kji_iterator i_it = Cijk.i_begin(j_it);
	   i_it != Cijk.i_end(j_it); ++i_it) {
	ordinal_type i = index(i_it);
	value_type c = value(i_it) / (pce_norms[j]*pce_norms[k]);
	val += pce[i]*c;
      }
      A(k,j) = val;
    }
  }

  // Each column i is given by projection of the i-th order monomial 
  // onto original basis
  vector_type u0 = Teuchos::getCol(Teuchos::View, K, 0);
  u0(0) = 1.0;
  for (ordinal_type i=1; i<pce_sz; i++)
    u0(i) = 0.0;
  for (ordinal_type k=1; k<pce_sz; k++) {
    vector_type u = Teuchos::getCol(Teuchos::View, K, k);
    vector_type up = Teuchos::getCol(Teuchos::View, K, k-1);
    u.multiply(Teuchos::NO_TRANS, Teuchos::NO_TRANS, 1.0, A, up, 0.0);
  }
  /*
  for (ordinal_type j=0; j<pce_sz; j++) {
    for (ordinal_type i=0; i<pce_sz; i++) {
      value_type val = 0.0;
      for (ordinal_type k=0; k<nqp; k++)
	val += weights[k]*std::pow(pce_vals[k],j)*basis_values[k][i];
      K(i,j) = val;
    }
  }
  */

  std::cout << K << std::endl << std::endl;

  // Compute QR factorization of K
  ordinal_type ws_size, info;
  value_type ws_size_query;
  Teuchos::Array<value_type> tau(pce_sz);
  Teuchos::LAPACK<ordinal_type,value_type> lapack;
  lapack.GEQRF(pce_sz, pce_sz, K.values(), K.stride(), &tau[0], 
	       &ws_size_query, -1, &info);
  TEUCHOS_TEST_FOR_EXCEPTION(info != 0, std::logic_error, 
		     "GEQRF returned value " << info);
  ws_size = static_cast<ordinal_type>(ws_size_query);
  Teuchos::Array<value_type> work(ws_size);
  lapack.GEQRF(pce_sz, pce_sz, K.values(), K.stride(), &tau[0], 
	       &work[0], ws_size, &info);
  TEUCHOS_TEST_FOR_EXCEPTION(info != 0, std::logic_error, 
		     "GEQRF returned value " << info);
  
//.........这里部分代码省略.........
开发者ID:00liujj,项目名称:trilinos,代码行数:101,代码来源:Stokhos_MonoProjPCEBasisImp.hpp


注:本文中的teuchos::LAPACK::GEQRF方法示例由纯净天空整理自Github/MSDocs等开源代码及文档管理平台,相关代码片段筛选自各路编程大神贡献的开源项目,源码版权归原作者所有,传播和使用请参考对应项目的License;未经允许,请勿转载。