int main(int argc, char *argv[])
{
    LibUtilities::SessionReaderSharedPtr session;
    string vDriverModule;
    DriverSharedPtr drv;

    try
    {
        // Create session reader.
        session = LibUtilities::SessionReader::CreateInstance(argc, argv);

        // Create driver
        session->LoadSolverInfo("Driver", vDriverModule, "Standard");
        drv = GetDriverFactory().CreateInstance(vDriverModule, session);

        // Execute driver
        drv->Execute();

        // Finalise session
        session->Finalise();
    }
    catch (const std::runtime_error& e)
    {
        return 1;
    }
    catch (const std::string& eStr)
    {
        cout << "Error: " << eStr << endl;
    }

    return 0;
}

        vector<ForcingSharedPtr> Forcing::Load(
                            const LibUtilities::SessionReaderSharedPtr& pSession,
                            const Array<OneD, MultiRegions::ExpListSharedPtr>& pFields,
                            const unsigned int& pNumForcingFields)
        {
            vector<ForcingSharedPtr> vForceList;

            if (!pSession->DefinesElement("Nektar/Forcing"))
            {
                return vForceList;
            }

            TiXmlElement* vForcing = pSession->GetElement("Nektar/Forcing");
            if (vForcing)
            {
                unsigned int vNumForcingFields = pNumForcingFields;
                if (!pNumForcingFields)
                {
                    vNumForcingFields = pFields.num_elements();
                }

                TiXmlElement* vForce = vForcing->FirstChildElement("FORCE");
                while (vForce)
                {
                    string vType = vForce->Attribute("TYPE");

                    vForceList.push_back(GetForcingFactory().CreateInstance(
                                            vType, pSession, pFields,
                                            vNumForcingFields, vForce));
                    vForce = vForce->NextSiblingElement("FORCE");
                }
            }
            return vForceList;
        }

/**
 * @param       pSession    Session reader for IO
 * @param       pParams     Parameters of filter
 */
FilterBenchmark::FilterBenchmark(
        const LibUtilities::SessionReaderSharedPtr &pSession,
        const std::map<std::string, std::string> &pParams)
    : Filter(pSession)
{
    ASSERTL0(pParams.find("ThresholdValue") != pParams.end(),
             "Missing parameter 'ThresholdValue'.");
    m_thresholdValue = atof(pParams.find("ThresholdValue")->second.c_str());
    ASSERTL0(pParams.find("InitialValue") != pParams.end(),
             "Missing parameter 'InitialValue'.");
    m_initialValue = atof(pParams.find("InitialValue")->second.c_str());
    ASSERTL0(!(pParams.find("OutputFile")->second.empty()),
             "Missing parameter 'OutputFile'.");
    m_outputFile = pParams.find("OutputFile")->second;

    m_startTime = 0.0;
    if (pParams.find("StartTime") != pParams.end())
    {
        m_startTime = atof(pParams.find("StartTime")->second.c_str());
    }

    m_fld = MemoryManager<LibUtilities::FieldIO>::AllocateSharedPtr(
                                                        pSession->GetComm());

}

        /**
         * @brief Set up filter with output file and frequency parameters.
         *
         * @param pSession  Current session.
         * @param pParams   Map of parameters defined in XML file.
         */
        FilterEnergy1D::FilterEnergy1D(
            const LibUtilities::SessionReaderSharedPtr &pSession,
            const std::map<std::string, std::string> &pParams) :
            Filter(pSession),
            m_index(0)
        {
            std::string outName;
            if (pParams.find("OutputFile") == pParams.end())
            {
                outName = m_session->GetSessionName();
            }
            else
            {
                ASSERTL0(!(pParams.find("OutputFile")->second.empty()),
                         "Missing parameter 'OutputFile'.");
                outName = pParams.find("OutputFile")->second;
            }

            if (pParams.find("OutputFrequency") == pParams.end())
            {
                m_outputFrequency = 1;
            }
            else
            {
                m_outputFrequency =
                    atoi(pParams.find("OutputFrequency")->second.c_str());
            }

            outName += ".eny";

            ASSERTL0(pSession->GetComm()->GetSize() == 1,
                     "The 1D energy filter currently only works in serial.");

            m_out.open(outName.c_str());
        }

        FilterEnergyBase::FilterEnergyBase(
            const LibUtilities::SessionReaderSharedPtr &pSession,
            const std::map<std::string, std::string> &pParams,
            const bool pConstDensity)
            : Filter        (pSession),
              m_index       (-1),
              m_homogeneous (false),
              m_planes      (),
              m_constDensity(pConstDensity)
        {
            std::string outName;
            if (pParams.find("OutputFile") == pParams.end())
            {
                outName = m_session->GetSessionName();
            }
            else
            {
                ASSERTL0(!(pParams.find("OutputFile")->second.empty()),
                         "Missing parameter 'OutputFile'.");
                outName = pParams.find("OutputFile")->second;
            }

            m_comm = pSession->GetComm();
            outName += ".eny";
            if (m_comm->GetRank() == 0)
            {
                m_outFile.open(outName.c_str());
                ASSERTL0(m_outFile.good(), "Unable to open: '" + outName + "'");
                m_outFile.setf(ios::scientific, ios::floatfield);
                m_outFile << "# Time                Kinetic energy        "
                          << "Enstrophy"
                          << endl
                          << "# ---------------------------------------------"
                          << "--------------"
                          << endl;
            }
            pSession->LoadParameter("LZ", m_homogeneousLength, 0.0);

            ASSERTL0(pParams.find("OutputFrequency") != pParams.end(),
                     "Missing parameter 'OutputFrequency'.");
            m_outputFrequency = atoi(
                pParams.find("OutputFrequency")->second.c_str());
        }

