C++ teuchos::Array类代码示例

  /** \brief Constructor.

      @param[in]     parlist is a parameter list specifying inputs

      parlist should contain sublists "SOL"->"Risk Measure"->"Convex Combination Risk Measure" and
      within the "Convex Combination Risk Measure" sublist should have the following parameters
      \li "Convex Combination Parameters" (greater than 0 and sum to 1)
      \li Sublists labeled 1 to n with risk measure definitions.
  */
  ConvexCombinationRiskMeasure(Teuchos::ParameterList &parlist)
    : RiskMeasure<Real>(), size_(0), firstReset_(true) {
    Teuchos::ParameterList &list
      = parlist.sublist("SOL").sublist("Risk Measure").sublist("Convex Combination Risk Measure");
    // Get convex combination parameters
    Teuchos::Array<Real> lambda
      = Teuchos::getArrayFromStringParameter<Real>(list,"Convex Combination Parameters");
    lambda_ = lambda.toVector();
    size_ = lambda_.size();
    // Build risk measures
    risk_.clear(); risk_.resize(size_,Teuchos::null);
    parlist_.clear(); parlist_.resize(size_);
    for (uint i = 0; i < size_; ++i) {
      std::ostringstream convert;
      convert << i;
      std::string si = convert.str();
      Teuchos::ParameterList &ilist = list.sublist(si);
      std::string name = ilist.get<std::string>("Name");
      parlist_[i].sublist("SOL").sublist("Risk Measure").set("Name",name);
      parlist_[i].sublist("SOL").sublist("Risk Measure").sublist(name) = ilist;
      risk_[i] = RiskMeasureFactory<Real>(parlist_[i]);
    }
    // Check inputs
    checkInputs();
  }

Adapter *
Adapter::initialize(const Teuchos::RCP<Teuchos::ParameterList>& catalystParams)
{
  // Validate parameters against list for this specific class
  catalystParams->validateParameters(*getValidAdapterParameters(),0);

  if (Adapter::instance)
    delete Adapter::instance;
  Adapter::instance = new Adapter();

  // Register our Grid class with Catalyst so that it can be used in a pipeline.
  if (vtkClientServerInterpreterInitializer *intInit =
      vtkClientServerInterpreterInitializer::GetInitializer()) {
    if (vtkClientServerInterpreter *interp = intInit->GetGlobalInterpreter()) {
      interp->AddNewInstanceFunction("Grid", Private::MakeGrid);
    }
  }

  // Load pipeline file
  Teuchos::Array<std::string> files =
      catalystParams->get<Teuchos::Array<std::string> >("Pipeline Files");
  typedef Teuchos::Array<std::string>::const_iterator FileIterT;
  for (FileIterT it = files.begin(), itEnd = files.end(); it != itEnd; ++it)
    Adapter::instance->addPythonScriptPipeline(*it);

  return Adapter::instance;
}

Teuchos::Array<int> getMyBlockLIDs(
    const Teuchos::ParameterList &blockingParams,
    const Albany::AbstractDiscretization &disc)
{
  Teuchos::Array<int> result;

  const std::string nodeSetLabel = blockingParams.get<std::string>("Node Set");
  const int dofRank = blockingParams.get<int>("Dof");

  const Albany::NodeSetList &nodeSets = disc.getNodeSets();
  const Albany::NodeSetList::const_iterator it = nodeSets.find(nodeSetLabel);
  TEUCHOS_TEST_FOR_EXCEPT(it == nodeSets.end());
  {
    typedef Albany::NodeSetList::mapped_type NodeSetEntryList;
    const NodeSetEntryList &nodeEntries = it->second;

    for (NodeSetEntryList::const_iterator jt = nodeEntries.begin(); jt != nodeEntries.end(); ++jt) {
      typedef NodeSetEntryList::value_type NodeEntryList;
      const NodeEntryList &entries = *jt;
      result.push_back(entries[dofRank]);
    }
  }

  return result;
}

