C++ Array::size方法代码示例

      //! Print eigenpair
      void print(std::ostream& os) const {
	os << eig_val << ", ";
	for (std::size_t i=0; i<eig_pairs.size()-1; i++) {
	  os << "(";
	  eig_pairs[i].eig_func->print(os);
	  os << ") * ";
	}
	os << "(";
	eig_pairs[eig_pairs.size()-1].eig_func->print(os);
	os << ")";
      }

//---------------------------------------------------------------------------//
TEUCHOS_UNIT_TEST( NanoflannTree, dim_3_factory_test )
{
    int dim = 3;
    int num_points = 10;
    int num_coords = dim*num_points;

    Teuchos::Array<double> coords(num_coords);
    for ( int i = 0; i < num_points; ++i )
    {
        coords[dim*i] = 1.0*i;
        coords[dim*i+1] = 1.0;
        coords[dim*i+2] = 1.0;
    }

    int max_leaf_size = 3;
    Teuchos::RCP<DataTransferKit::StaticSearchTree> tree =
        DataTransferKit::SearchTreeFactory::createStaticTree(
            3, coords(), max_leaf_size );

    Teuchos::Array<double> p1( dim );
    p1[0] = 4.9;
    p1[1] = 1.0;
    p1[2] = 1.0;
    Teuchos::Array<double> p2( dim );
    p2[0] = 11.4;
    p2[1] = 1.0;
    p2[2] = 1.0;

    int num_neighbors = 1;
    Teuchos::Array<unsigned> nnearest =
        tree->nnSearch( p1(), num_neighbors );
    TEST_EQUALITY( num_neighbors, nnearest.size() );
    TEST_EQUALITY( 5, nnearest[0] );

    nnearest = tree->nnSearch( p2(), num_neighbors );
    TEST_EQUALITY( num_neighbors, nnearest.size() );
    TEST_EQUALITY( 9, nnearest[0] );

    double radius = 1.1;
    nnearest = tree->radiusSearch( p1(), radius );
    TEST_EQUALITY( 3, nnearest.size() );
    TEST_EQUALITY( 5, nnearest[0] )
    TEST_EQUALITY( 4, nnearest[1] )
    TEST_EQUALITY( 6, nnearest[2] )

    nnearest = tree->radiusSearch( p2(), radius );
    TEST_EQUALITY( 0, nnearest.size() );

    radius = 2.41;
    nnearest = tree->radiusSearch( p2(), radius );
    TEST_EQUALITY( 1, nnearest.size() );
    TEST_EQUALITY( 9, nnearest[0] )
}

/*
 * Construct from an array of vectors and a specifier for the 
 * domain space. 
 */
template <class Scalar> inline
MultiVectorOperator<Scalar>
::MultiVectorOperator(const Teuchos::Array<Vector<Scalar> >& cols,
  const VectorSpace<Scalar>& domain)
  : LinearOpWithSpaces<Scalar>(domain, cols[0].space()),
    cols_(cols)
{
  TEUCHOS_TEST_FOR_EXCEPTION(cols.size() == 0, std::runtime_error,
    "empty multivector given to MultiVectorOperator ctor");
  for (int i=1; i<cols.size(); i++)
  {
    TEUCHOS_TEST_FOR_EXCEPTION(cols[i].space() != cols[0].space(), std::runtime_error,
      "inconsistent vector spaces in  MultiVectorOperator ctor");
  }
}

