Nektar++
Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Private Member Functions | List of all members
Nektar::DiffusionLDGNS Class Reference

#include <DiffusionLDGNS.h>

Inheritance diagram for Nektar::DiffusionLDGNS:
[legend]

Static Public Member Functions

static DiffusionSharedPtr create (std::string diffType)
 

Static Public Attributes

static std::string type
 

Protected Member Functions

 DiffusionLDGNS ()=default
 
void v_InitObject (LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields) override
 
void v_Diffuse (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
 Calculate weak DG Diffusion in the LDG form for the Navier-Stokes (NS) equations: More...
 
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
 
void v_DiffuseCalcDerivative (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd) override
 Diffusion Flux, calculate the physical derivatives. More...
 
void v_DiffuseVolumeFlux (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, TensorOfArray3D< NekDouble > &VolumeFlux, Array< OneD, int > &nonZeroIndex) override
 Diffusion Volume Flux. More...
 
void v_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, const Array< OneD, Array< OneD, NekDouble > > &pBwd, Array< OneD, int > &nonZeroIndex) override
 Diffusion term Trace Flux. More...
 
void v_SetHomoDerivs (Array< OneD, Array< OneD, NekDouble > > &deriv) override
 
TensorOfArray3D< NekDouble > & v_GetFluxTensor () override
 
- Protected Member Functions inherited from Nektar::SolverUtils::Diffusion
virtual SOLVER_UTILS_EXPORT ~Diffusion ()=default
 
virtual SOLVER_UTILS_EXPORT void v_InitObject (LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)=0
 
virtual SOLVER_UTILS_EXPORT void v_Diffuse (const std::size_t nConvectiveFields, 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)=0
 
virtual SOLVER_UTILS_EXPORT void v_DiffuseCoeffs (const std::size_t nConvectiveFields, 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)
 
virtual SOLVER_UTILS_EXPORT void v_DiffuseCalcDerivative (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
 Diffusion Flux, calculate the physical derivatives. More...
 
virtual SOLVER_UTILS_EXPORT void v_DiffuseVolumeFlux (const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, TensorOfArray3D< NekDouble > &VolumeFlux, Array< OneD, int > &nonZeroIndex)
 Diffusion Volume Flux. More...
 
virtual SOLVER_UTILS_EXPORT void v_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, const Array< OneD, Array< OneD, NekDouble > > &pBwd, Array< OneD, int > &nonZeroIndex)
 Diffusion term Trace Flux. More...
 
virtual void v_SetHomoDerivs (Array< OneD, Array< OneD, NekDouble > > &deriv)
 
virtual TensorOfArray3D< NekDouble > & v_GetFluxTensor ()
 
virtual SOLVER_UTILS_EXPORT const Array< OneD, const Array< OneD, NekDouble > > & v_GetTraceNormal ()
 

Protected Attributes

NekDouble m_C11
 Penalty coefficient for LDGNS. More...
 
Array< OneD, NekDoublem_traceOneOverH
 h scaling for penalty term More...
 
Array< OneD, Array< OneD, NekDouble > > m_traceVel
 
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
 
LibUtilities::SessionReaderSharedPtr m_session
 
NekDouble m_Twall
 
EquationOfStateSharedPtr m_eos
 Equation of system for computing temperature. More...
 
TensorOfArray3D< NekDoublem_viscTensor
 
Array< OneD, Array< OneD, NekDouble > > m_homoDerivs
 
std::size_t m_spaceDim
 
std::size_t m_diffDim
 
- Protected Attributes inherited from Nektar::SolverUtils::Diffusion
Array< OneD, NekDoublem_divVel
 Params for Ducros sensor. More...
 
Array< OneD, NekDoublem_divVelSquare
 
Array< OneD, NekDoublem_curlVelSquare
 
DiffusionFluxVecCB m_fluxVector
 
DiffusionFluxVecCBNS m_fluxVectorNS
 
DiffusionFluxPenaltyNS m_fluxPenaltyNS
 
DiffusionFluxCons m_FunctorDiffusionfluxCons
 
DiffusionFluxCons m_FunctorDiffusionfluxConsTrace
 
SpecialBndTreat m_SpecialBndTreat
 
DiffusionSymmFluxCons m_FunctorSymmetricfluxCons
 
Array< OneD, Array< OneD, NekDouble > > m_gridVelocityTrace
 