void FunctionTests::checkFunctionsAgree(FunctionPtr f1, FunctionPtr f2, BasisCachePtr basisCache) {
  ASSERT_EQ(f1->rank(), f2->rank())
    << "f1->rank() " << f1->rank() << " != f2->rank() " << f2->rank() << endl;

  int rank = f1->rank();
  int numCells = basisCache->getPhysicalCubaturePoints().dimension(0);
  int numPoints = basisCache->getPhysicalCubaturePoints().dimension(1);
  int spaceDim = basisCache->getPhysicalCubaturePoints().dimension(2);
  Teuchos::Array<int> dim;
  dim.append(numCells);
  dim.append(numPoints);
  for (int i=0; i<rank; i++) {
    dim.append(spaceDim);
  }
  FieldContainer<double> f1Values(dim);
  FieldContainer<double> f2Values(dim);
  f1->values(f1Values,basisCache);
  f2->values(f2Values,basisCache);
  
  double tol = 1e-14;
  double maxDiff;
  EXPECT_TRUE(fcsAgree(f1Values,f2Values,tol,maxDiff))
    << "Test failed: f1 and f2 disagree; maxDiff " << maxDiff << ".\n"
    << "f1Values: \n" << f1Values
    << "f2Values: \n" << f2Values;
}

// ============================================================================
Teuchos::RCP<Recti::Domain::Abstract>
Recti::Domain::Factory::
buildPolygon() const
{
    const Teuchos::ParameterList & pointsList =
            pList_.sublist("Vertices");
    
    Teuchos::Array<Point> vertices;
    
    Teuchos::ParameterList::ConstIterator k;
    for ( k=pointsList.begin(); k!=pointsList.end(); ++k )
    {
        Teuchos::ParameterEntry e = pointsList.entry(k);
        TEUCHOS_ASSERT( e.isList() );
        
        std::string label = pointsList.name(k);
        const Teuchos::ParameterList & coordinates =
                pointsList.sublist(label);
        
        Point X = Teuchos::tuple( coordinates.get<double>("x"),
                                  coordinates.get<double>("y"),
                                  0.0
                                );
        vertices.push_back( X );
    }
    
    return Teuchos::rcp( new Polygon( vertices ) );
}

  Teuchos::RCP<const Tpetra::Map<LO,GO,Node> >
  createMeshMap (const LO& blockSize, const Tpetra::Map<LO,GO,Node>& pointMap)
  {
    typedef Teuchos::OrdinalTraits<Tpetra::global_size_t> TOT;
    typedef Tpetra::Map<LO,GO,Node> map_type;

    //calculate mesh GIDs
    Teuchos::ArrayView<const GO> pointGids = pointMap.getNodeElementList();
    Teuchos::Array<GO> meshGids;
    GO indexBase = pointMap.getIndexBase();

    // Use hash table to track whether we've encountered this GID previously.  This will happen
    // when striding through the point DOFs in a block.  It should not happen otherwise.
    // I don't use sort/make unique because I don't want to change the ordering.
    meshGids.reserve(pointGids.size());
    Tpetra::Details::HashTable<GO,int> hashTable(pointGids.size());
    for (int i=0; i<pointGids.size(); ++i) {
      GO meshGid = (pointGids[i]-indexBase) / blockSize + indexBase;
      if (hashTable.get(meshGid) == -1) {
        hashTable.add(meshGid,1);   //(key,value)
        meshGids.push_back(meshGid);
      }
    }

    Teuchos::RCP<const map_type> meshMap = Teuchos::rcp( new map_type(TOT::invalid(), meshGids(), 0, pointMap.getComm()) );
    return meshMap;

  }

Teuchos::RCP<const Tpetra_Map>
Albany::SolutionResponseFunction::
buildCulledMapT(const Tpetra_Map& x_mapT,
	       const Teuchos::Array<int>& keepDOF) const
{
  int numKeepDOF = std::accumulate(keepDOF.begin(), keepDOF.end(), 0);
  int Neqns = keepDOF.size();
  int N = x_mapT.getNodeNumElements(); // x_mapT 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

  Teuchos::ArrayView<const GO> gidsT = x_mapT.getNodeElementList();
  Teuchos::Array<GO> 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++] = gidsT[(inode*Neqns)+ieqn];
  // end cull

  Teuchos::RCP<const Tpetra_Map> x_new_mapT = Tpetra::createNonContigMapWithNode<LO, GO, KokkosNode> (gids_new, x_mapT.getComm(), x_mapT.getNode());

  return x_new_mapT;

}

