Nektar::Utilities::ProcessEquiSpacedOutput Class Reference

This processing module interpolates one field to another. More...

#include <ProcessEquiSpacedOutput.h>

Inheritance diagram for Nektar::Utilities::ProcessEquiSpacedOutput:

Collaboration diagram for Nektar::Utilities::ProcessEquiSpacedOutput:

Public Member Functions
	ProcessEquiSpacedOutput (FieldSharedPtr f)
virtual	~ProcessEquiSpacedOutput ()
virtual void	Process (po::variables_map &vm)
	Write mesh to output file.
Public Member Functions inherited from Nektar::Utilities::ProcessModule
	ProcessModule ()
	ProcessModule (FieldSharedPtr p_f)
	ProcessModule (MeshSharedPtr p_m)
Public Member Functions inherited from Nektar::Utilities::Module
	Module (FieldSharedPtr p_f)
void	RegisterConfig (string key, string value)
	Register a configuration option with a module.
void	PrintConfig ()
	Print out all configuration options for a module.
void	SetDefaults ()
	Sets default configuration options for those which have not been set.
bool	GetRequireEquiSpaced (void)
void	SetRequireEquiSpaced (bool pVal)
void	EvaluateTriFieldAtEquiSpacedPts (LocalRegions::ExpansionSharedPtr &exp, const Array< OneD, const NekDouble > &infield, Array< OneD, NekDouble > &outfield)
	Module (MeshSharedPtr p_m)
virtual void	Process ()=0
void	RegisterConfig (string key, string value)
void	PrintConfig ()
void	SetDefaults ()
MeshSharedPtr	GetMesh ()
virtual void	ProcessVertices ()
	Extract element vertices.

Static Public Member Functions
static boost::shared_ptr< Module >	create (FieldSharedPtr f)
	Creates an instance of this class.

Static Public Attributes
static ModuleKey	className

Protected Member Functions
	ProcessEquiSpacedOutput ()
void	SetupEquiSpacedField (void)

Additional Inherited Members
Protected Attributes inherited from Nektar::Utilities::Module
FieldSharedPtr	m_f
	Field object.
map< string, ConfigOption >	m_config
	List of configuration values.
bool	m_requireEquiSpaced
MeshSharedPtr	m_mesh
	Mesh object.

Detailed Description

This processing module interpolates one field to another.

Definition at line 48 of file ProcessEquiSpacedOutput.h.

Constructor & Destructor Documentation

Nektar::Utilities::ProcessEquiSpacedOutput::ProcessEquiSpacedOutput

(

FieldSharedPtr

f

)

Definition at line 57 of file ProcessEquiSpacedOutput.cpp.

References Nektar::Utilities::Module::m_config.

    : ProcessModule(f)
{
    f->m_fieldPts = MemoryManager<FieldPts>::AllocateSharedPtr();
    f->m_setUpEquiSpacedFields = true;

    m_config["tetonly"] = ConfigOption(true, "NotSet",
                                "Only process tetrahedral elements");

}

Nektar::Utilities::ProcessEquiSpacedOutput::~ProcessEquiSpacedOutput ( )

virtual

Definition at line 68 of file ProcessEquiSpacedOutput.cpp.

{
}

Nektar::Utilities::ProcessEquiSpacedOutput::ProcessEquiSpacedOutput ( )

inlineprotected

Definition at line 64 of file ProcessEquiSpacedOutput.h.

{};

Member Function Documentation

static boost::shared_ptr<Module> Nektar::Utilities::ProcessEquiSpacedOutput::create ( FieldSharedPtr  f )

inlinestatic

Creates an instance of this class.

Reimplemented in Nektar::Utilities::ProcessIsoContour.

Definition at line 52 of file ProcessEquiSpacedOutput.h.

                                                                    {
                return MemoryManager<ProcessEquiSpacedOutput>::
                                                        AllocateSharedPtr(f);
            }

void Nektar::Utilities::ProcessEquiSpacedOutput::Process ( po::variables_map &  vm )

virtual

Write mesh to output file.

Implements Nektar::Utilities::Module.

Reimplemented in Nektar::Utilities::ProcessIsoContour.

Definition at line 72 of file ProcessEquiSpacedOutput.cpp.

References SetupEquiSpacedField().

{
    SetupEquiSpacedField();
}