NekDouble m_time = 0.0
 

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...
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::SolverUtils::Diffusion
SOLVER_UTILS_EXPORT void InitObject (LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
 
SOLVER_UTILS_EXPORT void Diffuse (const std::size_t nConvectiveFields, 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=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble > > &pBwd=NullNekDoubleArrayOfArray)
 
SOLVER_UTILS_EXPORT void Diffuse (const std::size_t nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, NekDouble time, const Array< OneD, Array< OneD, NekDouble > > &pFwd=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble > > &pBwd=NullNekDoubleArrayOfArray)
 
SOLVER_UTILS_EXPORT void DiffuseCoeffs (const std::size_t nConvectiveFields, 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=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble > > &pBwd=NullNekDoubleArrayOfArray)
 
SOLVER_UTILS_EXPORT void DiffuseCoeffs (const std::size_t nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, NekDouble time, const Array< OneD, Array< OneD, NekDouble > > &pFwd=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble > > &pBwd=NullNekDoubleArrayOfArray)
 
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)
 
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. More...
 
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. More...
 
SOLVER_UTILS_EXPORT void SetHomoDerivs (Array< OneD, Array< OneD, NekDouble > > &deriv)
 
SOLVER_UTILS_EXPORT TensorOfArray3D< NekDouble > & GetFluxTensor ()
 
SOLVER_UTILS_EXPORT const Array< OneD, const Array< OneD, NekDouble > > & GetTraceNormal ()
 Get trace normal. More...
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetFluxVector (FuncPointerT func, ObjectPointerT obj)
 
SOLVER_UTILS_EXPORT void SetFluxVector (DiffusionFluxVecCB fluxVector)
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetFluxVectorNS (FuncPointerT func, ObjectPointerT obj)
 
void SetFluxVectorNS (DiffusionFluxVecCBNS fluxVector)
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetFluxPenaltyNS (FuncPointerT func, ObjectPointerT obj)
 
void SetFluxPenaltyNS (DiffusionFluxPenaltyNS flux)
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetDiffusionFluxCons (FuncPointerT func, ObjectPointerT obj)
 
void SetDiffusionFluxCons (DiffusionFluxCons flux)
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetDiffusionFluxConsTrace (FuncPointerT func, ObjectPointerT obj)
 
void SetDiffusionFluxConsTrace (DiffusionFluxCons flux)
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetSpecialBndTreat (FuncPointerT func, ObjectPointerT obj)
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetDiffusionSymmFluxCons (FuncPointerT func, ObjectPointerT obj)
 
SOLVER_UTILS_EXPORT void SetGridVelocityTrace (Array< OneD, Array< OneD, NekDouble > > &gridVelocityTrace)
 

Detailed Description

Definition at line 47 of file DiffusionLDGNS.h.

Constructor & Destructor Documentation

◆ DiffusionLDGNS()

Nektar::DiffusionLDGNS::DiffusionLDGNS ( )
protecteddefault

Referenced by create().

Member Function Documentation

◆ ApplyBCsO1()

void Nektar::DiffusionLDGNS::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 
)
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);
524 Vmath::Vvtvp(nBndEdgePts, &m_traceNormals[k][id2],
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() ==
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() ==
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.
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Array< OneD, Array< OneD, NekDouble > > m_gridVelocityTrace
Definition: Diffusion.h:367
double NekDouble
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 Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273
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()

void Nektar::DiffusionLDGNS::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 
)
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() ==
791 !boost::iequals(
792 fields[var]->GetBndConditions()[i]->GetUserDefined(),
793 "WallAdiabatic") &&
794 !boost::iequals(
795 fields[var]->GetBndConditions()[i]->GetUserDefined(),
796 "WallRotational"))
797 {
798 Vmath::Vmul(nBndEdgePts, &m_traceNormals[dir][id2], 1,
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 {
807 ASSERTL0(false, "Neumann bcs not implemented for LDGNS");
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 {
816 if ((var == m_spaceDim + 1))
817 {
818 Vmath::Zero(nBndEdgePts, &penaltyflux[id2], 1);
819 }
820 else
821 {
822
823 Vmath::Vmul(nBndEdgePts, &m_traceNormals[dir][id2], 1,
824 &qFwd[id2], 1, &penaltyflux[id2], 1);
825 }
826 }
827 }
828 }
829}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208

