68{
69 int i, j;
70 int nTracePts =
m_fields[0]->GetTrace()->GetNpoints();
71 int nVariables = physarray.size();
73
74 const Array<OneD, const int> &traceBndMap =
m_fields[0]->GetTraceBndMap();
75
76
77
78 Array<OneD, NekDouble> Vn(nTracePts, 0.0);
79 Array<OneD, NekDouble> Vel(nTracePts, 0.0);
80 for (i = 0; i < nDimensions; ++i)
81 {
82 Vmath::Vdiv(nTracePts, Fwd[i + 1], 1, Fwd[0], 1, Vel, 1);
84 }
85
86
87
88 Array<OneD, NekDouble> absVel(nTracePts, 0.0);
89 m_varConv->GetAbsoluteVelocity(Fwd, absVel);
90
91
92 Array<OneD, NekDouble> soundSpeed(nTracePts);
93 m_varConv->GetSoundSpeed(Fwd, soundSpeed);
94
95
96 Array<OneD, NekDouble> Mach(nTracePts, 0.0);
97 Vmath::Vdiv(nTracePts, Vn, 1, soundSpeed, 1, Mach, 1);
99
100
101 int e, id1, id2, npts, pnt;
103
104
105 for (e = 0;
107 {
110 ->GetExp(e)
111 ->GetTotPoints();
112 id1 =
114 id2 =
116
117
119 Array<OneD, NekDouble> rho(npts, Fwd[0] + id2);
120 Array<OneD, NekDouble> Ei(npts);
122
123
124 for (i = 0; i < npts; i++)
125 {
126 pnt = id2 + i;
127
128
129 if (Mach[pnt] < 0.99)
130 {
131
133 for (j = 1; j < nVariables - 1; ++j)
134 {
135 Ek += 0.5 * (Fwd[j][pnt] * Fwd[j][pnt]) / Fwd[0][pnt];
136 }
137
138 rhoeb = Fwd[0][pnt] * Ei[i] + Ek;
139
140
141 for (j = 0; j < nVariables - 1; ++j)
142 {
145 ->UpdatePhys())[id1 + i] = Fwd[j][pnt];
146 }
147
150 ->UpdatePhys())[id1 + i] = rhoeb;
151 }
152
153 else
154 {
155 for (j = 0; j < nVariables; ++j)
156 {
157
160 ->UpdatePhys())[id1 + i] = Fwd[j][pnt];
161 }
162 }
163 }
164 }
165}
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.