void Nektar::Utilities::ProcessEquiSpacedOutput::SetupEquiSpacedField ( void  )

protected

Definition at line 78 of file ProcessEquiSpacedOutput.cpp.

References ASSERTL0, Nektar::StdRegions::eBwd, Nektar::StdRegions::eDir1BwdDir1_Dir2BwdDir2, Nektar::StdRegions::eDir1BwdDir1_Dir2FwdDir2, Nektar::StdRegions::eDir1BwdDir2_Dir2BwdDir1, Nektar::StdRegions::eDir1BwdDir2_Dir2FwdDir1, Nektar::StdRegions::eFwd, Nektar::LibUtilities::eHexahedron, Nektar::LibUtilities::ePrism, Nektar::Utilities::ePtsTetBlock, Nektar::Utilities::ePtsTriBlock, Nektar::LibUtilities::ePyramid, Nektar::LibUtilities::eQuadrilateral, Nektar::LibUtilities::eSegment, Nektar::LibUtilities::eTetrahedron, Nektar::LibUtilities::eTriangle, Nektar::LibUtilities::StdSegData::getNumberOfCoefficients(), Nektar::Utilities::Module::m_config, Nektar::Utilities::Module::m_f, npts, and Nektar::NullNekDouble1DArray.

Referenced by Process(), and Nektar::Utilities::ProcessIsoContour::Process().

