56 session = LibUtilities::SessionReader::CreateInstance(argc, argv);
63 string sessionName = session->GetSessionName();
64 string outFile = sessionName +
".fld";
65 unsigned int nSteps = session->GetParameter(
"NumSteps");
66 NekDouble delta_t = session->GetParameter(
"TimeStep");
67 NekDouble epsilon = session->GetParameter(
"epsilon" );
70 graph = SpatialDomains::MeshGraph::Read(session);
77 unsigned int nq = field->GetNpoints();
79 field->GetCoords(x0,x1,x2);
82 icond = session->GetFunction(
"InitialConditions",
"u");
83 icond->Evaluate(x0, x1, x2, 0.0, field->UpdatePhys());
89 Vmath::Zero(field->GetNcoeffs(), field->UpdateCoeffs(), 1);
92 for (
unsigned int n = 0; n < nSteps; ++n)
94 Vmath::Smul(nq, -1.0/delta_t/epsilon, field->GetPhys(), 1,
95 field->UpdatePhys(), 1);
97 field->HelmSolve(field->GetPhys(), field->UpdateCoeffs(),
100 field->BwdTrans(field->GetCoeffs(), field->UpdatePhys());
104 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef
105 = field->GetFieldDefinitions();
106 std::vector<std::vector<NekDouble> > FieldData(FieldDef.size());
107 for(
int i = 0; i < FieldDef.size(); ++i)
109 FieldDef[i]->m_fields.push_back(
"u");
110 field->AppendFieldData(FieldDef[i], FieldData[i]);
112 fld->Write(outFile, FieldDef, FieldData);
115 ex_sol = session->GetFunction(
"ExactSolution",0);
123 ex_sol->Evaluate(x0, x1, x2, (nSteps)*delta_t, exact);
127 cout <<
"L inf error: "
128 << field->Linf(field->GetPhys(), exact) << endl;
129 cout <<
"L 2 error: "
130 << field->L2(field->GetPhys(), exact) << endl;
131 cout <<
"H 1 error: "
132 << field->H1(field->GetPhys(), exact) << endl;
139 catch (
const std::runtime_error& e)
143 catch (
const std::string& eStr)
145 cout <<
"Error: " << eStr << endl;
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
boost::shared_ptr< ContField2D > ContField2DSharedPtr
std::map< ConstFactorType, NekDouble > ConstFactorMap
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
boost::shared_ptr< FieldIO > FieldIOSharedPtr
boost::shared_ptr< Equation > EquationSharedPtr
void Zero(int n, T *x, const int incx)
Zero vector.
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
static FlagList NullFlagList
An empty flag list.