C++ EquationSystems::get_mesh_data方法代码示例

//.........这里部分代码省略.........
            // of the element.
            QGauss qface(dim-1, SECOND);

            // Tell the finte element object to use our
            // quadrature rule.
            fe_face->attach_quadrature_rule (&qface);

            // The value of the shape functions at the quadrature
            // points.
            const std::vector<std::vector<Real> > & phi_face =
              fe_face->get_phi();

            // The Jacobian// Quadrature Weight at the quadrature
            // points on the face.
            const std::vector<Real> & JxW_face = fe_face->get_JxW();

            // Compute the shape function values on the element
            // face.
            fe_face->reinit(elem, side);

            // For the Robin BCs consider a normal admittance an=1
            // at some parts of the bounfdary
            const Real an_value = 1.;

            // Loop over the face quadrature points for integration.
            for (unsigned int qp=0; qp<qface.n_points(); qp++)
              {

                // Element matrix contributrion due to precribed
                // admittance boundary conditions.
                for (unsigned int i=0; i<phi_face.size(); i++)
                  for (unsigned int j=0; j<phi_face.size(); j++)
                    Ce(i,j) += rho*an_value*JxW_face[qp]*phi_face[i][qp]*phi_face[j][qp];
              }
          }

      // If this assembly program were to be used on an adaptive mesh,
      // we would have to apply any hanging node constraint equations
      // by uncommenting the following lines:
      // std::vector<unsigned int> dof_indicesC = dof_indices;
      // std::vector<unsigned int> dof_indicesM = dof_indices;
      // dof_map.constrain_element_matrix (Ke, dof_indices);
      // dof_map.constrain_element_matrix (Ce, dof_indicesC);
      // dof_map.constrain_element_matrix (Me, dof_indicesM);


      // Finally, simply add the contributions to the additional
      // matrices and vector.
      stiffness.add_matrix (Ke, dof_indices);
      damping.add_matrix (Ce, dof_indices);
      mass.add_matrix (Me, dof_indices);

      // For the overall matrix, explicitly zero the entries where
      // we added values in the other ones, so that we have
      // identical sparsity footprints.
      matrix.add_matrix(zero_matrix, dof_indices);

    } // loop el


  // It now remains to compute the rhs. Here, we simply
  // get the normal velocities values on the boundary from the
  // mesh data.
  {
    LOG_SCOPE("rhs", "assemble_helmholtz");

    // get a reference to the mesh data.
    const MeshData & mesh_data = es.get_mesh_data();

    // We will now loop over all nodes. In case nodal data
    // for a certain node is available in the MeshData, we simply
    // adopt the corresponding value for the rhs vector.
    // Note that normal data was given in the mesh data file,
    // i.e. one value per node
    libmesh_assert_equal_to (mesh_data.n_val_per_node(), 1);

    MeshBase::const_node_iterator       node_it  = mesh.nodes_begin();
    const MeshBase::const_node_iterator node_end = mesh.nodes_end();

    for ( ; node_it != node_end; ++node_it)
      {
        // the current node pointer
        const Node * node = *node_it;

        // check if the mesh data has data for the current node
        // and do for all components
        if (mesh_data.has_data(node))
          for (unsigned int comp=0; comp<node->n_comp(0, 0); comp++)
            {
              // the dof number
              unsigned int dn = node->dof_number(0, 0, comp);

              // set the nodal value
              freq_indep_rhs.set(dn, mesh_data.get_data(node)[0]);
            }
      }
  }

  // All done!
}