{

    if(m_f->m_verbose)
    {
        cout << "Interpolating fields to equispaced" << endl;
    }

    int coordim  = m_f->m_exp[0]->GetCoordim(0);
    int shapedim = m_f->m_exp[0]->GetExp(0)->GetShapeDimension();
    int npts     = m_f->m_exp[0]->GetTotPoints();
    Array<OneD, Array<OneD, NekDouble> > coords(3);

    int nel = m_f->m_exp[0]->GetExpSize();

    // set up the number of points in each element
    int newpoints;
    int newtotpoints = 0;

    Array<OneD,int> conn;
    int prevNcoeffs = 0;
    int prevNpoints = 0;
    int cnt = 0;

    // identify face 1 connectivity for prisms
    map<int,StdRegions::Orientation > face0orient;
    set<int> prismorient;
    LocalRegions::ExpansionSharedPtr e;

    for(int i = 0; i < nel; ++i)
    {
        e = m_f->m_exp[0]->GetExp(i);
        if(e->DetShapeType() == LibUtilities::ePrism)
        {
            StdRegions::Orientation forient = e->GetForient(0);
            int fid = e->GetGeom()->GetFid(0);
            if(face0orient.count(fid))
            { // face 1 meeting face 1 so reverse this id
                prismorient.insert(i);
            }
            else
            {
                // just store if Dir 1 is fwd or bwd
                if((forient == StdRegions::eDir1BwdDir1_Dir2FwdDir2) ||
                   (forient == StdRegions::eDir1BwdDir1_Dir2BwdDir2) ||
                   (forient == StdRegions::eDir1BwdDir2_Dir2FwdDir1) ||
                   (forient == StdRegions::eDir1BwdDir2_Dir2BwdDir1))
                {
                    face0orient[fid] = StdRegions::eBwd;
                }
                else
                {
                    face0orient[fid] = StdRegions::eFwd;
                }
            }
        }
    }

    for(int i = 0; i < nel; ++i)
    {
        e = m_f->m_exp[0]->GetExp(i);
        if(e->DetShapeType() == LibUtilities::ePrism)
        {
            int fid = e->GetGeom()->GetFid(2);
            // check to see if face 2 meets face 1
            if(face0orient.count(fid))
            {
                // check to see how face 2 is orientated
                StdRegions::Orientation forient2 = e->GetForient(2);
                StdRegions::Orientation forient0 = face0orient[fid];

                // If dir 1 or forient2 is bwd then check agains
                // face 1 value
                if((forient2 == StdRegions::eDir1BwdDir1_Dir2FwdDir2) ||
                   (forient2 == StdRegions::eDir1BwdDir1_Dir2BwdDir2) ||
                   (forient2 == StdRegions::eDir1BwdDir2_Dir2FwdDir1) ||
                   (forient2 == StdRegions::eDir1BwdDir2_Dir2BwdDir1))
                {
                    if(forient0 == StdRegions::eFwd)
                    {
                        prismorient.insert(i);
                    }
                }
                else
                {
                    if(forient0 == StdRegions::eBwd)
                    {
                        prismorient.insert(i);
                    }
                }
            }
        }
    }

    for(int i = 0; i < nel; ++i)
    {
        e = m_f->m_exp[0]->GetExp(i);
        if(m_config["tetonly"].m_beenSet)
        {
            if(m_f->m_exp[0]->GetExp(i)->DetShapeType() !=
                    LibUtilities::eTetrahedron)
            {
                if(i == nel-1)
                {
                    break;
                }
                else
                {
                    continue;
                }
            }
        }

        switch(e->DetShapeType())
        {
        case LibUtilities::eSegment:
            {
                int npoints0 = e->GetBasis(0)->GetNumPoints();

                newpoints = LibUtilities::StdSegData::
                                    getNumberOfCoefficients(npoints0);
                m_f->m_fieldPts->m_npts.push_back(newpoints);
                newtotpoints += newpoints;
            }
            break;
        case LibUtilities::eTriangle:
            {
                int np0 = e->GetBasis(0)->GetNumPoints();
                int np1 = e->GetBasis(1)->GetNumPoints();
                int np = max(np0,np1);
                newpoints     = LibUtilities::StdTriData::
                                    getNumberOfCoefficients(np,np);
                m_f->m_fieldPts->m_npts.push_back(newpoints);
                newtotpoints += newpoints;
            }
            break;
        case LibUtilities::eQuadrilateral:
            {
                int np0 = e->GetBasis(0)->GetNumPoints();
                int np1 = e->GetBasis(1)->GetNumPoints();
                int np = max(np0,np1);

                newpoints  = LibUtilities::StdQuadData::
                                    getNumberOfCoefficients(np,np);
                m_f->m_fieldPts->m_npts.push_back(newpoints);
                newtotpoints += newpoints;
            }
            break;
        case LibUtilities::eTetrahedron:
            {
                int np0 = e->GetBasis(0)->GetNumPoints();
                int np1 = e->GetBasis(1)->GetNumPoints();
                int np2 = e->GetBasis(2)->GetNumPoints();
                int np = max(np0,max(np1,np2));

                newpoints  = LibUtilities::StdTetData::
                                    getNumberOfCoefficients(np,np,np);
                m_f->m_fieldPts->m_npts.push_back(newpoints);
                newtotpoints += newpoints;
            }
            break;
        case LibUtilities::ePrism:
            {
                int np0 = e->GetBasis(0)->GetNumPoints();
                int np1 = e->GetBasis(1)->GetNumPoints();
                int np2 = e->GetBasis(2)->GetNumPoints();
                int np = max(np0,max(np1,np2));

                newpoints  = LibUtilities::StdPrismData::
                                    getNumberOfCoefficients(np,np,np);
                m_f->m_fieldPts->m_npts.push_back(newpoints);
                newtotpoints += newpoints;
            }
            break;
        case LibUtilities::ePyramid:
            {
                int np0 = e->GetBasis(0)->GetNumPoints();
                int np1 = e->GetBasis(1)->GetNumPoints();
                int np2 = e->GetBasis(2)->GetNumPoints();
                int np = max(np0,max(np1,np2));

                newpoints     = LibUtilities::StdPyrData::
                                    getNumberOfCoefficients(np,np,np);
                m_f->m_fieldPts->m_npts.push_back(newpoints);
                newtotpoints += newpoints;
            }
            break;
        case LibUtilities::eHexahedron:
            {
                int np0 = e->GetBasis(0)->GetNumPoints();
                int np1 = e->GetBasis(1)->GetNumPoints();
                int np2 = e->GetBasis(2)->GetNumPoints();
                int np = max(np0,max(np1,np2));

                newpoints     = LibUtilities::StdPyrData::
                                    getNumberOfCoefficients(np,np,np);
                m_f->m_fieldPts->m_npts.push_back(newpoints);
                newtotpoints += newpoints;
            }
            break;
        default:
            {
                ASSERTL0(false,"Points not known");
            }
        }

        if(e->DetShapeType() == LibUtilities::ePrism)
        {
            bool standard = true;

            if(prismorient.count(i))
            {
                standard = false; // reverse direction
            }

            e->GetSimplexEquiSpacedConnectivity(conn,standard);
        }
        else
        {

            if((prevNcoeffs != e->GetNcoeffs()) ||
               (prevNpoints != e->GetTotPoints()))
            {
                prevNcoeffs = e->GetNcoeffs();
                prevNpoints = e->GetTotPoints();

                e->GetSimplexEquiSpacedConnectivity(conn);
            }
        }
        Array<OneD, int> newconn(conn.num_elements());
        for(int j = 0; j < conn.num_elements(); ++j)
        {
            newconn[j] = conn[j] + cnt;
        }

        m_f->m_fieldPts->m_ptsConn.push_back(newconn);
        cnt += newpoints;
    }

    m_f->m_fieldPts->m_ptsDim  = coordim;
    if(m_f->m_fielddef.size())
    {
        m_f->m_fieldPts->m_nFields = m_f->m_exp.size();
    }
    else // just the mesh points
    {
        m_f->m_fieldPts->m_nFields = 0;
    }

    m_f->m_fieldPts->m_pts = Array<OneD, Array<OneD, NekDouble> >(
                                        m_f->m_fieldPts->m_nFields + coordim);

    for(int i = 0; i < m_f->m_fieldPts->m_nFields + coordim; ++i)
    {
        m_f->m_fieldPts->m_pts[i] = Array<OneD, NekDouble>(newtotpoints);
    }

    // Interpolate coordinates
    for(int i = 0; i < coordim; ++i)
    {
        coords[i] = Array<OneD, NekDouble>(npts);
    }

    for(int i = coordim; i < 3; ++i)
    {
        coords[i] = NullNekDouble1DArray;
    }

    m_f->m_exp[0]->GetCoords(coords[0],coords[1],coords[2]);

    int nq1 = m_f->m_exp[0]->GetTotPoints();

    Array<OneD, NekDouble> x1(nq1);
    Array<OneD, NekDouble> y1(nq1);
    Array<OneD, NekDouble> z1(nq1);

    m_f->m_exp[0]->GetCoords(x1, y1, z1);

    if (shapedim == 2)
    {
        m_f->m_fieldPts->m_ptype = ePtsTriBlock;
    }
    else if (shapedim == 3)
    {
        m_f->m_fieldPts->m_ptype = ePtsTetBlock;
    }

    Array<OneD, NekDouble> tmp;

    for(int n = 0; n < coordim; ++n)
    {
        cnt = 0;
        int cnt1 = 0;
        for(int i = 0; i < nel; ++i)
        {
            m_f->m_exp[0]->GetExp(i)->PhysInterpToSimplexEquiSpaced(
                                        coords[n] + cnt,
                                        tmp = m_f->m_fieldPts->m_pts[n] + cnt1);
            cnt1 += m_f->m_fieldPts->m_npts[i];
            cnt  += m_f->m_exp[0]->GetExp(i)->GetTotPoints();
        }
    }

    if(m_f->m_fielddef.size())
    {
        ASSERTL0(m_f->m_fielddef[0]->m_fields.size() == m_f->m_exp.size(),
                 "More expansion defined than fields");

        for(int n = 0; n < m_f->m_exp.size(); ++n)
        {
            cnt = 0;
            int cnt1 = 0;
            Array<OneD, const NekDouble> phys = m_f->m_exp[n]->GetPhys();
            for(int i = 0; i < nel; ++i)
            {
                m_f->m_exp[0]->GetExp(i)->PhysInterpToSimplexEquiSpaced(
                            phys + cnt,
                            tmp = m_f->m_fieldPts->m_pts[coordim + n] + cnt1);
                cnt1 += m_f->m_fieldPts->m_npts[i];
                cnt  += m_f->m_exp[0]->GetExp(i)->GetTotPoints();
            }

            // Set up Variable string.
            m_f->m_fieldPts->m_fields.push_back(
                                    m_f->m_fielddef[0]->m_fields[n]);
        }
    }
}

Member Data Documentation

ModuleKey Nektar::Utilities::ProcessEquiSpacedOutput::className

static

Initial value:

    GetModuleFactory().RegisterCreatorFunction(
           ModuleKey(eProcessModule, "equispacedoutput"),
           ProcessEquiSpacedOutput::create,
           "Write data as equi-spaced output using simplices to represent the data for connecting points")

Definition at line 56 of file ProcessEquiSpacedOutput.h.