57 #include <boost/format.hpp>
58 # include <boost/function.hpp>
88 Loki::ClassLevelLockable> Type;
89 return Type::Instance();
97 EquationSystem::EquationSystem(
99 : m_comm (pSession->GetComm()),
100 m_session (pSession),
105 const vector<std::string> filenames =
m_session->GetFilenames();
107 for(
int i = 0; i < filenames.size(); ++i)
109 string sessionname =
"SessionName";
110 sessionname += boost::lexical_cast<std::string>(i);
112 m_fieldMetaDataMap[
"ChkFileNum"] = boost::lexical_cast<std::string>(0);
148 if (
m_session->DefinesSolverInfo(
"HOMOGENEOUS"))
150 std::string HomoStr =
m_session->GetSolverInfo(
"HOMOGENEOUS");
153 if ((HomoStr ==
"HOMOGENEOUS1D") || (HomoStr ==
"Homogeneous1D")
154 || (HomoStr ==
"1D") || (HomoStr ==
"Homo1D"))
160 if(
m_session->DefinesSolverInfo(
"ModeType"))
162 m_session->MatchSolverInfo(
"ModeType",
"SingleMode",
164 m_session->MatchSolverInfo(
"ModeType",
"HalfMode",
166 m_session->MatchSolverInfo(
"ModeType",
"MultipleModes",
171 if (
m_session->DefinesSolverInfo(
"ModeType"))
187 ASSERTL0(
false,
"SolverInfo ModeType not valid");
196 if ((HomoStr ==
"HOMOGENEOUS2D") || (HomoStr ==
"Homogeneous2D")
197 || (HomoStr ==
"2D") || (HomoStr ==
"Homo2D"))
207 if ((HomoStr ==
"HOMOGENEOUS3D") || (HomoStr ==
"Homogeneous3D")
208 || (HomoStr ==
"3D") || (HomoStr ==
"Homo3D"))
220 m_session->MatchSolverInfo(
"DEALIASING",
"True",
224 m_session->MatchSolverInfo(
"DEALIASING",
"On",
235 m_session->MatchSolverInfo(
"SPECTRALHPDEALIASING",
"True",
239 m_session->MatchSolverInfo(
"SPECTRALHPDEALIASING",
"On",
245 if (
m_session->DefinesSolverInfo(
"PROJECTION"))
247 std::string ProjectStr =
m_session->GetSolverInfo(
"PROJECTION");
249 if ((ProjectStr ==
"Continuous") || (ProjectStr ==
"Galerkin") ||
250 (ProjectStr ==
"CONTINUOUS") || (ProjectStr ==
"GALERKIN"))
254 else if ((ProjectStr ==
"MixedCGDG") ||
255 (ProjectStr ==
"Mixed_CG_Discontinuous"))
259 else if(ProjectStr ==
"DisContinuous")
265 ASSERTL0(
false,
"PROJECTION value not recognised");
270 cerr <<
"Projection type not specified in SOLVERINFO,"
271 "defaulting to continuous Galerkin" << endl;
279 int nvariables =
m_session->GetVariables().size();
280 bool DeclareCoeffPhysArrays =
true;
307 for (i = 0; i <
m_fields.num_elements(); i++)
320 for (i = 0; i <
m_fields.num_elements(); i++)
347 for(i = 0; i <
m_fields.num_elements(); i++)
375 for (i = 0; i <
m_fields.num_elements(); i++)
377 if(
m_session->GetVariable(i).compare(
"w")
414 for (i = 0; i <
m_fields.num_elements(); i++)
435 DeclareCoeffPhysArrays,
438 for (i = 1; i < m_fields.num_elements(); i++)
441 SameExpansions(
m_session->GetVariable(0),
448 DeclareCoeffPhysArrays,
457 DeclareCoeffPhysArrays,
476 m_fields[0]->GetTrace()->
494 for (i = 1; i < m_fields.num_elements(); i++)
531 for(i = 1; i < m_fields.num_elements(); ++i)
533 m_fields[i]->GetTrace();
539 ASSERTL0(
false,
"Expansion dimension not recognised");
562 for (i = 0; i <
m_fields.num_elements(); i++)
575 for (i = 0; i <
m_fields.num_elements(); i++)
595 for (i = 0; i <
m_fields.num_elements(); i++)
607 for (i = 0; i <
m_fields.num_elements(); i++)
623 "3D fully periodic problems not implemented yet");
627 for (i = 0; i <
m_fields.num_elements(); i++)
638 ASSERTL0(
false,
"Expansion dimension not recognised");
662 m_session->LoadParameter(
"NumQuadPointsError",
690 std::string pFunctionName,
695 "Function '" + pFunctionName +
"' does not exist.");
697 std::vector<std::string> vFieldNames =
m_session->GetVariables();
699 for(
int i = 0 ; i < vFieldNames.size(); i++)
712 std::vector<std::string> pFieldNames,
714 const std::string &pFunctionName,
718 ASSERTL1(pFieldNames.size() == pFields.num_elements(),
719 "Function '" + pFunctionName +
720 "' variable list size mismatch with array storage.");
722 "Function '" + pFunctionName +
"' does not exist.");
724 for (
int i = 0; i < pFieldNames.size(); i++)
737 std::vector<std::string> pFieldNames,
739 const std::string &pFunctionName,
744 "Function '" + pFunctionName +
"' does not exist.");
745 ASSERTL0(pFieldNames.size() == pFields.num_elements(),
746 "Field list / name list size mismatch.");
748 for (
int i = 0; i < pFieldNames.size(); i++)
751 pFunctionName, pTime, domain);
752 pFields[i]->FwdTrans_IterPerExp(pFields[i]->GetPhys(),
753 pFields[i]->UpdateCoeffs());
759 const std::string &pFunctionName,
764 "Function '" + pFunctionName +
"' does not exist.");
767 vType =
m_session->GetFunctionType(pFunctionName, pFieldName, domain);
775 std::string filename =
776 m_session->GetFunctionFilename(pFunctionName, pFieldName, domain);
778 if (boost::filesystem::path(filename).extension() !=
".pts")
793 const string &pFunctionName,
798 "Function '" + pFunctionName +
"' does not exist.");
800 unsigned int nq =
m_fields[0]->GetNpoints();
801 if (pArray.num_elements() < nq)
807 vType =
m_session->GetFunctionType(pFunctionName, pFieldName, domain);
809 "vType not eFunctionTypeExpression");
819 m_session->GetFunction(pFunctionName, pFieldName, domain);
821 ffunc->Evaluate(x0, x1, x2, pTime, pArray);
826 const string &pFunctionName,
832 "Function '" + pFunctionName +
"' does not exist.");
834 unsigned int nq =
m_fields[0]->GetNpoints();
835 if (pArray.num_elements() < nq)
841 vType =
m_session->GetFunctionType(pFunctionName, pFieldName, domain);
844 "vType not eFunctionTypeFile or eFunctionTypeTransientFile");
846 std::string filename =
847 m_session->GetFunctionFilename(pFunctionName, pFieldName, domain);
848 std::string fileVar =
m_session->GetFunctionFilenameVariable(
849 pFunctionName, pFieldName, domain);
851 if (fileVar.length() == 0)
853 fileVar = pFieldName;
862 #if (defined _WIN32 && _MSC_VER < 1900)
865 unsigned int old_exponent_format;
866 old_exponent_format = _set_output_format(_TWO_DIGIT_EXPONENT);
868 _set_output_format(old_exponent_format);
875 ASSERTL0(
false,
"Invalid Filename in function \"" + pFunctionName +
876 "\", variable \"" + fileVar +
"\"")
880 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
881 std::vector<std::vector<NekDouble> > FieldData;
884 int numexp =
m_fields[0]->GetExpSize();
886 for (
int i = 0; i < numexp; ++i)
888 ElementGIDs[i] =
m_fields[0]->GetExp(i)->GetGeom()->GetGlobalID();
894 std::string funcFilename =
895 m_session->GetFunctionFilename(pFunctionName, pFieldName, domain);
909 pts_fld->Import(filename, FieldDef, FieldData,
918 pts_fld->Import(filename, FieldDef, FieldData,
932 for (
int i = 0; i < FieldDef.size(); ++i)
936 for (
int j = 0; j < FieldDef[i]->m_fields.size(); ++j)
938 if (FieldDef[i]->
m_fields[j] == fileVar)
947 FieldDef[i], FieldData[i], FieldDef[i]->
m_fields[idx], vCoeffs);
951 cout <<
"Field " + fileVar +
" not found." << endl;
955 m_fields[0]->BwdTrans_IterPerExp(vCoeffs, pArray);
960 const string &pFunctionName,
966 "Function '" + pFunctionName +
"' does not exist.");
968 unsigned int nq =
m_fields[0]->GetNpoints();
969 if (pArray.num_elements() < nq)
975 vType =
m_session->GetFunctionType(pFunctionName, pFieldName, domain);
978 "vType not eFunctionTypeFile or eFunctionTypeTransientFile");
980 std::string filename =
981 m_session->GetFunctionFilename(pFunctionName, pFieldName, domain);
982 std::string fileVar =
m_session->GetFunctionFilenameVariable(
983 pFunctionName, pFieldName, domain);
985 if (fileVar.length() == 0)
987 fileVar = pFieldName;
996 #if (defined _WIN32 && _MSC_VER < 1900)
999 unsigned int old_exponent_format;
1000 old_exponent_format = _set_output_format(_TWO_DIGIT_EXPONENT);
1002 _set_output_format(old_exponent_format);
1009 ASSERTL0(
false,
"Invalid Filename in function \"" + pFunctionName +
1010 "\", variable \"" + fileVar +
"\"")
1018 std::string funcFilename =
1019 m_session->GetFunctionFilename(pFunctionName, pFieldName, domain);
1031 LoadPts(funcFilename, filename, outPts);
1036 LoadPts(funcFilename, filename, outPts);
1045 vector<string> fieldNames = outPts->GetFieldNames();
1046 for (fieldInd = 0; fieldInd < fieldNames.size(); ++fieldInd)
1048 if (outPts->GetFieldName(fieldInd) == pFieldName)
1053 ASSERTL0(fieldInd != fieldNames.size(),
"field not found");
1055 pArray = outPts->GetPts(fieldInd + outPts->GetDim());
1059 std::string filename,
1062 unsigned int nq =
m_fields[0]->GetNpoints();
1066 ptsIO.
Import(filename, inPts);
1069 inPts->GetNFields());
1070 for (
int i = 0; i < inPts->GetDim() + inPts->GetNFields(); ++i)
1074 if (inPts->GetDim() == 1)
1078 else if (inPts->GetDim() == 2)
1080 m_fields[0]->GetCoords(pts[0], pts[1]);
1082 else if (inPts->GetDim() == 3)
1084 m_fields[0]->GetCoords(pts[0], pts[1], pts[2]);
1087 inPts->GetDim(), inPts->GetFieldNames(), pts);
1094 if (
m_comm->GetRank() == 0)
1100 if (
m_comm->GetRank() == 0)
1103 if (
GetSession()->DefinesCmdLineArgument(
"verbose"))
1120 std::string pFieldName,
1121 const std::string &pFunctionName,
1125 "Function '" + pFunctionName +
"' does not exist.");
1130 vType =
m_session->GetFunctionType(pFunctionName, pFieldName);
1134 =
m_session->GetFunction(pFunctionName, pFieldName,domain);
1135 retVal = ffunc->GetExpression();
1139 std::string filename
1140 =
m_session->GetFunctionFilename(pFunctionName, pFieldName,domain);
1141 retVal =
"from file " + filename;
1154 std::string varName;
1155 int nvariables =
m_fields.num_elements();
1156 for (
int i = 0; i < nvariables; ++i)
1159 m_fields[i]->EvaluateBoundaryConditions(time, varName);
1179 if (
m_fields[field]->GetPhysState() ==
false)
1185 if (exactsoln.num_elements())
1189 else if (
m_session->DefinesFunction(
"ExactSolution"))
1195 "ExactSolution",
m_time);
1204 if (Normalised ==
true)
1212 L2error = sqrt(L2error*L2error/Vol);
1238 if (
m_fields[field]->GetPhysState() ==
false)
1244 if (exactsoln.num_elements())
1248 else if (
m_session->DefinesFunction(
"ExactSolution"))
1254 "ExactSolution",
m_time);
1267 Linferror = L2INF[1];
1306 int ErrorCoordim = ErrorExp->GetCoordim(0);
1307 int ErrorNq = ErrorExp->GetTotPoints();
1313 switch(ErrorCoordim)
1316 ErrorExp->GetCoords(ErrorXc0);
1319 ErrorExp->GetCoords(ErrorXc0, ErrorXc1);
1322 ErrorExp->GetCoords(ErrorXc0, ErrorXc1, ErrorXc2);
1326 m_session->GetFunction(
"ExactSolution", field);
1331 exSol->Evaluate(ErrorXc0,ErrorXc1,ErrorXc2,
m_time,ErrorSol);
1335 ErrorExp->BwdTrans_IterPerExp(
m_fields[field]->GetCoeffs(),
1336 ErrorExp->UpdatePhys());
1338 L2INF[0] = ErrorExp->L2 (ErrorExp->GetPhys(), ErrorSol);
1339 L2INF[1] = ErrorExp->Linf(ErrorExp->GetPhys(), ErrorSol);
1352 bool dumpInitialConditions,
1355 if (
m_session->GetComm()->GetRank() == 0)
1357 cout <<
"Initial Conditions:" << endl;
1360 if (
m_session->DefinesFunction(
"InitialConditions"))
1363 "InitialConditions",
m_time, domain);
1365 if (
m_session->GetComm()->GetRank() == 0)
1368 for (
int i = 0; i <
m_fields.num_elements(); ++i)
1370 std::string varName =
m_session->GetVariable(i);
1371 cout <<
" - Field " << varName <<
": "
1379 int nq =
m_fields[0]->GetNpoints();
1380 for (
int i = 0; i <
m_fields.num_elements(); i++)
1386 if (
m_session->GetComm()->GetRank() == 0)
1388 cout <<
" - Field " <<
m_session->GetVariable(i)
1389 <<
": 0 (default)" << endl;
1408 "ExactSolution array size mismatch.");
1409 Vmath::Zero(outfield.num_elements(), outfield, 1);
1410 if (
m_session->DefinesFunction(
"ExactSolution"))
1413 "ExactSolution", time);
1431 std::vector<std::string> vel;
1432 vel.push_back(
"Vx");
1433 vel.push_back(
"Vy");
1434 vel.push_back(
"Vz");
1447 int nvariables =
m_session->GetVariables().size();
1456 for(i = 0; i <
m_base.num_elements(); i++)
1476 for (i = 0 ; i <
m_base.num_elements(); i++)
1480 ::AllocateSharedPtr(
1485 m_base[i]->SetWaveSpace(
true);
1497 for (i = 0 ; i <
m_base.num_elements(); i++)
1501 ::AllocateSharedPtr(
1506 m_base[i]->SetWaveSpace(
true);
1517 for (i = 0 ; i <
m_base.num_elements(); i++)
1521 ::AllocateSharedPtr(
1526 m_base[i]->SetWaveSpace(
false);
1539 for (i = 1 ; i < m_base.num_elements(); i++)
1556 for (i = 1 ; i < m_base.num_elements(); i++)
1565 ASSERTL0(
false,
"Expansion dimension not recognised");
1577 for(i = 0 ; i <
m_base.num_elements(); i++)
1588 for(i = 0 ; i <
m_base.num_elements(); i++)
1600 ASSERTL0(
false,
"Expansion dimension not recognised");
1608 std::string pInfile,
1611 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
1612 std::vector<std::vector<NekDouble> > FieldData;
1617 base_fld->Import(pInfile,FieldDef,FieldData);
1619 int nvar =
m_session->GetVariables().size();
1620 if (
m_session->DefinesSolverInfo(
"HOMOGENEOUS"))
1622 std::string HomoStr =
m_session->GetSolverInfo(
"HOMOGENEOUS");
1625 for (
int j = 0; j < nvar; ++j)
1627 for(
int i = 0; i < FieldDef.size(); ++i)
1629 bool flag = FieldDef[i]->m_fields[j] ==
1631 ASSERTL0(flag, (std::string(
"Order of ") + pInfile
1632 + std::string(
" data and that defined in "
1633 "the session differs")).c_str());
1635 m_base[j]->ExtractDataToCoeffs(FieldDef[i], FieldData[i],
1637 m_base[j]->UpdateCoeffs());
1698 for (
int i = 0; i <
m_fields.num_elements(); i++)
1710 for (
int i = 0; i <
m_fields.num_elements(); i++)
1732 m_fields[0]->IProductWRTDerivBase(F,outarray);
1746 int ndim = F.num_elements();
1747 int nPointsTot =
m_fields[0]->GetNpoints();
1752 for (
int i = 0; i < ndim; ++i)
1759 m_fields[0]->IProductWRTBase(div, outarray);
1776 int ndim = V.num_elements();
1778 int nPointsTot =
m_fields[0]->GetNpoints();
1786 m_fields[0]->IProductWRTBase(tmp, outarray,
1791 m_fields[0]->IProductWRTBase_IterPerExp(tmp, outarray);
1809 int ndim = V.num_elements();
1814 int nPointsTot =
m_fields[0]->GetNpoints();
1818 if (wk.num_elements())
1832 Vmath::Vmul(nPointsTot, grad0, 1, V[0], 1, outarray,1);
1836 m_fields[0]->PhysDeriv(u, grad0, grad1);
1837 Vmath::Vmul (nPointsTot, grad0, 1, V[0], 1, outarray, 1);
1839 outarray, 1, outarray, 1);
1844 m_fields[0]->PhysDeriv(u,grad0,grad1,grad2);
1845 Vmath::Vmul (nPointsTot, grad0, 1, V[0], 1, outarray, 1);
1847 outarray, 1, outarray, 1);
1849 outarray, 1, outarray, 1);
1852 ASSERTL0(
false,
"dimension unknown");
1869 bool NumericalFluxIncludesNormal,
1870 bool InFieldIsInPhysSpace,
1881 nvariables =
m_fields.num_elements();
1887 for(i = 0; i < nVelDim; ++i)
1894 if (InFieldIsInPhysSpace ==
true)
1896 for (i = 0; i < nvariables; ++i)
1898 physfield[i] = InField[i];
1904 for(i = 0; i < nvariables; ++i)
1908 m_fields[i]->BwdTrans(InField[i],physfield[i]);
1913 for (i = 0; i < nvariables; ++i)
1925 if (NumericalFluxIncludesNormal ==
true)
1929 for (i = 0; i < nvariables; ++i)
1939 for (i = 0; i < nvariables; ++i)
1942 m_fields[i]->AddTraceIntegral(numflux[i],OutField[i]);
1953 for (i = 0; i < nvariables; ++i)
1964 for(i = 0; i < nvariables; ++i)
1967 m_fields[i]->AddTraceIntegral(numfluxX[i], numfluxY[i],
1983 bool NumericalFluxIncludesNormal,
1984 bool InFieldIsInPhysSpace)
1990 int nvariables =
m_fields.num_elements();
2002 for (j = 0; j < nqvar; ++j)
2007 for (i = 0; i< nvariables; ++i)
2021 if (InFieldIsInPhysSpace ==
true)
2023 for (i = 0; i < nvariables; ++i)
2025 ufield[i] = InField[i];
2031 for (i = 0; i < nvariables; ++i)
2035 m_fields[i]->BwdTrans(InField[i],ufield[i]);
2045 for (j = 0; j < nqvar; ++j)
2047 for (i = 0; i < nvariables; ++i)
2057 ufield[i], 1, fluxvector[k], 1);
2069 m_fields[i]->AddTraceIntegral(flux[j][i], qcoeffs);
2084 m_fields[i]->MultiplyByElmtInvMass(qcoeffs, qcoeffs);
2087 m_fields[i]->BwdTrans(qcoeffs, qfield[j][i]);
2099 for (i = 0; i < nvariables; ++i)
2107 for (j = 0; j < nqvar; ++j)
2112 qfield[j][i], 1, fluxvector[k], 1);
2133 m_fields[i]->AddTraceIntegral(flux[0][i], OutField[i]);
2145 boost::lexical_cast<std::string>(n);
2157 std::vector<std::string> &variables)
2160 boost::lexical_cast<std::string>(n);
2161 WriteFld(outname, field, fieldcoeffs, variables);
2171 boost::lexical_cast<std::string>(n);
2182 std::vector<Array<OneD, NekDouble> > fieldcoeffs(
2184 std::vector<std::string> variables(
m_fields.num_elements());
2186 for (
int i = 0; i <
m_fields.num_elements(); ++i)
2190 fieldcoeffs[i] =
m_fields[i]->UpdateCoeffs();
2218 const std::string &outname,
2221 std::vector<std::string> &variables)
2223 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef
2224 = field->GetFieldDefinitions();
2225 std::vector<std::vector<NekDouble> > FieldData(FieldDef.size());
2228 for(
int j = 0; j < fieldcoeffs.size(); ++j)
2230 for(
int i = 0; i < FieldDef.size(); ++i)
2233 FieldDef[i]->m_fields.push_back(variables[j]);
2234 field->AppendFieldData(FieldDef[i], FieldData[i],
2258 mapping->Output( fieldMetaDataMap, outname);
2260 m_fld->Write(outname, FieldDef, FieldData, fieldMetaDataMap,
true);
2271 const std::string &infile,
2274 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
2275 std::vector<std::vector<NekDouble> > FieldData;
2278 field_fld->Import(infile,FieldDef,FieldData);
2281 for(
int j = 0; j < pFields.num_elements(); ++j)
2284 pFields[j]->UpdateCoeffs(),1);
2286 for(
int i = 0; i < FieldDef.size(); ++i)
2290 std::string(
"Order of ") + infile
2291 + std::string(
" data and that defined in "
2292 "m_boundaryconditions differs"));
2294 pFields[j]->ExtractDataToCoeffs(FieldDef[i], FieldData[i],
2296 pFields[j]->UpdateCoeffs());
2298 pFields[j]->BwdTrans(pFields[j]->GetCoeffs(),
2299 pFields[j]->UpdatePhys());
2313 const std::string &infile,
2317 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
2318 std::vector<std::vector<NekDouble> > FieldData;
2324 ASSERTL0(ndomains*nvariables == pFields.num_elements(),
"Number of fields does not match the number of variables and domains");
2327 for(
int j = 0; j < ndomains; ++j)
2329 for(
int i = 0; i < nvariables; ++i)
2333 for(
int n = 0; n < FieldDef.size(); ++n)
2336 std::string(
"Order of ") + infile
2337 + std::string(
" data and that defined in "
2338 "m_boundaryconditions differs"));
2340 pFields[j*nvariables+i]->ExtractDataToCoeffs(FieldDef[n], FieldData[n],
2342 pFields[j*nvariables+i]->UpdateCoeffs());
2344 pFields[j*nvariables+i]->BwdTrans(pFields[j*nvariables+i]->GetCoeffs(),
2345 pFields[j*nvariables+i]->UpdatePhys());
2356 const std::string &infile,
2358 std::string &pFieldName)
2360 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
2361 std::vector<std::vector<NekDouble> > FieldData;
2365 field_fld->Import(infile,FieldDef,FieldData);
2368 Vmath::Zero(pField->GetNcoeffs(),pField->UpdateCoeffs(),1);
2370 for(
int i = 0; i < FieldDef.size(); ++i)
2373 for(
int j = 0; j < FieldData.size(); ++j)
2375 if (FieldDef[i]->
m_fields[j] == pFieldName)
2380 ASSERTL1(idx >= 0,
"Field " + pFieldName +
" not found.");
2382 pField->ExtractDataToCoeffs(FieldDef[i], FieldData[i],
2384 pField->UpdateCoeffs());
2386 pField->BwdTrans(pField->GetCoeffs(), pField->UpdatePhys());
2397 const std::string &infile,
2398 std::vector< std::string> &fieldStr,
2402 ASSERTL0(fieldStr.size() <= coeffs.num_elements(),
2403 "length of fieldstr should be the same as pFields");
2405 std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef;
2406 std::vector<std::vector<NekDouble> > FieldData;
2410 field_fld->Import(infile,FieldDef,FieldData);
2413 for(
int j = 0; j < fieldStr.size(); ++j)
2415 Vmath::Zero(coeffs[j].num_elements(),coeffs[j],1);
2416 for(
int i = 0; i < FieldDef.size(); ++i)
2418 m_fields[0]->ExtractDataToCoeffs(FieldDef[i], FieldData[i],
2419 fieldStr[j], coeffs[j]);
2435 if (
m_session->GetComm()->GetSize() > 1)
2476 if (
m_session->DefinesSolverInfo(
"UpwindType"))
2479 m_session->GetSolverInfo(
"UpwindType"));
2482 if (
m_session->DefinesSolverInfo(
"AdvectionType"))
2484 std::string AdvectionType;
2485 AdvectionType =
m_session->GetSolverInfo(
"AdvectionType");
2487 GetClassDescription(AdvectionType));
2501 "Mixed Continuous Galerkin and Discontinuous");
2504 if (
m_session->DefinesSolverInfo(
"DiffusionType"))
2506 std::string DiffusionType;
2507 DiffusionType =
m_session->GetSolverInfo(
"DiffusionType");
2509 GetClassDescription(DiffusionType));
2526 string::const_iterator it1=s1.begin();
2527 string::const_iterator it2=s2.begin();
2530 while ( (it1!=s1.end()) && (it2!=s2.end()) )
2532 if(::toupper(*it1) != ::toupper(*it2))
2535 return (::toupper(*it1) < ::toupper(*it2)) ? -1 : 1;
2543 size_t size1=s1.size();
2544 size_t size2=s2.size();
2552 return (size1 < size2) ? -1 : 1;
2565 ASSERTL0(
false,
"v_GetFluxVector: This function is not valid "
2566 "for the Base class");
2570 const int i,
const int j,
2574 ASSERTL0(
false,
"v_GetqFluxVector: This function is not valid "
2575 "for the Base class");
2584 ASSERTL0(
false,
"v_GetFluxVector: This function is not valid "
2585 "for the Base class");
2592 ASSERTL0(
false,
"v_NumericalFlux: This function is not valid "
2593 "for the Base class");
2601 ASSERTL0(
false,
"v_NumericalFlux: This function is not valid "
2602 "for the Base class");
2609 ASSERTL0(
false,
"v_NumFluxforScalar: This function is not valid "
2610 "for the Base class");
2618 ASSERTL0(
false,
"v_NumFluxforVector: This function is not valid "
2619 "for the Base class");
2624 ASSERTL0(
false,
"This function is not valid for the Base class");
2631 std::vector<std::string> &variables)
SOLVER_UTILS_EXPORT void EvaluateFunctionPts(std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, const NekDouble &pTime=0.0, const int domain=0)
virtual SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &l)
Virtual function for generating summary information.
bool m_singleMode
Flag to determine if single homogeneous mode is used.
SOLVER_UTILS_EXPORT void NumFluxforScalar(const Array< OneD, Array< OneD, NekDouble > > &ufield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &uflux)
f SOLVER_UTILS_EXPORT void AdvectionNonConservativeForm(const Array< OneD, Array< OneD, NekDouble > > &V, const Array< OneD, const NekDouble > &u, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wk=NullNekDouble1DArray)
Compute the non-conservative advection.
#define ASSERTL0(condition, msg)
SOLVER_UTILS_EXPORT void InitialiseBaseFlow(Array< OneD, Array< OneD, NekDouble > > &base)
Perform initialisation of the base flow.
A class that contains algorithms for interpolation between pts fields, expansions and different meshe...
LibUtilities::NekFactory< std::string, EquationSystem, const LibUtilities::SessionReaderSharedPtr & > EquationSystemFactory
Datatype of the NekFactory used to instantiate classes derived from the EquationSystem class...
SOLVER_UTILS_EXPORT std::string DescribeFunction(std::string pFieldName, const std::string &pFunctionName, const int domain)
Provide a description of a function for a given field name.
static boost::shared_ptr< MeshGraph > Read(const LibUtilities::SessionReaderSharedPtr &pSession, DomainRangeShPtr &rng=NullDomainRangeShPtr)
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
SOLVER_UTILS_EXPORT void LoadPts(std::string funcFilename, std::string filename, Nektar::LibUtilities::PtsFieldSharedPtr &outPts)
NekDouble m_time
Current time of simulation.
SOLVER_UTILS_EXPORT void WeakAdvectionNonConservativeForm(const Array< OneD, Array< OneD, NekDouble > > &V, const Array< OneD, const NekDouble > &u, Array< OneD, NekDouble > &outarray, bool UseContCoeffs=false)
Compute the inner product .
std::map< std::string, std::pair< std::string, loadedFldField > > m_loadedFldFields
void Import(const std::string &infilename, std::vector< FieldDefinitionsSharedPtr > &fielddefs, std::vector< std::vector< NekDouble > > &fielddata, FieldMetaDataMap &fieldinfomap, const Array< OneD, int > &ElementIDs)
This function allows for data to be imported from an FLD file when a session and/or communicator is n...
NekDouble m_timestep
Time step size.
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
SOLVER_UTILS_EXPORT int NoCaseStringCompare(const std::string &s1, const std::string &s2)
Perform a case-insensitive string comparison.
std::vector< std::pair< std::string, std::string > > SummaryList
SOLVER_UTILS_EXPORT int GetNvariables()
bool m_halfMode
Flag to determine if half homogeneous mode is used.
NekDouble m_LhomZ
physical length in Z direction (if homogeneous)
virtual SOLVER_UTILS_EXPORT NekDouble v_L2Error(unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray, bool Normalised=false)
Virtual function for the L_2 error computation between fields and a given exact solution.
Array< OneD, bool > m_checkIfSystemSingular
Flag to indicate if the fields should be checked for singularity.
int m_expdim
Expansion dimension.
SOLVER_UTILS_EXPORT void Checkpoint_Output(const int n)
Write checkpoint file of m_fields.
DiffusionFactory & GetDiffusionFactory()
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
enum MultiRegions::ProjectionType m_projectionType
Type of projection; e.g continuous or discontinuous.
Principle Modified Functions .
std::map< std::string, std::pair< std::string, LibUtilities::PtsFieldSharedPtr > > m_loadedPtsFields
pts fields we already read from disk: {funcFilename: (filename, ptsfield)}
SOLVER_UTILS_EXPORT Array< OneD, NekDouble > ErrorExtraPoints(unsigned int field)
Compute error (L2 and L_inf) over an larger set of quadrature points return [L2 Linf].
SOLVER_UTILS_EXPORT void WeakAdvectionGreensDivergenceForm(const Array< OneD, Array< OneD, NekDouble > > &F, Array< OneD, NekDouble > &outarray)
Compute the inner product .
NekDouble m_checktime
Time between checkpoints.
bool m_useFFT
Flag to determine if FFT is used for homogeneous transform.
NekDouble m_LhomY
physical length in Y direction (if homogeneous)
boost::shared_ptr< ContField2D > ContField2DSharedPtr
bool m_specHP_dealiasing
Flag to determine if dealisising is usde for the Spectral/hp element discretisation.
SOLVER_UTILS_EXPORT int GetNumExpModes()
virtual SOLVER_UTILS_EXPORT void v_ExtraFldOutput(std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
SOLVER_UTILS_EXPORT void PrintProgressbar(const int position, const int goal) const
virtual SOLVER_UTILS_EXPORT NekDouble v_LinfError(unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray)
Virtual function for the L_inf error computation between fields and a given exact solution...
int m_npointsZ
number of points in Z direction (if homogeneous)
virtual SOLVER_UTILS_EXPORT ~EquationSystem()
Destructor.
std::string m_sessionName
Name of the session.
int m_nchk
Number of checkpoints written so far.
SOLVER_UTILS_EXPORT void WeakAdvectionDivergenceForm(const Array< OneD, Array< OneD, NekDouble > > &F, Array< OneD, NekDouble > &outarray)
Compute the inner product .
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
std::map< std::string, std::string > FieldMetaDataMap
int m_checksteps
Number of steps between checkpoints.
SOLVER_UTILS_EXPORT void EvaluateFunctionExp(std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, const NekDouble &pTime=0.0, const int domain=0)
boost::shared_ptr< PtsField > PtsFieldSharedPtr
LibUtilities::CommSharedPtr m_comm
Communicator.
Gauss Radau pinned at x=-1, .
NekDouble m_fintime
Finish time of the simulation.
virtual SOLVER_UTILS_EXPORT void v_InitObject()
Initialisation object for EquationSystem.
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_tanbasis
2 x m_spacedim x nq
int m_steps
Number of steps to take.
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Array holding trace normals for DG simulations in the forwards direction.
int m_HomoDirec
number of homogenous directions
SOLVER_UTILS_EXPORT void FwdTransFields()
SOLVER_UTILS_EXPORT void WeakDGDiffusion(const Array< OneD, Array< OneD, NekDouble > > &InField, Array< OneD, Array< OneD, NekDouble > > &OutField, bool NumericalFluxIncludesNormal=true, bool InFieldIsInPhysSpace=false)
Calculate weak DG Diffusion in the LDG form.
SOLVER_UTILS_EXPORT void EvaluateFunctionFld(std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, const NekDouble &pTime=0.0, const int domain=0)
1D Evenly-spaced points using Fourier Fit
bool m_multipleModes
Flag to determine if use multiple homogenenous modes are used.
Fourier Modified expansions with just the real part of the first mode .
virtual SOLVER_UTILS_EXPORT Array< OneD, bool > v_GetSystemSingularChecks()
SOLVER_UTILS_EXPORT void NumericalFlux(Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &numflux)
SOLVER_UTILS_EXPORT void SetUpBaseFields(SpatialDomains::MeshGraphSharedPtr &mesh)
This class is the abstraction of a global continuous two- dimensional spectral/hp element expansion w...
int m_npointsY
number of points in Y direction (if homogeneous)
Abstraction of a global continuous one-dimensional spectral/hp element expansion which approximates t...
NekMatrix< NekMatrix< NekMatrix< NekDouble, StandardMatrixTag >, ScaledMatrixTag >, BlockMatrixTag > DNekScalBlkMat
virtual SOLVER_UTILS_EXPORT void v_TransPhysToCoeff()
Virtual function for transformation to coefficient space.
virtual SOLVER_UTILS_EXPORT void v_GetFluxVector(const int i, Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &flux)
Principle Modified Functions .
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
virtual SOLVER_UTILS_EXPORT void v_Output(void)
SpatialDomains::BoundaryConditionsSharedPtr m_boundaryConditions
Pointer to boundary conditions object.
boost::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
int m_spacedim
Spatial dimension (>= expansion dim).
This class is the abstraction of a global discontinuous two- dimensional spectral/hp element expansio...
LibUtilities::FieldMetaDataMap m_fieldMetaDataMap
Map to identify relevant solver info to dump in output fields.
AdvectionFactory & GetAdvectionFactory()
Gets the factory for initialising advection objects.
Defines a specification for a set of points.
void Neg(int n, T *x, const int incx)
Negate x = -x.
boost::shared_ptr< FieldIO > FieldIOSharedPtr
virtual SOLVER_UTILS_EXPORT void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
SOLVER_UTILS_EXPORT LibUtilities::SessionReaderSharedPtr GetSession()
Get Session name.
SOLVER_UTILS_EXPORT void GetFluxVector(const int i, Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &flux)
boost::shared_ptr< ExpList2D > ExpList2DSharedPtr
Shared pointer to an ExpList2D object.
SOLVER_UTILS_EXPORT void EvaluateFunction(Array< OneD, Array< OneD, NekDouble > > &pArray, std::string pFunctionName, const NekDouble pTime=0.0, const int domain=0)
Evaluates a function as specified in the session file.
GLOBAL_MAPPING_EXPORT typedef boost::shared_ptr< Mapping > MappingSharedPtr
A shared pointer to a Mapping object.
void Import(const string &inFile, PtsFieldSharedPtr &ptsField, FieldMetaDataMap &fieldmetadatamap=NullFieldMetaDataMap)
Import a pts field from file.
boost::shared_ptr< Equation > EquationSharedPtr
virtual SOLVER_UTILS_EXPORT void v_DoInitialise()
Virtual function for initialisation implementation.
bool m_homogen_dealiasing
Flag to determine if dealiasing is used for homogeneous simulations.
EquationSystemFactory & GetEquationSystemFactory()
Fourier Modified expansions with just the imaginary part of the first mode .
SOLVER_UTILS_EXPORT void ImportFldToMultiDomains(const std::string &infile, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int ndomains)
Input field data from the given file to multiple domains.
SOLVER_UTILS_EXPORT void SetBoundaryConditions(NekDouble time)
Evaluates the boundary conditions at the given time.
MultiRegions::Direction const DirCartesianMap[]
static boost::shared_ptr< FieldIO > CreateDefault(const LibUtilities::SessionReaderSharedPtr session)
Returns an object for the default FieldIO method.
SOLVER_UTILS_EXPORT void WriteFld(const std::string &outname)
Write field data to the given filename.
static boost::shared_ptr< FieldIO > CreateForFile(const LibUtilities::SessionReaderSharedPtr session, const std::string &filename)
Construct a FieldIO object for a given input filename.
SOLVER_UTILS_EXPORT int GetNpoints()
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
Fourier ModifiedExpansion with just the first mode .
Array< OneD, MultiRegions::ExpListSharedPtr > m_base
Base fields.
SOLVER_UTILS_EXPORT int GetTraceNpoints()
SOLVER_UTILS_EXPORT int GetNcoeffs()
1D Non Evenly-spaced points for Single Mode analysis
LibUtilities::FieldIOSharedPtr m_fld
Field input/output.
virtual SOLVER_UTILS_EXPORT void v_NumFluxforVector(const Array< OneD, Array< OneD, NekDouble > > &ufield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &qfield, Array< OneD, Array< OneD, NekDouble > > &qflux)
virtual SOLVER_UTILS_EXPORT bool v_NegatedOp()
Virtual function to identify if operator is negated in DoSolve.
virtual SOLVER_UTILS_EXPORT MultiRegions::ExpListSharedPtr v_GetPressure(void)
SOLVER_UTILS_EXPORT void NumFluxforVector(const Array< OneD, Array< OneD, NekDouble > > &ufield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &qfield, Array< OneD, Array< OneD, NekDouble > > &qflux)
SOLVER_UTILS_EXPORT void ZeroPhysFields()
SOLVER_UTILS_EXPORT void SessionSummary(SummaryList &vSummary)
Write out a session summary.
SpatialDomains::MeshGraphSharedPtr m_graph
Pointer to graph defining mesh.
virtual SOLVER_UTILS_EXPORT void v_NumFluxforScalar(const Array< OneD, Array< OneD, NekDouble > > &ufield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &uflux)
SOLVER_UTILS_EXPORT void WeakDGAdvection(const Array< OneD, Array< OneD, NekDouble > > &InField, Array< OneD, Array< OneD, NekDouble > > &OutField, bool NumericalFluxIncludesNormal=true, bool InFieldIsInPhysSpace=false, int nvariables=0)
Calculate the weak discontinuous Galerkin advection.
Used to lookup the create function in NekManager.
SOLVER_UTILS_EXPORT void Checkpoint_BaseFlow(const int n)
Write base flow file of m_fields.
SOLVER_UTILS_EXPORT void ImportFldBase(std::string pInfile, SpatialDomains::MeshGraphSharedPtr pGraph)
virtual SOLVER_UTILS_EXPORT void v_DoSolve()
Virtual function for solve implementation.
std::map< std::string, FieldUtils::Interpolator > m_interpolators
Map of interpolator objects.
boost::shared_ptr< ContField3D > ContField3DSharedPtr
void Zero(int n, T *x, const int incx)
Zero vector.
int m_NumQuadPointsError
Number of Quadrature points used to work out the error.
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
static FieldMetaDataMap NullFieldMetaDataMap
NekMatrix< NekMatrix< NekDouble, StandardMatrixTag >, ScaledMatrixTag > DNekScalMat
static GLOBAL_MAPPING_EXPORT MappingSharedPtr Load(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
Return a pointer to the mapping, creating it on first call.
virtual SOLVER_UTILS_EXPORT void v_EvaluateExactSolution(unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
virtual SOLVER_UTILS_EXPORT void v_TransCoeffToPhys()
Virtual function for transformation to physical space.
Describes the specification for a Basis.
virtual SOLVER_UTILS_EXPORT void v_NumericalFlux(Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &numflux)
enum HomogeneousType m_HomogeneousType
1D Gauss-Lobatto-Legendre quadrature points
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Provides a generic Factory class.
SOLVER_UTILS_EXPORT void ImportFld(const std::string &infile, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
Input field data from the given file.