70{
71 boost::ignore_unused(time);
72
73 int i, j;
74 int nTracePts =
m_fields[0]->GetTrace()->GetNpoints();
75 int nVariables = physarray.size();
77
78 const Array<OneD, const int> &traceBndMap =
m_fields[0]->GetTraceBndMap();
79
80
81
82 Array<OneD, NekDouble> Vn(nTracePts, 0.0);
83 Array<OneD, NekDouble> Vel(nTracePts, 0.0);
84 for (i = 0; i < nDimensions; ++i)
85 {
86 Vmath::Vdiv(nTracePts, Fwd[i + 1], 1, Fwd[0], 1, Vel, 1);
88 }
89
90
91
92 Array<OneD, NekDouble> absVel(nTracePts, 0.0);
93 m_varConv->GetAbsoluteVelocity(Fwd, absVel);
94
95
96 Array<OneD, NekDouble> soundSpeed(nTracePts);
97 m_varConv->GetSoundSpeed(Fwd, soundSpeed);
98
99
100 Array<OneD, NekDouble> Mach(nTracePts, 0.0);
101 Vmath::Vdiv(nTracePts, Vn, 1, soundSpeed, 1, Mach, 1);
103
104
105 int e, id1, id2, npts, pnt;
107
108
109 for (e = 0;
111 {
114 ->GetExp(e)
115 ->GetTotPoints();
116 id1 =
118 id2 =
120
121
123 Array<OneD, NekDouble> rho(npts, Fwd[0] + id2);
124 Array<OneD, NekDouble> Ei(npts);
126
127
128 for (i = 0; i < npts; i++)
129 {
130 pnt = id2 + i;
131
132
133 if (Mach[pnt] < 0.99)
134 {
135
137 for (j = 1; j < nVariables - 1; ++j)
138 {
139 Ek += 0.5 * (Fwd[j][pnt] * Fwd[j][pnt]) / Fwd[0][pnt];
140 }
141
142 rhoeb = Fwd[0][pnt] * Ei[i] + Ek;
143
144
145 for (j = 0; j < nVariables - 1; ++j)
146 {
149 ->UpdatePhys())[id1 + i] = Fwd[j][pnt];
150 }
151
154 ->UpdatePhys())[id1 + i] = rhoeb;
155 }
156
157 else
158 {
159 for (j = 0; j < nVariables; ++j)
160 {
161
164 ->UpdatePhys())[id1 + i] = Fwd[j][pnt];
165 }
166 }
167 }
168 }
169}
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.