#include "tutorials/scl-12/src/electrostatics.hpp"
|  | 
|  | OpEnergyDensity (const std::string &field_name, boost::shared_ptr< MatrixDouble > grad_u_negative, boost::shared_ptr< VectorDouble > energy_densiity_ptr, boost::shared_ptr< std::map< int, BlockData > > perm_block_sets_ptr, boost::shared_ptr< DataAtIntegrationPts > common_data_ptr) | 
|  | 
| MoFEMErrorCode | doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data) | 
|  | 
◆ OpEnergyDensity()
  
  | 
        
          | OpEnergyDensity::OpEnergyDensity | ( | const std::string & | field_name, |  
          |  |  | boost::shared_ptr< MatrixDouble > | grad_u_negative, |  
          |  |  | boost::shared_ptr< VectorDouble > | energy_densiity_ptr, |  
          |  |  | boost::shared_ptr< std::map< int, BlockData > > | perm_block_sets_ptr, |  
          |  |  | boost::shared_ptr< DataAtIntegrationPts > | common_data_ptr |  
          |  | ) |  |  |  | inline | 
 
Definition at line 197 of file electrostatics.hpp.
  206    std::fill(&doEntities[MBVERTEX], &doEntities[MBMAXTYPE], false);
  207    doEntities[MBVERTEX] = true;
  208  }
DomainEle::UserDataOperator DomainEleOp
constexpr auto field_name
boost::shared_ptr< MatrixDouble > gradUNegative
boost::shared_ptr< DataAtIntegrationPts > commonDataPtr
boost::shared_ptr< std::map< int, BlockData > > permBlockSetsPtr
boost::shared_ptr< VectorDouble > energyDensity
 
 
 
◆ doWork()
Definition at line 210 of file electrostatics.hpp.
  213                                                            {
  215 
  216    size_t nb_gauss_pts = getGaussPts().size2();
  220    auto t_negative_grad_u = getFTensor1FromMat<SPACE_DIM>(*
gradUNegative);
 
  222    double blockPermittivity = 0.0;
  224      if (
n.second.domainEnts.find(getFEEntityHandle()) !=
 
  225          n.second.domainEnts.end()) {
 
  226        blockPermittivity = 
n.second.epsPermit;
 
  227      }
  228    }
  229    
  230    for (int gg = 0; gg != nb_gauss_pts; gg++) {
  231      t_energy_density =
  232          0.5 * t_negative_grad_u(
I) * t_negative_grad_u(
I) * blockPermittivity;
 
  233      ++t_negative_grad_u;
  234      ++t_energy_density;
  235    }
  236 
  238  }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
const double n
refractive index of diffusive medium
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
constexpr IntegrationType I
 
 
◆ commonDataPtr
◆ energyDensity
  
  | 
        
          | boost::shared_ptr<VectorDouble> OpEnergyDensity::energyDensity |  | private | 
 
 
◆ gradUNegative
  
  | 
        
          | boost::shared_ptr<MatrixDouble> OpEnergyDensity::gradUNegative |  | private | 
 
 
◆ permBlockSetsPtr
  
  | 
        
          | boost::shared_ptr<std::map<int, BlockData> > OpEnergyDensity::permBlockSetsPtr |  | private | 
 
 
The documentation for this struct was generated from the following file: