73{
74 int i, j;
75 int nTracePts =
m_fields[0]->GetTrace()->GetNpoints();
76 int nVariables = physarray.size();
78
79 const Array<OneD, const int> &traceBndMap =
m_fields[0]->GetTraceBndMap();
80
81
82
83 Array<OneD, NekDouble> Vn(nTracePts, 0.0);
84 Array<OneD, NekDouble> Vel(nTracePts, 0.0);
85 for (i = 0; i < nDimensions; ++i)
86 {
87 Vmath::Vdiv(nTracePts, Fwd[i + 1], 1, Fwd[0], 1, Vel, 1);
89 }
90
91
92
93 Array<OneD, NekDouble> absVel(nTracePts, 0.0);
94 m_varConv->GetAbsoluteVelocity(Fwd, absVel);
95
96
97 Array<OneD, NekDouble>
pressure(nTracePts);
98 Array<OneD, NekDouble> soundSpeed(nTracePts);
99
101 m_varConv->GetSoundSpeed(Fwd, soundSpeed);
102
103
104 Array<OneD, NekDouble> Mach(nTracePts, 0.0);
105 Vmath::Vdiv(nTracePts, Vn, 1, soundSpeed, 1, Mach, 1);
107
108
109 int e, id1, id2, npts, pnt;
111
112
113 for (e = 0;
115 {
118 ->GetExp(e)
119 ->GetTotPoints();
120 id1 =
122 id2 =
124
125
126 Array<OneD, NekDouble> tmpPressure(npts,
pressure + id2);
127 Array<OneD, NekDouble> rho(npts, Fwd[0] + id2);
128 Array<OneD, NekDouble> Ei(npts);
129 m_varConv->GetEFromRhoP(rho, tmpPressure, Ei);
130
131
132 for (i = 0; i < npts; i++)
133 {
134 pnt = id2 + i;
135
136
137 if (Mach[pnt] < 0.99)
138 {
139
140 for (j = 0; j < nVariables - 1; ++j)
141 {
145 }
146
147
149 for (j = 1; j < nVariables - 1; ++j)
150 {
151 Ek += 0.5 *
155 }
156
157 rhoeb = Fwd[0][pnt] * Ei[i] + Ek;
158
161 ->UpdatePhys())[id1 + i] = rhoeb;
162 }
163
164 else
165 {
166 for (j = 0; j < nVariables; ++j)
167 {
168
171 ->UpdatePhys())[id1 + i] = Fwd[j][pnt];
172 }
173 }
174 }
175 }
176}
int m_spacedim
Space dimension.
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Trace normals.
VariableConverterSharedPtr m_varConv
Auxiliary object to convert variables.
void Vabs(int n, const T *x, const int incx, T *y, const int incy)
vabs: y = |x|
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
void Vdiv(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x/y.