void
Albany::SamplingBasedScalarResponseFunction::
evaluateSGResponse(
  const double curr_time,
  const Stokhos::EpetraVectorOrthogPoly* sg_xdot,
  const Stokhos::EpetraVectorOrthogPoly* sg_xdotdot,
  const Stokhos::EpetraVectorOrthogPoly& sg_x,
  const Teuchos::Array<ParamVec>& p,
  const Teuchos::Array<int>& sg_p_index,
  const Teuchos::Array< Teuchos::Array<SGType> >& sg_p_vals,
  Stokhos::EpetraVectorOrthogPoly& sg_g)
{
  RCP<const Epetra_BlockMap> x_map = sg_x.coefficientMap();
  RCP<Epetra_Vector> xdot;
  if (sg_xdot != NULL)
    xdot = rcp(new Epetra_Vector(*x_map));
  RCP<Epetra_Vector> xdotdot;
  if (sg_xdotdot != NULL)
    xdotdot = rcp(new Epetra_Vector(*x_map));
  Epetra_Vector x(*x_map);
  Teuchos::Array<ParamVec> pp = p;
  
  RCP<const Epetra_BlockMap> g_map = sg_g.coefficientMap();
  Epetra_Vector g(*g_map);

  // Get quadrature data
  const Teuchos::Array<double>& norms = basis->norm_squared();
  const Teuchos::Array< Teuchos::Array<double> >& points = 
    quad->getQuadPoints();
  const Teuchos::Array<double>& weights = quad->getQuadWeights();
  const Teuchos::Array< Teuchos::Array<double> >& vals = 
    quad->getBasisAtQuadPoints();
  int nqp = points.size();

  // Compute sg_g via quadrature
  sg_g.init(0.0);
  SGConverter c(this, commT);
  for (int qp=0; qp<nqp; qp++) {

    // Evaluate sg_x, sg_xdot at quadrature point
    sg_x.evaluate(vals[qp], x);
    if (sg_xdot != NULL)
      sg_xdot->evaluate(vals[qp], *xdot);
    if (sg_xdotdot != NULL)
      sg_xdotdot->evaluate(vals[qp], *xdotdot);

    // Evaluate parameters at quadrature point
    for (int i=0; i<sg_p_index.size(); i++) {
      int ii = sg_p_index[i];
      for (unsigned int j=0; j<pp[ii].size(); j++)
	pp[ii][j].baseValue = sg_p_vals[ii][j].evaluate(points[qp], vals[qp]);
    }

    // Compute response at quadrature point
    c.evaluateResponse(curr_time, xdot.get(), xdotdot.get(), x, pp, g);

    // Add result into integral
    sg_g.sumIntoAllTerms(weights[qp], vals[qp], norms, g);
  }
}

void
Albany::SamplingBasedScalarResponseFunction::
evaluateMPResponse(
  const double curr_time,
  const Stokhos::ProductEpetraVector* mp_xdot,
  const Stokhos::ProductEpetraVector* mp_xdotdot,
  const Stokhos::ProductEpetraVector& mp_x,
  const Teuchos::Array<ParamVec>& p,
  const Teuchos::Array<int>& mp_p_index,
  const Teuchos::Array< Teuchos::Array<MPType> >& mp_p_vals,
  Stokhos::ProductEpetraVector& mp_g)
{
  Teuchos::Array<ParamVec> pp = p;
  const Epetra_Vector* xdot = NULL;
  const Epetra_Vector* xdotdot = NULL;

  SGConverter c(this, commT);
  for (int i=0; i<mp_x.size(); i++) {

    for (int k=0; k<mp_p_index.size(); k++) {
      int kk = mp_p_index[k];
      for (unsigned int j=0; j<pp[kk].size(); j++)
	pp[kk][j].baseValue = mp_p_vals[kk][j].coeff(i);
    }

    if (mp_xdot != NULL)
      xdot = mp_xdot->getCoeffPtr(i).get();
    if (mp_xdotdot != NULL)
      xdotdot = mp_xdotdot->getCoeffPtr(i).get();
    
    // Evaluate response function
    c.evaluateResponse(curr_time, xdot, xdotdot, mp_x[i], pp, mp_g[i]);
  }
}

void getD2_ij_FromMultiplicities(int &i, int &j, const Teuchos::Array<int> &multiplicities)
{
  bool i_assigned = false;
  for (int k=0; k<multiplicities.size(); k++)
  {
    if (multiplicities[k] == 2)
    {
      i = k;
      j = k;
      return;
    }
    else if (multiplicities[k] == 1)
    {
      if (! i_assigned)
      {
        i = k;
        i_assigned = true;
      }
      else
      {
        j = k;
        return;
      }
    }
  }
  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, "i,j not identified from multiplicities...")
}

 RiskVector( Teuchos::ParameterList &parlist,
       const Teuchos::RCP<Vector<Real> > &vec,
       const Real stat = 1. )
   : vec_(vec), augmented_(false), nStat_(0) {
   stat_.clear();
   dual_vec1_ = vec->dual().clone();
   std::string type = parlist.sublist("SOL").sublist("Risk Measure").get("Name","CVaR");
   if ( type == "CVaR"                       ||
        type == "HMCR"                       ||
        type == "KL Divergence"              ||
        type == "Moreau-Yosida CVaR"         ||
        type == "Log-Exponential Quadrangle" ||
        type == "Log-Quantile Quadrangle"    ||
        type == "Quantile-Based Quadrangle"  ||
        type == "Truncated Mean Quadrangle" ) {
     augmented_ = true;
     nStat_     = 1;
     stat_.resize(nStat_,stat);
   }
   else if ( type == "Mixed-Quantile Quadrangle" ) {
     Teuchos::ParameterList &list
       = parlist.sublist("SOL").sublist("Risk Measure").sublist("Mixed-Quantile Quadrangle");
     Teuchos::Array<Real> prob
       = Teuchos::getArrayFromStringParameter<Real>(list,"Probability Array");
     augmented_ = true;
     nStat_     = prob.size();
     stat_.resize(nStat_,stat);
   }
   else if ( type == "Quantile-Radius Quadrangle" ) {
     augmented_ = true;
     nStat_     = 2;
     stat_.resize(nStat_,stat);
   }
 }

