87 session = LibUtilities::SessionReader::CreateInstance(argc, argv);
94 sessionName = session->GetSessionName();
95 scheme = session->GetSolverInfo(
"TimeIntegrationMethod");
96 nSteps = session->GetParameter(
"NumSteps");
97 delta_t = session->GetParameter(
"TimeStep");
98 epsilon = session->GetParameter(
"epsilon");
99 lambda = 1.0/delta_t/epsilon;
102 graph = SpatialDomains::MeshGraph::Read(session);
109 fields[0] = field->UpdatePhys();
112 unsigned int nq = field->GetNpoints();
114 field->GetCoords(x0,x1,x2);
118 = session->GetFunction(
"InitialConditions",
"u");
119 icond->Evaluate(x0,x1,x2,0.0,field->UpdatePhys());
130 CreateInstance(scheme);
135 u = IntScheme->InitializeScheme(delta_t, fields, 0.0, ode);
138 Vmath::Zero(field->GetNcoeffs(), field->UpdateCoeffs(), 1);
140 for (
int n = 0; n < nSteps; ++n)
142 fields = IntScheme->TimeIntegrate(n, delta_t, u, ode);
144 Vmath::Vcopy(field->GetNpoints(), fields[0], 1, field->UpdatePhys(), 1);
161 for (
int i = 0; i < inarray.num_elements(); ++i)
169 field->HelmSolve(outarray[i],
170 field->UpdateCoeffs(),
174 field->BwdTrans (field->GetCoeffs(), outarray[i]);
181 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef
182 = field->GetFieldDefinitions();
183 std::vector<std::vector<NekDouble> > FieldData(FieldDef.size());
184 for(
int i = 0; i < FieldDef.size(); ++i)
186 FieldDef[i]->m_fields.push_back(
"u");
187 field->AppendFieldData(FieldDef[i], FieldData[i]);
189 fld->Write(session->GetSessionName() +
".fld", FieldDef, FieldData);
196 unsigned int nq = field->GetNpoints();
198 field->GetCoords(x0,x1,x2);
201 session->GetFunction(
"ExactSolution",0);
207 ex_sol->Evaluate(x0, x1, x2, (nSteps)*delta_t, exact);
210 cout <<
"L inf error: "
211 << field->Linf(field->GetPhys(), exact) << endl;
212 cout <<
"L 2 error: "
213 << field->L2(field->GetPhys(), exact) << endl;
214 cout <<
"H 1 error: "
215 << field->H1(field->GetPhys(), exact) << endl;
220 int main(
int argc,
char *argv[])
227 catch (
const std::runtime_error& e)
231 catch (
const std::string& eStr)
233 cout <<
"Error: " << eStr << endl;
LibUtilities::TimeIntegrationWrapperSharedPtr IntScheme
MultiRegions::ContField2DSharedPtr field
Diffusion(int argc, char *argv[])
int main(int argc, char *argv[])
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
boost::shared_ptr< TimeIntegrationWrapper > TimeIntegrationWrapperSharedPtr
LibUtilities::TimeIntegrationSolutionSharedPtr u
LibUtilities::FieldIOSharedPtr fld
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.
Array< OneD, Array< OneD, NekDouble > > fields
boost::shared_ptr< FieldIO > FieldIOSharedPtr
TimeIntegrationWrapperFactory & GetTimeIntegrationWrapperFactory()
boost::shared_ptr< Equation > EquationSharedPtr
LibUtilities::SessionReaderSharedPtr session
void DoImplicitSolve(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time, const NekDouble lambda)
SpatialDomains::MeshGraphSharedPtr graph
LibUtilities::TimeIntegrationSchemeOperators ode
boost::shared_ptr< TimeIntegrationSolution > TimeIntegrationSolutionSharedPtr
void Zero(int n, T *x, const int incx)
Zero vector.
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
static FlagList NullFlagList
An empty flag list.