int main(int argc, char *argv[])
{
    if(argc != 2)
    {
        fprintf(stderr,"Usage: ./Aliasing file.xml \n");
        fprintf(stderr,"\t Method will read intiial conditions section of .xml file for input \n");
        exit(1);
    }

    LibUtilities::SessionReaderSharedPtr session;
    string vDriverModule;
    DriverSharedPtr drv;  
    try
    {
        // Create session reader.
        session = LibUtilities::SessionReader::CreateInstance(argc, argv);
        
        // Create driver
        session->LoadSolverInfo("Driver", vDriverModule, "Standard");
        drv = GetDriverFactory().CreateInstance(vDriverModule, session);


        EquationSystemSharedPtr EqSys = drv->GetEqu()[0];
        IncNavierStokesSharedPtr IncNav = EqSys->as<IncNavierStokes>();
        
        IncNav->SetInitialConditions(0.0,false);
        Array<OneD, MultiRegions::ExpListSharedPtr> fields = IncNav->UpdateFields();
        
        int i;
        int nConvectiveFields = IncNav->GetNConvectiveFields();
        int nphys = fields[0]->GetTotPoints();
        Array<OneD, Array<OneD, NekDouble> > VelFields(nConvectiveFields);
        Array<OneD, Array<OneD, NekDouble> > NonLinear(nConvectiveFields);
        Array<OneD, Array<OneD, NekDouble> > NonLinearDealiased(nConvectiveFields);
        
        for(i = 0; i < nConvectiveFields; ++i)
        {
            VelFields[i] = fields[i]->UpdatePhys();
            NonLinear[i] = Array<OneD, NekDouble> (nphys);
            NonLinearDealiased[i] = Array<OneD, NekDouble> (nphys);
        }

        boost::shared_ptr<NavierStokesAdvection> A
            = boost::dynamic_pointer_cast<NavierStokesAdvection>(IncNav->GetAdvObject());

        if (!A)
        {
            cout << "Must use non-linear Navier-Stokes advection" << endl;
            exit(-1);
        }

        // calculate non-linear terms without dealiasing
        A->SetSpecHPDealiasing(false);
        A->Advect(nConvectiveFields, fields,
                                            VelFields, VelFields, 
                                            NonLinear, 0.0);


        // calculate non-linear terms with dealiasing
        A->SetSpecHPDealiasing(true);
        A->Advect(nConvectiveFields, fields,
                                            VelFields, VelFields, 
                                            NonLinearDealiased, 0.0);

        // Evaulate Difference and put into fields;
        for(i = 0; i < nConvectiveFields; ++i)
        {
            Vmath::Vsub(nphys,NonLinearDealiased[i],1,NonLinear[i],1,NonLinear[i],1);
            fields[i]->FwdTrans_IterPerExp(NonLinear[i],fields[i]->UpdateCoeffs());
            // Need to reset varibale name for output
            string name = "NL_Aliasing_"+session->GetVariable(i);
            session->SetVariable(i,name.c_str());
        }


        // Reset session name for output file
        std::string outname = IncNav->GetSessionName();
        
        outname += "_NonLinear_Aliasing";
        IncNav->ResetSessionName(outname);
        IncNav->Output();

    }
    catch (const std::runtime_error&)
    {
        return 1;
    }
    catch (const std::string& eStr)
    {
        cout << "Error: " << eStr << endl;
    }


    return 0;
}

