// ========================================================================
void process_surface_entity(Ioss::SideSet *sset, stk::mesh::MetaData &meta,
                            stk::mesh::EntityRank sset_rank)
{
    assert(sset->type() == Ioss::SIDESET);
    Ioss::SideSet *fs = dynamic_cast<Ioss::SideSet *>(sset);
    assert(fs != NULL);
    const Ioss::SideBlockContainer& blocks = fs->get_side_blocks();
    stk::io::default_part_processing(blocks, meta, sset_rank);

    stk::mesh::Part* const fs_part = meta.get_part(sset->name());
    STKIORequire(fs_part != NULL);

    stk::mesh::Field<double, stk::mesh::ElementNode> *distribution_factors_field = NULL;
    bool surface_df_defined = false; // Has the surface df field been defined yet?


    int block_count = sset->block_count();
    for (int i=0; i < block_count; i++) {
        Ioss::SideBlock *side_block = sset->get_block(i);
        if (stk::io::include_entity(side_block)) {
            stk::mesh::Part * const side_block_part = meta.get_part(side_block->name());
            STKIORequire(side_block_part != NULL);
            meta.declare_part_subset(*fs_part, *side_block_part);

            const stk::mesh::EntityRank part_rank = side_block_part->primary_entity_rank();

            if (side_block->field_exists("distribution_factors")) {
                if (!surface_df_defined) {
                    std::string field_name = sset->name() + "_distribution_factors";
                    distribution_factors_field =
                        &meta.declare_field<stk::mesh::Field<double, stk::mesh::ElementNode> >(static_cast<stk::topology::rank_t>(part_rank), field_name);
                    stk::io::set_distribution_factor_field(*fs_part, *distribution_factors_field);
                    surface_df_defined = true;
                }
                stk::io::set_distribution_factor_field(*side_block_part, *distribution_factors_field);
                int side_node_count = side_block->topology()->number_nodes();
                stk::mesh::put_field(*distribution_factors_field,
                                     *side_block_part, side_node_count);
            }

            /** \todo IMPLEMENT truly handle fields... For this case we
             * are just defining a field for each transient field that is
             * present in the mesh...
             */
            stk::io::define_io_fields(side_block, Ioss::Field::TRANSIENT,
                                      *side_block_part,
                                      part_rank);
        }
    }
}
Esempio n. 2
0
    inline void StkMeshIoBroker::add_heartbeat_global(size_t index,
						      const std::string &name,
						      const boost::any *value,
						      stk::util::ParameterType::Type type)
    {
      STKIORequire(index < m_heartbeat.size());
      m_heartbeat[index]->add_global_ref(name, value, type);
    }
