#include <ForcingMovingReferenceFrame.h>

Public Member Functions
	Newmark_BetaSolver ()

	~Newmark_BetaSolver ()

void	SetNewmarkBeta (NekDouble beta, NekDouble gamma, NekDouble dt, DNekMatSharedPtr M, DNekMatSharedPtr C, DNekMatSharedPtr K, std::set< int > DirDoFs)

void	Solve (Array< OneD, Array< OneD, NekDouble > > u, Array< OneD, NekDouble > force, std::map< int, NekDouble > motionPrescribed)

Public Attributes
int	m_rows

Array< OneD, NekDouble >	m_coeffs

DNekMatSharedPtr	m_M

DNekMatSharedPtr	m_C

DNekMatSharedPtr	m_K

DNekMatSharedPtr	m_coeffMatrix

DNekMatSharedPtr	m_inverseMatrix

Detailed Description

Definition at line 65 of file ForcingMovingReferenceFrame.h.

Constructor & Destructor Documentation

◆ Newmark_BetaSolver()

Nektar::SolverUtils::Newmark_BetaSolver::Newmark_BetaSolver ( )

inline

Definition at line 68 of file ForcingMovingReferenceFrame.h.

68{};

◆ ~Newmark_BetaSolver()

Nektar::SolverUtils::Newmark_BetaSolver::~Newmark_BetaSolver ( )

inline

Definition at line 69 of file ForcingMovingReferenceFrame.h.

69{};

Member Function Documentation

◆ SetNewmarkBeta()

void Nektar::SolverUtils::Newmark_BetaSolver::SetNewmarkBeta	(	NekDouble	beta,
		NekDouble	gamma,
		NekDouble	dt,
		DNekMatSharedPtr	M,
		DNekMatSharedPtr	C,
		DNekMatSharedPtr	K,
		std::set< int >	DirDoFs
	)

Definition at line 1011 of file ForcingMovingReferenceFrame.cpp.

{
    m_coeffs    = Array<OneD, NekDouble>(5, 0.);
    m_coeffs[0] = 1. / (gamma * dt);
    m_coeffs[1] = 1. / gamma - 1.;
    m_coeffs[2] = beta * dt / gamma;
    m_coeffs[3] = dt * (1. - beta / gamma);
    m_coeffs[4] = (0.5 - beta / gamma) * dt * dt;
 
    m_rows = M->GetRows();
    m_coeffMatrix =
        MemoryManager<DNekMat>::AllocateSharedPtr(m_rows, m_rows, 0.0, eFULL);
    m_inverseMatrix =
        MemoryManager<DNekMat>::AllocateSharedPtr(m_rows, m_rows, 0.0, eFULL);
    for (int i = 0; i < m_rows; ++i)
    {
        for (int j = 0; j < m_rows; ++j)
        {
            NekDouble value = m_coeffs[0] * M->GetValue(i, j) +
                              C->GetValue(i, j) +
                              m_coeffs[2] * K->GetValue(i, j);
            m_coeffMatrix->SetValue(i, j, value);
            if (DirDoFs.find(i) != DirDoFs.end() ||
                DirDoFs.find(j) != DirDoFs.end())
            {
                value = (i == j) ? 1. : 0.;
            }
            m_inverseMatrix->SetValue(i, j, value);
        }
    }
    m_inverseMatrix->Invert();
    m_M = M;
    m_C = C;
    m_K = K;
}

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::LibUtilities::beta, Nektar::eFULL, m_C, m_coeffMatrix, m_coeffs, m_inverseMatrix, m_K, m_M, and m_rows.

Referenced by Nektar::SolverUtils::ForcingMovingReferenceFrame::InitBodySolver().

◆ Solve()

void Nektar::SolverUtils::Newmark_BetaSolver::Solve	(	Array< OneD, Array< OneD, NekDouble > >	u,
		Array< OneD, NekDouble >	force,
		std::map< int, NekDouble >	motionPrescribed
	)

Definition at line 1050 of file ForcingMovingReferenceFrame.cpp.

{
    NekVector<NekDouble> bm(m_rows);
    for (int i = 0; i < m_rows; ++i)
    {
        bm[i] = m_coeffs[0] * u[1][i] + m_coeffs[1] * u[2][i];
    }
    NekVector<NekDouble> fbm(m_rows);
    Multiply(fbm, *m_M, bm);
    NekVector<NekDouble> bk(m_rows);
    for (int i = 0; i < m_rows; ++i)
    {
        bk[i] = u[0][i] + m_coeffs[3] * u[1][i] + m_coeffs[4] * u[2][i];
    }
    NekVector<NekDouble> fbk(m_rows);
    Multiply(fbk, *m_K, bk);
    NekVector<NekDouble> rhs(m_rows);
    for (int i = 0; i < m_rows; ++i)
    {
        rhs[i] = force[i] - fbk[i] + fbm[i];
    }
    // apply Dirichelet DoFs
    for (int i = 0; i < 3; ++i)
    {
        if (motionPrescribed.find(i) != motionPrescribed.end())
        {
            rhs[i] = motionPrescribed[i];
        }
        else
        {
            for (auto it : motionPrescribed)
            {
                rhs[i] -= it.second * m_coeffMatrix->GetValue(it.first, i);
            }
        }
    }
    // solve
    NekVector<NekDouble> b(m_rows);
    Multiply(b, *m_inverseMatrix, rhs);
    for (int i = 0; i < m_rows; ++i)
    {
        u[1][i] = b[i];
        u[0][i] = m_coeffs[2] * u[1][i] + bk[i];
        u[2][i] = m_coeffs[0] * u[1][i] - bm[i];
    }
}