int main(int argc, char *argv[])
{
    Array<OneD,NekDouble>  fce;
    Array<OneD,NekDouble>  xc0,xc1,xc2;

    if(argc < 2)
    {
        fprintf(stderr,"Usage: XmlToTecplot meshfile\n");
        exit(1);
    }

    LibUtilities::SessionReader::RegisterCmdLineFlag(
        "multi-zone", "m", "Output multi-zone format (one element per zone).");

    LibUtilities::SessionReaderSharedPtr vSession
            = LibUtilities::SessionReader::CreateInstance(argc, argv);

    //----------------------------------------------
    // Read in mesh from input file
    string meshfile(argv[argc-1]);
    SpatialDomains::MeshGraphSharedPtr graphShPt
        = SpatialDomains::MeshGraph::Read(vSession);
    //----------------------------------------------

    //----------------------------------------------
    // Set up Expansion information
    SpatialDomains::ExpansionMap emap = graphShPt->GetExpansions();
    SpatialDomains::ExpansionMapIter it;

    for (it = emap.begin(); it != emap.end(); ++it)
    {
        for (int i = 0; i < it->second->m_basisKeyVector.size(); ++i)
        {
            LibUtilities::BasisKey  tmp1 = it->second->m_basisKeyVector[i];
            LibUtilities::PointsKey tmp2 = tmp1.GetPointsKey();
            it->second->m_basisKeyVector[i] = LibUtilities::BasisKey(
                tmp1.GetBasisType(), tmp1.GetNumModes(),
                LibUtilities::PointsKey(tmp1.GetNumModes(),
                                        LibUtilities::ePolyEvenlySpaced));
        }
    }
    //----------------------------------------------

    //----------------------------------------------
    // Define Expansion
    int expdim  = graphShPt->GetMeshDimension();
    Array<OneD, MultiRegions::ExpListSharedPtr> Exp(1);

    switch(expdim)
    {
        case 1:
        {
            MultiRegions::ExpList1DSharedPtr Exp1D;
            Exp1D = MemoryManager<MultiRegions::ExpList1D>::AllocateSharedPtr(vSession,graphShPt);
            Exp[0] = Exp1D;
            break;
        }
        case 2:
        {
            if(vSession->DefinesSolverInfo("HOMOGENEOUS"))
            {
                std::string HomoStr = vSession->GetSolverInfo("HOMOGENEOUS");
                MultiRegions::ExpList3DHomogeneous1DSharedPtr Exp3DH1;

                ASSERTL0(
                    HomoStr == "HOMOGENEOUS1D" || HomoStr == "Homogeneous1D" ||
                    HomoStr == "1D"            || HomoStr == "Homo1D",
                    "Only 3DH1D supported for XML output currently.");

                int nplanes;
                vSession->LoadParameter("HomModesZ", nplanes);

                // choose points to be at evenly spaced points at nplanes + 1
                // points
                const LibUtilities::PointsKey Pkey(
                    nplanes + 1, LibUtilities::ePolyEvenlySpaced);
                const LibUtilities::BasisKey  Bkey(
                    LibUtilities::eFourier, nplanes, Pkey);
                NekDouble lz = vSession->GetParameter("LZ");

                Exp3DH1 = MemoryManager<MultiRegions::ExpList3DHomogeneous1D>
                    ::AllocateSharedPtr(
                        vSession, Bkey, lz, false, false, graphShPt);
                Exp[0] = Exp3DH1;
            }
            else
            {
                MultiRegions::ExpList2DSharedPtr Exp2D;
                Exp2D = MemoryManager<MultiRegions::ExpList2D>::AllocateSharedPtr(vSession,graphShPt);
                Exp[0] =  Exp2D;
            }
            break;
        }
        case 3:
        {
            MultiRegions::ExpList3DSharedPtr Exp3D;
            Exp3D = MemoryManager<MultiRegions::ExpList3D>::AllocateSharedPtr(vSession,graphShPt);
            Exp[0] =  Exp3D;
            break;
        }
//.........这里部分代码省略.........

int main(int argc, char *argv[])
{
    LibUtilities::SessionReaderSharedPtr vSession
            = LibUtilities::SessionReader::CreateInstance(argc, argv);

    LibUtilities::CommSharedPtr vComm = vSession->GetComm();

    MultiRegions::DisContField3DSharedPtr Exp,Fce;
    int     i, nq, coordim;
    Array<OneD,NekDouble>  fce;
    Array<OneD,NekDouble>  xc0,xc1,xc2;
    StdRegions::ConstFactorMap factors;

    if(argc < 2)
    {
        fprintf(stderr,"Usage: PostProcHDG3D  meshfile [solntype]\n");
        exit(1);
    }

    //----------------------------------------------
    // Read in mesh from input file
    SpatialDomains::MeshGraphSharedPtr graph3D = MemoryManager<SpatialDomains::MeshGraph3D>::AllocateSharedPtr(vSession);
    //----------------------------------------------

    //----------------------------------------------
    // Print summary of solution details
    factors[StdRegions::eFactorLambda] = vSession->GetParameter("Lambda");
    factors[StdRegions::eFactorTau] = 1.0;
    const SpatialDomains::ExpansionMap &expansions = graph3D->GetExpansions();
    LibUtilities::BasisKey bkey0
                            = expansions.begin()->second->m_basisKeyVector[0];

	//MAY NEED ADJUSTMENT FOR VARIOUS ELEMENT TYPES
	int num_modes = bkey0.GetNumModes();
	int num_points = bkey0.GetNumPoints();

    if (vComm->GetRank() == 0)
    {
        cout << "Solving 3D Helmholtz:"  << endl;
        cout << "         Lambda     : " << factors[StdRegions::eFactorLambda] << endl;
        cout << "         No. modes  : " << num_modes << endl;
        cout << "         No. points : " << num_points << endl;
        cout << endl;
    }

    //----------------------------------------------
    // Define Expansion
    //----------------------------------------------
    Exp = MemoryManager<MultiRegions::DisContField3D>::
        AllocateSharedPtr(vSession,graph3D,vSession->GetVariable(0));
    //----------------------------------------------
    Timing("Read files and define exp ..");

    //----------------------------------------------
    // Set up coordinates of mesh for Forcing function evaluation
    coordim = Exp->GetCoordim(0);
    nq      = Exp->GetTotPoints();

    xc0 = Array<OneD,NekDouble>(nq,0.0);
    xc1 = Array<OneD,NekDouble>(nq,0.0);
    xc2 = Array<OneD,NekDouble>(nq,0.0);

    switch(coordim)
    {
    case 1:
        Exp->GetCoords(xc0);
        break;
    case 2:
        Exp->GetCoords(xc0,xc1);
        break;
    case 3:
        Exp->GetCoords(xc0,xc1,xc2);
        break;
    }
    //----------------------------------------------

    //----------------------------------------------
    // Define forcing function for first variable defined in file
    fce = Array<OneD,NekDouble>(nq);
    LibUtilities::EquationSharedPtr ffunc
                                    = vSession->GetFunction("Forcing", 0);

    ffunc->Evaluate(xc0, xc1, xc2, fce);

    //----------------------------------------------


    //----------------------------------------------
    // Setup expansion containing the  forcing function
    Fce = MemoryManager<MultiRegions::DisContField3D>::AllocateSharedPtr(*Exp);
    Fce->SetPhys(fce);
    //----------------------------------------------
    Timing("Define forcing ..");

    //----------------------------------------------
    // Helmholtz solution taking physical forcing
    Exp->HelmSolve(Fce->GetPhys(), Exp->UpdateCoeffs(), NullFlagList, factors);
    //----------------------------------------------

    Timing("Helmholtz Solve ..");
//.........这里部分代码省略.........

        void Forcing::EvaluateFunction(
                Array<OneD, MultiRegions::ExpListSharedPtr>       pFields,
                LibUtilities::SessionReaderSharedPtr              pSession,
                std::string                                       pFieldName,
                Array<OneD, NekDouble>&                           pArray,
                const std::string&                                pFunctionName,
                NekDouble                                         pTime)
        {
            ASSERTL0(pSession->DefinesFunction(pFunctionName),
                     "Function '" + pFunctionName + "' does not exist.");

            unsigned int nq = pFields[0]->GetNpoints();
            if (pArray.num_elements() != nq)
            {
                pArray = Array<OneD, NekDouble> (nq);
            }

            LibUtilities::FunctionType vType;
            vType = pSession->GetFunctionType(pFunctionName, pFieldName);
            if (vType == LibUtilities::eFunctionTypeExpression)
            {
                Array<OneD, NekDouble> x0(nq);
                Array<OneD, NekDouble> x1(nq);
                Array<OneD, NekDouble> x2(nq);

                pFields[0]->GetCoords(x0, x1, x2);
                LibUtilities::EquationSharedPtr ffunc =
                        pSession->GetFunction(pFunctionName, pFieldName);

                ffunc->Evaluate(x0, x1, x2, pTime, pArray);
            }
            else if (vType == LibUtilities::eFunctionTypeFile)
            {
                std::string filename = pSession->GetFunctionFilename(
                                                    pFunctionName,
                                                    pFieldName);

                std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
                std::vector<std::vector<NekDouble> > FieldData;
                Array<OneD, NekDouble> vCoeffs(pFields[0]->GetNcoeffs());
                Vmath::Zero(vCoeffs.num_elements(), vCoeffs, 1);

                LibUtilities::FieldIOSharedPtr fld =
                    MemoryManager<LibUtilities::FieldIO>::AllocateSharedPtr(m_session->GetComm());
                fld->Import(filename, FieldDef, FieldData);

                int idx = -1;
                for (int i = 0; i < FieldDef.size(); ++i)
                {
                    for (int j = 0; j < FieldDef[i]->m_fields.size(); ++j)
                    {
                        if (FieldDef[i]->m_fields[j] == pFieldName)
                        {
                            idx = j;
                        }
                    }

                    if (idx >= 0)
                    {
                        pFields[0]->ExtractDataToCoeffs(
                                                    FieldDef[i],
                                                    FieldData[i],
                                                    FieldDef[i]->m_fields[idx],
                                                    vCoeffs);
                    }
                    else
                    {
                        cout << "Field " + pFieldName + " not found." << endl;
                    }
                }
                pFields[0]->BwdTrans_IterPerExp(vCoeffs, pArray);
            }
        }

        /**
         * Read global optimisation parameters from a file and set up flags.
         * @param   fileName    File to read parameters from.
         */
        GlobalOptParam::GlobalOptParam(const LibUtilities::SessionReaderSharedPtr& pSession, const int dim,
                                         const Array<OneD, const int> &NumShapeElements):
            m_doGlobalMatOp(SIZE_OptimizeOperationType,false)
        {
            int i;
            int numShapes = 0;
            TiXmlDocument& doc = pSession->GetDocument();

            m_shapeNumElements = NumShapeElements;

            switch (dim)
            {
            case 1:
                numShapes = 1;
                ASSERTL0(false,"Needs setting up for dimension 1");
                break;
            case 2:
                numShapes = 2;
                m_shapeList = Array<OneD, LibUtilities::ShapeType>(numShapes);
                m_shapeList[0] = LibUtilities::eTriangle;
                m_shapeList[1] = LibUtilities::eQuadrilateral;
                break;
            case 3:
                numShapes = 4;
                m_shapeList = Array<OneD, LibUtilities::ShapeType>(numShapes);
                m_shapeList[0] = LibUtilities::eTetrahedron;
                m_shapeList[1] = LibUtilities::ePyramid;
                m_shapeList[2] = LibUtilities::ePrism;
                m_shapeList[3] = LibUtilities::eHexahedron;
                break;
            }

            m_doBlockMatOp = Array<OneD, Array<OneD,bool> > (SIZE_OptimizeOperationType);
            for(i = 0; i < SIZE_OptimizeOperationType; ++i)
            {
                m_doBlockMatOp[i] = Array<OneD, bool> (numShapes,false);
            }

            TiXmlHandle docHandle(&doc);
            TiXmlElement* master
                        = docHandle.FirstChildElement("NEKTAR").Element();
            ASSERTL0(master   , "Unable to find NEKTAR tag in file.");

            TiXmlElement* paramList
                        = docHandle.FirstChildElement("NEKTAR")
                            .FirstChildElement("GLOBALOPTIMIZATIONPARAMETERS")
                            .Element();

            // If no global optimisation parameters set, we're done
            if (!paramList)
            {
                return;
            }

            // Check if there is a reference an external file and if so, load it
            TiXmlElement* source
                        = paramList->FirstChildElement("SOURCE");
            if (source)
            {
                std::string sourceFile = source->Attribute("FILE");
                TiXmlDocument docSource;
                bool loadOkay = docSource.LoadFile(sourceFile);
                ASSERTL0(loadOkay, (std::string("Unable to load file: ") +
                                                sourceFile).c_str());
                TiXmlHandle docSourceHandle(&docSource);
                master = docHandle.FirstChildElement("NEKTAR").Element();
                ASSERTL0(master   , "Unable to find NEKTAR tag in file.");

                paramList = docHandle.FirstChildElement("NEKTAR")
                                .FirstChildElement("GLOBALOPTIMIZATIONPARAMETERS")
                                .Element();
                ASSERTL0(paramList, std::string("Specified source file '"
                        + sourceFile + "' is missing an "
                        "GLOBALOPTIMIZATIONPARAMETERS tag.").c_str());
            }

            int n;
            for(n = 0; n < SIZE_OptimizeOperationType; n++)
            {
                TiXmlElement* operationType = paramList->FirstChildElement(
                                      OptimizationOperationTypeMap[n]);

                if(operationType)
                {
                    TiXmlElement* arrayElement = operationType
                                        ->FirstChildElement("DO_GLOBAL_MAT_OP");
                    if(arrayElement)
                    {
                        int value;
                        int err;

                        err = arrayElement->QueryIntAttribute("VALUE", &value);
                        ASSERTL0(err == TIXML_SUCCESS,(
                           std::string("Unable to read DO_GLOBAL_MAT_OP "
                                       "attribute VALUE for ")
                         + std::string(OptimizationOperationTypeMap[n])
//.........这里部分代码省略.........

int main(int argc, char *argv[])
{
    string fname = std::string(argv[2]);
    int fdot = fname.find_last_of('.');
    if (fdot != std::string::npos)
    {
        string ending = fname.substr(fdot);

        // If .chk or .fld we exchange the extension in the output file.
        // For all other files (e.g. .bse) we append the extension to avoid
        // conflicts.
        if (ending == ".chk" || ending == ".fld")
        {
            fname = fname.substr(0,fdot);
        }
    }

    fname = fname + ".txt";

    int cnt;
    int id1, id2;
    int i, j, n, e, b;
    Array<OneD, NekDouble> auxArray;

    int nBndEdgePts, nBndEdges, nBndRegions;

    if (argc < 3)
    {
        fprintf(stderr,
                "Usage: ExtractSurface3DCFS meshfile fieldFile\n");
        fprintf(stderr,
                "Extracts a surface from a 3D fld file"
                "(only for CompressibleFlowSolver and purely 3D .fld files)\n");
        exit(1);
    }

    LibUtilities::SessionReaderSharedPtr vSession
        = LibUtilities::SessionReader::CreateInstance(3, argv);

    std::string                         m_ViscosityType;

    NekDouble                           m_gamma;
    NekDouble                           m_pInf;
    NekDouble                           m_rhoInf;
    NekDouble                           m_uInf;
    NekDouble                           m_vInf;
    NekDouble                           m_wInf;
    NekDouble                           m_gasConstant;
    NekDouble                           m_Twall;
    NekDouble                           m_mu;
    NekDouble                           m_thermalConductivity;

    int m_spacedim = 3;
    int nDimensions = m_spacedim;
    int phys_offset;

    // Get gamma parameter from session file.
    ASSERTL0(vSession->DefinesParameter("Gamma"),
             "Compressible flow sessions must define a Gamma parameter.");
    vSession->LoadParameter("Gamma", m_gamma, 1.4);

    // Get E0 parameter from session file.
    ASSERTL0(vSession->DefinesParameter("pInf"),
             "Compressible flow sessions must define a pInf parameter.");
    vSession->LoadParameter("pInf", m_pInf, 101325);

    // Get rhoInf parameter from session file.
    ASSERTL0(vSession->DefinesParameter("rhoInf"),
             "Compressible flow sessions must define a rhoInf parameter.");
    vSession->LoadParameter("rhoInf", m_rhoInf, 1.225);

    // Get uInf parameter from session file.
    ASSERTL0(vSession->DefinesParameter("uInf"),
             "Compressible flow sessions must define a uInf parameter.");
    vSession->LoadParameter("uInf", m_uInf, 0.1);

    // Get vInf parameter from session file.
    if (m_spacedim == 2 || m_spacedim == 3)
    {
        ASSERTL0(vSession->DefinesParameter("vInf"),
                 "Compressible flow sessions must define a vInf parameter"
                 "for 2D/3D problems.");
        vSession->LoadParameter("vInf", m_vInf, 0.0);
    }

    // Get wInf parameter from session file.
    if (m_spacedim == 3)
    {
        ASSERTL0(vSession->DefinesParameter("wInf"),
                 "Compressible flow sessions must define a wInf parameter"
                 "for 3D problems.");
        vSession->LoadParameter("wInf", m_wInf, 0.0);
    }

    vSession->LoadParameter ("GasConstant",   m_gasConstant,   287.058);
    vSession->LoadParameter ("Twall",         m_Twall,         300.15);
    vSession->LoadSolverInfo("ViscosityType", m_ViscosityType, "Constant");
    vSession->LoadParameter ("mu",            m_mu,            1.78e-05);
    vSession->LoadParameter ("thermalConductivity",
                             m_thermalConductivity, 0.0257);
//.........这里部分代码省略.........

int main(int argc, char *argv[])
{
    unsigned int i,j;

    if(argc < 3)
    {
        fprintf(stderr,"Usage: FldToPts  meshfile  fieldfile(s)\n");
        exit(1);
    }


    int nExtraPoints;
    LibUtilities::SessionReaderSharedPtr vSession
            = LibUtilities::SessionReader::CreateInstance(argc, argv);

    vSession->LoadParameter("OutputExtraPoints",nExtraPoints,0);

    //----------------------------------------------
    // Read in mesh from input file
    SpatialDomains::MeshGraphSharedPtr graphShPt = SpatialDomains::MeshGraph::Read(vSession);//->GetFilename(), false);
    //----------------------------------------------

    for (int n = 2; n < argc; ++n)
    {
        string fname = std::string(argv[n]);
        int fdot = fname.find_last_of('.');
        if (fdot != std::string::npos)
        {
            string ending = fname.substr(fdot);
            if (ending == ".chk" || ending == ".fld")
            {
                fname = fname.substr(0,fdot);
            }
        }
        fname = fname + ".pts";

        if (argc > 3)
        {
            if (fexist(fname.c_str()))
            {
                cout << "Skipping converted file: " << argv[n] << endl;
                continue;
            }
            cout << "Processing " << argv[n] << endl;
        }

        //----------------------------------------------
        // Import field file.
        string fieldfile(argv[n]);
        vector<SpatialDomains::FieldDefinitionsSharedPtr> fielddef;
        vector<vector<NekDouble> > fielddata;
        graphShPt->Import(fieldfile,fielddef,fielddata);
        bool useFFT = false;
		bool dealiasing = false;
        //----------------------------------------------

        //----------------------------------------------
        // Set up Expansion information
        for(i = 0; i < fielddef.size(); ++i)
        {
            vector<LibUtilities::PointsType> ptype;
            for(j = 0; j < 3; ++j)
            {
                ptype.push_back(LibUtilities::ePolyEvenlySpaced);
            }
            
            fielddef[i]->m_pointsDef = true;
            fielddef[i]->m_points    = ptype; 
            
            vector<unsigned int> porder;
            if(fielddef[i]->m_numPointsDef == false)
            {
                for(j = 0; j < fielddef[i]->m_numModes.size(); ++j)
                {
                    porder.push_back(fielddef[i]->m_numModes[j]+nExtraPoints);
                }
                
                fielddef[i]->m_numPointsDef = true;
            }
            else
            {
                for(j = 0; j < fielddef[i]->m_numPoints.size(); ++j)
                {
                    porder.push_back(fielddef[i]->m_numPoints[j]+nExtraPoints);
                }
            }
            fielddef[i]->m_numPoints = porder;
            
        }
        graphShPt->SetExpansions(fielddef);
        //----------------------------------------------


        //----------------------------------------------
        // Define Expansion
        int expdim  = graphShPt->GetMeshDimension();
        int nfields = fielddef[0]->m_fields.size();
        Array<OneD, MultiRegions::ExpListSharedPtr> Exp(nfields);

        switch(expdim)
//.........这里部分代码省略.........

/**
 * Main function.
 *
 * Usage: VtkToFld session.xml input.vtk output.fld [options]
 */
int main(int argc, char* argv[])
{
    // Set up available options
    po::options_description desc("Available options");
    desc.add_options()
        ("help,h",         "Produce this help message.")
        ("name,n", po::value<string>()->default_value("Intensity"),
                "Name of field in VTK file to use for intensity.")
        ("outname,m", po::value<string>()->default_value("intensity"),
                "Name of field in output FLD file.")
        ("precision,p",  po::value<double>()->default_value(1),
             "Precision of vertex matching.");

    po::options_description hidden("Hidden options");
    hidden.add_options()
        ("file",   po::value<vector<string> >(), "Input filename");

    po::options_description cmdline_options;
    cmdline_options.add(desc).add(hidden);

    po::positional_options_description p;
    p.add("file", -1);

    po::variables_map vm;

    // Parse command-line options
    try
    {
        po::store(po::command_line_parser(argc, argv).
                  options(cmdline_options).positional(p).run(), vm);
        po::notify(vm);
    }
    catch (const std::exception& e)
    {
        cerr << e.what() << endl;
        cerr << desc;
        return 1;
    }

    if ( vm.count("help") || vm.count("file") == 0 ||
                             vm["file"].as<vector<string> >().size() != 3) {
        cerr << "Usage: VtkToFld session.xml intensity.vtk output.fld [options]"
             << endl;
        cerr << desc;
        return 1;
    }

    // Extract command-line argument values
    std::vector<std::string> vFiles = vm["file"].as<vector<string> >();
    const string infile  = vFiles[1];
    const string outfile = vFiles[2];
    const double factor  = vm["precision"].as<double>();
    const string name    = vm["name"].as<string>();
    const string outname = vm["outname"].as<string>();

    std::vector<std::string> vFilenames;
    LibUtilities::SessionReaderSharedPtr vSession;
    SpatialDomains::MeshGraphSharedPtr graph2D;
    MultiRegions::ExpList2DSharedPtr Exp;

    vFilenames.push_back(vFiles[0]);
    vSession = LibUtilities::SessionReader::CreateInstance(2, argv, vFilenames);

    try
    {
        //----------------------------------------------
        // Read in mesh from input file
        graph2D = MemoryManager<SpatialDomains::MeshGraph2D>::
                    AllocateSharedPtr(vSession);
        //----------------------------------------------

        //----------------------------------------------
        // Define Expansion
        Exp = MemoryManager<MultiRegions::ExpList2D>::
                    AllocateSharedPtr(vSession,graph2D);
        //----------------------------------------------

        //----------------------------------------------
        // Set up coordinates of mesh
        int coordim = Exp->GetCoordim(0);
        int nq      = Exp->GetNpoints();

        Array<OneD, NekDouble> xc0(nq,0.0);
        Array<OneD, NekDouble> xc1(nq,0.0);
        Array<OneD, NekDouble> xc2(nq,0.0);

        switch(coordim)
        {
        case 2:
            Exp->GetCoords(xc0,xc1);
            break;
        case 3:
            Exp->GetCoords(xc0,xc1,xc2);
            break;
        default:
//.........这里部分代码省略.........

int main(int argc, char *argv[])
{
    LibUtilities::SessionReaderSharedPtr vSession
            = LibUtilities::SessionReader::CreateInstance(argc, argv);

    MultiRegions::ContField2DSharedPtr Exp,Fce;
    int     nq,  coordim;
    Array<OneD,NekDouble>  fce; 
    Array<OneD,NekDouble>  xc0,xc1,xc2; 
    NekDouble  lambda;
    NekDouble  ax,ay;

    if((argc != 2)&&(argc != 3))
    {
        fprintf(stderr,"Usage: SteadyLinearAdvectionReaction2D  meshfile [SysSolnType]\n");
        exit(1);
    }

    //----------------------------------------------
    // Read in mesh from input file
    SpatialDomains::MeshGraphSharedPtr graph2D = MemoryManager<SpatialDomains::MeshGraph2D>::AllocateSharedPtr(vSession);
    //----------------------------------------------

    //----------------------------------------------
    // Get Advection Velocity
    ax = vSession->GetParameter("Advection_x");
    ay = vSession->GetParameter("Advection_y");
    //----------------------------------------------

    //----------------------------------------------
    // Print summary of solution details
    lambda = vSession->GetParameter("Lambda");
    cout << "            Lambda         : " << lambda << endl;
    const SpatialDomains::ExpansionMap &expansions = graph2D->GetExpansions();
    LibUtilities::BasisKey bkey0
                            = expansions.begin()->second->m_basisKeyVector[0];
    LibUtilities::BasisKey bkey1
                            = expansions.begin()->second->m_basisKeyVector[1];
    cout << "Solving Steady 2D LinearAdvection :"  << endl; 
    cout << "            Advection_x    : " << ax << endl; 
    cout << "            Advection_y    : " << ay << endl; 
    cout << "            Expansion      : (" << LibUtilities::BasisTypeMap[bkey0.GetBasisType()] <<","<< LibUtilities::BasisTypeMap[bkey1.GetBasisType()]  << ")" << endl;
    cout << "            No. modes      : " << bkey0.GetNumModes() << endl;
    cout << endl;
    //----------------------------------------------
   
    //----------------------------------------------
    // Define Expansion 
    Exp = MemoryManager<MultiRegions::ContField2D>::
        AllocateSharedPtr(vSession,graph2D,vSession->GetVariable(0));
    //----------------------------------------------

    Timing("Read files and define exp ..");
    
    //----------------------------------------------
    // Set up coordinates of mesh for Forcing function evaluation
    coordim = Exp->GetCoordim(0);
    nq      = Exp->GetTotPoints();
    
    xc0 = Array<OneD,NekDouble>(nq,0.0);
    xc1 = Array<OneD,NekDouble>(nq,0.0);
    xc2 = Array<OneD,NekDouble>(nq,0.0);
    
    switch(coordim)
    {
    case 1:
        Exp->GetCoords(xc0);
        break;
    case 2:
        Exp->GetCoords(xc0,xc1);
        break;
    case 3:
        Exp->GetCoords(xc0,xc1,xc2);
        break;
    }

    Array<OneD, Array< OneD, NekDouble> > Vel(2);
    Vel[0] = Array<OneD, NekDouble> (nq,ax);
    Vel[1] = Array<OneD, NekDouble> (nq,ay);
    //----------------------------------------------
    
    //----------------------------------------------
    // Define forcing function for first variable defined in file 
    fce = Array<OneD,NekDouble>(nq);
    LibUtilities::EquationSharedPtr ffunc = vSession->GetFunction("Forcing",0);

    ffunc->Evaluate(xc0,xc1,xc2,fce);
    //----------------------------------------------

    //----------------------------------------------
    // Setup expansion containing the  forcing function
    Fce = MemoryManager<MultiRegions::ContField2D>::AllocateSharedPtr(*Exp);
    Fce->SetPhys(fce);
    //----------------------------------------------
    Timing("Define forcing ..");
  
    //----------------------------------------------
    // Helmholtz solution taking physical forcing 
    Exp->LinearAdvectionReactionSolve(Vel, Fce->GetPhys(), Exp->UpdateCoeffs(), lambda, MultiRegions::eGlobal);
    //----------------------------------------------
//.........这里部分代码省略.........

int main(int argc, char *argv[])
{
    int cnt;
    int id1, id2;
    int i, j, e, b;
    
    int nBndEdgePts, nBndEdges, nBndRegions;
        
    if (argc < 3)
    {
        fprintf(stderr,
                "Usage: ExtractSurface2DCFS meshfile fieldFile\n");
        fprintf(stderr,
                "Extracts a surface from a 2D fld file" 
                "(only for CompressibleFlowSolver and purely 2D .fld files)\n");
        exit(1);
    }

    LibUtilities::SessionReaderSharedPtr vSession
            = LibUtilities::SessionReader::CreateInstance(3, argv);

    //--------------------------------------------------------------------------
    // Read in mesh from input file
    string meshfile(argv[1]);
    SpatialDomains::MeshGraphSharedPtr graphShPt = 
        SpatialDomains::MeshGraph::Read(vSession);
    //--------------------------------------------------------------------------

    //--------------------------------------------------------------------------
    // Import field file
    string                                          fieldFile(argv[2]);
    vector<LibUtilities::FieldDefinitionsSharedPtr> fieldDef;
    vector<vector<NekDouble> >                      fieldData;
    
    LibUtilities::Import(fieldFile, fieldDef, fieldData);
    //----------------------------------