//---------------------------------------------------------------------------//
// Get the parameteric center of an entity.
void MoabEntityLocalMap::parametricCenter(
    const Entity& entity,
    Teuchos::Array<double>& center ) const
{
    moab::EntityType moab_type = d_moab_mesh->get_moab()->type_from_handle(
	MoabHelpers::extractEntity(entity) );

    switch( moab_type )
    {
	case moab::MBTRI:
	    center.assign( 2, 1.0 / 3.0 );
	    break;
	case moab::MBQUAD:
	    center.assign( 2, 0.0 );
	    break;
	case moab::MBTET:
	    center.assign( 3, 1.0 / 6.0 );
	    break;
	case moab::MBHEX:
	    center.assign( 3, 0.0 );
	    break;
	default:
	    center.resize( 0 );
	    break;
    }
}

// =============================================================================
Teuchos::RCP<Tpetra::Map<ORD> >
VIO::TpetraMesh::Mesh::
getElemsToNodesMap_() const
{
  // create list of elements that need to be accessible from this process


  // Make sure that *all entries that belong to any of the elements in
  // this core are accessible.
  // First mark all the nodes that need to be accessible:
  TEUCHOS_ASSERT( !elems_.is_null() );
  TEUCHOS_ASSERT( !nodes_.is_null() );
  int numNodes = nodes_.size();
  Teuchos::Array<bool> mustBeAccessible( numNodes );
  for ( int k=0; k<elems_.size(); k++ )
      for ( int l=0; l<elems_[k].size(); l++ )
          mustBeAccessible[ elems_[k][l] ] = true;
  // now create the list
  Teuchos::Array<ORD> entryList;
  for ( int k=0; k<numNodes; k++ )
      if ( mustBeAccessible[k] )
          entryList.append( k );

  Teuchos::RCP<Tpetra::Map<ORD> > map =
      Teuchos::rcp( new Tpetra::Map<ORD>( Teuchos::OrdinalTraits<ORD>::invalid(), entryList(), 0, comm_ )
                  );

  return map;
}

/*!
 * \brief Find the set of entities a point neighbors.
 */
void CoarseLocalSearch::search( const Teuchos::ArrayView<const double>& point,
				const Teuchos::ParameterList& parameters,
				Teuchos::Array<Entity>& neighbors ) const
{
    // Find the leaf of nearest neighbors.
    int num_neighbors = 100;
    if ( parameters.isParameter("Coarse Local Search kNN") )
    {
	num_neighbors = parameters.get<int>("Coarse Local Search kNN");
    }
    num_neighbors = 
	std::min( num_neighbors, Teuchos::as<int>(d_entity_map.size()) );
    Teuchos::Array<unsigned> local_neighbors = 
	d_tree->nnSearch( point, num_neighbors );

    // Extract the neighbors.
    neighbors.resize( local_neighbors.size() );
    Teuchos::Array<unsigned>::const_iterator local_it;
    Teuchos::Array<Entity>::iterator entity_it;
    for ( local_it = local_neighbors.begin(),
	 entity_it = neighbors.begin();
	  local_it != local_neighbors.end();
	  ++local_it, ++entity_it )
    {
	DTK_CHECK( d_entity_map.count(*local_it) );
	*entity_it = d_entity_map.find(*local_it)->second;
    }
}

Teuchos::RCP<const Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node> > create_tridiag_graph(LocalOrdinal num_rows_per_proc)
{
  Teuchos::RCP<const Teuchos::Comm<int> > comm = getDefaultComm();

  const global_size_t INVALID = Teuchos::OrdinalTraits<global_size_t>::invalid();
  const LocalOrdinal indexBase = 0;

  Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rowmap = Teuchos::rcp(new Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node>(INVALID, num_rows_per_proc, indexBase, comm));

  Teuchos::RCP<Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node> > crsgraph = Teuchos::rcp(new Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>(rowmap, 0));

  Teuchos::Array<GlobalOrdinal> cols;
 for(LocalOrdinal l_row = 0; (size_t) l_row<rowmap->getNodeNumElements(); l_row++) {
    GlobalOrdinal g_row = rowmap->getGlobalElement(l_row);
    if (g_row == rowmap->getMinAllGlobalIndex() ||
        g_row == rowmap->getMaxAllGlobalIndex()) cols.resize(2);
    else cols.resize(3);

    size_t coloffset = 0;
    if (g_row > rowmap->getMinAllGlobalIndex()) cols[coloffset++] = g_row-1;
    cols[coloffset++] = g_row;
    if (g_row < rowmap->getMaxAllGlobalIndex()) cols[coloffset++] = g_row+1;

    crsgraph->insertGlobalIndices(g_row, cols());
  }

  crsgraph->fillComplete();

  return crsgraph;
}