Teuchos::RCP<Epetra_Map>
Albany::SolutionResponseFunction::
buildCulledMap(const Epetra_Map& x_map,
	       const Teuchos::Array<int>& keepDOF) const
{
  int numKeepDOF = std::accumulate(keepDOF.begin(), keepDOF.end(), 0);
  int Neqns = keepDOF.size();
  int N = x_map.NumMyElements(); // x_map is map for solution vector

  TEUCHOS_ASSERT( !(N % Neqns) ); // Assume that all the equations for
                                  // a given node are on the assigned
                                  // processor. I.e. need to ensure
                                  // that N is exactly Neqns-divisible

  int nnodes = N / Neqns;          // number of fem nodes
  int N_new = nnodes * numKeepDOF; // length of local x_new

  int *gids = x_map.MyGlobalElements(); // Fill local x_map into gids array
  Teuchos::Array<int> gids_new(N_new);
  int idx = 0;
  for ( int inode = 0; inode < N/Neqns ; ++inode) // For every node
    for ( int ieqn = 0; ieqn < Neqns; ++ieqn )  // Check every dof on the node
      if ( keepDOF[ieqn] == 1 )  // then want to keep this dof
	gids_new[idx++] = gids[(inode*Neqns)+ieqn];
  // end cull

  Teuchos::RCP<Epetra_Map> x_new_map =
    Teuchos::rcp( new Epetra_Map( -1, N_new, &gids_new[0], 0, x_map.Comm() ) );

  return x_new_map;
}

 void values(FieldContainer<double> &values, BasisCachePtr sliceBasisCache) {
   vector<GlobalIndexType> sliceCellIDs = sliceBasisCache->cellIDs();
   
   Teuchos::Array<int> dim;
   values.dimensions(dim);
   dim[0] = 1; // one cell
   Teuchos::Array<int> offset(dim.size());
   
   for (int cellOrdinal = 0; cellOrdinal < sliceCellIDs.size(); cellOrdinal++) {
     offset[0] = cellOrdinal;
     int enumeration = values.getEnumeration(offset);
     FieldContainer<double>valuesForCell(dim,&values[enumeration]);
     GlobalIndexType sliceCellID = sliceCellIDs[cellOrdinal];
     int numPoints = sliceBasisCache->getPhysicalCubaturePoints().dimension(1);
     int spaceDim = sliceBasisCache->getPhysicalCubaturePoints().dimension(2);
     FieldContainer<double> spaceTimePhysicalPoints(1,numPoints,spaceDim+1);
     for (int ptOrdinal=0; ptOrdinal<numPoints; ptOrdinal++) {
       for (int d=0; d<spaceDim; d++) {
         spaceTimePhysicalPoints(0,ptOrdinal,d) = sliceBasisCache->getPhysicalCubaturePoints()(cellOrdinal,ptOrdinal,d);
       }
       spaceTimePhysicalPoints(0,ptOrdinal,spaceDim) = _t;
     }
     
     GlobalIndexType cellID = _cellIDMap[sliceCellID];
     BasisCachePtr spaceTimeBasisCache = BasisCache::basisCacheForCell(_spaceTimeMesh, cellID);
     
     FieldContainer<double> spaceTimeRefPoints(1,numPoints,spaceDim+1);
     CamelliaCellTools::mapToReferenceFrame(spaceTimeRefPoints, spaceTimePhysicalPoints, _spaceTimeMesh, cellID);
     spaceTimeRefPoints.resize(numPoints,spaceDim+1);
     spaceTimeBasisCache->setRefCellPoints(spaceTimeRefPoints);
     _spaceTimeFunction->values(valuesForCell, spaceTimeBasisCache);
   }
 }