References ASSERTL0, Nektar::SpatialDomains::eDirichlet, Nektar::SpatialDomains::eNeumann, Nektar::SpatialDomains::ePeriodic, m_spaceDim, m_traceNormals, Vmath::Vmul(), and Vmath::Zero().

Referenced by NumericalFluxO2().

◆ create()

static DiffusionSharedPtr Nektar::DiffusionLDGNS::create ( std::string  diffType)
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()

void Nektar::DiffusionLDGNS::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 
)
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
414 for (std::size_t j = 0; j < m_spaceDim; ++j)
415 {
416 for (std::size_t i = 0; i < nScalars; ++i)
417 {
418 Vmath::Vmul(nTracePts, m_traceNormals[j], 1, numflux[i], 1,
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.

References ApplyBCsO1(), m_spaceDim, m_traceNormals, and Vmath::Vmul().

Referenced by v_DiffuseCalcDerivative().

◆ NumericalFluxO2()

void Nektar::DiffusionLDGNS::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 
)
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
728 Vmath::Vmul(nTracePts, m_traceNormals[j], 1, qfluxtemp, 1,
729 qfluxtemp, 1);
730
731 // Add penalty term
732 Vmath::Vvtvp(nTracePts, m_traceOneOverH, 1, fluxPen[i - 1], 1,
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.
Array< OneD, NekDouble > m_traceOneOverH
h scaling for penalty term
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()

void Nektar::DiffusionLDGNS::v_Diffuse ( const std::size_t  nConvectiveFields,
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 
)
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

References v_DiffuseCoeffs().

◆ v_DiffuseCalcDerivative()

void Nektar::DiffusionLDGNS::v_DiffuseCalcDerivative ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  fields,
const Array< OneD, Array< OneD, NekDouble > > &  inarray,
TensorOfArray3D< NekDouble > &  qfields,
const Array< OneD, Array< OneD, NekDouble > > &  pFwd,
const Array< OneD, Array< OneD, NekDouble > > &  pBwd 
)
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
330 for (std::size_t j = 0; j < m_spaceDim; ++j)
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
355 if (m_diffDim == 1)
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
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 Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.hpp:292

References m_diffDim, m_homoDerivs, m_spaceDim, Vmath::Neg(), and NumericalFluxO1().

◆ v_DiffuseCoeffs()

void Nektar::DiffusionLDGNS::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 
)
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
252 for (std::size_t j = 0; j < m_spaceDim; ++j)
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
265 m_viscTensor = TensorOfArray3D<NekDouble>{m_spaceDim};
266 Array<OneD, Array<OneD, NekDouble>> viscousFlux{nConvectiveFields};
267
268 for (std::size_t j = 0; j < m_spaceDim; ++j)
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
282 DiffuseVolumeFlux(fields, inarray, derivativesO1, m_viscTensor);
283
284 // Compute u from q_{\eta} and q_{\xi}
285 // Obtain numerical fluxes
286 // Note: derivativesO1 is not used in DiffuseTraceFlux
287 DiffuseTraceFlux(fields, inarray, m_viscTensor, derivativesO1, viscousFlux,
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
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()

void Nektar::DiffusionLDGNS::v_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,
const Array< OneD, Array< OneD, NekDouble > > &  pBwd,
Array< OneD, int > &  nonZeroIndex 
)
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.

References NumericalFluxO2().

◆ v_DiffuseVolumeFlux()

void Nektar::DiffusionLDGNS::v_DiffuseVolumeFlux ( const Array< OneD, MultiRegions::ExpListSharedPtr > &  fields,
const Array< OneD, Array< OneD, NekDouble > > &  inarray,
TensorOfArray3D< NekDouble > &  qfields,
TensorOfArray3D< NekDouble > &  VolumeFlux,
Array< OneD, int > &  nonZeroIndex 
)
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()

TensorOfArray3D< NekDouble > & Nektar::DiffusionLDGNS::v_GetFluxTensor ( )
inlineoverrideprotectedvirtual

Reimplemented from Nektar::SolverUtils::Diffusion.

Definition at line 129 of file DiffusionLDGNS.h.

130 {
131 return m_viscTensor;
132 }

References m_viscTensor.

◆ v_InitObject()

void Nektar::DiffusionLDGNS::v_InitObject ( LibUtilities::SessionReaderSharedPtr  pSession,
Array< OneD, MultiRegions::ExpListSharedPtr pFields 
)
overrideprotectedvirtual

Implements Nektar::SolverUtils::Diffusion.

Definition at line 49 of file DiffusionLDGNS.cpp.

52{
53 m_session = pSession;
54 m_session->LoadParameter("Twall", m_Twall, 300.15);
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
66 m_diffDim = m_spaceDim - nDim;
67
68 m_traceVel = Array<OneD, Array<OneD, NekDouble>>{
69 m_spaceDim}; // @TODO: What is this trace vel used for???
70 m_traceNormals = Array<OneD, Array<OneD, NekDouble>>{m_spaceDim};
71 m_gridVelocityTrace = Array<OneD, Array<OneD, NekDouble>>{m_spaceDim};
72 for (std::size_t i = 0; i < m_spaceDim; ++i)
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;
82 m_session->LoadSolverInfo("EquationOfState", eosType, "IdealGas");
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 {
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 {
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 {
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 {
135 ASSERTL0(false, "Dimension out of bound.")
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
191 m_session->LoadParameter("LDGNSc11", m_C11, 1.0);
192 Vmath::Sdiv(nTracePts, -m_C11, traceH, 1, m_traceOneOverH, 1);
193}
LibUtilities::SessionReaderSharedPtr m_session
Array< OneD, Array< OneD, NekDouble > > m_traceVel
NekDouble m_C11
Penalty coefficient for LDGNS.
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
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.
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()

void Nektar::DiffusionLDGNS::v_SetHomoDerivs ( Array< OneD, Array< OneD, NekDouble > > &  deriv)
inlineoverrideprotectedvirtual

Reimplemented from Nektar::SolverUtils::Diffusion.

Definition at line 124 of file DiffusionLDGNS.h.

125 {
126 m_homoDerivs = deriv;
127 }

References m_homoDerivs.

Member Data Documentation

◆ m_C11

NekDouble Nektar::DiffusionLDGNS::m_C11
protected

Penalty coefficient for LDGNS.

Definition at line 59 of file DiffusionLDGNS.h.

Referenced by v_InitObject().

◆ m_diffDim

std::size_t Nektar::DiffusionLDGNS::m_diffDim
protected

Definition at line 77 of file DiffusionLDGNS.h.

Referenced by v_DiffuseCalcDerivative(), and v_InitObject().

◆ m_eos

EquationOfStateSharedPtr Nektar::DiffusionLDGNS::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

Array<OneD, Array<OneD, NekDouble> > Nektar::DiffusionLDGNS::m_homoDerivs
protected

Definition at line 74 of file DiffusionLDGNS.h.

Referenced by v_DiffuseCalcDerivative(), and v_SetHomoDerivs().

◆ m_session

LibUtilities::SessionReaderSharedPtr Nektar::DiffusionLDGNS::m_session
protected

Definition at line 66 of file DiffusionLDGNS.h.

Referenced by v_InitObject().

◆ m_spaceDim

std::size_t Nektar::DiffusionLDGNS::m_spaceDim
protected

◆ m_traceNormals

Array<OneD, Array<OneD, NekDouble> > Nektar::DiffusionLDGNS::m_traceNormals
protected

◆ m_traceOneOverH

Array<OneD, NekDouble> Nektar::DiffusionLDGNS::m_traceOneOverH
protected

h scaling for penalty term

Definition at line 62 of file DiffusionLDGNS.h.

Referenced by NumericalFluxO2(), and v_InitObject().

◆ m_traceVel

Array<OneD, Array<OneD, NekDouble> > Nektar::DiffusionLDGNS::m_traceVel
protected

Definition at line 64 of file DiffusionLDGNS.h.

Referenced by v_InitObject().

◆ m_Twall

NekDouble Nektar::DiffusionLDGNS::m_Twall
protected

Definition at line 67 of file DiffusionLDGNS.h.

Referenced by ApplyBCsO1(), and v_InitObject().

◆ m_viscTensor

TensorOfArray3D<NekDouble> Nektar::DiffusionLDGNS::m_viscTensor
protected

Definition at line 72 of file DiffusionLDGNS.h.

Referenced by v_DiffuseCoeffs(), and v_GetFluxTensor().

◆ type

std::string Nektar::DiffusionLDGNS::type
static
Initial value:
=
"LDGNS", DiffusionLDGNS::create, "LDG for Navier-Stokes")
static DiffusionSharedPtr create(std::string diffType)
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
DiffusionFactory & GetDiffusionFactory()
Definition: Diffusion.cpp:39

Definition at line 55 of file DiffusionLDGNS.h.