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