44#include <boost/algorithm/string.hpp>
60 const std::weak_ptr<EquationSystem> &pEquation,
const ParamMap &pParams)
61 :
Filter(pSession, pEquation)
64 auto it = pParams.find(
"OutputFile");
65 if (it == pParams.end())
71 ASSERTL0(it->second.length() > 0,
"Missing parameter 'OutputFile'.");
82 it = pParams.find(
"OutputFrequency");
83 if (it == pParams.end())
94 it = pParams.find(
"StartTime");
95 if (it == pParams.end())
109 it = pParams.find(
"OutputAllPlanes");
110 if (it == pParams.end())
116 std::string sOption = it->second.c_str();
118 (boost::iequals(sOption,
"yes"));
127 it = pParams.find(
"Boundary");
128 ASSERTL0(it != pParams.end(),
"Missing parameter 'Boundary");
129 ASSERTL0(it->second.length() > 0,
"Empty parameter 'Boundary'.");
144 for (
int i = 0; i < 3; ++i)
152 std::stringstream directionStream;
153 std::string directionString;
155 for (
int i = 0; i < 3; ++i)
157 std::stringstream tmp;
159 std::string dir =
"Direction" + tmp.str();
160 it = pParams.find(dir);
161 if (it != pParams.end())
163 ASSERTL0(!(it->second.empty()),
"Missing parameter '" + dir +
"'.");
164 directionStream.str(it->second);
167 directionStream.clear();
170 for (
int j = 0; j < 3; j++)
172 directionStream >> directionString;
173 if (!directionString.empty())
182 for (
int j = 0; j < 3; j++)
189 for (
int i = 0; i < 3; ++i)
191 for (
int j = 0; j < 3; ++j)
199 it = pParams.find(
"PivotPoint");
200 if (it == pParams.end())
202 it = pParams.find(
"MomentPoint");
205 if (it != pParams.end())
207 ASSERTL0(!(it->second.empty()),
"Missing parameter 'PivotPoint'.");
208 std::stringstream pivotPointStream;
209 std::string pivotPointString;
210 pivotPointStream.str(it->second);
212 for (
int j = 0; j < 3; ++j)
214 pivotPointStream >> pivotPointString;
215 if (!pivotPointString.empty())
247 (std::string(
"Error reading boundary region definition:") +
251 std::string IndString =
256 (std::string(
"Unable to read boundary regions index "
257 "range for FilterAeroForces: ") +
263 unsigned int numBoundaryRegions = pFields[0]->GetBndConditions().size();
265 numBoundaryRegions, 0);
272 for (
auto &it : bregions)
286 pFields[0]->GetPlane(0)->GetBoundaryToElmtMap(
m_BCtoElmtID,
300 m_nPlanes = pFields[0]->GetHomogeneousBasis()->GetZ().size();
317 for (j = 0; j < IDs.size(); ++j)
341 int expdim = pFields[0]->GetGraph()->GetMeshDimension();
342 int momdim = (expdim == 2) ? 1 : 3;
343 if (vComm->GetRank() == 0)
347 m_session->MatchSolverInfo(
"Driver",
"Adaptive", adaptive,
false);
356 m_outputStream <<
"# Forces and moments acting on bodies" << endl;
357 for (
int i = 0; i < expdim; i++)
360 <<
" Direction" << i + 1 <<
" = (";
381 m_outputStream <<
"# Boundary regions: " << IndString.c_str() << endl;
385 for (
int i = 1; i <= expdim; i++)
406 for (
int i = 1; i <= momdim; i++)
421 if (
m_session->DefinesSolverInfo(
"HomoStrip"))
425 "Output forces on all planes not compatible with HomoStrips");
459 int expdim = pFields[0]->GetGraph()->GetMeshDimension();
460 int momdim = (expdim == 2) ? 1 : 3;
464 for (
int i = 0; i < expdim; i++)
468 Ft[i] = Fp[i] + Fv[i];
474 for (
int i = 0; i < momdim; i++)
478 Mt[i] = Mp[i] + Mv[i];
490 if (
m_session->DefinesSolverInfo(
"HomoStrip"))
492 m_session->LoadParameter(
"Strip_Z", nstrips);
493 colSize = vComm->GetColumnComm()->GetSize();
494 idOffset = colSize / nstrips;
499 m_session->LoadParameter(
"HomModesZ", nZ);
503 if (vComm->GetRank() == 0)
508 for (
int plane = 0; plane <
m_nPlanes; plane++)
514 for (
int i = 0; i < expdim; i++)
524 for (
int i = 0; i < momdim; i++)
541 for (
int i = 0; i < expdim; i++)
550 for (
int i = 0; i < momdim; i++)
565 if (
m_session->DefinesSolverInfo(
"HomoStrip"))
574 for (
int i = 1; i < nstrips; i++)
576 int id = i * idOffset;
577 vComm->GetColumnComm()->Recv(
id, Fp);
578 vComm->GetColumnComm()->Recv(
id, Fv);
579 vComm->GetColumnComm()->Recv(
id, Ft);
580 vComm->GetColumnComm()->Recv(
id, Mp);
581 vComm->GetColumnComm()->Recv(
id, Mv);
582 vComm->GetColumnComm()->Recv(
id, Mt);
586 for (
int j = 0; j < expdim; j++)
595 for (
int j = 0; j < momdim; j++)
617 if (
m_session->DefinesSolverInfo(
"HomoStrip") &&
618 (vComm->GetRowComm()->GetRank() == 0))
622 for (
int i = 1; i < nstrips; ++i)
624 int rid = i * idOffset;
625 int sid = vComm->GetColumnComm()->GetRank();
628 vComm->GetColumnComm()->Send(0, Fp);
629 vComm->GetColumnComm()->Send(0, Fv);
630 vComm->GetColumnComm()->Send(0, Ft);
631 vComm->GetColumnComm()->Send(0, Mp);
632 vComm->GetColumnComm()->Send(0, Mv);
633 vComm->GetColumnComm()->Send(0, Mt);
647 if (pFields[0]->GetComm()->GetRank() == 0)
680 int local_planes = ZIDs.size();
681 int expdim = pFields[0]->GetGraph()->GetMeshDimension();
686 for (
int plane = 0; plane < local_planes; plane++)
688 for (
int dir = 0; dir < expdim; dir++)
690 Aeroforces[plane + dir * local_planes] =
697 for (
int plane = 0; plane < local_planes; plane++)
699 for (
int dir = 0; dir < expdim; dir++)
701 Aeroforces[plane + dir * local_planes] =
m_Ft[dir];
724 int local_planes = ZIDs.size();
725 int expdim = pFields[0]->GetGraph()->GetMeshDimension();
726 int momdim = (expdim == 2) ? 1 : 3;
731 for (
int plane = 0; plane < local_planes; plane++)
733 for (
int dir = 0; dir < momdim; ++dir)
735 moments[plane + dir * local_planes] =
742 for (
int plane = 0; plane < local_planes; ++plane)
744 for (
int dir = 0; dir < momdim; dir++)
746 moments[plane + dir * local_planes] =
m_Mt[dir];
773 auto equ =
m_equ.lock();
774 ASSERTL0(equ,
"Weak pointer expired");
776 auto fluidEqu = std::dynamic_pointer_cast<FluidInterface>(equ);
777 ASSERTL0(fluidEqu,
"Aero forces filter is incompatible with this solver.");
782 bnu::matrix<NekDouble> projMat = bnu::identity_matrix<NekDouble>(3, 3);
786 fluidEqu->GetMovingFrameProjectionMat(projMat);
789 for (
int idir = 0; idir < 3; ++idir)
791 bnu::vector<NekDouble> v0 = bnu::zero_vector<NekDouble>(3);
792 bnu::vector<NekDouble> v1 = bnu::zero_vector<NekDouble>(3);
794 for (
int j = 0; j < 3; ++j)
799 v1 = bnu::prec_prod(projMat, v0);
801 for (
int j = 0; j < 3; ++j)
808 int cnt, elmtid, nq, offset, boundary;
810 int physTot = pFields[0]->GetNpoints();
811 int expdim = pFields[0]->GetGraph()->GetMeshDimension();
812 int momdim = (expdim == 2) ? 1 : 3;
848 m_session->DefinesSolverInfo(
"HomoStrip")
849 ? vComm->GetColumnComm()->GetColumnComm()
850 : vComm->GetColumnComm();
863 for (
int i = 0; i < expdim; i++)
874 for (
int i = 0; i < momdim; ++i)
891 if (
m_session->DefinesParameter(
"Kinvis"))
893 rho = (
m_session->DefinesParameter(
"rho"))
896 mu = rho *
m_session->GetParameter(
"Kinvis");
906 for (
int i = 0; i < pFields.size(); ++i)
910 pFields[i]->SetWaveSpace(
false);
912 pFields[i]->BwdTrans(pFields[i]->GetCoeffs(), pFields[i]->UpdatePhys());
913 pFields[i]->SetPhysState(
true);
921 BndConds = pFields[0]->GetPlane(0)->GetBndConditions();
922 BndExp = pFields[0]->GetPlane(0)->GetBndCondExpansions();
926 BndConds = pFields[0]->GetBndConditions();
927 BndExp = pFields[0]->GetBndCondExpansions();
932 for (
int plane = 0; plane <
m_nPlanes; plane++)
941 for (
int n = 0; n < pFields.size(); n++)
943 fields[n] = pFields[n]->GetPlane(
m_planesID[plane]);
948 for (
int n = 0; n < pFields.size(); n++)
950 fields[n] = pFields[n];
955 for (
int n = 0; n < physfields.size(); ++n)
957 physfields[n] = fields[n]->GetPhys();
959 for (
int n = 0; n < nVel; ++n)
964 fluidEqu->GetVelocity(physfields,
velocity);
965 fluidEqu->GetPressure(physfields,
pressure);
968 for (
int n = cnt = 0; n < BndConds.size(); n++)
972 for (
int i = 0; i < BndExp[n]->GetExpSize(); ++i, cnt++)
975 elmt = fields[0]->GetExp(elmtid);
976 nq = elmt->GetTotPoints();
977 offset = fields[0]->GetPhys_Offset(elmtid);
980 for (
int j = 0; j < expdim; j++)
988 for (
int j = 0; j < expdim; j++)
990 for (
int k = 0; k < expdim; k++)
992 grad[j * expdim + k] =
994 elmt->PhysDeriv(k, velElmt[j],
995 grad[j * expdim + k]);
1008 for (
int j = 0; j < 3; ++j)
1012 elmt->GetCoords(coords[0], coords[1], coords[2]);
1020 elmt->GetTraceNormal(boundary);
1024 BndExp[n]->GetExp(i);
1027 int nbc = bc->GetTotPoints();
1031 elmt->GetTracePhysVals(boundary, bc, pElmt, Pb);
1033 for (
int j = 0; j < expdim * expdim; ++j)
1036 elmt->GetTracePhysVals(boundary, bc, grad[j],
1040 for (
int j = 0; j < 3; ++j)
1043 elmt->GetTracePhysVals(boundary, bc, coords[j],
1054 for (
int j = 0; j < expdim; j++)
1063 for (
int j = 0; j < expdim; j++)
1066 for (
int k = 0; k < expdim; k++)
1069 normals[k], 1, fv[j], 1, fv[j], 1);
1071 normals[k], 1, fv[j], 1, fv[j], 1);
1073 if (!fluidEqu->HasConstantDensity())
1096 1, coordsb[1], 1, mp[0], 1);
1099 1, coordsb[1], 1, mv[0], 1);
1113 1, coordsb[2], 1, mp[0], 1);
1115 1, coordsb[2], 1, mv[0], 1);
1123 1, coordsb[0], 1, mp[1], 1);
1125 1, coordsb[0], 1, mv[1], 1);
1133 1, coordsb[1], 1, mp[2], 1);
1135 1, coordsb[1], 1, mv[2], 1);
1139 for (
int j = 0; j < expdim; j++)
1142 BndExp[n]->GetExp(i)->Integral(fp[j]);
1144 BndExp[n]->GetExp(i)->Integral(fv[j]);
1146 for (
int j = 0; j < momdim; ++j)
1149 BndExp[n]->GetExp(i)->Integral(mp[j]);
1151 BndExp[n]->GetExp(i)->Integral(mv[j]);
1157 cnt += BndExp[n]->GetExpSize();
1166 for (
int i = 0; i < expdim; i++)
1174 for (
int i = 0; i < momdim; ++i)
1186 for (
size_t i = 0; i <
m_Ft.size(); ++i)
1192 for (
size_t i = 0; i <
m_Mt.size(); ++i)
1194 int j =
m_Mt.size() - 1 - i;
1199 fluidEqu->SetAeroForce(aeroforces);
1202 for (
int plane = 0; plane <
m_nPlanes; plane++)
1206 for (
int i = 0; i < expdim; i++)
1208 for (
int j = 0; j < expdim; j++)
1215 for (
int i = 0; i < expdim; i++)
1224 for (
int i = 0; i < 3; ++i)
1228 for (
int j = 0; j < 3; ++j)
1235 for (
int i = 0; i < 3; ++i)
1244 for (
int plane = 0; plane <
m_nPlanes; plane++)
1246 for (
int i = 0; i < expdim; i++)
1253 for (
int i = 0; i < 3; ++i)
1265 for (
int i = 0; i < expdim; i++)
1271 for (
int i = 0; i < momdim; i++)
1281 for (
int i = 0; i < pFields.size(); ++i)
1283 pFields[i]->SetWaveSpace(
true);
1284 pFields[i]->HomogeneousFwdTrans(pFields[i]->GetTotPoints(),
1285 pFields[i]->GetPhys(),
1286 pFields[i]->UpdatePhys());
1299 int cnt, elmtid, offset, boundary;
1301 int physTot = pFields[0]->GetNpoints();
1302 int expdim = pFields[0]->GetGraph()->GetMeshDimension();
1303 int momdim = (expdim == 2) ? 1 : 3;
1346 m_session->DefinesSolverInfo(
"HomoStrip")
1347 ? vComm->GetColumnComm()->GetColumnComm()
1348 : vComm->GetColumnComm();
1362 for (
int i = 0; i < expdim; i++)
1373 for (
int i = 0; i < momdim; ++i)
1390 if (
m_session->DefinesParameter(
"Kinvis"))
1392 rho = (
m_session->DefinesParameter(
"rho"))
1395 mu = rho *
m_session->GetParameter(
"Kinvis");
1405 for (
int i = 0; i < pFields.size(); ++i)
1409 pFields[i]->SetWaveSpace(
false);
1411 pFields[i]->BwdTrans(pFields[i]->GetCoeffs(), pFields[i]->UpdatePhys());
1412 pFields[i]->SetPhysState(
true);
1420 BndConds = pFields[0]->GetPlane(0)->GetBndConditions();
1421 BndExp = pFields[0]->GetPlane(0)->GetBndCondExpansions();
1425 BndConds = pFields[0]->GetBndConditions();
1426 BndExp = pFields[0]->GetBndCondExpansions();
1441 Exp3DH1 = std::dynamic_pointer_cast<
1443 for (
int i = 0; i < nVel * nVel; i++)
1456 for (
int i = 0; i < 3; ++i)
1468 Exp2D = std::dynamic_pointer_cast<MultiRegions::ExpList>(
1470 for (
int i = 0; i < nVel * nVel; i++)
1484 for (
int i = 0; i < 3; ++i)
1499 std::dynamic_pointer_cast<MultiRegions::ExpList>(pFields[0]);
1500 for (
int i = 0; i < nVel * nVel; i++)
1512 for (
int i = 0; i < 3; ++i)
1525 "Expansion dimension not supported by FilterAeroForces");
1534 m_mapping->GetCartesianCoordinates(coords[0]->UpdatePhys(),
1535 coords[1]->UpdatePhys(),
1536 coords[2]->UpdatePhys());
1539 for (
int i = 0; i < nVel * nVel; ++i)
1541 Vmath::Vmul(physTot, tmp[i], 1, pFields[nVel]->GetPhys(), 1,
1542 P[i]->UpdatePhys(), 1);
1547 for (
int i = 0; i < nVel; ++i)
1550 Vmath::Vcopy(physTot, pFields[i]->GetPhys(), 1, wk[i], 1);
1552 m_mapping->ApplyChristoffelContravar(wk, tmp);
1553 for (
int i = 0; i < nVel; ++i)
1555 for (
int k = 0; k < nVel; ++k)
1558 grad[i * nVel + k]->UpdatePhys());
1561 grad[i * nVel + k]->UpdatePhys(), 1,
1562 grad[i * nVel + k]->UpdatePhys(), 1);
1566 for (
int i = 0; i < nVel; ++i)
1568 for (
int k = 0; k < nVel; ++k)
1570 Vmath::Vcopy(physTot, grad[i * nVel + k]->GetPhys(), 1, wk[k], 1);
1573 for (
int j = 0; j < nVel; ++j)
1575 Vmath::Vcopy(physTot, wk[j], 1, grad[i * nVel + j]->UpdatePhys(),
1581 m_mapping->GetJacobian(Jac->UpdatePhys());
1584 for (
int i = 0; i < nVel; ++i)
1587 for (
int k = 0; k < nVel; ++k)
1594 m_mapping->ContravarToCartesian(wk, wk);
1596 for (
int k = 0; k < nVel; k++)
1598 Vmath::Vcopy(physTot, wk[k], 1, C[k * nVel + i]->UpdatePhys(), 1);
1608 for (
int plane = 0; plane <
m_nPlanes; ++plane)
1617 for (
int n = 0; n < nVel * nVel; ++n)
1620 gradPlane[n] = grad[n]->GetPlane(
m_planesID[plane]);
1621 CPlane[n] = C[n]->GetPlane(
m_planesID[plane]);
1623 for (
int n = 0; n < 3; ++n)
1625 coordsPlane[n] = coords[n]->GetPlane(
m_planesID[plane]);
1631 for (
int n = 0; n < nVel * nVel; ++n)
1634 gradPlane[n] = grad[n];
1637 for (
int n = 0; n < 3; ++n)
1639 coordsPlane[n] = coords[n];
1645 for (
int n = cnt = 0; n < BndConds.size(); n++)
1649 for (
int i = 0; i < BndExp[n]->GetExpSize(); ++i, cnt++)
1652 elmt = PPlane[0]->GetExp(elmtid);
1653 offset = PPlane[0]->GetPhys_Offset(elmtid);
1656 for (
int j = 0; j < nVel * nVel; j++)
1658 PElmt[j] = PPlane[j]->GetPhys() + offset;
1659 gradElmt[j] = gradPlane[j]->GetPhys() + offset;
1660 CElmt[j] = CPlane[j]->GetPhys() + offset;
1662 for (
int j = 0; j < 3; ++j)
1664 coordsElmt[j] = coordsPlane[j]->GetPhys() + offset;
1666 JacElmt = JacPlane->GetPhys() + offset;
1675 normals = elmt->GetTraceNormal(boundary);
1679 BndExp[n]->GetExp(i);
1682 int nbc = bc->GetTotPoints();
1685 for (
int j = 0; j < nVel * nVel; ++j)
1688 elmt->GetTracePhysVals(boundary, bc, gradElmt[j],
1692 elmt->GetTracePhysVals(boundary, bc, PElmt[j],
1696 elmt->GetTracePhysVals(boundary, bc, CElmt[j],
1699 for (
int j = 0; j < 3; ++j)
1702 elmt->GetTracePhysVals(boundary, bc, coordsElmt[j],
1707 coordsBnd[j], 1, coordsBnd[j], 1);
1710 elmt->GetTracePhysVals(boundary, bc, JacElmt, JacBnd);
1715 for (
int j = 0; j < nVel; j++)
1719 for (
int k = 0; k < expdim; k++)
1722 normals[k], 1, fp[j], 1, fp[j], 1);
1729 for (
int j = 0; j < nVel; j++)
1732 for (
int k = 0; k < expdim; k++)
1735 normals[k], 1, fv[j], 1, fv[j], 1);
1737 normals[k], 1, fv[j], 1, fv[j], 1);
1743 for (
int k = 0; k < nVel; k++)
1750 for (
int k = 0; k < nVel; k++)
1753 for (
int j = 0; j < nVel; j++)
1756 1, wk[k], 1, wk[k], 1);
1759 for (
int k = 0; k < nVel; k++)
1764 for (
int k = 0; k < nVel; k++)
1767 for (
int j = 0; j < nVel; j++)
1770 1, wk[k], 1, wk[k], 1);
1773 for (
int k = 0; k < nVel; k++)
1791 fptmp, 1, coordsBnd[1], 1, mp[0], 1);
1793 fvtmp, 1, coordsBnd[1], 1, mv[0], 1);
1807 fptmp, 1, coordsBnd[2], 1, mp[0], 1);
1809 fvtmp, 1, coordsBnd[2], 1, mv[0], 1);
1817 fptmp, 1, coordsBnd[0], 1, mp[1], 1);
1819 fvtmp, 1, coordsBnd[0], 1, mv[1], 1);
1827 fptmp, 1, coordsBnd[1], 1, mp[2], 1);
1829 fvtmp, 1, coordsBnd[1], 1, mv[2], 1);
1833 for (
int j = 0; j < expdim; j++)
1836 BndExp[n]->GetExp(i)->Integral(fp[j]);
1838 BndExp[n]->GetExp(i)->Integral(fv[j]);
1842 for (
int j = 0; j < momdim; ++j)
1845 BndExp[n]->GetExp(i)->Integral(mp[j]);
1847 BndExp[n]->GetExp(i)->Integral(mv[j]);
1853 cnt += BndExp[n]->GetExpSize();
1862 for (
int i = 0; i < expdim; ++i)
1870 for (
int i = 0; i < momdim; ++i)
1880 for (
int plane = 0; plane <
m_nPlanes; plane++)
1884 for (
int i = 0; i < expdim; i++)
1886 for (
int j = 0; j < expdim; j++)
1893 for (
int i = 0; i < expdim; i++)
1903 for (
int i = 0; i < 3; ++i)
1907 for (
int j = 0; j < 3; ++j)
1915 for (
int i = 0; i < 3; ++i)
1924 for (
int plane = 0; plane <
m_nPlanes; plane++)
1926 for (
int i = 0; i < expdim; i++)
1933 for (
int i = 0; i < 3; ++i)
1945 for (
int i = 0; i < expdim; i++)
1951 for (
int i = 0; i < momdim; i++)
1961 for (
int i = 0; i < pFields.size(); ++i)
1963 pFields[i]->SetWaveSpace(
true);
1964 pFields[i]->HomogeneousFwdTrans(pFields[i]->GetTotPoints(),
1965 pFields[i]->GetPhys(),
1966 pFields[i]->UpdatePhys());
#define ASSERTL0(condition, msg)
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.
NekDouble Evaluate() const
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expa...
static bool GenerateSeqVector(const std::string &str, std::vector< unsigned int > &out)
Takes a comma-separated compressed string and converts it to entries in a vector.
void v_Finalise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time) override
Array< OneD, int > m_planesID
Array< OneD, Array< OneD, NekDouble > > m_Fvplane
void v_Update(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time) override
void CalculateForces(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
Array< OneD, int > m_BCtoElmtID
int m_nPlanes
number of planes for homogeneous1D expansion
SOLVER_UTILS_EXPORT ~FilterAeroForces() override
Array< OneD, NekDouble > m_pivotPoint
Array< OneD, NekDouble > m_Mp
bool m_outputAllPlanes
if using a homogeneous1D expansion, determine if should output all planes or just the average
void CalculateForcesMapping(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
static FilterSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const std::weak_ptr< EquationSystem > &pEquation, const std::map< std::string, std::string > &pParams)
Creates an instance of this class.
std::vector< unsigned int > m_boundaryRegionsIdList
ID's of boundary regions where we want the forces.
Array< OneD, Array< OneD, NekDouble > > m_directions
Array< OneD, NekDouble > m_Fv
std::ofstream m_outputStream
std::vector< bool > m_boundaryRegionIsInList
Determines if a given Boundary Region is in m_boundaryRegionsIdList.
Array< OneD, NekDouble > m_Mv
std::string m_BoundaryString
static std::string className
Name of the class.
SOLVER_UTILS_EXPORT FilterAeroForces(const LibUtilities::SessionReaderSharedPtr &pSession, const std::weak_ptr< EquationSystem > &pEquation, const std::map< std::string, std::string > &pParams)
Array< OneD, int > m_BCtoTraceID
Array< OneD, Array< OneD, NekDouble > > m_Mpplane
bool v_IsTimeDependent() override
SOLVER_UTILS_EXPORT void GetMoments(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, Array< OneD, NekDouble > &moments, const NekDouble &time)
Array< OneD, Array< OneD, NekDouble > > m_directions0
GlobalMapping::MappingSharedPtr m_mapping
Array< OneD, Array< OneD, NekDouble > > m_Ftplane
unsigned int m_outputFrequency
void v_Initialise(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time) override
Array< OneD, NekDouble > m_Ft
SOLVER_UTILS_EXPORT void GetForces(const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, Array< OneD, NekDouble > &Aeroforces, const NekDouble &time)
Array< OneD, NekDouble > m_Fp
Array< OneD, Array< OneD, NekDouble > > m_Mtplane
Array< OneD, NekDouble > m_Mt
Array< OneD, Array< OneD, NekDouble > > m_Fpplane
Array< OneD, Array< OneD, NekDouble > > m_Mvplane
LibUtilities::SessionReaderSharedPtr m_session
const std::weak_ptr< EquationSystem > m_equ
std::map< std::string, std::string > ParamMap
const BoundaryRegionCollection & GetBoundaryRegions(void) const
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
@ P
Monomial polynomials .
std::shared_ptr< Expansion > ExpansionSharedPtr
const std::vector< NekDouble > velocity
MultiRegions::Direction const DirCartesianMap[]
std::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
std::shared_ptr< ExpList3DHomogeneous1D > ExpList3DHomogeneous1DSharedPtr
Shared pointer to an ExpList3DHomogeneous1D object.
FilterFactory & GetFilterFactory()
std::map< int, BoundaryRegionShPtr > BoundaryRegionCollection
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
static Array< OneD, NekDouble > NullNekDouble1DArray
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.
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
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
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 Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
void Vvtvvtp(int n, const T *v, int incv, const T *w, int incw, const T *x, int incx, const T *y, int incy, T *z, int incz)
vvtvvtp (vector times vector plus vector times vector):
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
scalarT< T > sqrt(scalarT< T > in)