#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()
|
protecteddefault |
Referenced by create().
Member Function Documentation
◆ ApplyBCsO1()
|
private |
Imposes appropriate bcs for the 1st order derivatives.
Definition at line 428 of file DiffusionLDGNS.cpp.
434{
435 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
436 std::size_t nScalars = inarray.size();
437
438 Array<OneD, NekDouble> tmp1{nTracePts, 0.0};
439 Array<OneD, NekDouble> tmp2{nTracePts, 0.0};
440 Array<OneD, NekDouble> Tw{nTracePts, m_Twall};
441
442 Array<OneD, Array<OneD, NekDouble>> scalarVariables{nScalars};
443
444 // Extract internal values of the scalar variables for Neumann bcs
445 for (std::size_t i = 0; i < nScalars; ++i)
446 {
447 scalarVariables[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
448 }
449
450 // Compute boundary conditions for velocities
451 for (std::size_t i = 0; i < nScalars - 1; ++i)
452 {
453 // Note that cnt has to loop on nBndRegions and nBndEdges
454 // and has to be reset to zero per each equation
455 std::size_t cnt = 0;
456 std::size_t nBndRegions = fields[i + 1]->GetBndCondExpansions().size();
457 for (std::size_t j = 0; j < nBndRegions; ++j)
458 {
459 if (fields[i + 1]
460 ->GetBndConditions()[j]
461 ->GetBoundaryConditionType() == SpatialDomains::ePeriodic)
462 {
463 continue;
464 }
465
466 std::size_t nBndEdges =
467 fields[i + 1]->GetBndCondExpansions()[j]->GetExpSize();
468 for (std::size_t e = 0; e < nBndEdges; ++e)
469 {
470 std::size_t nBndEdgePts = fields[i + 1]
471 ->GetBndCondExpansions()[j]
472 ->GetExp(e)
473 ->GetTotPoints();
474
475 std::size_t id1 =
476 fields[i + 1]->GetBndCondExpansions()[j]->GetPhys_Offset(e);
477
478 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
479 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(
480 cnt++));
481
482 if (boost::iequals(
483 fields[i]->GetBndConditions()[j]->GetUserDefined(),
484 "WallViscous") ||
485 boost::iequals(
486 fields[i]->GetBndConditions()[j]->GetUserDefined(),
487 "WallAdiabatic") ||
488 boost::iequals(
489 fields[i]->GetBndConditions()[j]->GetUserDefined(),
490 "WallRotational"))
491 {
492 // Reinforcing bcs for velocity in case of Wall bcs
493 for (int pt = 0; pt < nBndEdgePts; ++pt)
494 {
495 scalarVariables[i][id2 + pt] = m_gridVelocityTrace
496 [i][id2 + pt]; // If movement this equals trace grid
497 // velocity otherwise 0
498 }
499 }
500 else if (boost::iequals(
501 fields[i]->GetBndConditions()[j]->GetUserDefined(),
502 "Wall") ||
503 boost::iequals(
504 fields[i]->GetBndConditions()[j]->GetUserDefined(),
505 "Symmetry"))
506 {
507 // Symmetry bc: normal velocity is zero
508 // get all velocities at once because we need u.n
509 if (i == 0)
510 {
511 // tmp1 = -(u.n)
512 Vmath::Zero(nBndEdgePts, tmp1, 1);
513 for (std::size_t k = 0; k < nScalars - 1; ++k)
514 {
515 Vmath::Vdiv(nBndEdgePts,
516 &(fields[k + 1]
517 ->GetBndCondExpansions()[j]
518 ->UpdatePhys())[id1],
519 1,
520 &(fields[0]
521 ->GetBndCondExpansions()[j]
522 ->UpdatePhys())[id1],
523 1, &scalarVariables[k][id2], 1);
525 1, &scalarVariables[k][id2], 1,
526 &tmp1[0], 1, &tmp1[0], 1);
527 }
528 Vmath::Smul(nBndEdgePts, -1.0, &tmp1[0], 1, &tmp1[0],
529 1);
530
531 // u_i - (u.n)n_i
532 for (std::size_t k = 0; k < nScalars - 1; ++k)
533 {
534 Vmath::Vvtvp(nBndEdgePts, &tmp1[0], 1,
535 &m_traceNormals[k][id2], 1,
536 &scalarVariables[k][id2], 1,
537 &scalarVariables[k][id2], 1);
538 }
539 }
540 }
541 else if (fields[i]
542 ->GetBndConditions()[j]
543 ->GetBoundaryConditionType() ==
545 {
546 // Imposing velocity bcs if not Wall
547 Vmath::Vdiv(nBndEdgePts,
548 &(fields[i + 1]
549 ->GetBndCondExpansions()[j]
550 ->UpdatePhys())[id1],
551 1,
552 &(fields[0]
553 ->GetBndCondExpansions()[j]
554 ->UpdatePhys())[id1],
555 1, &scalarVariables[i][id2], 1);
556 }
557
558 // For Dirichlet boundary condition: uflux = u_bcs
559 if (fields[i]
560 ->GetBndConditions()[j]
561 ->GetBoundaryConditionType() ==
563 {
564 Vmath::Vcopy(nBndEdgePts, &scalarVariables[i][id2], 1,
565 &fluxO1[i][id2], 1);
566 }
567
568 // For Neumann boundary condition: uflux = u_+
569 else if ((fields[i]->GetBndConditions()[j])
570 ->GetBoundaryConditionType() ==
572 {
573 Vmath::Vcopy(nBndEdgePts, &pFwd[i][id2], 1, &fluxO1[i][id2],
574 1);
575 }
576
577 // Building kinetic energy to be used for T bcs
578 Vmath::Vmul(nBndEdgePts, &scalarVariables[i][id2], 1,
579 &scalarVariables[i][id2], 1, &tmp1[id2], 1);
580
581 Vmath::Smul(nBndEdgePts, 0.5, &tmp1[id2], 1, &tmp1[id2], 1);
582
583 Vmath::Vadd(nBndEdgePts, &tmp2[id2], 1, &tmp1[id2], 1,
584 &tmp2[id2], 1);
585 }
586 }
587 }
588
589 // Compute boundary conditions for temperature
590 std::size_t cnt = 0;
591 std::size_t nBndRegions = fields[nScalars]->GetBndCondExpansions().size();
592 for (std::size_t j = 0; j < nBndRegions; ++j)
593 {
594 if (fields[nScalars]
595 ->GetBndConditions()[j]
596 ->GetBoundaryConditionType() == SpatialDomains::ePeriodic)
597 {
598 continue;
599 }
600
601 std::size_t nBndEdges =
602 fields[nScalars]->GetBndCondExpansions()[j]->GetExpSize();
603 for (std::size_t e = 0; e < nBndEdges; ++e)
604 {
605 std::size_t nBndEdgePts = fields[nScalars]
606 ->GetBndCondExpansions()[j]
607 ->GetExp(e)
608 ->GetTotPoints();
609
610 std::size_t id1 =
611 fields[nScalars]->GetBndCondExpansions()[j]->GetPhys_Offset(e);
612
613 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
614 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
615
616 // Imposing Temperature Twall at the wall
617 if (boost::iequals(
618 fields[nScalars]->GetBndConditions()[j]->GetUserDefined(),
619 "WallViscous"))
620 {
621 Vmath::Vcopy(nBndEdgePts, &Tw[0], 1,
622 &scalarVariables[nScalars - 1][id2], 1);
623 }
624 // Imposing Temperature through condition on the Energy
625 // for no wall boundaries (e.g. farfield)
626 else if (fields[nScalars]
627 ->GetBndConditions()[j]
628 ->GetBoundaryConditionType() ==
630 {
631 // Use equation of state to evaluate temperature
632 NekDouble rho, ene;
633 for (std::size_t n = 0; n < nBndEdgePts; ++n)
634 {
635 rho = fields[0]
636 ->GetBndCondExpansions()[j]
637 ->GetPhys()[id1 + n];
638 ene = fields[nScalars]
639 ->GetBndCondExpansions()[j]
640 ->GetPhys()[id1 + n] /
641 rho -
642 tmp2[id2 + n];
643 scalarVariables[nScalars - 1][id2 + n] =
644 m_eos->GetTemperature(rho, ene);
645 }
646 }
647
648 // For Dirichlet boundary condition: uflux = u_bcs
649 if (fields[nScalars]
650 ->GetBndConditions()[j]
651 ->GetBoundaryConditionType() ==
652 SpatialDomains::eDirichlet &&
653 (!boost::iequals(
654 fields[nScalars]->GetBndConditions()[j]->GetUserDefined(),
655 "WallAdiabatic") ||
656 boost::iequals(
657 fields[nScalars]->GetBndConditions()[j]->GetUserDefined(),
658 "WallRotational")))
659 {
660 Vmath::Vcopy(nBndEdgePts, &scalarVariables[nScalars - 1][id2],
661 1, &fluxO1[nScalars - 1][id2], 1);
662 }
663
664 // For Neumann boundary condition: uflux = u_+
665 else if (((fields[nScalars]->GetBndConditions()[j])
666 ->GetBoundaryConditionType() ==
667 SpatialDomains::eNeumann) ||
668 boost::iequals(fields[nScalars]
669 ->GetBndConditions()[j]
670 ->GetUserDefined(),
671 "WallAdiabatic") ||
672 boost::iequals(fields[nScalars]
673 ->GetBndConditions()[j]
674 ->GetUserDefined(),
675 "WallRotational"))
676 {
677 Vmath::Vcopy(nBndEdgePts, &pFwd[nScalars - 1][id2], 1,
678 &fluxO1[nScalars - 1][id2], 1);
679 }
680 }
681 }
682}
EquationOfStateSharedPtr m_eos
Equation of system for computing temperature.
Definition: DiffusionLDGNS.h:70
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Definition: DiffusionLDGNS.h:65
Array< OneD, Array< OneD, NekDouble > > m_gridVelocityTrace
Definition: Diffusion.h:367
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, Nektar::SolverUtils::Diffusion::m_gridVelocityTrace, 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 751 of file DiffusionLDGNS.cpp.
758{
759 std::size_t cnt = 0;
760 std::size_t nBndRegions = fields[var]->GetBndCondExpansions().size();
761 // Loop on the boundary regions to apply appropriate bcs
762 for (std::size_t i = 0; i < nBndRegions; ++i)
763 {
764 // Number of boundary regions related to region 'i'
765 std::size_t nBndEdges =
766 fields[var]->GetBndCondExpansions()[i]->GetExpSize();
767
768 if (fields[var]->GetBndConditions()[i]->GetBoundaryConditionType() ==
770 {
771 continue;
772 }
773
774 // Weakly impose bcs by modifying flux values
775 for (std::size_t e = 0; e < nBndEdges; ++e)
776 {
777 std::size_t nBndEdgePts = fields[var]
778 ->GetBndCondExpansions()[i]
779 ->GetExp(e)
780 ->GetTotPoints();
781
782 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
783 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
784
785 // In case of Dirichlet bcs:
786 // uflux = gD
787 if (fields[var]
788 ->GetBndConditions()[i]
789 ->GetBoundaryConditionType() ==
790 SpatialDomains::eDirichlet &&
791 !boost::iequals(
792 fields[var]->GetBndConditions()[i]->GetUserDefined(),
793 "WallAdiabatic") &&
794 !boost::iequals(
795 fields[var]->GetBndConditions()[i]->GetUserDefined(),
796 "WallRotational"))
797 {
799 &qFwd[id2], 1, &penaltyflux[id2], 1);
800 }
801 // 3.4) In case of Neumann bcs:
802 // uflux = u+
803 else if ((fields[var]->GetBndConditions()[i])
804 ->GetBoundaryConditionType() ==
806 {
808 }
809 else if (boost::iequals(
810 fields[var]->GetBndConditions()[i]->GetUserDefined(),
811 "WallAdiabatic") ||
812 boost::iequals(
813 fields[var]->GetBndConditions()[i]->GetUserDefined(),
814 "WallRotational"))
815 {
817 {
818 Vmath::Zero(nBndEdgePts, &penaltyflux[id2], 1);
819 }
820 else
821 {
822
824 &qFwd[id2], 1, &penaltyflux[id2], 1);
825 }
826 }
827 }
828 }
829}
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 }
DiffusionLDGNS()=default
std::shared_ptr< Diffusion > DiffusionSharedPtr
A shared pointer to an EquationSystem object.
Definition: Diffusion.h:55
References DiffusionLDGNS().
◆ NumericalFluxO1()
|
private |
Builds the numerical flux for the 1st order derivatives.
Definition at line 390 of file DiffusionLDGNS.cpp.
396{
397 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
398 std::size_t nScalars = inarray.size();
399
400 // Upwind
401 Array<OneD, Array<OneD, NekDouble>> numflux{nScalars};
402 for (std::size_t i = 0; i < nScalars; ++i)
403 {
404 numflux[i] = {pFwd[i]};
405 }
406
407 // Modify the values in case of boundary interfaces
408 if (fields[0]->GetBndCondExpansions().size())
409 {
410 ApplyBCsO1(fields, inarray, pFwd, pBwd, numflux);
411 }
412
413 // Splitting the numerical flux into the dimensions
415 {
416 for (std::size_t i = 0; i < nScalars; ++i)
417 {
419 numericalFluxO1[j][i], 1);
420 }
421 }
422}
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:428
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 688 of file DiffusionLDGNS.cpp.
694{
695 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
696 std::size_t nVariables = fields.size();
697
698 // Initialize penalty flux
699 Array<OneD, Array<OneD, NekDouble>> fluxPen{nVariables - 1};
700 for (std::size_t i = 0; i < nVariables - 1; ++i)
701 {
702 fluxPen[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
703 }
704
705 // Get penalty flux
706 m_fluxPenaltyNS(uFwd, uBwd, fluxPen);
707
708 // Evaluate Riemann flux
709 // qflux = \hat{q} \cdot u = q \cdot n
710 // Notice: i = 1 (first row of the viscous tensor is zero)
711
712 Array<OneD, NekDouble> qFwd{nTracePts};
713 Array<OneD, NekDouble> qBwd{nTracePts};
714 Array<OneD, NekDouble> qtemp{nTracePts};
715 Array<OneD, NekDouble> qfluxtemp{nTracePts, 0.0};
716 std::size_t nDim = fields[0]->GetCoordim(0);
717 for (std::size_t i = 1; i < nVariables; ++i)
718 {
719 for (std::size_t j = 0; j < nDim; ++j)
720 {
721 // Compute qFwd and qBwd value of qfield in position 'ji'
722 fields[i]->GetFwdBwdTracePhys(qfield[j][i], qFwd, qBwd);
723
724 // Downwind
725 Vmath::Vcopy(nTracePts, qBwd, 1, qfluxtemp, 1);
726
727 // Multiply the Riemann flux by the trace normals
729 qfluxtemp, 1);
730
731 // Add penalty term
733 qfluxtemp, 1, qfluxtemp, 1);
734
735 // Impose weak boundary condition with flux
736 if (fields[0]->GetBndCondExpansions().size())
737 {
738 ApplyBCsO2(fields, i, j, qfield[j][i], qFwd, qBwd, qfluxtemp);
739 }
740
741 // Store the final flux into qflux
742 Vmath::Vadd(nTracePts, qfluxtemp, 1, qflux[i], 1, qflux[i], 1);
743 }
744 }
745}
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:751
DiffusionFluxPenaltyNS m_fluxPenaltyNS
Definition: Diffusion.h:362
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 208 of file DiffusionLDGNS.cpp.
215{
216 std::size_t nCoeffs = fields[0]->GetNcoeffs();
217 Array<OneD, Array<OneD, NekDouble>> tmp2{nConvectiveFields};
218 for (std::size_t i = 0; i < nConvectiveFields; ++i)
219 {
220 tmp2[i] = Array<OneD, NekDouble>{nCoeffs, 0.0};
221 }
222 v_DiffuseCoeffs(nConvectiveFields, fields, inarray, tmp2, pFwd, pBwd);
223 for (std::size_t i = 0; i < nConvectiveFields; ++i)
224 {
225 fields[i]->BwdTrans(tmp2[i], outarray[i]);
226 }
227}
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:229
References v_DiffuseCoeffs().
◆ v_DiffuseCalcDerivative()
|
overrideprotectedvirtual |
Diffusion Flux, calculate the physical derivatives.
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 313 of file DiffusionLDGNS.cpp.
319{
320 std::size_t nDim = fields[0]->GetCoordim(0);
321 std::size_t nCoeffs = fields[0]->GetNcoeffs();
322 std::size_t nScalars = inarray.size();
323 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
324 std::size_t nConvectiveFields = fields.size();
325
326 Array<OneD, NekDouble> tmp1{nCoeffs};
327 Array<OneD, Array<OneD, NekDouble>> tmp2{nConvectiveFields};
328 TensorOfArray3D<NekDouble> numericalFluxO1{m_spaceDim};
329
331 {
332 numericalFluxO1[j] = Array<OneD, Array<OneD, NekDouble>>{nScalars};
333
334 for (std::size_t i = 0; i < nScalars; ++i)
335 {
336 numericalFluxO1[j][i] = Array<OneD, NekDouble>{nTracePts, 0.0};
337 }
338 }
339
340 NumericalFluxO1(fields, inarray, numericalFluxO1, pFwd, pBwd);
341
342 for (std::size_t j = 0; j < nDim; ++j)
343 {
344 for (std::size_t i = 0; i < nScalars; ++i)
345 {
346 fields[i]->IProductWRTDerivBase(j, inarray[i], tmp1);
347 Vmath::Neg(nCoeffs, tmp1, 1);
348 fields[i]->AddTraceIntegral(numericalFluxO1[j][i], tmp1);
349 fields[i]->SetPhysState(false);
350 fields[i]->MultiplyByElmtInvMass(tmp1, tmp1);
351 fields[i]->BwdTrans(tmp1, qfields[j][i]);
352 }
353 }
354 // For 3D Homogeneous 1D only take derivatives in 3rd direction
356 {
357 for (std::size_t i = 0; i < nScalars; ++i)
358 {
359 qfields[2][i] = m_homoDerivs[i];
360 }
361 }
362}
Array< OneD, Array< OneD, NekDouble > > m_homoDerivs
Definition: DiffusionLDGNS.h:74
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:390
References m_diffDim, m_homoDerivs, m_spaceDim, Vmath::Neg(), and NumericalFluxO1().
◆ v_DiffuseCoeffs()
|
overrideprotectedvirtual |
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 229 of file DiffusionLDGNS.cpp.
236{
237
238 if (fields[0]->GetGraph()->GetMovement()->GetMoveFlag()) // i.e. if
239 // m_ALESolver
240 {
241 fields[0]->GetTrace()->GetNormals(m_traceNormals);
242 }
243
244 std::size_t nDim = fields[0]->GetCoordim(0);
245 std::size_t nPts = fields[0]->GetTotPoints();
246 std::size_t nCoeffs = fields[0]->GetNcoeffs();
247 std::size_t nScalars = inarray.size();
248 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
249
250 TensorOfArray3D<NekDouble> derivativesO1{m_spaceDim};
251
253 {
254 derivativesO1[j] = Array<OneD, Array<OneD, NekDouble>>{nScalars};
255
256 for (std::size_t i = 0; i < nScalars; ++i)
257 {
258 derivativesO1[j][i] = Array<OneD, NekDouble>{nPts, 0.0};
259 }
260 }
261
262 DiffuseCalcDerivative(fields, inarray, derivativesO1, pFwd, pBwd);
263
264 // Initialisation viscous tensor
266 Array<OneD, Array<OneD, NekDouble>> viscousFlux{nConvectiveFields};
267
269 {
270 m_viscTensor[j] = Array<OneD, Array<OneD, NekDouble>>{nScalars + 1};
271 for (std::size_t i = 0; i < nScalars + 1; ++i)
272 {
273 m_viscTensor[j][i] = Array<OneD, NekDouble>{nPts, 0.0};
274 }
275 }
276
277 for (std::size_t i = 0; i < nConvectiveFields; ++i)
278 {
279 viscousFlux[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
280 }
281
283
284 // Compute u from q_{\eta} and q_{\xi}
285 // Obtain numerical fluxes
286 // Note: derivativesO1 is not used in DiffuseTraceFlux
288 pFwd, pBwd);
289 Array<OneD, NekDouble> tmpOut{nCoeffs};
290 Array<OneD, Array<OneD, NekDouble>> tmpIn{nDim};
291
292 for (std::size_t i = 0; i < nConvectiveFields; ++i)
293 {
294 // Temporary fix to call collection op
295 // we can reorder m_viscTensor later
296 for (std::size_t j = 0; j < nDim; ++j)
297 {
298 tmpIn[j] = m_viscTensor[j][i];
299 }
300 fields[i]->IProductWRTDerivBase(tmpIn, tmpOut);
301
302 // Evaulate <\phi, \hat{F}\cdot n> - outarray[i]
303 Vmath::Neg(nCoeffs, tmpOut, 1);
304 fields[i]->AddTraceIntegral(viscousFlux[i], tmpOut);
305 fields[i]->SetPhysState(false);
306 if (!fields[0]->GetGraph()->GetMovement()->GetMoveFlag())
307 {
308 fields[i]->MultiplyByElmtInvMass(tmpOut, outarray[i]);
309 }
310 }
311}
TensorOfArray3D< NekDouble > m_viscTensor
Definition: DiffusionLDGNS.h:72
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:209
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:222
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:233
References Nektar::SolverUtils::Diffusion::DiffuseCalcDerivative(), Nektar::SolverUtils::Diffusion::DiffuseTraceFlux(), Nektar::SolverUtils::Diffusion::DiffuseVolumeFlux(), m_spaceDim, m_traceNormals, m_viscTensor, and Vmath::Neg().
Referenced by v_Diffuse().
◆ v_DiffuseTraceFlux()
|
overrideprotectedvirtual |
Diffusion term Trace Flux.
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 373 of file DiffusionLDGNS.cpp.
382{
383 NumericalFluxO2(fields, qfields, TraceFlux, pFwd, pBwd);
384}
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:688
References NumericalFluxO2().
◆ v_DiffuseVolumeFlux()
|
overrideprotectedvirtual |
Diffusion Volume Flux.
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 364 of file DiffusionLDGNS.cpp.
369{
370 m_fluxVectorNS(inarray, qfields, VolumeFlux);
371}
DiffusionFluxVecCBNS m_fluxVectorNS
Definition: Diffusion.h:361
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 49 of file DiffusionLDGNS.cpp.
52{
53 m_session = pSession;
55
56 // Setting up the normals
57 std::size_t nDim = pFields[0]->GetCoordim(0);
58 std::size_t nTracePts = pFields[0]->GetTrace()->GetTotPoints();
59
60 m_spaceDim = nDim;
61 if (pSession->DefinesSolverInfo("HOMOGENEOUS"))
62 {
63 m_spaceDim = 3;
64 }
65
67
68 m_traceVel = Array<OneD, Array<OneD, NekDouble>>{
73 {
74 m_traceVel[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
75 m_traceNormals[i] = Array<OneD, NekDouble>{nTracePts};
76 m_gridVelocityTrace[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
77 }
78 pFields[0]->GetTrace()->GetNormals(m_traceNormals);
79
80 // Create equation of state object
81 std::string eosType;
84
85 // Set up {h} reference on the trace for penalty term
86 //
87 // Note, this shold be replaced with something smarter when merging
88 // LDG with IP
89
90 // Get min h per element
91 std::size_t nElements = pFields[0]->GetExpSize();
92 Array<OneD, NekDouble> hEle{nElements, 1.0};
93 for (std::size_t e = 0; e < nElements; e++)
94 {
95 NekDouble h{1.0e+10};
96 std::size_t expDim = pFields[0]->GetShapeDimension();
97 switch (expDim)
98 {
99 case 3:
100 {
101 LocalRegions::Expansion3DSharedPtr exp3D;
102 exp3D = pFields[0]->GetExp(e)->as<LocalRegions::Expansion3D>();
103 for (std::size_t i = 0; i < exp3D->GetNtraces(); ++i)
104 {
105 h = std::min(
106 h, exp3D->GetGeom3D()->GetEdge(i)->GetVertex(0)->dist(*(
107 exp3D->GetGeom3D()->GetEdge(i)->GetVertex(1))));
108 }
109 break;
110 }
111
112 case 2:
113 {
114 LocalRegions::Expansion2DSharedPtr exp2D;
115 exp2D = pFields[0]->GetExp(e)->as<LocalRegions::Expansion2D>();
116 for (std::size_t i = 0; i < exp2D->GetNtraces(); ++i)
117 {
118 h = std::min(
119 h, exp2D->GetGeom2D()->GetEdge(i)->GetVertex(0)->dist(*(
120 exp2D->GetGeom2D()->GetEdge(i)->GetVertex(1))));
121 }
122 break;
123 }
124 case 1:
125 {
126 LocalRegions::Expansion1DSharedPtr exp1D;
127 exp1D = pFields[0]->GetExp(e)->as<LocalRegions::Expansion1D>();
128
129 h = std::min(h, exp1D->GetGeom1D()->GetVertex(0)->dist(
130 *(exp1D->GetGeom1D()->GetVertex(1))));
131 break;
132 }
133 default:
134 {
136 }
137 }
138
139 // Store scaling
140 hEle[e] = h;
141 }
142 // Expand h from elements to points
143 std::size_t nPts = pFields[0]->GetTotPoints();
144 Array<OneD, NekDouble> hElePts{nPts, 0.0};
145 Array<OneD, NekDouble> tmp;
146 for (std::size_t e = 0; e < pFields[0]->GetExpSize(); e++)
147 {
148 std::size_t nElmtPoints = pFields[0]->GetExp(e)->GetTotPoints();
149 std::size_t physOffset = pFields[0]->GetPhys_Offset(e);
150 Vmath::Fill(nElmtPoints, hEle[e], tmp = hElePts + physOffset, 1);
151 }
152 // Get Fwd and Bwd traces
153 Array<OneD, NekDouble> Fwd{nTracePts, 0.0};
154 Array<OneD, NekDouble> Bwd{nTracePts, 0.0};
155 pFields[0]->GetFwdBwdTracePhys(hElePts, Fwd, Bwd);
156 // Fix boundaries
157 std::size_t cnt = 0;
158 std::size_t nBndRegions = pFields[0]->GetBndCondExpansions().size();
159 // Loop on the boundary regions
160 for (std::size_t i = 0; i < nBndRegions; ++i)
161 {
162 // Number of boundary regions related to region 'i'
163 std::size_t nBndEdges =
164 pFields[0]->GetBndCondExpansions()[i]->GetExpSize();
165
166 if (pFields[0]->GetBndConditions()[i]->GetBoundaryConditionType() ==
168 {
169 continue;
170 }
171
172 // Get value from interior
173 for (std::size_t e = 0; e < nBndEdges; ++e)
174 {
175 std::size_t nBndEdgePts = pFields[0]
176 ->GetBndCondExpansions()[i]
177 ->GetExp(e)
178 ->GetTotPoints();
179
180 std::size_t id2 = pFields[0]->GetTrace()->GetPhys_Offset(
181 pFields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
182
183 Vmath::Vcopy(nBndEdgePts, &Fwd[id2], 1, &Bwd[id2], 1);
184 }
185 }
186 // Get average of traces
187 Array<OneD, NekDouble> traceH{nTracePts, 1.0};
188 m_traceOneOverH = Array<OneD, NekDouble>{nTracePts, 1.0};
189 Vmath::Svtsvtp(nTracePts, 0.5, Fwd, 1, 0.5, Bwd, 1, traceH, 1);
190 // Multiply by coefficient = - C11 / h
193}
LibUtilities::SessionReaderSharedPtr m_session
Definition: DiffusionLDGNS.h:66
Array< OneD, Array< OneD, NekDouble > > m_traceVel
Definition: DiffusionLDGNS.h:64
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
Definition: BasicUtils/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:39
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, Nektar::SolverUtils::Diffusion::m_gridVelocityTrace, 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 59 of file DiffusionLDGNS.h.
Referenced by v_InitObject().
◆ m_diffDim
|
protected |
Definition at line 77 of file DiffusionLDGNS.h.
Referenced by v_DiffuseCalcDerivative(), and v_InitObject().
◆ m_eos
|
protected |
Equation of system for computing temperature.
Definition at line 70 of file DiffusionLDGNS.h.
Referenced by ApplyBCsO1(), and v_InitObject().
◆ m_homoDerivs
Definition at line 74 of file DiffusionLDGNS.h.
Referenced by v_DiffuseCalcDerivative(), and v_SetHomoDerivs().
◆ m_session
|
protected |
Definition at line 66 of file DiffusionLDGNS.h.
Referenced by v_InitObject().
◆ m_spaceDim
|
protected |
Definition at line 76 of file DiffusionLDGNS.h.
Referenced by ApplyBCsO2(), NumericalFluxO1(), v_DiffuseCalcDerivative(), v_DiffuseCoeffs(), and v_InitObject().
◆ m_traceNormals
Definition at line 65 of file DiffusionLDGNS.h.
Referenced by ApplyBCsO1(), ApplyBCsO2(), NumericalFluxO1(), NumericalFluxO2(), v_DiffuseCoeffs(), and v_InitObject().
◆ m_traceOneOverH
h scaling for penalty term
Definition at line 62 of file DiffusionLDGNS.h.
Referenced by NumericalFluxO2(), and v_InitObject().
◆ m_traceVel
Definition at line 64 of file DiffusionLDGNS.h.
Referenced by v_InitObject().
◆ m_Twall
|
protected |
Definition at line 67 of file DiffusionLDGNS.h.
Referenced by ApplyBCsO1(), and v_InitObject().
◆ m_viscTensor
|
protected |
Definition at line 72 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: BasicUtils/NekFactory.hpp:197
DiffusionFactory & GetDiffusionFactory()
Definition: Diffusion.cpp:39
Definition at line 55 of file DiffusionLDGNS.h.
