Skip to content

Commit d01db7d

Browse files
bigfootedbigfooted
authored
Merge branch 'develop' into fix_symmetry_logic_with_deform_marker
2 parents 9ce7cec + f33b248 commit d01db7d

11 files changed

Lines changed: 595 additions & 639 deletions

File tree

‎SU2_CFD/src/numerics/NEMO/CNEMONumerics.cpp‎

Lines changed: 68 additions & 77 deletions
Original file line numberDiff line numberDiff line change
@@ -84,36 +84,35 @@ void CNEMONumerics::GetInviscidProjFlux(const su2double *val_U,
8484
const su2double *val_V,
8585
const su2double *val_normal,
8686
su2double *val_Proj_Flux) {
87-
unsigned short iSpecies, iVar;
88-
su2double rho, u, v, w, rhoEve, P, H;
87+
8988
const su2double *rhos;
9089

9190
/*--- Initialize vectors ---*/
92-
for (iVar = 0; iVar < nVar; iVar++)
91+
for (auto iVar = 0ul; iVar < nVar; iVar++)
9392
val_Proj_Flux[iVar] = 0.0;
9493

9594
/*--- Rename for convienience ---*/
96-
rho = val_V[RHO_INDEX];
97-
u = val_V[VEL_INDEX];
98-
v = val_V[VEL_INDEX+1];
99-
w = val_V[VEL_INDEX+2];
100-
P = val_V[P_INDEX];
101-
H = val_V[H_INDEX];
102-
rhoEve = val_U[nSpecies+nDim+1];
95+
const su2double rho = val_V[RHO_INDEX];
96+
const su2double u = val_V[VEL_INDEX];
97+
const su2double v = val_V[VEL_INDEX+1];
98+
const su2double w = val_V[VEL_INDEX+2];
99+
const su2double P = val_V[P_INDEX];
100+
const su2double H = val_V[H_INDEX];
101+
const su2double rhoEve = val_U[nSpecies+nDim+1];
103102
rhos = &val_V[RHOS_INDEX];
104103

105104
if (nDim == 2) {
106105

107106
/*--- iDim = 0 (x-direction) ---*/
108-
for (iSpecies = 0; iSpecies < nSpecies; iSpecies++)
107+
for (auto iSpecies = 0ul; iSpecies < nSpecies; iSpecies++)
109108
val_Proj_Flux[iSpecies] = (rhos[iSpecies]*u) * val_normal[0];
110109
val_Proj_Flux[nSpecies] = (rho*u*u + P) * val_normal[0];
111110
val_Proj_Flux[nSpecies+1] = (rho*u*v) * val_normal[0];
112111
val_Proj_Flux[nSpecies+2] = (rho*u*H) * val_normal[0];
113112
val_Proj_Flux[nSpecies+3] = (rhoEve*u) * val_normal[0];
114113

115114
/*---- iDim = 1 (y-direction) ---*/
116-
for (iSpecies = 0; iSpecies < nSpecies; iSpecies++)
115+
for (auto iSpecies = 0ul; iSpecies < nSpecies; iSpecies++)
117116
val_Proj_Flux[iSpecies] += (rhos[iSpecies]*v) * val_normal[1];
118117
val_Proj_Flux[nSpecies] += (rho*v*u) * val_normal[1];
119118
val_Proj_Flux[nSpecies+1] += (rho*v*v + P) * val_normal[1];
@@ -123,7 +122,7 @@ void CNEMONumerics::GetInviscidProjFlux(const su2double *val_U,
123122
else {
124123

125124
/*--- iDim = 0 (x-direction) ---*/
126-
for (iSpecies = 0; iSpecies < nSpecies; iSpecies++)
125+
for (auto iSpecies = 0ul; iSpecies < nSpecies; iSpecies++)
127126
val_Proj_Flux[iSpecies] = (rhos[iSpecies]*u) * val_normal[0];
128127
val_Proj_Flux[nSpecies] = (rho*u*u + P) * val_normal[0];
129128
val_Proj_Flux[nSpecies+1] = (rho*u*v) * val_normal[0];
@@ -132,7 +131,7 @@ void CNEMONumerics::GetInviscidProjFlux(const su2double *val_U,
132131
val_Proj_Flux[nSpecies+4] = (rhoEve*u) * val_normal[0];
133132

134133
/*--- iDim = 0 (y-direction) ---*/
135-
for (iSpecies = 0; iSpecies < nSpecies; iSpecies++)
134+
for (auto iSpecies = 0ul; iSpecies < nSpecies; iSpecies++)
136135
val_Proj_Flux[iSpecies] += (rhos[iSpecies]*v) * val_normal[1];
137136
val_Proj_Flux[nSpecies] += (rho*v*u) * val_normal[1];
138137
val_Proj_Flux[nSpecies+1] += (rho*v*v + P) * val_normal[1];
@@ -141,7 +140,7 @@ void CNEMONumerics::GetInviscidProjFlux(const su2double *val_U,
141140
val_Proj_Flux[nSpecies+4] += (rhoEve*v) * val_normal[1];
142141

143142
/*--- iDim = 0 (z-direction) ---*/
144-
for (iSpecies = 0; iSpecies < nSpecies; iSpecies++)
143+
for (auto iSpecies = 0ul; iSpecies < nSpecies; iSpecies++)
145144
val_Proj_Flux[iSpecies] += (rhos[iSpecies]*w) * val_normal[2];
146145
val_Proj_Flux[nSpecies] += (rho*w*u) * val_normal[2];
147146
val_Proj_Flux[nSpecies+1] += (rho*w*v) * val_normal[2];
@@ -157,33 +156,32 @@ void CNEMONumerics::GetInviscidProjJac(const su2double *val_U, const su2doubl
157156
su2double **val_Proj_Jac_Tensor) {
158157

159158
const su2double *rhos;
160-
161159
rhos = &val_V[RHOS_INDEX];
162160

163161
/*--- Initialize the Jacobian tensor ---*/
164-
for (unsigned short iVar = 0; iVar < nVar; iVar++)
165-
for (unsigned short jVar = 0; jVar < nVar; jVar++)
162+
for (auto iVar = 0ul; iVar < nVar; iVar++)
163+
for (auto jVar = 0ul; jVar < nVar; jVar++)
166164
val_Proj_Jac_Tensor[iVar][jVar] = 0.0;
167165

168166
/*--- Rename for convenience ---*/
169-
su2double rho = val_V[RHO_INDEX];
170-
su2double H = val_V[H_INDEX];
171-
su2double rhoEve = val_U[nSpecies+nDim+1];
167+
const su2double rho = val_V[RHO_INDEX];
168+
const su2double H = val_V[H_INDEX];
169+
const su2double rhoEve = val_U[nSpecies+nDim+1];
172170

173171
su2double u[MAXNDIM];
174-
for (unsigned short iDim = 0; iDim < nDim; iDim++)
172+
for (auto iDim = 0ul; iDim < nDim; iDim++)
175173
u[iDim] = val_V[VEL_INDEX+iDim];
176174

177175
/*--- Calculate projected velocity ---*/
178176
su2double proj_vel = GeometryToolbox::DotProduct(nDim, u, val_normal);
179177

180178
/*--- Species density rows ---*/
181-
for (unsigned short iSpecies = 0; iSpecies < nSpecies; iSpecies++) {
182-
for (unsigned short jSpecies = 0; jSpecies < nSpecies; jSpecies++) {
179+
for (auto iSpecies = 0ul; iSpecies < nSpecies; iSpecies++) {
180+
for (auto jSpecies = 0ul; jSpecies < nSpecies; jSpecies++) {
183181
val_Proj_Jac_Tensor[iSpecies][jSpecies] += -(rhos[iSpecies]/rho) * proj_vel;
184182
}
185183
val_Proj_Jac_Tensor[iSpecies][iSpecies] += proj_vel;
186-
for (unsigned short iDim = 0; iDim < nDim; iDim++) {
184+
for (auto iDim = 0ul; iDim < nDim; iDim++) {
187185
val_Proj_Jac_Tensor[iSpecies][nSpecies+iDim] += (rhos[iSpecies]/rho) * val_normal[iDim];
188186
val_Proj_Jac_Tensor[nSpecies+iDim][iSpecies] += val_dPdU[iSpecies]*val_normal[iDim] - proj_vel*u[iDim];
189187
}
@@ -192,8 +190,8 @@ void CNEMONumerics::GetInviscidProjJac(const su2double *val_U, const su2doubl
192190
}
193191

194192
/*--- Momentum rows ---*/
195-
for (unsigned short iDim = 0; iDim < nDim; iDim++) {
196-
for (unsigned short jDim = 0; jDim < nDim; jDim++) {
193+
for (auto iDim = 0ul; iDim < nDim; iDim++) {
194+
for (auto jDim = 0ul; jDim < nDim; jDim++) {
197195
val_Proj_Jac_Tensor[nSpecies+iDim][nSpecies+jDim] += val_dPdU[nSpecies+jDim]*val_normal[iDim] + u[iDim]*val_normal[jDim];
198196
}
199197
val_Proj_Jac_Tensor[nSpecies+iDim][nSpecies+iDim] += proj_vel;
@@ -202,18 +200,18 @@ void CNEMONumerics::GetInviscidProjJac(const su2double *val_U, const su2doubl
202200
}
203201

204202
/*--- Total energy row ---*/
205-
for (unsigned short iDim = 0; iDim < nDim; iDim++)
203+
for (auto iDim = 0ul; iDim < nDim; iDim++)
206204
val_Proj_Jac_Tensor[nSpecies+nDim][nSpecies+iDim] += val_dPdU[nSpecies+iDim]*proj_vel + H*val_normal[iDim];
207205
val_Proj_Jac_Tensor[nSpecies+nDim][nSpecies+nDim] += (1+val_dPdU[nSpecies+nDim])*proj_vel;
208206
val_Proj_Jac_Tensor[nSpecies+nDim][nSpecies+nDim+1] += val_dPdU[nSpecies+nDim+1] *proj_vel;
209207

210208
/*--- Vib.-el. energy row ---*/
211-
for (unsigned short iDim = 0; iDim < nDim; iDim++)
209+
for (auto iDim = 0ul; iDim < nDim; iDim++)
212210
val_Proj_Jac_Tensor[nSpecies+nDim+1][nSpecies+iDim] = rhoEve/rho*val_normal[iDim];
213211
val_Proj_Jac_Tensor[nSpecies+nDim+1][nSpecies+nDim+1] = proj_vel;
214212

215-
for (unsigned short iVar = 0; iVar < nVar; iVar++)
216-
for (unsigned short jVar = 0; jVar < nVar; jVar++)
213+
for (auto iVar = 0ul; iVar < nVar; iVar++)
214+
for (auto jVar = 0ul; jVar < nVar; jVar++)
217215
val_Proj_Jac_Tensor[iVar][jVar] = val_scale * val_Proj_Jac_Tensor[iVar][jVar];
218216
}
219217

@@ -240,17 +238,17 @@ void CNEMONumerics::GetViscousProjFlux(const su2double *val_primvar,
240238
su2activematrix Flux_Tensor(nVar,nDim);
241239

242240
/*--- Initialize ---*/
243-
for (auto iVar = 0; iVar < nVar; iVar++) {
241+
for (auto iVar = 0ul; iVar < nVar; iVar++) {
244242
Proj_Flux_Tensor[iVar] = 0.0;
245-
for (auto iDim = 0; iDim < nDim; iDim++)
243+
for (auto iDim = 0ul; iDim < nDim; iDim++)
246244
Flux_Tensor[iVar][iDim] = 0.0;
247245
}
248246

249247
/*--- Rename variables for convenience ---*/
250248
const auto& Ds = val_diffusioncoeff;
251249
const su2double mu = val_lam_viscosity+val_eddy_viscosity;
252-
su2double ktr = val_therm_conductivity;
253-
su2double kve = val_therm_conductivity_ve;
250+
const su2double ktr = val_therm_conductivity;
251+
const su2double kve = val_therm_conductivity_ve;
254252
const su2double rho = val_primvar[RHO_INDEX];
255253
const su2double T = val_primvar[T_INDEX];
256254
const su2double Tve = val_primvar[TVE_INDEX];
@@ -260,8 +258,8 @@ void CNEMONumerics::GetViscousProjFlux(const su2double *val_primvar,
260258

261259
/*--- Pre-compute mixture quantities ---*/ //TODO
262260
su2double Vector[MAXNDIM] = {0.0};
263-
for (auto iDim = 0; iDim < nDim; iDim++) {
264-
for (auto iSpecies = 0; iSpecies < nHeavy; iSpecies++) {
261+
for (auto iDim = 0ul; iDim < nDim; iDim++) {
262+
for (auto iSpecies = 0ul; iSpecies < nHeavy; iSpecies++) {
265263
Vector[iDim] += rho*Ds[iSpecies]*GV[RHOS_INDEX+iSpecies][iDim];
266264
}
267265
}
@@ -270,23 +268,23 @@ void CNEMONumerics::GetViscousProjFlux(const su2double *val_primvar,
270268
ComputeStressTensor(nDim,tau,val_gradprimvar+VEL_INDEX, mu);
271269

272270
/*--- Populate entries in the viscous flux vector ---*/
273-
for (auto iDim = 0; iDim < nDim; iDim++) {
271+
for (auto iDim = 0ul; iDim < nDim; iDim++) {
274272

275273
/*--- Species diffusion velocity ---*/
276-
for (auto iSpecies = 0; iSpecies < nHeavy; iSpecies++) {
274+
for (auto iSpecies = 0ul; iSpecies < nHeavy; iSpecies++) {
277275
Flux_Tensor[iSpecies][iDim] = rho*Ds[iSpecies]*GV[RHOS_INDEX+iSpecies][iDim]
278276
- V[RHOS_INDEX+iSpecies]*Vector[iDim];
279277
}
280278

281279
/*--- Shear-stress/momentum related terms ---*/
282280
Flux_Tensor[nSpecies+nDim][iDim] = 0.0;
283-
for (auto jDim = 0; jDim < nDim; jDim++) {
281+
for (auto jDim = 0ul; jDim < nDim; jDim++) {
284282
Flux_Tensor[nSpecies+jDim][iDim] = tau[iDim][jDim];
285283
Flux_Tensor[nSpecies+nDim][iDim] += tau[iDim][jDim]*val_primvar[VEL_INDEX+jDim];
286284
}
287285

288286
/*--- Diffusion terms ---*/
289-
for (auto iSpecies = 0; iSpecies < nHeavy; iSpecies++) {
287+
for (auto iSpecies = 0ul; iSpecies < nHeavy; iSpecies++) {
290288
Flux_Tensor[nSpecies+nDim][iDim] += Flux_Tensor[iSpecies][iDim] * hs[iSpecies];
291289
Flux_Tensor[nSpecies+nDim+1][iDim] += Flux_Tensor[iSpecies][iDim] * val_eve[iSpecies];
292290
}
@@ -296,8 +294,8 @@ void CNEMONumerics::GetViscousProjFlux(const su2double *val_primvar,
296294
Flux_Tensor[nSpecies+nDim+1][iDim] += kve*GV[TVE_INDEX][iDim];
297295
}
298296

299-
for (auto iVar = 0; iVar < nVar; iVar++) {
300-
for (auto iDim = 0; iDim < nDim; iDim++) {
297+
for (auto iVar = 0ul; iVar < nVar; iVar++) {
298+
for (auto iDim = 0ul; iDim < nDim; iDim++) {
301299
Proj_Flux_Tensor[iVar] += Flux_Tensor[iVar][iDim]*val_normal[iDim];
302300
}
303301
}
@@ -356,32 +354,28 @@ void CNEMONumerics::GetPMatrix(const su2double *U, const su2double *V, const su2
356354
const su2double *val_normal, const su2double *l, const su2double *m,
357355
su2double **val_p_tensor) const {
358356

359-
// P matrix is equivalent to the L matrix in Gnoffo
360-
unsigned short iSpecies, iDim, iVar, jVar;
361-
su2double sqvel, rho, a, a2, eve;
362-
su2double vU, vV, vW;
363-
364357
/*--- Initialize the P matrix to zero ---*/
365-
for (iVar = 0; iVar < nVar; iVar++)
366-
for (jVar = 0; jVar < nVar; jVar++)
358+
for (auto iVar = 0ul; iVar < nVar; iVar++)
359+
for (auto jVar = 0ul; jVar < nVar; jVar++)
367360
val_p_tensor[iVar][jVar] = 0.0;
368361

369362
/*--- Pre-compute useful quantities ---*/
370-
sqvel = 0.0;
371-
rho = V[RHO_INDEX];
372-
eve = U[nSpecies+nDim+1]/rho;
373-
vU = 0.0; vV = 0.0; vW = 0.0;
374-
for (iDim = 0; iDim < nDim; iDim++) {
363+
su2double sqvel = 0.0;
364+
const su2double rho = V[RHO_INDEX];
365+
const su2double eve = U[nSpecies+nDim+1]/rho;
366+
su2double vU, vV, vW;
367+
vU = vV = vW = 0.0;
368+
for (auto iDim = 0ul; iDim < nDim; iDim++) {
375369
vU += V[VEL_INDEX+iDim] * val_normal[iDim];
376370
vV += V[VEL_INDEX+iDim] * l[iDim];
377371
vW += V[VEL_INDEX+iDim] * m[iDim];
378372
sqvel += V[VEL_INDEX+iDim] * V[VEL_INDEX+iDim];
379373
}
380-
a = V[A_INDEX];
381-
a2 = V[A_INDEX]*V[A_INDEX];
374+
const su2double a = V[A_INDEX];
375+
const su2double a2 = V[A_INDEX]*V[A_INDEX];
382376

383377
if(nDim == 2) {
384-
for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) {
378+
for (auto iSpecies = 0ul; iSpecies < nSpecies; iSpecies++) {
385379
val_p_tensor[iSpecies][iSpecies] += 1.0/a2;
386380
val_p_tensor[iSpecies][nSpecies] += 0.0;
387381
val_p_tensor[iSpecies][nSpecies+1] += V[RHOS_INDEX+iSpecies] / (2.0*rho*a2);
@@ -395,7 +389,7 @@ void CNEMONumerics::GetPMatrix(const su2double *U, const su2double *V, const su2
395389
val_p_tensor[nSpecies+3][iSpecies] += 0.0;
396390
}
397391

398-
for (iDim = 0; iDim < nDim; iDim++){
392+
for (auto iDim = 0ul; iDim < nDim; iDim++){
399393
val_p_tensor[nSpecies+iDim][nSpecies] += l[iDim];
400394
val_p_tensor[nSpecies+iDim][nSpecies+1] += (V[VEL_INDEX+iDim]+a*val_normal[iDim]) / (2.0*a2);
401395
val_p_tensor[nSpecies+iDim][nSpecies+2] += (V[VEL_INDEX+iDim]-a*val_normal[iDim]) / (2.0*a2);
@@ -414,7 +408,7 @@ void CNEMONumerics::GetPMatrix(const su2double *U, const su2double *V, const su2
414408

415409
} else {
416410

417-
for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) {
411+
for (auto iSpecies = 0ul; iSpecies < nSpecies; iSpecies++) {
418412
val_p_tensor[iSpecies][iSpecies] = 1.0/a2;
419413
val_p_tensor[iSpecies][nSpecies] = 0.0;
420414
val_p_tensor[iSpecies][nSpecies+1] = 0.0;
@@ -430,7 +424,7 @@ void CNEMONumerics::GetPMatrix(const su2double *U, const su2double *V, const su2
430424
val_p_tensor[nSpecies+4][iSpecies] = 0.0;
431425
}
432426

433-
for (iDim = 0; iDim < nDim; iDim++){
427+
for (auto iDim = 0ul; iDim < nDim; iDim++){
434428
val_p_tensor[nSpecies+iDim][nSpecies] = l[iDim];
435429
val_p_tensor[nSpecies+iDim][nSpecies+1] = m[iDim];
436430
val_p_tensor[nSpecies+iDim][nSpecies+2] = (V[VEL_INDEX+iDim]+a*val_normal[iDim]) / (2.0*a2);
@@ -456,30 +450,27 @@ void CNEMONumerics::GetPMatrix_inv(const su2double *U, const su2double *V, const
456450
const su2double *val_normal, const su2double *l, const su2double *m,
457451
su2double **val_invp_tensor) const {
458452

459-
unsigned short iSpecies, jSpecies, iDim, iVar, jVar;
460-
su2double rho, a, a2, eve;
461-
su2double vU, vV, vW;
462-
463-
for (iVar = 0; iVar < nVar; iVar++)
464-
for (jVar = 0; jVar < nVar; jVar++)
453+
for (auto iVar = 0ul; iVar < nVar; iVar++)
454+
for (auto jVar = 0ul; jVar < nVar; jVar++)
465455
val_invp_tensor[iVar][jVar] = 0.0;
466456

467457
/*--- Pre-compute useful quantities ---*/
468-
rho = V[RHO_INDEX];
469-
eve = U[nSpecies+nDim+1]/rho;
470-
vU = 0.0; vV = 0.0; vW = 0.0;
471-
for (iDim = 0; iDim < nDim; iDim++) {
458+
const su2double rho = V[RHO_INDEX];
459+
const su2double eve = U[nSpecies+nDim+1]/rho;
460+
su2double vU, vV, vW;
461+
vU = vV = vW = 0.0;
462+
for (auto iDim = 0ul; iDim < nDim; iDim++) {
472463
vU += V[VEL_INDEX+iDim] * val_normal[iDim];
473464
vV += V[VEL_INDEX+iDim] * l[iDim];
474465
vW += V[VEL_INDEX+iDim] * m[iDim];
475466
}
476-
a = V[A_INDEX];
477-
a2 = V[A_INDEX]*V[A_INDEX];
467+
const su2double a = V[A_INDEX];
468+
const su2double a2 = V[A_INDEX]*V[A_INDEX];
478469

479470
if (nDim == 2) {
480471

481-
for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) {
482-
for (jSpecies = 0; jSpecies < nSpecies; jSpecies++) {
472+
for (auto iSpecies = 0ul; iSpecies < nSpecies; iSpecies++) {
473+
for (auto jSpecies = 0ul; jSpecies < nSpecies; jSpecies++) {
483474
val_invp_tensor[iSpecies][jSpecies] += -(V[RHOS_INDEX+iSpecies]/rho) * val_dPdU[jSpecies];
484475
}
485476
val_invp_tensor[iSpecies][iSpecies] += a2;
@@ -516,8 +507,8 @@ void CNEMONumerics::GetPMatrix_inv(const su2double *U, const su2double *V, const
516507

517508
} else {
518509

519-
for (iSpecies = 0; iSpecies < nSpecies; iSpecies++) {
520-
for (jSpecies = 0; jSpecies < nSpecies; jSpecies++) {
510+
for (auto iSpecies = 0ul; iSpecies < nSpecies; iSpecies++) {
511+
for (auto jSpecies = 0ul; jSpecies < nSpecies; jSpecies++) {
521512
val_invp_tensor[iSpecies][jSpecies] += -(V[RHOS_INDEX+iSpecies]/rho) * val_dPdU[jSpecies];
522513
}
523514
val_invp_tensor[iSpecies][iSpecies] += a2;

0 commit comments

Comments
 (0)