// ========================================================================
void process_elementblocks(Ioss::Region &region, stk::mesh::BulkData &bulk)
{
    const Ioss::ElementBlockContainer& elem_blocks = region.get_element_blocks();

    for(Ioss::ElementBlockContainer::const_iterator it = elem_blocks.begin();
            it != elem_blocks.end(); ++it) {
        Ioss::ElementBlock *entity = *it;

        if (stk::io::include_entity(entity)) {
            const std::string &name = entity->name();
            const stk::mesh::MetaData& meta = stk::mesh::MetaData::get(bulk);
            stk::mesh::Part* const part = meta.get_part(name);
            STKIORequire(part != NULL);

            const stk::topology topo = part->topology();
            if (topo == stk::topology::INVALID_TOPOLOGY) {
                std::ostringstream msg ;
                msg << " INTERNAL_ERROR: Part " << part->name() << " returned INVALID from get_topology()";
                throw std::runtime_error( msg.str() );
            }

            std::vector<int> elem_ids ;
            std::vector<int> connectivity ;

            entity->get_field_data("ids", elem_ids);
            entity->get_field_data("connectivity", connectivity);

            size_t element_count = elem_ids.size();
            int nodes_per_elem = topo.num_nodes();

            stk::mesh::EntityIdVector connectivity2(nodes_per_elem);

            std::vector<int>::const_iterator connBegin = connectivity.begin();
            std::vector<stk::mesh::Entity> elements(element_count);
            for(size_t i=0; i<element_count; ++i, connBegin += nodes_per_elem) {
                std::copy(connBegin, connBegin + nodes_per_elem, connectivity2.begin());
                elements[i] = stk::mesh::declare_element(bulk, *part, elem_ids[i], connectivity2);
            }

            // For this example, we are just taking all attribute fields
            // found on the io database and populating fields on the
            // corresponding mesh part.  In practice, would probably be
            // selective about which attributes to use...
            Ioss::NameList names;
            entity->field_describe(Ioss::Field::ATTRIBUTE, &names);
            for (Ioss::NameList::const_iterator I = names.begin(); I != names.end(); ++I) {
                if (*I == "attribute" && names.size() > 1)
                    continue;
                stk::mesh::FieldBase *field = meta.get_field<stk::mesh::FieldBase>(stk::topology::ELEMENT_RANK, *I);
                stk::io::field_data_from_ioss(bulk, field, elements, entity, *I);

            }
        }
    }
}
// ========================================================================
void process_nodesets(Ioss::Region &region, stk::mesh::MetaData &meta)
{
    const Ioss::NodeSetContainer& node_sets = region.get_nodesets();
    stk::io::default_part_processing(node_sets, meta, stk::topology::NODE_RANK);

    /** \todo REFACTOR should "distribution_factor" be a default field
     * that is automatically declared on all objects that it exists
     * on as is done in current framework?
     */
    stk::mesh::Field<double> & distribution_factors_field =
        meta.declare_field<stk::mesh::Field<double> >(stk::topology::NODE_RANK, "distribution_factors");

    /** \todo REFACTOR How to associate distribution_factors field
     * with the nodeset part if a node is a member of multiple
     * nodesets
     */

    for(Ioss::NodeSetContainer::const_iterator it = node_sets.begin();
            it != node_sets.end(); ++it) {
        Ioss::NodeSet *entity = *it;

        if (stk::io::include_entity(entity)) {
            stk::mesh::Part* const part = meta.get_part(entity->name());
            STKIORequire(part != NULL);
            STKIORequire(entity->field_exists("distribution_factors"));

            stk::mesh::put_field(distribution_factors_field, *part);

            /** \todo IMPLEMENT truly handle fields... For this case we
             * are just defining a field for each transient field that is
             * present in the mesh...
             */
            stk::io::define_io_fields(entity, Ioss::Field::TRANSIENT,
                                      *part,
                                      part->primary_entity_rank() ) ;
        }
    }
}
// ========================================================================
void process_nodesets(Ioss::Region &region, stk::mesh::BulkData &bulk)
{
    // Should only process nodes that have already been defined via the element
    // blocks connectivity lists.
    const Ioss::NodeSetContainer& node_sets = region.get_nodesets();

    for(Ioss::NodeSetContainer::const_iterator it = node_sets.begin();
            it != node_sets.end(); ++it) {
        Ioss::NodeSet *entity = *it;

        if (stk::io::include_entity(entity)) {
            const std::string & name = entity->name();
            const stk::mesh::MetaData& meta = stk::mesh::MetaData::get(bulk);
            stk::mesh::Part* const part = meta.get_part(name);
            STKIORequire(part != NULL);
            stk::mesh::PartVector add_parts( 1 , part );

            std::vector<int> node_ids ;
            int node_count = entity->get_field_data("ids", node_ids);

            std::vector<stk::mesh::Entity> nodes(node_count);
            for(int i=0; i<node_count; ++i) {
                nodes[i] = bulk.get_entity( stk::topology::NODE_RANK, node_ids[i] );
                if (bulk.is_valid(nodes[i])) {
                    bulk.declare_entity(stk::topology::NODE_RANK, node_ids[i], add_parts );
                }
            }

            /** \todo REFACTOR Application would probably store this field
             * (and others) somewhere after the declaration instead of
             * looking it up each time it is needed.
             */
            stk::mesh::Field<double> *df_field =
                meta.get_field<stk::mesh::Field<double> >(stk::topology::NODE_RANK, "distribution_factors");

            if (df_field != NULL) {
                stk::io::field_data_from_ioss(bulk, df_field, nodes, entity, "distribution_factors");
            }
        }
    }
}
// ========================================================================
void process_elementblocks(Ioss::Region &region, stk::mesh::MetaData &meta)
{
    const stk::mesh::EntityRank element_rank = stk::topology::ELEMENT_RANK;

    const Ioss::ElementBlockContainer& elem_blocks = region.get_element_blocks();
    stk::io::default_part_processing(elem_blocks, meta, element_rank);

    // Parts were created above, now handle element block specific
    // information (topology, attributes, ...);
    for(Ioss::ElementBlockContainer::const_iterator it = elem_blocks.begin();
            it != elem_blocks.end(); ++it) {
        Ioss::ElementBlock *entity = *it;

        if (stk::io::include_entity(entity)) {
            stk::mesh::Part* const part = meta.get_part(entity->name());
            STKIORequire(part != NULL);

            const stk::mesh::EntityRank part_rank = part->primary_entity_rank();

            // Element Block attributes (if any)...
            /** \todo IMPLEMENT truly handle attribute fields... For this
             * case we are just defining a field for each attribute field
             * that is present in the mesh...
             */
            stk::io::define_io_fields(entity, Ioss::Field::ATTRIBUTE,
                                      *part,
                                      part_rank);

            /** \todo IMPLEMENT truly handle fields... For this case we
             * are just defining a field for each transient field that is
             * present in the mesh...
             */
            stk::io::define_io_fields(entity, Ioss::Field::TRANSIENT,
                                      *part,
                                      part_rank);
        }
    }
}
Esempio n. 7
0
 inline void StkMeshIoBroker::process_heartbeat_output(size_t index, int step, double time)
 {
   STKIORequire(index < m_heartbeat.size());
   m_heartbeat[index]->process_output(step, time);
 }