#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. | |
| 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. | |
| 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. | |
| 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. | |
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 430 of file DiffusionLDGNS.cpp.
436{
437 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
438 std::size_t nScalars = inarray.size();
439
440 Array<OneD, NekDouble> tmp1{nTracePts, 0.0};
441 Array<OneD, NekDouble> tmp2{nTracePts, 0.0};
442 Array<OneD, NekDouble> Tw{nTracePts, m_Twall};
443
444 Array<OneD, Array<OneD, NekDouble>> scalarVariables{nScalars};
445
446 // Extract internal values of the scalar variables for Neumann bcs
447 for (std::size_t i = 0; i < nScalars; ++i)
448 {
449 scalarVariables[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
450 }
451
452 // Compute boundary conditions for velocities
453 for (std::size_t i = 0; i < nScalars - 1; ++i)
454 {
455 // Note that cnt has to loop on nBndRegions and nBndEdges
456 // and has to be reset to zero per each equation
457 std::size_t cnt = 0;
458 std::size_t nBndRegions = fields[i + 1]->GetBndCondExpansions().size();
459 for (std::size_t j = 0; j < nBndRegions; ++j)
460 {
461 if (fields[i + 1]
462 ->GetBndConditions()[j]
463 ->GetBoundaryConditionType() == SpatialDomains::ePeriodic)
464 {
465 continue;
466 }
467
468 std::size_t nBndEdges =
469 fields[i + 1]->GetBndCondExpansions()[j]->GetExpSize();
470 for (std::size_t e = 0; e < nBndEdges; ++e)
471 {
472 std::size_t nBndEdgePts = fields[i + 1]
473 ->GetBndCondExpansions()[j]
474 ->GetExp(e)
475 ->GetTotPoints();
476
477 std::size_t id1 =
478 fields[i + 1]->GetBndCondExpansions()[j]->GetPhys_Offset(e);
479
480 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
481 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(
482 cnt++));
483
484 if (boost::iequals(
485 fields[i]->GetBndConditions()[j]->GetUserDefined(),
486 "WallViscous") ||
487 boost::iequals(
488 fields[i]->GetBndConditions()[j]->GetUserDefined(),
489 "WallAdiabatic") ||
490 boost::iequals(
491 fields[i]->GetBndConditions()[j]->GetUserDefined(),
492 "WallRotational"))
493 {
494 // Reinforcing bcs for velocity in case of Wall bcs
495 for (int pt = 0; pt < nBndEdgePts; ++pt)
496 {
497 scalarVariables[i][id2 + pt] =
498 m_gridVelocityTrace[i][id2 + pt];
499 // If movement this equals trace grid
500 // velocity otherwise 0
501 }
502 }
503 else if (boost::iequals(
504 fields[i]->GetBndConditions()[j]->GetUserDefined(),
505 "Wall") ||
506 boost::iequals(
507 fields[i]->GetBndConditions()[j]->GetUserDefined(),
508 "Symmetry"))
509 {
510 // Symmetry bc: normal velocity is zero
511 // get all velocities at once because we need u.n
512 if (i == 0)
513 {
514 // tmp1 = -(u.n)
515 Vmath::Zero(nBndEdgePts, tmp1, 1);
516 for (std::size_t k = 0; k < nScalars - 1; ++k)
517 {
518 Vmath::Vdiv(nBndEdgePts,
519 &(fields[k + 1]
520 ->GetBndCondExpansions()[j]
521 ->UpdatePhys())[id1],
522 1,
523 &(fields[0]
524 ->GetBndCondExpansions()[j]
525 ->UpdatePhys())[id1],
526 1, &scalarVariables[k][id2], 1);
528 1, &scalarVariables[k][id2], 1,
529 &tmp1[0], 1, &tmp1[0], 1);
530 }
531 Vmath::Smul(nBndEdgePts, -1.0, &tmp1[0], 1, &tmp1[0],
532 1);
533
534 // u_i - (u.n)n_i
535 for (std::size_t k = 0; k < nScalars - 1; ++k)
536 {
537 Vmath::Vvtvp(nBndEdgePts, &tmp1[0], 1,
538 &m_traceNormals[k][id2], 1,
539 &scalarVariables[k][id2], 1,
540 &scalarVariables[k][id2], 1);
541 }
542 }
543 }
544 else if (fields[i]
545 ->GetBndConditions()[j]
546 ->GetBoundaryConditionType() ==
548 {
549 // Imposing velocity bcs if not Wall
550 Vmath::Vdiv(nBndEdgePts,
551 &(fields[i + 1]
552 ->GetBndCondExpansions()[j]
553 ->UpdatePhys())[id1],
554 1,
555 &(fields[0]
556 ->GetBndCondExpansions()[j]
557 ->UpdatePhys())[id1],
558 1, &scalarVariables[i][id2], 1);
559 }
560
561 // For Dirichlet boundary condition: uflux = u_bcs
562 if (fields[i]
563 ->GetBndConditions()[j]
564 ->GetBoundaryConditionType() ==
566 {
567 Vmath::Vcopy(nBndEdgePts, &scalarVariables[i][id2], 1,
568 &fluxO1[i][id2], 1);
569 }
570
571 // For Neumann boundary condition: uflux = u_+
572 else if ((fields[i]->GetBndConditions()[j])
573 ->GetBoundaryConditionType() ==
575 {
576 Vmath::Vcopy(nBndEdgePts, &pFwd[i][id2], 1, &fluxO1[i][id2],
577 1);
578 }
579
580 // Building kinetic energy to be used for T bcs
581 Vmath::Vmul(nBndEdgePts, &scalarVariables[i][id2], 1,
582 &scalarVariables[i][id2], 1, &tmp1[id2], 1);
583
584 Vmath::Smul(nBndEdgePts, 0.5, &tmp1[id2], 1, &tmp1[id2], 1);
585
586 Vmath::Vadd(nBndEdgePts, &tmp2[id2], 1, &tmp1[id2], 1,
587 &tmp2[id2], 1);
588 }
589 }
590 }
591
592 // Compute boundary conditions for temperature
593 std::size_t cnt = 0;
594 std::size_t nBndRegions = fields[nScalars]->GetBndCondExpansions().size();
595 for (std::size_t j = 0; j < nBndRegions; ++j)
596 {
597 if (fields[nScalars]
598 ->GetBndConditions()[j]
599 ->GetBoundaryConditionType() == SpatialDomains::ePeriodic)
600 {
601 continue;
602 }
603
604 std::size_t nBndEdges =
605 fields[nScalars]->GetBndCondExpansions()[j]->GetExpSize();
606 for (std::size_t e = 0; e < nBndEdges; ++e)
607 {
608 std::size_t nBndEdgePts = fields[nScalars]
609 ->GetBndCondExpansions()[j]
610 ->GetExp(e)
611 ->GetTotPoints();
612
613 std::size_t id1 =
614 fields[nScalars]->GetBndCondExpansions()[j]->GetPhys_Offset(e);
615
616 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
617 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
618
619 // Imposing Temperature Twall at the wall
620 if (boost::iequals(
621 fields[nScalars]->GetBndConditions()[j]->GetUserDefined(),
622 "WallViscous"))
623 {
624 Vmath::Vcopy(nBndEdgePts, &Tw[0], 1,
625 &scalarVariables[nScalars - 1][id2], 1);
626 }
627 // Imposing Temperature through condition on the Energy
628 // for no wall boundaries (e.g. farfield)
629 else if (fields[nScalars]
630 ->GetBndConditions()[j]
631 ->GetBoundaryConditionType() ==
633 {
634 // Use equation of state to evaluate temperature
635 NekDouble rho, ene;
636 for (std::size_t n = 0; n < nBndEdgePts; ++n)
637 {
638 rho = fields[0]
639 ->GetBndCondExpansions()[j]
640 ->GetPhys()[id1 + n];
641 ene = fields[nScalars]
642 ->GetBndCondExpansions()[j]
643 ->GetPhys()[id1 + n] /
644 rho -
645 tmp2[id2 + n];
646 scalarVariables[nScalars - 1][id2 + n] =
647 m_eos->GetTemperature(rho, ene);
648 }
649 }
650
651 // For Dirichlet boundary condition: uflux = u_bcs
652 if (fields[nScalars]
653 ->GetBndConditions()[j]
654 ->GetBoundaryConditionType() ==
655 SpatialDomains::eDirichlet &&
656 (!boost::iequals(
657 fields[nScalars]->GetBndConditions()[j]->GetUserDefined(),
658 "WallAdiabatic") ||
659 boost::iequals(
660 fields[nScalars]->GetBndConditions()[j]->GetUserDefined(),
661 "WallRotational")))
662 {
663 Vmath::Vcopy(nBndEdgePts, &scalarVariables[nScalars - 1][id2],
664 1, &fluxO1[nScalars - 1][id2], 1);
665 }
666
667 // For Neumann boundary condition: uflux = u_+
668 else if (((fields[nScalars]->GetBndConditions()[j])
669 ->GetBoundaryConditionType() ==
670 SpatialDomains::eNeumann) ||
671 boost::iequals(fields[nScalars]
672 ->GetBndConditions()[j]
673 ->GetUserDefined(),
674 "WallAdiabatic") ||
675 boost::iequals(fields[nScalars]
676 ->GetBndConditions()[j]
677 ->GetUserDefined(),
678 "WallRotational"))
679 {
680 Vmath::Vcopy(nBndEdgePts, &pFwd[nScalars - 1][id2], 1,
681 &fluxO1[nScalars - 1][id2], 1);
682 }
683 }
684 }
685}
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 754 of file DiffusionLDGNS.cpp.
761{
762 std::size_t cnt = 0;
763 std::size_t nBndRegions = fields[var]->GetBndCondExpansions().size();
764 // Loop on the boundary regions to apply appropriate bcs
765 for (std::size_t i = 0; i < nBndRegions; ++i)
766 {
767 // Number of boundary regions related to region 'i'
768 std::size_t nBndEdges =
769 fields[var]->GetBndCondExpansions()[i]->GetExpSize();
770
771 if (fields[var]->GetBndConditions()[i]->GetBoundaryConditionType() ==
773 {
774 continue;
775 }
776
777 // Weakly impose bcs by modifying flux values
778 for (std::size_t e = 0; e < nBndEdges; ++e)
779 {
780 std::size_t nBndEdgePts = fields[var]
781 ->GetBndCondExpansions()[i]
782 ->GetExp(e)
783 ->GetTotPoints();
784
785 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
786 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
787
788 // In case of Dirichlet bcs:
789 // uflux = gD
790 if (fields[var]
791 ->GetBndConditions()[i]
792 ->GetBoundaryConditionType() ==
793 SpatialDomains::eDirichlet &&
794 !boost::iequals(
795 fields[var]->GetBndConditions()[i]->GetUserDefined(),
796 "WallAdiabatic") &&
797 !boost::iequals(
798 fields[var]->GetBndConditions()[i]->GetUserDefined(),
799 "WallRotational"))
800 {
802 &qFwd[id2], 1, &penaltyflux[id2], 1);
803 }
804 // 3.4) In case of Neumann bcs:
805 // uflux = u+
806 else if ((fields[var]->GetBndConditions()[i])
807 ->GetBoundaryConditionType() ==
809 {
811 }
812 else if (boost::iequals(
813 fields[var]->GetBndConditions()[i]->GetUserDefined(),
814 "WallAdiabatic") ||
815 boost::iequals(
816 fields[var]->GetBndConditions()[i]->GetUserDefined(),
817 "WallRotational"))
818 {
820 {
821 Vmath::Zero(nBndEdgePts, &penaltyflux[id2], 1);
822 }
823 else
824 {
825
827 &qFwd[id2], 1, &penaltyflux[id2], 1);
828 }
829 }
830 }
831 }
832}
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 392 of file DiffusionLDGNS.cpp.
398{
399 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
400 std::size_t nScalars = inarray.size();
401
402 // Upwind
403 Array<OneD, Array<OneD, NekDouble>> numflux{nScalars};
404 for (std::size_t i = 0; i < nScalars; ++i)
405 {
406 numflux[i] = {pFwd[i]};
407 }
408
409 // Modify the values in case of boundary interfaces
410 if (fields[0]->GetBndCondExpansions().size())
411 {
412 ApplyBCsO1(fields, inarray, pFwd, pBwd, numflux);
413 }
414
415 // Splitting the numerical flux into the dimensions
417 {
418 for (std::size_t i = 0; i < nScalars; ++i)
419 {
421 numericalFluxO1[j][i], 1);
422 }
423 }
424}
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:430
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 691 of file DiffusionLDGNS.cpp.
697{
698 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
699 std::size_t nVariables = fields.size();
700
701 // Initialize penalty flux
702 Array<OneD, Array<OneD, NekDouble>> fluxPen{nVariables - 1};
703 for (std::size_t i = 0; i < nVariables - 1; ++i)
704 {
705 fluxPen[i] = Array<OneD, NekDouble>{nTracePts, 0.0};
706 }
707
708 // Get penalty flux
709 m_fluxPenaltyNS(uFwd, uBwd, fluxPen);
710
711 // Evaluate Riemann flux
712 // qflux = \hat{q} \cdot u = q \cdot n
713 // Notice: i = 1 (first row of the viscous tensor is zero)
714
715 Array<OneD, NekDouble> qFwd{nTracePts};
716 Array<OneD, NekDouble> qBwd{nTracePts};
717 Array<OneD, NekDouble> qtemp{nTracePts};
718 Array<OneD, NekDouble> qfluxtemp{nTracePts, 0.0};
719 std::size_t nDim = fields[0]->GetCoordim(0);
720 for (std::size_t i = 1; i < nVariables; ++i)
721 {
722 for (std::size_t j = 0; j < nDim; ++j)
723 {
724 // Compute qFwd and qBwd value of qfield in position 'ji'
725 fields[i]->GetFwdBwdTracePhys(qfield[j][i], qFwd, qBwd);
726
727 // Downwind
728 Vmath::Vcopy(nTracePts, qBwd, 1, qfluxtemp, 1);
729
730 // Multiply the Riemann flux by the trace normals
732 qfluxtemp, 1);
733
734 // Add penalty term
736 qfluxtemp, 1, qfluxtemp, 1);
737
738 // Impose weak boundary condition with flux
739 if (fields[0]->GetBndCondExpansions().size())
740 {
741 ApplyBCsO2(fields, i, j, qfield[j][i], qFwd, qBwd, qfluxtemp);
742 }
743
744 // Store the final flux into qflux
745 Vmath::Vadd(nTracePts, qfluxtemp, 1, qflux[i], 1, qflux[i], 1);
746 }
747 }
748}
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:754
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 315 of file DiffusionLDGNS.cpp.
321{
322 std::size_t nDim = fields[0]->GetCoordim(0);
323 std::size_t nCoeffs = fields[0]->GetNcoeffs();
324 std::size_t nScalars = inarray.size();
325 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
326 std::size_t nConvectiveFields = fields.size();
327
328 Array<OneD, NekDouble> tmp1{nCoeffs};
329 Array<OneD, Array<OneD, NekDouble>> tmp2{nConvectiveFields};
330 TensorOfArray3D<NekDouble> numericalFluxO1{m_spaceDim};
331
333 {
334 numericalFluxO1[j] = Array<OneD, Array<OneD, NekDouble>>{nScalars};
335
336 for (std::size_t i = 0; i < nScalars; ++i)
337 {
338 numericalFluxO1[j][i] = Array<OneD, NekDouble>{nTracePts, 0.0};
339 }
340 }
341
342 NumericalFluxO1(fields, inarray, numericalFluxO1, pFwd, pBwd);
343
344 for (std::size_t j = 0; j < nDim; ++j)
345 {
346 for (std::size_t i = 0; i < nScalars; ++i)
347 {
348 fields[i]->IProductWRTDerivBase(j, inarray[i], tmp1);
349 Vmath::Neg(nCoeffs, tmp1, 1);
350 fields[i]->AddTraceIntegral(numericalFluxO1[j][i], tmp1);
351 fields[i]->SetPhysState(false);
352 fields[i]->MultiplyByElmtInvMass(tmp1, tmp1);
353 fields[i]->BwdTrans(tmp1, qfields[j][i]);
354 }
355 }
356 // For 3D Homogeneous 1D only take derivatives in 3rd direction
358 {
359 for (std::size_t i = 0; i < nScalars; ++i)
360 {
361 qfields[2][i] = m_homoDerivs[i];
362 }
363 }
364}
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:392
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
307 // If mesh is not distorted
308 if (!fields[0]->GetGraph()->GetMovement()->GetMeshDistortedFlag())
309 {
310 fields[i]->MultiplyByElmtInvMass(tmpOut, outarray[i]);
311 }
312 }
313}
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 375 of file DiffusionLDGNS.cpp.
384{
385 NumericalFluxO2(fields, qfields, TraceFlux, pFwd, pBwd);
386}
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:691
References NumericalFluxO2().
◆ v_DiffuseVolumeFlux()
|
overrideprotectedvirtual |
Diffusion Volume Flux.
Reimplemented from Nektar::SolverUtils::Diffusion.
Definition at line 366 of file DiffusionLDGNS.cpp.
371{
372 m_fluxVectorNS(inarray, qfields, VolumeFlux);
373}
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:142
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::PointGeom::dist(), Nektar::SpatialDomains::ePeriodic, Vmath::Fill(), Nektar::GetEquationOfStateFactory(), Nektar::LocalRegions::Expansion1D::GetGeom1D(), Nektar::SpatialDomains::Geometry::GetVertex(), 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:196
DiffusionFactory & GetDiffusionFactory()
Definition Diffusion.cpp:39
Definition at line 55 of file DiffusionLDGNS.h.
