Evaluation of the velocity gradient in the cartesian directions Du_x: traceFieldsAdded[10] Du_y: traceFieldsAdded[11] Dv_x: traceFieldsAdded[12] Dv_y: traceFieldsAdded[13]
71{
72 string fname;
73 po::options_description desc("Available options");
74 desc.add_options()("help,h", "Produce this help message.");
75
76 po::options_description hidden("Hidden options");
77 hidden.add_options()("input-file", po::value<vector<string>>(),
78 "input filename");
79
80 po::options_description cmdline_options;
81 cmdline_options.add(desc).add(hidden);
82
83 po::options_description visible("Allowed options");
84 visible.add(desc);
85
86 po::positional_options_description
p;
87 p.add(
"input-file", -1);
88
89 po::variables_map vm;
90
91 try
92 {
93 po::store(po::command_line_parser(argc, argv)
94 .options(cmdline_options)
97 vm);
98 po::notify(vm);
99 }
100 catch (const std::exception &e)
101 {
102 cerr << e.what() << endl;
103 cerr << desc;
104 return 1;
105 }
106
107 std::vector<std::string> filenames;
108 if (vm.count("input-file"))
109 {
110 filenames = vm["input-file"].as<std::vector<std::string>>();
111 }
112
113 if (vm.count("help") || vm.count("input-file") != 1)
114 {
115 cerr << "Description: Extracts a surface from a 2D .fld file created "
116 "by the CompressibleFlowSolver and places it into a .cfs file"
117 << endl;
118 cerr << "Usage: ExtractSurface2DCFS [options] meshFile fieldFile"
119 << endl;
120 cout << desc;
121 return 1;
122 }
123
125 LibUtilities::SessionReader::CreateInstance(argc, argv);
127 SpatialDomains::MeshGraph::Read(vSession);
128
129 fname = vSession->GetSessionName() + ".cfs";
130
131
132 std::string fieldFile;
133 for (auto &file : filenames)
134 {
135 if (file.size() > 4 && (file.substr(file.size() - 4, 4) == ".fld" ||
136 file.substr(file.size() - 4, 4) == ".chk"))
137 {
138 fieldFile = file;
139 break;
140 }
141 }
142
143 int cnt;
144 int id1 = 0;
145 int id2 = 0;
146 int i, j, n, e, b;
148
149 int nBndEdgePts, nBndEdges, nBndRegions;
150
151 std::string m_ViscosityType;
161
162 int m_spacedim = 2;
163 int nDimensions = m_spacedim;
164 int phys_offset;
165
166
167 ASSERTL0(vSession->DefinesParameter(
"Gamma"),
168 "Compressible flow sessions must define a Gamma parameter.");
169 vSession->LoadParameter("Gamma", m_gamma, 1.4);
170
171
172 ASSERTL0(vSession->DefinesParameter(
"pInf"),
173 "Compressible flow sessions must define a pInf parameter.");
174 vSession->LoadParameter("pInf", m_pInf, 101325);
175
176
177 ASSERTL0(vSession->DefinesParameter(
"rhoInf"),
178 "Compressible flow sessions must define a rhoInf parameter.");
179 vSession->LoadParameter(
"rhoInf",
m_rhoInf, 1.225);
180
181
182 ASSERTL0(vSession->DefinesParameter(
"uInf"),
183 "Compressible flow sessions must define a uInf parameter.");
184 vSession->LoadParameter(
"uInf",
m_uInf, 0.1);
185
186
187 if (m_spacedim == 2 || m_spacedim == 3)
188 {
189 ASSERTL0(vSession->DefinesParameter(
"vInf"),
190 "Compressible flow sessions must define a vInf parameter"
191 "for 2D/3D problems.");
192 vSession->LoadParameter(
"vInf",
m_vInf, 0.0);
193 }
194
195
196 if (m_spacedim == 3)
197 {
198 ASSERTL0(vSession->DefinesParameter(
"wInf"),
199 "Compressible flow sessions must define a wInf parameter"
200 "for 3D problems.");
201 vSession->LoadParameter("wInf", m_wInf, 0.0);
202 }
203
204 vSession->LoadParameter("GasConstant", m_gasConstant, 287.058);
205 vSession->LoadParameter(
"Twall",
m_Twall, 300.15);
206 vSession->LoadSolverInfo("ViscosityType", m_ViscosityType, "Constant");
207 vSession->LoadParameter(
"mu",
m_mu, 1.78e-05);
208
209
210
211 vector<LibUtilities::FieldDefinitionsSharedPtr> fieldDef;
212 vector<vector<NekDouble>> fieldData;
213
216
217
218
219
220 vector<vector<LibUtilities::PointsType>> pointsType;
221 for (i = 0; i < fieldDef.size(); ++i)
222 {
223 vector<LibUtilities::PointsType> ptype;
224 for (j = 0; j < 2; ++j)
225 {
227 }
228 pointsType.push_back(ptype);
229 }
230 graphShPt->SetExpansionInfo(fieldDef, pointsType);
231
232
233
234
235
236 int nfields = vSession->GetVariables().size();
239
240 for (i = 0; i < pFields.size(); i++)
241 {
242 pFields[i] =
244 vSession, graphShPt, vSession->GetVariable(i));
245 }
246
249 graphShPt);
250
251 Exp[0] = Exp2D;
252
253 for (i = 1; i < nfields; ++i)
254 {
255 Exp[i] =
257 }
258
259 int nSolutionPts = pFields[0]->GetNpoints();
260 int nTracePts = pFields[0]->GetTrace()->GetTotPoints();
261 int nElements = pFields[0]->GetExpSize();
262
264
268
272
276
277 pFields[0]->GetCoords(x, y,
z);
278
279 pFields[0]->ExtractTracePhys(x, traceX);
280 pFields[0]->ExtractTracePhys(y, traceY);
281 pFields[0]->ExtractTracePhys(
z, traceZ);
282
283
284
285
289
290
291 for (j = 0; j < nfields; ++j)
292 {
296
297 for (i = 0; i < fieldData.size(); ++i)
298 {
299 Exp[j]->ExtractDataToCoeffs(fieldDef[i], fieldData[i],
300 fieldDef[i]->m_fields[j],
301 Exp[j]->UpdateCoeffs());
302 }
303 Exp[j]->BwdTrans(Exp[j]->GetCoeffs(), Exp[j]->UpdatePhys());
304 Vmath::Vcopy(nSolutionPts, Exp[j]->GetPhys(), 1, uFields[j], 1);
305 pFields[0]->ExtractTracePhys(uFields[j], traceFields[j]);
306 }
307
308
309
310 int nfieldsAdded = 20;
313
314 for (j = 0; j < nfieldsAdded; ++j)
315 {
318 }
319
320
321
322
323
324
325
326
328 for (i = 0; i < nDimensions; ++i)
329 {
331 }
332 pFields[0]->GetTrace()->GetNormals(m_traceNormals);
333
335 for (i = 0; i < nDimensions; ++i)
336 {
338 }
339
340
341 Vmath::Vcopy(nTracePts, &m_traceNormals[0][0], 1, &traceFieldsAdded[0][0],
342 1);
343
344
345 Vmath::Vcopy(nTracePts, &m_traceNormals[1][0], 1, &traceFieldsAdded[1][0],
346 1);
347
348
349 Vmath::Vcopy(nTracePts, &m_traceNormals[1][0], 1, &m_traceTangents[0][0],
350 1);
351 Vmath::Neg(nTracePts, &m_traceTangents[0][0], 1);
352
353 Vmath::Vcopy(nTracePts, &m_traceTangents[0][0], 1, &traceFieldsAdded[2][0],
354 1);
355
356
357 Vmath::Vcopy(nTracePts, &m_traceNormals[0][0], 1, &m_traceTangents[1][0],
358 1);
359
360 Vmath::Vcopy(nTracePts, &m_traceTangents[1][0], 1, &traceFieldsAdded[3][0],
361 1);
362
363
364
365
366
367
370
371 for (i = 0; i < m_spacedim; i++)
372 {
373 Vmath::Vmul(nSolutionPts, &uFields[i + 1][0], 1, &uFields[i + 1][0], 1,
374 &tmp[0], 1);
375
376 Vmath::Smul(nSolutionPts, 0.5, &tmp[0], 1, &tmp[0], 1);
377
379 }
380
382 1);
383
386
388
389
390 pFields[0]->ExtractTracePhys(
pressure, traceFieldsAdded[4]);
391
392
393
394
395
396
398
400 &temperature[0], 1);
401
402 NekDouble GasConstantInv = 1.0 / m_gasConstant;
403 Vmath::Smul(nSolutionPts, GasConstantInv, &temperature[0], 1,
404 &temperature[0], 1);
405
406
407 pFields[0]->ExtractTracePhys(temperature, traceFieldsAdded[5]);
408
409
410
411
412
415
416 for (i = 0; i < nDimensions; ++i)
417 {
420 }
421
422 for (i = 0; i < nDimensions; ++i)
423 {
424 for (n = 0; n < nElements; n++)
425 {
426 phys_offset = pFields[0]->GetPhys_Offset(n);
427
428 pFields[i]->GetExp(n)->PhysDeriv(i, temperature + phys_offset,
429 auxArray =
430 Dtemperature[i] + phys_offset);
431 }
432
433 pFields[0]->ExtractTracePhys(Dtemperature[i], traceDtemperature[i]);
434 }
435
436 for (i = 0; i < nDimensions; ++i)
437 {
439 &traceDtemperature[i][0], 1, &tmp[0], 1);
440
441 Vmath::Vadd(nTracePts, &traceFieldsAdded[6][0], 1, &tmp[0], 1,
442 &traceFieldsAdded[6][0], 1);
443 }
444
445
446
447
448
449
450
453
454 for (i = 0; i < nDimensions; ++i)
455 {
458 }
459
460 for (i = 0; i < nDimensions; ++i)
461 {
462 for (n = 0; n < nElements; n++)
463 {
464 phys_offset = pFields[0]->GetPhys_Offset(n);
465
466 pFields[i]->GetExp(n)->PhysDeriv(i,
pressure + phys_offset,
467 auxArray =
468 Dpressure[i] + phys_offset);
469 }
470
471 pFields[0]->ExtractTracePhys(Dpressure[i], traceDpressure[i]);
472 }
473
474
475 for (i = 0; i < nDimensions; ++i)
476 {
477 Vmath::Vmul(nTracePts, &m_traceTangents[i][0], 1, &traceDpressure[i][0],
478 1, &tmp[0], 1);
479
480 Vmath::Vadd(nTracePts, &traceFieldsAdded[7][0], 1, &tmp[0], 1,
481 &traceFieldsAdded[7][0], 1);
482 }
483
484
485 Vmath::Vcopy(nTracePts, &traceDpressure[0][0], 1, &traceFieldsAdded[8][0],
486 1);
487
488
489 Vmath::Vcopy(nTracePts, &traceDpressure[1][0], 1, &traceFieldsAdded[9][0],
490 1);
491
492
493
494
495
496
497
500 nDimensions);
502
503 for (i = 0; i < nDimensions; ++i)
504 {
508
510 1);
511
512 for (j = 0; j < nDimensions; ++j)
513 {
516 }
517 }
518
519 for (i = 0; i < nDimensions; ++i)
520 {
521 for (j = 0; j < nDimensions; ++j)
522 {
523 for (n = 0; n < nElements; n++)
524 {
525 phys_offset = pFields[0]->GetPhys_Offset(n);
526
527 pFields[i]->GetExp(n)->PhysDeriv(j,
velocity[i] + phys_offset,
528 auxArray = Dvelocity[i][j] +
529 phys_offset);
530 }
531
532
533 pFields[0]->ExtractTracePhys(Dvelocity[i][j], traceDvelocity[i][j]);
534 }
535 }
536
538 &traceFieldsAdded[10][0], 1);
540 &traceFieldsAdded[11][0], 1);
542 &traceFieldsAdded[12][0], 1);
544 &traceFieldsAdded[13][0], 1);
545
546
547
548
549
550
551
552
554
555
559
560
561 if (m_ViscosityType == "Variable")
562 {
566
567 for (int i = 0; i < nSolutionPts; ++i)
568 {
569 ratio = temperature[i] / T_star;
570 mu[i] = mu_star * ratio *
sqrt(ratio) * (T_star + 110.0) /
571 (temperature[i] + 110.0);
572 }
573 }
574 else
575 {
577 }
578
579
582
583
584 Vmath::Smul(nSolutionPts, 2.0, &mu[0], 1, &mu2[0], 1);
585
586
587 Vmath::Vadd(nSolutionPts, &divVel[0], 1, &Dvelocity[0][0][0], 1, &divVel[0],
588 1);
589 Vmath::Vadd(nSolutionPts, &divVel[0], 1, &Dvelocity[1][1][0], 1, &divVel[0],
590 1);
591
592
593 Vmath::Smul(nSolutionPts, lambda, &divVel[0], 1, &divVel[0], 1);
594 Vmath::Vmul(nSolutionPts, &mu[0], 1, &divVel[0], 1, &divVel[0], 1);
595
596
597
598 for (j = 0; j < m_spacedim; ++j)
599 {
602
603 Vmath::Vmul(nSolutionPts, &mu2[0], 1, &Dvelocity[j][j][0], 1,
604 &temp[j][0], 1);
605
606 Vmath::Vadd(nSolutionPts, &temp[j][0], 1, &divVel[0], 1, &Sgg[j][0], 1);
607 }
608
609
611
612
613 Vmath::Vadd(nSolutionPts, &Dvelocity[0][1][0], 1, &Dvelocity[1][0][0], 1,
614 &Sxy[0], 1);
615
616
617 Vmath::Vmul(nSolutionPts, &mu[0], 1, &Sxy[0], 1, &Sxy[0], 1);
618
619 pFields[0]->ExtractTracePhys(Sgg[0], traceFieldsAdded[14]);
620 pFields[0]->ExtractTracePhys(Sgg[1], traceFieldsAdded[15]);
621 pFields[0]->ExtractTracePhys(Sxy, traceFieldsAdded[16]);
622
623
624
625
632
634
635
636 Vmath::Vcopy(nTracePts, &m_traceNormals[0][0], 1, &cosTeta[0], 1);
637
638
639 Vmath::Vcopy(nTracePts, &m_traceNormals[1][0], 1, &sinTeta[0], 1);
640
641
642 Vmath::Vsub(nTracePts, &traceFieldsAdded[14][0], 1,
643 &traceFieldsAdded[15][0], 1, &sigma_diff[0], 1);
644
645
646 Vmath::Vmul(nTracePts, &cosTeta[0], 1, &sinTeta[0], 1, &tmpTeta[0], 1);
647 Vmath::Smul(nTracePts, 2.0, &tmpTeta[0], 1, &sin2Teta[0], 1);
648
649
650 Vmath::Vmul(nTracePts, &cosTeta[0], 1, &cosTeta[0], 1, &cos2Teta[0], 1);
651 Vmath::Vmul(nTracePts, &sinTeta[0], 1, &sinTeta[0], 1, &tmpTeta[0], 1);
652 Vmath::Vsub(nTracePts, &cos2Teta[0], 1, &tmpTeta[0], 1, &cos2Teta[0], 1);
653
654
655 Vmath::Smul(nTracePts, -0.5, &sigma_diff[0], 1, &sigma_diff[0], 1);
656 Vmath::Vmul(nTracePts, &sigma_diff[0], 1, &sin2Teta[0], 1, &tau_t[0], 1);
657 Vmath::Vmul(nTracePts, &traceFieldsAdded[16][0], 1, &cos2Teta[0], 1,
658 &tmpTeta[0], 1);
659 Vmath::Vadd(nTracePts, &tau_t[0], 1, &tmpTeta[0], 1, &tau_t[0], 1);
660
661 Vmath::Vcopy(nTracePts, &tau_t[0], 1, &traceFieldsAdded[17][0], 1);
662
663
664
665
666
667 pFields[0]->ExtractTracePhys(mu, traceFieldsAdded[18]);
668
669
670
671
673
674
676
678 Vmath::Smul(nSolutionPts, gamma, soundspeed, 1, soundspeed, 1);
679 Vmath::Vsqrt(nSolutionPts, soundspeed, 1, soundspeed, 1);
680
681
683
684 for (int i = 0; i < m_spacedim; ++i)
685 {
686 Vmath::Vvtvp(nSolutionPts, uFields[i + 1], 1, uFields[i + 1], 1, mach,
687 1, mach, 1);
688 }
689
690 Vmath::Vdiv(nSolutionPts, mach, 1, uFields[0], 1, mach, 1);
691 Vmath::Vdiv(nSolutionPts, mach, 1, uFields[0], 1, mach, 1);
693 Vmath::Vdiv(nSolutionPts, mach, 1, soundspeed, 1, mach, 1);
694
695 pFields[0]->ExtractTracePhys(mach, traceFieldsAdded[19]);
696
697
698
699
700 if (pFields[0]->GetBndCondExpansions().size())
701 {
702 id1 = 0;
703 cnt = 0;
704 nBndRegions = pFields[0]->GetBndCondExpansions().size();
705 for (b = 0; b < nBndRegions; ++b)
706 {
707 nBndEdges = pFields[0]->GetBndCondExpansions()[b]->GetExpSize();
708 for (e = 0; e < nBndEdges; ++e)
709 {
710 nBndEdgePts = pFields[0]
711 ->GetBndCondExpansions()[b]
712 ->GetExp(e)
713 ->GetNumPoints(0);
714
715 id2 = pFields[0]->GetTrace()->GetPhys_Offset(
716 pFields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(
717 cnt++));
718
719 if (pFields[0]->GetBndConditions()[b]->GetUserDefined() ==
720 "WallViscous" ||
721 pFields[0]->GetBndConditions()[b]->GetUserDefined() ==
722 "WallAdiabatic" ||
723 pFields[0]->GetBndConditions()[b]->GetUserDefined() ==
724 "Wall")
725 {
726 Vmath::Vcopy(nBndEdgePts, &traceX[id2], 1, &surfaceX[id1],
727 1);
728
729 Vmath::Vcopy(nBndEdgePts, &traceY[id2], 1, &surfaceY[id1],
730 1);
731
732 Vmath::Vcopy(nBndEdgePts, &traceZ[id2], 1, &surfaceZ[id1],
733 1);
734
735 id1 += nBndEdgePts;
736 }
737 }
738 }
739 }
740
741
742 if (pFields[0]->GetBndCondExpansions().size())
743 {
744 for (j = 0; j < nfields; ++j)
745 {
746 id1 = 0;
747 cnt = 0;
748 nBndRegions = pFields[j]->GetBndCondExpansions().size();
749 for (b = 0; b < nBndRegions; ++b)
750 {
751 nBndEdges = pFields[j]->GetBndCondExpansions()[b]->GetExpSize();
752 for (e = 0; e < nBndEdges; ++e)
753 {
754 nBndEdgePts = pFields[j]
755 ->GetBndCondExpansions()[b]
756 ->GetExp(e)
757 ->GetNumPoints(0);
758
759 id2 = pFields[j]->GetTrace()->GetPhys_Offset(
760 pFields[j]->GetTraceMap()->GetBndCondIDToGlobalTraceID(
761 cnt++));
762
763 if (pFields[j]->GetBndConditions()[b]->GetUserDefined() ==
764 "WallViscous" ||
765 pFields[j]->GetBndConditions()[b]->GetUserDefined() ==
766 "WallAdiabatic" ||
767 pFields[j]->GetBndConditions()[b]->GetUserDefined() ==
768 "Wall")
769 {
771 &surfaceFields[j][id1], 1);
772
773 id1 += nBndEdgePts;
774 }
775 }
776 }
777 }
778 }
779
780
781 if (pFields[0]->GetBndCondExpansions().size())
782 {
783 for (j = 0; j < nfieldsAdded; ++j)
784 {
785 id1 = 0;
786 cnt = 0;
787 nBndRegions = pFields[0]->GetBndCondExpansions().size();
788 for (b = 0; b < nBndRegions; ++b)
789 {
790 nBndEdges = pFields[0]->GetBndCondExpansions()[b]->GetExpSize();
791 for (e = 0; e < nBndEdges; ++e)
792 {
793 nBndEdgePts = pFields[0]
794 ->GetBndCondExpansions()[b]
795 ->GetExp(e)
796 ->GetNumPoints(0);
797
798 id2 = pFields[0]->GetTrace()->GetPhys_Offset(
799 pFields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(
800 cnt++));
801
802 if (pFields[0]->GetBndConditions()[b]->GetUserDefined() ==
803 "WallViscous" ||
804 pFields[0]->GetBndConditions()[b]->GetUserDefined() ==
805 "WallAdiabatic" ||
806 pFields[0]->GetBndConditions()[b]->GetUserDefined() ==
807 "Wall")
808 {
810 &surfaceFieldsAdded[j][id1], 1);
811
812 id1 += nBndEdgePts;
813 }
814 }
815 }
816 }
817 }
818
819
820
821 std::string vEquation = vSession->GetSolverInfo("EQType");
822
824 BndExp = pFields[0]->GetBndCondExpansions();
826
832
833 int GlobalIndex(0);
834
835 for (int i = 0; i < BndExp[0]->GetExpSize(); ++i)
836 {
838
840
845
846
847
851
855
857
858 for (int j = 0; j < nbc; ++j)
859 {
860
861 nxOnBnd[j] = surfaceFieldsAdded[0][GlobalIndex];
862 nyOnBnd[j] = surfaceFieldsAdded[1][GlobalIndex];
863 txOnBnd[j] = surfaceFieldsAdded[2][GlobalIndex];
864 tyOnBnd[j] = surfaceFieldsAdded[3][GlobalIndex];
865
866 PressurOnBnd[j] = surfaceFieldsAdded[4][GlobalIndex];
867
868 if (vEquation == "NavierStokesCFE")
869 {
870 ShearStressOnBnd[j] = surfaceFieldsAdded[17][GlobalIndex];
871 }
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887 GlobalIndex++;
888 }
889
890 Vmath::Vmul(nbc, PressurOnBnd, 1, nxOnBnd, 1, drag_p, 1);
891 Vmath::Vmul(nbc, PressurOnBnd, 1, nyOnBnd, 1, lift_p, 1);
892
893
894
895
896 Fxp += bc->Integral(drag_p);
897 Fyp += bc->Integral(lift_p);
898
899 if (vEquation == "NavierStokesCFE")
900 {
901 Vmath::Vmul(nbc, ShearStressOnBnd, 1, txOnBnd, 1, drag_v, 1);
902 Vmath::Vmul(nbc, ShearStressOnBnd, 1, tyOnBnd, 1, lift_v, 1);
903
904
905
906
907 Fxv += bc->Integral(drag_v);
908 Fyv += bc->Integral(lift_v);
909 }
910
911 Sref += bc->Integral(Unity);
912 }
913
914 cout << "\n Sref = " << Sref << endl;
915 Fxp = Fxp / Sref;
916 Fyp = Fyp / Sref;
917 Fxv = Fxv / Sref;
918 Fyv = Fyv / Sref;
919 cout << " Pressure drag (Fxp) = " << Fxp << endl;
920 cout << " Pressure lift (Fyp) = " << Fyp << endl;
921 cout << " Viscous drag (Fxv) = " << Fxv << endl;
922 cout << " Viscous lift (Fyv) = " << Fyv << endl;
923 cout << "\n ==> Total drag = " << Fxp + Fxv << endl;
924 cout << " ==> Total lift = " << Fyp + Fyv << "\n" << endl;
925
926
927
928
929
930
931 ofstream outfile;
932 outfile.open(fname.c_str());
933 outfile << "% x[m] "
934 << " \t"
935 << "y[m] "
936 << " \t"
937 << "z[m] "
938 << " \t"
939 << "nx[] "
940 << " \t"
941 << "ny[] "
942 << " \t"
943 << "tx[] "
944 << " \t"
945 << "ty[] "
946 << " \t"
947 << "rho[kg/m^3] "
948 << " \t"
949 << "rhou[kg/(m^2 s)] "
950 << " \t"
951 << "rhov[kg/(m^2 s)] "
952 << " \t"
953 << "E[Pa] "
954 << " \t"
955 << "p[Pa] "
956 << " \t"
957 << "T[k] "
958 << " \t"
959 << "dT/dn[k/m] "
960 << " \t"
961 << "dp/dT[Pa/m] "
962 << " \t"
963 << "dp/dx[Pa/m] "
964 << " \t"
965 << "dp/dy[Pa/m] "
966 << " \t"
967 << "du/dx[s^-1] "
968 << " \t"
969 << "du/dy[s^-1] "
970 << " \t"
971 << "dv/dx[s^-1] "
972 << " \t"
973 << "dv/dy[s^-1] "
974 << " \t"
975 << "tau_xx[Pa] "
976 << " \t"
977 << "tau_yy[Pa] "
978 << " \t"
979 << "tau_xy[Pa] "
980 << " \t"
981 << "tau_t[Pa] "
982 << " \t"
983 << "mu[Pa s] "
984 << " \t"
985 << "M[] "
986 << " \t"
987
988 << endl;
989 for (i = 0; i < id1; ++i)
990 {
991 outfile << scientific << setw(17) << setprecision(16) << surfaceX[i]
992 << " \t " << surfaceY[i] << " \t " << surfaceZ[i] << " \t "
993 << surfaceFieldsAdded[0][i] << " \t "
994 << surfaceFieldsAdded[1][i] << " \t "
995 << surfaceFieldsAdded[2][i] << " \t "
996 << surfaceFieldsAdded[3][i] << " \t " << surfaceFields[0][i]
997 << " \t " << surfaceFields[1][i] << " \t "
998 << surfaceFields[2][i] << " \t " << surfaceFields[3][i]
999 << " \t " << surfaceFieldsAdded[4][i] << " \t "
1000 << surfaceFieldsAdded[5][i] << " \t "
1001 << surfaceFieldsAdded[6][i] << " \t "
1002 << surfaceFieldsAdded[7][i] << " \t "
1003 << surfaceFieldsAdded[8][i] << " \t "
1004 << surfaceFieldsAdded[9][i] << " \t "
1005 << surfaceFieldsAdded[10][i] << " \t "
1006 << surfaceFieldsAdded[11][i] << " \t "
1007 << surfaceFieldsAdded[12][i] << " \t "
1008 << surfaceFieldsAdded[13][i] << " \t "
1009 << surfaceFieldsAdded[14][i] << " \t "
1010 << surfaceFieldsAdded[15][i] << " \t "
1011 << surfaceFieldsAdded[16][i] << " \t "
1012 << surfaceFieldsAdded[17][i] << " \t "
1013 << surfaceFieldsAdded[18][i] << " \t "
1014 << surfaceFieldsAdded[19][i]
1015 << " \t "
1016
1017 << endl;
1018 }
1019 outfile << endl << endl;
1020 outfile.close();
1021
1022 return 0;
1023}
#define ASSERTL0(condition, msg)
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
void Import(const std::string &infilename, std::vector< FieldDefinitionsSharedPtr > &fielddefs, std::vector< std::vector< NekDouble > > &fielddata, FieldMetaDataMap &fieldinfomap, const Array< OneD, int > &ElementIDs)
This function allows for data to be imported from an FLD file when a session and/or communicator is n...
std::map< std::string, std::string > FieldMetaDataMap
std::shared_ptr< SessionReader > SessionReaderSharedPtr
@ ePolyEvenlySpaced
1D Evenly-spaced points using Lagrange polynomial
const std::vector< NekDouble > velocity
std::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
std::shared_ptr< MeshGraph > MeshGraphSharedPtr
std::shared_ptr< StdExpansion1D > StdExpansion1DSharedPtr
std::vector< double > z(NPUPPER)
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
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 Neg(int n, T *x, const int incx)
Negate x = -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 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 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.
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
scalarT< T > sqrt(scalarT< T > in)