Teuchos::Array<EpetraGlobalIndex> getMyLIDs(
    const Epetra_BlockMap &map,
    const Teuchos::ArrayView<const EpetraGlobalIndex> &selectedGIDs)
{
  Teuchos::Array<EpetraGlobalIndex> sortedMyGIDs(map.MyGlobalElements(), map.MyGlobalElements() + map.NumMyElements());
  std::sort(sortedMyGIDs.begin(), sortedMyGIDs.end());

  Teuchos::Array<EpetraGlobalIndex> sortedSelectedGIDs(selectedGIDs);
  std::sort(sortedSelectedGIDs.begin(), sortedSelectedGIDs.end());

  Teuchos::Array<EpetraGlobalIndex> mySelectedGIDs;
  std::set_intersection(sortedMyGIDs.begin(), sortedMyGIDs.end(),
                        sortedSelectedGIDs.begin(), sortedSelectedGIDs.end(),
                        std::back_inserter(mySelectedGIDs));

  Teuchos::Array<EpetraGlobalIndex> result;
  result.reserve(mySelectedGIDs.size());

  std::transform(
      mySelectedGIDs.begin(), mySelectedGIDs.end(),
      std::back_inserter(result),
      std::bind1st(std::mem_fun_ref(static_cast<int(Epetra_BlockMap::*)(EpetraGlobalIndex) const>(&Epetra_BlockMap::LID)), map));

  return result;
}

Teuchos::Array<int>
Albany::NodeSetSolutionCullingStrategy::
selectedGIDs(const Epetra_BlockMap &sourceMap) const
{
  Teuchos::Array<int> result;
  {
    Teuchos::Array<int> mySelectedGIDs;
    {
      const NodeSetList &nodeSets = disc_->getNodeSets();

      const NodeSetList::const_iterator it = nodeSets.find(nodeSetLabel_);
      if (it != nodeSets.end()) {
        typedef NodeSetList::mapped_type NodeSetEntryList;
        const NodeSetEntryList &sampleNodeEntries = it->second;

        for (NodeSetEntryList::const_iterator jt = sampleNodeEntries.begin(); jt != sampleNodeEntries.end(); ++jt) {
          typedef NodeSetEntryList::value_type NodeEntryList;
          const NodeEntryList &sampleEntries = *jt;
          for (NodeEntryList::const_iterator kt = sampleEntries.begin(); kt != sampleEntries.end(); ++kt) {
            mySelectedGIDs.push_back(sourceMap.GID(*kt));
          }
        }
      }
    }

    const Epetra_Comm &comm = sourceMap.Comm();

    {
      int selectedGIDCount;
      {
        int mySelectedGIDCount = mySelectedGIDs.size();
        comm.SumAll(&mySelectedGIDCount, &selectedGIDCount, 1);
      }
      result.resize(selectedGIDCount);
    }

    const int ierr = Epetra::GatherAllV(
        comm,
        mySelectedGIDs.getRawPtr(), mySelectedGIDs.size(),
        result.getRawPtr(), result.size());
    TEUCHOS_ASSERT(ierr == 0);
  }

  std::sort(result.begin(), result.end());

  return result;
}

