#include <DiffusionLDGNS.h>
Inheritance diagram for Nektar::DiffusionLDGNS:
Static Public Member Functions | |
| static DiffusionSharedPtr | create (std::string diffType) |
Static Public Attributes | |
| static std::string | type |
Protected Member Functions inherited from Private Member Functions | |
| void | NumericalFluxO1 (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &numericalFluxO1, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd) |
| Builds the numerical flux for the 1st order derivatives. More... | |
| void | ApplyBCsO1 (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd, Array< OneD, Array< OneD, NekDouble > > &flux01) |
| Imposes appropriate bcs for the 1st order derivatives. More... | |
| void | NumericalFluxO2 (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, TensorOfArray3D< NekDouble > &qfield, Array< OneD, Array< OneD, NekDouble > > &qflux, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd) |
| Build the numerical flux for the 2nd order derivatives. More... | |
| void | ApplyBCsO2 (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const std::size_t var, const std::size_t dir, const Array< OneD, const NekDouble > &qfield, const Array< OneD, const NekDouble > &qFwd, const Array< OneD, const NekDouble > &qBwd, Array< OneD, NekDouble > &penaltyflux) |
| Imposes appropriate bcs for the 2nd order derivatives. More... | |
Detailed Description
Definition at line 47 of file DiffusionLDGNS.h.
Constructor & Destructor Documentation
◆ DiffusionLDGNS()
|
protected |
Member Function Documentation
◆ ApplyBCsO1()
|
private |
Imposes appropriate bcs for the 1st order derivatives.
Definition at line 420 of file DiffusionLDGNS.cpp.
426{
427 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
428 std::size_t nScalars = inarray.size();
429
430 Array<OneD, NekDouble> tmp1{nTracePts, 0.0};
431 Array<OneD, NekDouble> tmp2{nTracePts, 0.0};
432 Array<OneD, NekDouble> Tw{nTracePts, m_Twall};
433
434 Array<OneD, Array<OneD, NekDouble>> scalarVariables{nScalars};
435
436 // Extract internal values of the scalar variables for Neumann bcs
437 for (std::size_t i = 0; i < nScalars; ++i)
438 {
439 scalarVariables[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
440 }
441
442 // Compute boundary conditions for velocities
443 for (std::size_t i = 0; i < nScalars - 1; ++i)
444 {
445 // Note that cnt has to loop on nBndRegions and nBndEdges
446 // and has to be reset to zero per each equation
447 std::size_t cnt = 0;
448 std::size_t nBndRegions = fields[i + 1]->GetBndCondExpansions().size();
449 for (std::size_t j = 0; j < nBndRegions; ++j)
450 {
451 if (fields[i + 1]
452 ->GetBndConditions()[j]
453 ->GetBoundaryConditionType() == SpatialDomains::ePeriodic)
454 {
455 continue;
456 }
457
458 std::size_t nBndEdges =
459 fields[i + 1]->GetBndCondExpansions()[j]->GetExpSize();
460 for (std::size_t e = 0; e < nBndEdges; ++e)
461 {
462 std::size_t nBndEdgePts = fields[i + 1]
463 ->GetBndCondExpansions()[j]
464 ->GetExp(e)
465 ->GetTotPoints();
466
467 std::size_t id1 =
468 fields[i + 1]->GetBndCondExpansions()[j]->GetPhys_Offset(e);
469
470 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
471 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(
472 cnt++));
473
474 if (boost::iequals(
475 fields[i]->GetBndConditions()[j]->GetUserDefined(),
476 "WallViscous") ||
477 boost::iequals(
478 fields[i]->GetBndConditions()[j]->GetUserDefined(),
479 "WallAdiabatic"))
480 {
481 // Reinforcing bcs for velocity in case of Wall bcs
482 Vmath::Zero(nBndEdgePts, &scalarVariables[i][id2], 1);
483 }
484 else if (boost::iequals(
485 fields[i]->GetBndConditions()[j]->GetUserDefined(),
486 "Wall") ||
487 boost::iequals(
488 fields[i]->GetBndConditions()[j]->GetUserDefined(),
489 "Symmetry"))
490 {
491 // Symmetry bc: normal velocity is zero
492 // get all velocities at once because we need u.n
493 if (i == 0)
494 {
495 // tmp1 = -(u.n)
496 Vmath::Zero(nBndEdgePts, tmp1, 1);
497 for (std::size_t k = 0; k < nScalars - 1; ++k)
498 {
499 Vmath::Vdiv(nBndEdgePts,
500 &(fields[k + 1]
501 ->GetBndCondExpansions()[j]
502 ->UpdatePhys())[id1],
503 1,
504 &(fields[0]
505 ->GetBndCondExpansions()[j]
506 ->UpdatePhys())[id1],
507 1, &scalarVariables[k][id2], 1);
509 1, &scalarVariables[k][id2], 1,
510 &tmp1[0], 1, &tmp1[0], 1);
511 }
512 Vmath::Smul(nBndEdgePts, -1.0, &tmp1[0], 1, &tmp1[0],
513 1);
514
515 // u_i - (u.n)n_i
516 for (std::size_t k = 0; k < nScalars - 1; ++k)
517 {
518 Vmath::Vvtvp(nBndEdgePts, &tmp1[0], 1,
519 &m_traceNormals[k][id2], 1,
520 &scalarVariables[k][id2], 1,
521 &scalarVariables[k][id2], 1);
522 }
523 }
524 }
525 else if (fields[i]
526 ->GetBndConditions()[j]
527 ->GetBoundaryConditionType() ==
529 {
530 // Imposing velocity bcs if not Wall
531 Vmath::Vdiv(nBndEdgePts,
532 &(fields[i + 1]
533 ->GetBndCondExpansions()[j]
534 ->UpdatePhys())[id1],
535 1,
536 &(fields[0]
537 ->GetBndCondExpansions()[j]
538 ->UpdatePhys())[id1],
539 1, &scalarVariables[i][id2], 1);
540 }
541
542 // For Dirichlet boundary condition: uflux = u_bcs
543 if (fields[i]
544 ->GetBndConditions()[j]
545 ->GetBoundaryConditionType() ==
547 {
548 Vmath::Vcopy(nBndEdgePts, &scalarVariables[i][id2], 1,
549 &fluxO1[i][id2], 1);
550 }
551
552 // For Neumann boundary condition: uflux = u_+
553 else if ((fields[i]->GetBndConditions()[j])
554 ->GetBoundaryConditionType() ==
556 {
557 Vmath::Vcopy(nBndEdgePts, &pFwd[i][id2], 1, &fluxO1[i][id2],
558 1);
559 }
560
561 // Building kinetic energy to be used for T bcs
562 Vmath::Vmul(nBndEdgePts, &scalarVariables[i][id2], 1,
563 &scalarVariables[i][id2], 1, &tmp1[id2], 1);
564
565 Vmath::Smul(nBndEdgePts, 0.5, &tmp1[id2], 1, &tmp1[id2], 1);
566
567 Vmath::Vadd(nBndEdgePts, &tmp2[id2], 1, &tmp1[id2], 1,
568 &tmp2[id2], 1);
569 }
570 }
571 }
572
573 // Compute boundary conditions for temperature
574 std::size_t cnt = 0;
575 std::size_t nBndRegions = fields[nScalars]->GetBndCondExpansions().size();
576 for (std::size_t j = 0; j < nBndRegions; ++j)
577 {
578 if (fields[nScalars]
579 ->GetBndConditions()[j]
580 ->GetBoundaryConditionType() == SpatialDomains::ePeriodic)
581 {
582 continue;
583 }
584
585 std::size_t nBndEdges =
586 fields[nScalars]->GetBndCondExpansions()[j]->GetExpSize();
587 for (std::size_t e = 0; e < nBndEdges; ++e)
588 {
589 std::size_t nBndEdgePts = fields[nScalars]
590 ->GetBndCondExpansions()[j]
591 ->GetExp(e)
592 ->GetTotPoints();
593
594 std::size_t id1 =
595 fields[nScalars]->GetBndCondExpansions()[j]->GetPhys_Offset(e);
596
597 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
598 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
599
600 // Imposing Temperature Twall at the wall
601 if (boost::iequals(
602 fields[nScalars]->GetBndConditions()[j]->GetUserDefined(),
603 "WallViscous"))
604 {
605 Vmath::Vcopy(nBndEdgePts, &Tw[0], 1,
606 &scalarVariables[nScalars - 1][id2], 1);
607 }
608 // Imposing Temperature through condition on the Energy
609 // for no wall boundaries (e.g. farfield)
610 else if (fields[nScalars]
611 ->GetBndConditions()[j]
612 ->GetBoundaryConditionType() ==
614 {
615 // Use equation of state to evaluate temperature
616 NekDouble rho, ene;
617 for (std::size_t n = 0; n < nBndEdgePts; ++n)
618 {
619 rho = fields[0]
620 ->GetBndCondExpansions()[j]
621 ->GetPhys()[id1 + n];
622 ene = fields[nScalars]
623 ->GetBndCondExpansions()[j]
624 ->GetPhys()[id1 + n] /
625 rho -
626 tmp2[id2 + n];
627 scalarVariables[nScalars - 1][id2 + n] =
628 m_eos->GetTemperature(rho, ene);
629 }
630 }
631
632 // For Dirichlet boundary condition: uflux = u_bcs
633 if (fields[nScalars]
634 ->GetBndConditions()[j]
635 ->GetBoundaryConditionType() ==
636 SpatialDomains::eDirichlet &&
637 !boost::iequals(
638 fields[nScalars]->GetBndConditions()[j]->GetUserDefined(),
639 "WallAdiabatic"))
640 {
641 Vmath::Vcopy(nBndEdgePts, &scalarVariables[nScalars - 1][id2],
642 1, &fluxO1[nScalars - 1][id2], 1);
643 }
644
645 // For Neumann boundary condition: uflux = u_+
646 else if (((fields[nScalars]->GetBndConditions()[j])
647 ->GetBoundaryConditionType() ==
648 SpatialDomains::eNeumann) ||
649 boost::iequals(fields[nScalars]
650 ->GetBndConditions()[j]
651 ->GetUserDefined(),
652 "WallAdiabatic"))
653 {
654 Vmath::Vcopy(nBndEdgePts, &pFwd[nScalars - 1][id2], 1,
655 &fluxO1[nScalars - 1][id2], 1);
656 }
657 }
658 }
659}
EquationOfStateSharedPtr m_eos
Equation of system for computing temperature.
Definition: DiffusionLDGNS.h:72
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Definition: DiffusionLDGNS.h:67
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.
Definition: Vmath.hpp:72
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
Definition: Vmath.hpp:366
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.
Definition: Vmath.hpp:180
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.hpp:100
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.
Definition: Vmath.hpp:126
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825
References Nektar::SpatialDomains::eDirichlet, Nektar::SpatialDomains::eNeumann, Nektar::SpatialDomains::ePeriodic, m_eos, m_traceNormals, m_Twall, Vmath::Smul(), Vmath::Vadd(), Vmath::Vcopy(), Vmath::Vdiv(), Vmath::Vmul(), Vmath::Vvtvp(), and Vmath::Zero().
Referenced by NumericalFluxO1().
◆ ApplyBCsO2()
|
private |
Imposes appropriate bcs for the 2nd order derivatives.
Definition at line 728 of file DiffusionLDGNS.cpp.
735{
736 std::size_t cnt = 0;
737 std::size_t nBndRegions = fields[var]->GetBndCondExpansions().size();
738 // Loop on the boundary regions to apply appropriate bcs
739 for (std::size_t i = 0; i < nBndRegions; ++i)
740 {
741 // Number of boundary regions related to region 'i'
742 std::size_t nBndEdges =
743 fields[var]->GetBndCondExpansions()[i]->GetExpSize();
744
745 if (fields[var]->GetBndConditions()[i]->GetBoundaryConditionType() ==
747 {
748 continue;
749 }
750
751 // Weakly impose bcs by modifying flux values
752 for (std::size_t e = 0; e < nBndEdges; ++e)
753 {
754 std::size_t nBndEdgePts = fields[var]
755 ->GetBndCondExpansions()[i]
756 ->GetExp(e)
757 ->GetTotPoints();
758
759 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
760 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
761
762 // In case of Dirichlet bcs:
763 // uflux = gD
764 if (fields[var]
765 ->GetBndConditions()[i]
766 ->GetBoundaryConditionType() ==
767 SpatialDomains::eDirichlet &&
768 !boost::iequals(
769 fields[var]->GetBndConditions()[i]->GetUserDefined(),
770 "WallAdiabatic"))
771 {
773 &qFwd[id2], 1, &penaltyflux[id2], 1);
774 }
775 // 3.4) In case of Neumann bcs:
776 // uflux = u+
777 else if ((fields[var]->GetBndConditions()[i])
778 ->GetBoundaryConditionType() ==
780 {
782 }
783 else if (boost::iequals(
784 fields[var]->GetBndConditions()[i]->GetUserDefined(),
785 "WallAdiabatic"))
786 {
788 {
789 Vmath::Zero(nBndEdgePts, &penaltyflux[id2], 1);
790 }
791 else
792 {
793
795 &qFwd[id2], 1, &penaltyflux[id2], 1);
796 }
797 }
798 }
799 }
800}
References ASSERTL0, Nektar::SpatialDomains::eDirichlet, Nektar::SpatialDomains::eNeumann, Nektar::SpatialDomains::ePeriodic, m_spaceDim, m_traceNormals, Vmath::Vmul(), and Vmath::Zero().
Referenced by NumericalFluxO2().
◆ create()
|
inlinestatic |
Definition at line 50 of file DiffusionLDGNS.h.
51 {
53 }
std::shared_ptr< Diffusion > DiffusionSharedPtr
A shared pointer to an EquationSystem object.
Definition: Diffusion.h:54
References DiffusionLDGNS().
◆ NumericalFluxO1()
|
private |
Builds the numerical flux for the 1st order derivatives.
Definition at line 382 of file DiffusionLDGNS.cpp.
388{
389 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
390 std::size_t nScalars = inarray.size();
391
392 // Upwind
393 Array<OneD, Array<OneD, NekDouble>> numflux{nScalars};
394 for (std::size_t i = 0; i < nScalars; ++i)
395 {
396 numflux[i] = {pFwd[i]};
397 }
398
399 // Modify the values in case of boundary interfaces
400 if (fields[0]->GetBndCondExpansions().size())
401 {
402 ApplyBCsO1(fields, inarray, pFwd, pBwd, numflux);
403 }
404
405 // Splitting the numerical flux into the dimensions
407 {
408 for (std::size_t i = 0; i < nScalars; ++i)
409 {
411 numericalFluxO1[j][i], 1);
412 }
413 }
414}
void ApplyBCsO1(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd, Array< OneD, Array< OneD, NekDouble > > &flux01)
Imposes appropriate bcs for the 1st order derivatives.
Definition: DiffusionLDGNS.cpp:420
References ApplyBCsO1(), m_spaceDim, m_traceNormals, and Vmath::Vmul().
Referenced by v_DiffuseCalcDerivative().
◆ NumericalFluxO2()
|
private |
Build the numerical flux for the 2nd order derivatives.
Definition at line 665 of file DiffusionLDGNS.cpp.
671{
672 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
673 std::size_t nVariables = fields.size();
674
675 // Initialize penalty flux
676 Array<OneD, Array<OneD, NekDouble>> fluxPen{nVariables - 1};
677 for (std::size_t i = 0; i < nVariables - 1; ++i)
678 {
679 fluxPen[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
680 }
681
682 // Get penalty flux
683 m_fluxPenaltyNS(uFwd, uBwd, fluxPen);
684
685 // Evaluate Riemann flux
686 // qflux = \hat{q} \cdot u = q \cdot n
687 // Notice: i = 1 (first row of the viscous tensor is zero)
688
689 Array<OneD, NekDouble> qFwd{nTracePts};
690 Array<OneD, NekDouble> qBwd{nTracePts};
691 Array<OneD, NekDouble> qtemp{nTracePts};
692 Array<OneD, NekDouble> qfluxtemp{nTracePts, 0.0};
693 std::size_t nDim = fields[0]->GetCoordim(0);
694 for (std::size_t i = 1; i < nVariables; ++i)
695 {
696 for (std::size_t j = 0; j < nDim; ++j)
697 {
698 // Compute qFwd and qBwd value of qfield in position 'ji'
699 fields[i]->GetFwdBwdTracePhys(qfield[j][i], qFwd, qBwd);
700
701 // Downwind
702 Vmath::Vcopy(nTracePts, qBwd, 1, qfluxtemp, 1);
703
704 // Multiply the Riemann flux by the trace normals
706 qfluxtemp, 1);
707
708 // Add penalty term
710 qfluxtemp, 1, qfluxtemp, 1);
711
712 // Impose weak boundary condition with flux
713 if (fields[0]->GetBndCondExpansions().size())
714 {
715 ApplyBCsO2(fields, i, j, qfield[j][i], qFwd, qBwd, qfluxtemp);
716 }
717
718 // Store the final flux into qflux
719 Vmath::Vadd(nTracePts, qfluxtemp, 1, qflux[i], 1, qflux[i], 1);
720 }
721 }
722}
void ApplyBCsO2(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const std::size_t var, const std::size_t dir, const Array< OneD, const NekDouble > &qfield, const Array< OneD, const NekDouble > &qFwd, const Array< OneD, const NekDouble > &qBwd, Array< OneD, NekDouble > &penaltyflux)
Imposes appropriate bcs for the 2nd order derivatives.
Definition: DiffusionLDGNS.cpp:728
DiffusionFluxPenaltyNS m_fluxPenaltyNS
Definition: Diffusion.h:348
References ApplyBCsO2(), Nektar::SolverUtils::Diffusion::m_fluxPenaltyNS, m_traceNormals, m_traceOneOverH, Vmath::Vadd(), Vmath::Vcopy(), Vmath::Vmul(), and Vmath::Vvtvp().
Referenced by v_DiffuseTraceFlux().
◆ v_Diffuse()
|
overrideprotectedvirtual |
Calculate weak DG Diffusion in the LDG form for the Navier-Stokes (NS) equations:
\( \langle\psi, \hat{u}\cdot n\rangle - \langle\nabla\psi \cdot u\rangle \langle\phi, \hat{q}\cdot n\rangle - (\nabla \phi \cdot q) \rangle \)
The equations that need a diffusion operator are those related with the velocities and with the energy.
Implements Nektar::SolverUtils::Diffusion.
Definition at line 209 of file DiffusionLDGNS.cpp.
216{
217 std::size_t nCoeffs = fields[0]->GetNcoeffs();
218 Array<OneD, Array<OneD, NekDouble>> tmp2{nConvectiveFields};
219 for (std::size_t i = 0; i < nConvectiveFields; ++i)
220 {
221 tmp2[i] = Array<OneD, NekDouble>{nCoeffs, 0.0};
222 }
223 v_DiffuseCoeffs(nConvectiveFields, fields, inarray, tmp2, pFwd, pBwd);
224 for (std::size_t i = 0; i < nConvectiveFields; ++i)
225 {
226 fields[i]->BwdTrans(tmp2[i], outarray[i]);
227 }
228}
void v_DiffuseCoeffs(const std::size_t nConvective, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd) override
Definition: DiffusionLDGNS.cpp:230
References v_DiffuseCoeffs().
◆ v_DiffuseCalcDerivative()
|
overrideprotectedvirtual |
Diffusion Flux, calculate the physical derivatives.
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 305 of file DiffusionLDGNS.cpp.
311{
312 std::size_t nDim = fields[0]->GetCoordim(0);
313 std::size_t nCoeffs = fields[0]->GetNcoeffs();
314 std::size_t nScalars = inarray.size();
315 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
316 std::size_t nConvectiveFields = fields.size();
317
318 Array<OneD, NekDouble> tmp1{nCoeffs};
319 Array<OneD, Array<OneD, NekDouble>> tmp2{nConvectiveFields};
320 TensorOfArray3D<NekDouble> numericalFluxO1{m_spaceDim};
321
323 {
324 numericalFluxO1[j] = Array<OneD, Array<OneD, NekDouble>>{nScalars};
325
326 for (std::size_t i = 0; i < nScalars; ++i)
327 {
328 numericalFluxO1[j][i] = Array<OneD, NekDouble>{nTracePts, 0.0};
329 }
330 }
331
332 NumericalFluxO1(fields, inarray, numericalFluxO1, pFwd, pBwd);
333
334 for (std::size_t j = 0; j < nDim; ++j)
335 {
336 for (std::size_t i = 0; i < nScalars; ++i)
337 {
338 fields[i]->IProductWRTDerivBase(j, inarray[i], tmp1);
339 Vmath::Neg(nCoeffs, tmp1, 1);
340 fields[i]->AddTraceIntegral(numericalFluxO1[j][i], tmp1);
341 fields[i]->SetPhysState(false);
342 fields[i]->MultiplyByElmtInvMass(tmp1, tmp1);
343 fields[i]->BwdTrans(tmp1, qfields[j][i]);
344 }
345 }
346 // For 3D Homogeneous 1D only take derivatives in 3rd direction
348 {
349 for (std::size_t i = 0; i < nScalars; ++i)
350 {
351 qfields[2][i] = m_homoDerivs[i];
352 }
353 }
354}
Array< OneD, Array< OneD, NekDouble > > m_homoDerivs
Definition: DiffusionLDGNS.h:76
void NumericalFluxO1(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &numericalFluxO1, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
Builds the numerical flux for the 1st order derivatives.
Definition: DiffusionLDGNS.cpp:382
References m_diffDim, m_homoDerivs, m_spaceDim, Vmath::Neg(), and NumericalFluxO1().
◆ v_DiffuseCoeffs()
|
overrideprotectedvirtual |
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 230 of file DiffusionLDGNS.cpp.
237{
238 std::size_t nDim = fields[0]->GetCoordim(0);
239 std::size_t nPts = fields[0]->GetTotPoints();
240 std::size_t nCoeffs = fields[0]->GetNcoeffs();
241 std::size_t nScalars = inarray.size();
242 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
243
244 TensorOfArray3D<NekDouble> derivativesO1{m_spaceDim};
245
247 {
248 derivativesO1[j] = Array<OneD, Array<OneD, NekDouble>>{nScalars};
249
250 for (std::size_t i = 0; i < nScalars; ++i)
251 {
252 derivativesO1[j][i] = Array<OneD, NekDouble>{nPts, 0.0};
253 }
254 }
255
256 DiffuseCalcDerivative(fields, inarray, derivativesO1, pFwd, pBwd);
257
258 // Initialisation viscous tensor
260 Array<OneD, Array<OneD, NekDouble>> viscousFlux{nConvectiveFields};
261
263 {
264 m_viscTensor[j] = Array<OneD, Array<OneD, NekDouble>>{nScalars + 1};
265 for (std::size_t i = 0; i < nScalars + 1; ++i)
266 {
267 m_viscTensor[j][i] = Array<OneD, NekDouble>{nPts, 0.0};
268 }
269 }
270
271 for (std::size_t i = 0; i < nConvectiveFields; ++i)
272 {
273 viscousFlux[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
274 }
275
277
278 // Compute u from q_{\eta} and q_{\xi}
279 // Obtain numerical fluxes
280 // Note: derivativesO1 is not used in DiffuseTraceFlux
282 pFwd, pBwd);
283
284 Array<OneD, NekDouble> tmpOut{nCoeffs};
285 Array<OneD, Array<OneD, NekDouble>> tmpIn{nDim};
286
287 for (std::size_t i = 0; i < nConvectiveFields; ++i)
288 {
289 // Temporary fix to call collection op
290 // we can reorder m_viscTensor later
291 for (std::size_t j = 0; j < nDim; ++j)
292 {
293 tmpIn[j] = m_viscTensor[j][i];
294 }
295 fields[i]->IProductWRTDerivBase(tmpIn, tmpOut);
296
297 // Evaulate <\phi, \hat{F}\cdot n> - outarray[i]
298 Vmath::Neg(nCoeffs, tmpOut, 1);
299 fields[i]->AddTraceIntegral(viscousFlux[i], tmpOut);
300 fields[i]->SetPhysState(false);
301 fields[i]->MultiplyByElmtInvMass(tmpOut, outarray[i]);
302 }
303}
TensorOfArray3D< NekDouble > m_viscTensor
Definition: DiffusionLDGNS.h:74
SOLVER_UTILS_EXPORT void DiffuseCalcDerivative(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, const Array< OneD, Array< OneD, NekDouble > > &pFwd=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble > > &pBwd=NullNekDoubleArrayOfArray)
Definition: Diffusion.h:201
SOLVER_UTILS_EXPORT void DiffuseVolumeFlux(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, TensorOfArray3D< NekDouble > &VolumeFlux, Array< OneD, int > &nonZeroIndex=NullInt1DArray)
Diffusion Volume FLux.
Definition: Diffusion.h:214
SOLVER_UTILS_EXPORT void DiffuseTraceFlux(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, TensorOfArray3D< NekDouble > &VolumeFlux, Array< OneD, Array< OneD, NekDouble > > &TraceFlux, const Array< OneD, Array< OneD, NekDouble > > &pFwd=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble > > &pBwd=NullNekDoubleArrayOfArray, Array< OneD, int > &nonZeroIndex=NullInt1DArray)
Diffusion term Trace Flux.
Definition: Diffusion.h:225
References Nektar::SolverUtils::Diffusion::DiffuseCalcDerivative(), Nektar::SolverUtils::Diffusion::DiffuseTraceFlux(), Nektar::SolverUtils::Diffusion::DiffuseVolumeFlux(), m_spaceDim, m_viscTensor, and Vmath::Neg().
Referenced by v_Diffuse().
◆ v_DiffuseTraceFlux()
|
overrideprotectedvirtual |
Diffusion term Trace Flux.
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 365 of file DiffusionLDGNS.cpp.
374{
375 NumericalFluxO2(fields, qfields, TraceFlux, pFwd, pBwd);
376}
void NumericalFluxO2(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, TensorOfArray3D< NekDouble > &qfield, Array< OneD, Array< OneD, NekDouble > > &qflux, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
Build the numerical flux for the 2nd order derivatives.
Definition: DiffusionLDGNS.cpp:665
References NumericalFluxO2().
◆ v_DiffuseVolumeFlux()
|
overrideprotectedvirtual |
Diffusion Volume Flux.
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 356 of file DiffusionLDGNS.cpp.
361{
362 m_fluxVectorNS(inarray, qfields, VolumeFlux);
363}
DiffusionFluxVecCBNS m_fluxVectorNS
Definition: Diffusion.h:347
References Nektar::SolverUtils::Diffusion::m_fluxVectorNS.
◆ v_GetFluxTensor()
|
inlineoverrideprotectedvirtual |
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 129 of file DiffusionLDGNS.h.
References m_viscTensor.
◆ v_InitObject()
|
overrideprotectedvirtual |
Implements Nektar::SolverUtils::Diffusion.
Definition at line 53 of file DiffusionLDGNS.cpp.
56{
57 m_session = pSession;
59
60 // Setting up the normals
61 std::size_t nDim = pFields[0]->GetCoordim(0);
62 std::size_t nTracePts = pFields[0]->GetTrace()->GetTotPoints();
63
64 m_spaceDim = nDim;
65 if (pSession->DefinesSolverInfo("HOMOGENEOUS"))
66 {
67 m_spaceDim = 3;
68 }
69
71
75 {
76 m_traceVel[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
77 m_traceNormals[i] = Array<OneD, NekDouble>{nTracePts};
78 }
79 pFields[0]->GetTrace()->GetNormals(m_traceNormals);
80
81 // Create equation of state object
82 std::string eosType;
85
86 // Set up {h} reference on the trace for penalty term
87 //
88 // Note, this shold be replaced with something smarter when merging
89 // LDG with IP
90
91 // Get min h per element
92 std::size_t nElements = pFields[0]->GetExpSize();
93 Array<OneD, NekDouble> hEle{nElements, 1.0};
94 for (std::size_t e = 0; e < nElements; e++)
95 {
96 NekDouble h{1.0e+10};
97 std::size_t expDim = pFields[0]->GetShapeDimension();
98 switch (expDim)
99 {
100 case 3:
101 {
102 LocalRegions::Expansion3DSharedPtr exp3D;
103 exp3D = pFields[0]->GetExp(e)->as<LocalRegions::Expansion3D>();
104 for (std::size_t i = 0; i < exp3D->GetNtraces(); ++i)
105 {
106 h = std::min(
107 h, exp3D->GetGeom3D()->GetEdge(i)->GetVertex(0)->dist(*(
108 exp3D->GetGeom3D()->GetEdge(i)->GetVertex(1))));
109 }
110 break;
111 }
112
113 case 2:
114 {
115 LocalRegions::Expansion2DSharedPtr exp2D;
116 exp2D = pFields[0]->GetExp(e)->as<LocalRegions::Expansion2D>();
117 for (std::size_t i = 0; i < exp2D->GetNtraces(); ++i)
118 {
119 h = std::min(
120 h, exp2D->GetGeom2D()->GetEdge(i)->GetVertex(0)->dist(*(
121 exp2D->GetGeom2D()->GetEdge(i)->GetVertex(1))));
122 }
123 break;
124 }
125 case 1:
126 {
127 LocalRegions::Expansion1DSharedPtr exp1D;
128 exp1D = pFields[0]->GetExp(e)->as<LocalRegions::Expansion1D>();
129
130 h = std::min(h, exp1D->GetGeom1D()->GetVertex(0)->dist(
131 *(exp1D->GetGeom1D()->GetVertex(1))));
132 break;
133 }
134 default:
135 {
137 }
138 }
139
140 // Store scaling
141 hEle[e] = h;
142 }
143 // Expand h from elements to points
144 std::size_t nPts = pFields[0]->GetTotPoints();
145 Array<OneD, NekDouble> hElePts{nPts, 0.0};
146 Array<OneD, NekDouble> tmp;
147 for (std::size_t e = 0; e < pFields[0]->GetExpSize(); e++)
148 {
149 std::size_t nElmtPoints = pFields[0]->GetExp(e)->GetTotPoints();
150 std::size_t physOffset = pFields[0]->GetPhys_Offset(e);
151 Vmath::Fill(nElmtPoints, hEle[e], tmp = hElePts + physOffset, 1);
152 }
153 // Get Fwd and Bwd traces
154 Array<OneD, NekDouble> Fwd{nTracePts, 0.0};
155 Array<OneD, NekDouble> Bwd{nTracePts, 0.0};
156 pFields[0]->GetFwdBwdTracePhys(hElePts, Fwd, Bwd);
157 // Fix boundaries
158 std::size_t cnt = 0;
159 std::size_t nBndRegions = pFields[0]->GetBndCondExpansions().size();
160 // Loop on the boundary regions
161 for (std::size_t i = 0; i < nBndRegions; ++i)
162 {
163 // Number of boundary regions related to region 'i'
164 std::size_t nBndEdges =
165 pFields[0]->GetBndCondExpansions()[i]->GetExpSize();
166
167 if (pFields[0]->GetBndConditions()[i]->GetBoundaryConditionType() ==
169 {
170 continue;
171 }
172
173 // Get value from interior
174 for (std::size_t e = 0; e < nBndEdges; ++e)
175 {
176 std::size_t nBndEdgePts = pFields[0]
177 ->GetBndCondExpansions()[i]
178 ->GetExp(e)
179 ->GetTotPoints();
180
181 std::size_t id2 = pFields[0]->GetTrace()->GetPhys_Offset(
182 pFields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
183
184 Vmath::Vcopy(nBndEdgePts, &Fwd[id2], 1, &Bwd[id2], 1);
185 }
186 }
187 // Get average of traces
188 Array<OneD, NekDouble> traceH{nTracePts, 1.0};
189 m_traceOneOverH = Array<OneD, NekDouble>{nTracePts, 1.0};
190 Vmath::Svtsvtp(nTracePts, 0.5, Fwd, 1, 0.5, Bwd, 1, traceH, 1);
191 // Multiply by coefficient = - C11 / h
194}
LibUtilities::SessionReaderSharedPtr m_session
Definition: DiffusionLDGNS.h:68
Array< OneD, Array< OneD, NekDouble > > m_traceVel
Definition: DiffusionLDGNS.h:66
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:143
std::shared_ptr< Expansion2D > Expansion2DSharedPtr
Definition: Expansion1D.h:46
std::shared_ptr< Expansion1D > Expansion1DSharedPtr
Definition: Expansion1D.h:50
std::shared_ptr< Expansion3D > Expansion3DSharedPtr
Definition: Expansion2D.h:47
EquationOfStateFactory & GetEquationOfStateFactory()
Declaration of the equation of state factory singleton.
Definition: EquationOfState.cpp:41
void Svtsvtp(int n, const T alpha, const T *x, int incx, const T beta, const T *y, int incy, T *z, int incz)
Svtsvtp (scalar times vector plus scalar times vector):
Definition: Vmath.hpp:473
void Sdiv(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha/x.
Definition: Vmath.hpp:154
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.hpp:54
References ASSERTL0, Nektar::LibUtilities::NekFactory< tKey, tBase, tParam >::CreateInstance(), Nektar::SpatialDomains::ePeriodic, Vmath::Fill(), Nektar::GetEquationOfStateFactory(), Nektar::LocalRegions::Expansion1D::GetGeom1D(), m_C11, m_diffDim, m_eos, m_session, m_spaceDim, m_traceNormals, m_traceOneOverH, m_traceVel, m_Twall, Vmath::Sdiv(), Vmath::Svtsvtp(), and Vmath::Vcopy().
◆ v_SetHomoDerivs()
|
inlineoverrideprotectedvirtual |
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 124 of file DiffusionLDGNS.h.
References m_homoDerivs.
Member Data Documentation
◆ m_C11
|
protected |
Penalty coefficient for LDGNS.
Definition at line 61 of file DiffusionLDGNS.h.
Referenced by v_InitObject().
◆ m_diffDim
|
protected |
Definition at line 79 of file DiffusionLDGNS.h.
Referenced by v_DiffuseCalcDerivative(), and v_InitObject().
◆ m_eos
|
protected |
Equation of system for computing temperature.
Definition at line 72 of file DiffusionLDGNS.h.
Referenced by ApplyBCsO1(), and v_InitObject().
◆ m_homoDerivs
Definition at line 76 of file DiffusionLDGNS.h.
Referenced by v_DiffuseCalcDerivative(), and v_SetHomoDerivs().
◆ m_session
|
protected |
Definition at line 68 of file DiffusionLDGNS.h.
Referenced by v_InitObject().
◆ m_spaceDim
|
protected |
Definition at line 78 of file DiffusionLDGNS.h.
Referenced by ApplyBCsO2(), NumericalFluxO1(), v_DiffuseCalcDerivative(), v_DiffuseCoeffs(), and v_InitObject().
◆ m_traceNormals
Definition at line 67 of file DiffusionLDGNS.h.
Referenced by ApplyBCsO1(), ApplyBCsO2(), NumericalFluxO1(), NumericalFluxO2(), and v_InitObject().
◆ m_traceOneOverH
h scaling for penalty term
Definition at line 64 of file DiffusionLDGNS.h.
Referenced by NumericalFluxO2(), and v_InitObject().
◆ m_traceVel
Definition at line 66 of file DiffusionLDGNS.h.
Referenced by v_InitObject().
◆ m_Twall
|
protected |
Definition at line 69 of file DiffusionLDGNS.h.
Referenced by ApplyBCsO1(), and v_InitObject().
◆ m_viscTensor
|
protected |
Definition at line 74 of file DiffusionLDGNS.h.
Referenced by v_DiffuseCoeffs(), and v_GetFluxTensor().
◆ type
|
static |
Initial value:
=
static DiffusionSharedPtr create(std::string diffType)
Definition: DiffusionLDGNS.h:50
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:197
DiffusionFactory & GetDiffusionFactory()
Definition: Diffusion.cpp:39
Definition at line 55 of file DiffusionLDGNS.h.