void FELIX::Hydrology::constructNeumannEvaluators (const Teuchos::RCP<Albany::MeshSpecsStruct>& meshSpecs)
{
    // Note: we only enter this function if sidesets are defined in the mesh file
    // i.e. meshSpecs.ssNames.size() > 0

    Albany::BCUtils<Albany::NeumannTraits> nbcUtils;

    // Check to make sure that Neumann BCs are given in the input file
    if (!nbcUtils.haveBCSpecified(this->params))
    {
        return;
    }

    // Construct BC evaluators for all side sets and names
    // Note that the string index sets up the equation offset, so ordering is important

    std::vector<std::string> neumannNames(neq + 1);
    Teuchos::Array<Teuchos::Array<int> > offsets;
    offsets.resize(1);

    neumannNames[0] = "Hydraulic Potential";
    neumannNames[1] = "all";
    offsets[0].resize(1);
    offsets[0][0] = 0;

    // Construct BC evaluators for all possible names of conditions
    std::vector<std::string> condNames(1);
    condNames[0] = "neumann";

    nfm.resize(1); // FELIX problem only has one element block

    nfm[0] = nbcUtils.constructBCEvaluators(meshSpecs, neumannNames, dof_names, false, 0,
                                            condNames, offsets, dl,
                                            this->params, this->paramLib);
}

 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);
   }
 }

void BilinearFormUtility::weightCellBasisValues(FieldContainer<double> &basisValues, const FieldContainer<double> &weights, int offset) {
  // weights are (numCells, offset+numFields)
  // basisValues are (numCells, numFields, ...)
  int numCells = basisValues.dimension(0);
  int numFields = basisValues.dimension(1);
  
  Teuchos::Array<int> dimensions;
  basisValues.dimensions(dimensions);

  int numAffectedValues = 1;
  for (int dimIndex=2; dimIndex<dimensions.size(); dimIndex++) {
    numAffectedValues *= dimensions[dimIndex];
  }
  
  Teuchos::Array<int> index(dimensions.size(),0);
  
  for (int cellIndex=0; cellIndex < numCells; cellIndex++) {
    index[0] = cellIndex;
    for (int fieldIndex=0; fieldIndex < numFields; fieldIndex++) {
      index[1] = fieldIndex;
      int enumIndex = basisValues.getEnumeration(index);
      for (int valIndex=enumIndex; valIndex < numAffectedValues + enumIndex; valIndex++) {
        basisValues[valIndex] *= weights(cellIndex,fieldIndex+offset);
      }
    }
  }
  
}

void
Stokhos::PecosOneDOrthogPolyBasis<ordinal_type,value_type>::
getQuadPoints(ordinal_type quad_order,
	      Teuchos::Array<value_type>& quad_points,
	      Teuchos::Array<value_type>& quad_weights,
	      Teuchos::Array< Teuchos::Array<value_type> >& quad_values) const
{
  // This assumes Gaussian quadrature
  ordinal_type num_points = 
    static_cast<ordinal_type>(std::ceil((quad_order+1)/2.0));
  const Pecos::RealArray& gp = pecosPoly->collocation_points(num_points);
  const Pecos::RealArray& gw = pecosPoly->type1_collocation_weights(num_points); 
  quad_points.resize(num_points);
  quad_weights.resize(num_points);
  for (ordinal_type i=0; i<num_points; i++) {
    quad_points[i] = gp[i];
    quad_weights[i] = gw[i];
  }
  
  // Evalute basis at gauss points
  quad_values.resize(num_points);
  for (ordinal_type i=0; i<num_points; i++) {
    quad_values[i].resize(p+1);
    evaluateBases(quad_points[i], quad_values[i]);
  }
}