void
Albany::SamplingBasedScalarResponseFunction::
evaluateMPGradient(
  const double current_time,
  const Stokhos::ProductEpetraVector* mp_xdot,
  const Stokhos::ProductEpetraVector* mp_xdotdot,
  const Stokhos::ProductEpetraVector& mp_x,
  const Teuchos::Array<ParamVec>& p,
  const Teuchos::Array<int>& mp_p_index,
  const Teuchos::Array< Teuchos::Array<MPType> >& mp_p_vals,
  ParamVec* deriv_p,
  Stokhos::ProductEpetraVector* mp_g,
  Stokhos::ProductEpetraMultiVector* mp_dg_dx,
  Stokhos::ProductEpetraMultiVector* mp_dg_dxdot,
  Stokhos::ProductEpetraMultiVector* mp_dg_dxdotdot,
  Stokhos::ProductEpetraMultiVector* mp_dg_dp)
{
  Teuchos::Array<ParamVec> pp = p;
  const Epetra_Vector* xdot = NULL;
  const Epetra_Vector* xdotdot = NULL;
  Epetra_Vector* g = NULL;
  Epetra_MultiVector* dg_dx = NULL;
  Epetra_MultiVector* dg_dxdot = NULL;
  Epetra_MultiVector* dg_dxdotdot = NULL;
  Epetra_MultiVector* dg_dp = NULL;
  for (int i=0; i<mp_x.size(); i++) {

    for (int k=0; k<mp_p_index.size(); k++) {
      int kk = mp_p_index[k];
      for (unsigned int j=0; j<pp[kk].size(); j++) {
	pp[kk][j].baseValue = mp_p_vals[kk][j].coeff(i);
	if (deriv_p != NULL) {
	  for (unsigned int l=0; l<deriv_p->size(); l++)
	    if ((*deriv_p)[l].family->getName() == pp[kk][j].family->getName())
	      (*deriv_p)[l].baseValue = pp[kk][j].baseValue;
	}
      }
    }

    if (mp_xdot != NULL)
      xdot = mp_xdot->getCoeffPtr(i).get();
    if (mp_xdotdot != NULL)
      xdotdot = mp_xdotdot->getCoeffPtr(i).get();
    if (mp_g != NULL)
      g = mp_g->getCoeffPtr(i).get();
    if (mp_dg_dx != NULL)
      dg_dx = mp_dg_dx->getCoeffPtr(i).get();
    if(mp_dg_dxdot != NULL)
      dg_dxdot = mp_dg_dxdot->getCoeffPtr(i).get();
    if(mp_dg_dxdotdot != NULL)
      dg_dxdotdot = mp_dg_dxdotdot->getCoeffPtr(i).get();
    if (mp_dg_dp != NULL)
      dg_dp = mp_dg_dp->getCoeffPtr(i).get();
    
    // Evaluate response function
    evaluateGradient(current_time, xdot, xdotdot,  mp_x[i], pp, deriv_p, 
		     g, dg_dx, dg_dxdot, dg_dxdotdot, dg_dp);
  }
}

 void checkInputs(void) const {
   int pSize = prob_.size(), cSize = coeff_.size();
   TEUCHOS_TEST_FOR_EXCEPTION((pSize!=cSize),std::invalid_argument,
     ">>> ERROR (ROL::MixedQuantileQuadrangle): Probability and coefficient arrays have different sizes!");
   Real sum(0), zero(0), one(1);
   for (int i = 0; i < pSize; i++) {
     TEUCHOS_TEST_FOR_EXCEPTION((prob_[i]>one || prob_[i]<zero), std::invalid_argument,
       ">>> ERROR (ROL::MixedQuantileQuadrangle): Element of probability array out of range!");
     TEUCHOS_TEST_FOR_EXCEPTION((coeff_[i]>one || coeff_[i]<zero), std::invalid_argument,
       ">>> ERROR (ROL::MixedQuantileQuadrangle): Element of coefficient array out of range!");
     sum += coeff_[i];
   }
   TEUCHOS_TEST_FOR_EXCEPTION((std::abs(sum-one) > std::sqrt(ROL_EPSILON<Real>())),std::invalid_argument,
     ">>> ERROR (ROL::MixedQuantileQuadrangle): Coefficients do not sum to one!");
   TEUCHOS_TEST_FOR_EXCEPTION(plusFunction_ == Teuchos::null, std::invalid_argument,
     ">>> ERROR (ROL::MixedQuantileQuadrangle): PlusFunction pointer is null!");
 }

 void initialize(void) {
   size_ = prob_.size();
   // Initialize temporary storage
   Real zero(0);
   xvar_.clear(); xvar_.resize(size_,zero);
   vvar_.clear(); vvar_.resize(size_,zero);
   vec_.clear();  vec_.resize(size_,zero);
 }

 void interpolate(const std::vector<double>& positions,
   std::vector<double>& results) const 
   {
     Teuchos::Array<double> in(positions.size());
     for (int i=0; i<in.size(); i++) in[i] = positions[i];
     Teuchos::Array<double> out;
     interpolate(in, out);
     for (int i=0; i<out.size(); i++) results[i] = out[i];
   }