Commit 80a0cde7 authored by Gianmarco Mengaldo's avatar Gianmarco Mengaldo
Browse files

Diffusion and Advection classes cleaned up. Removed DiffusionLFR

parent 4f5b0266
......@@ -53,28 +53,7 @@ namespace Nektar
{
v_InitObject(pSession, pFields);
}
void Advection::SetupMetrics(
const LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields)
{
v_SetupMetrics(pSession, pFields);
}
void Advection::SetupCFunctions(
const LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields)
{
v_SetupCFunctions(pSession, pFields);
}
void Advection::SetupInterpolationMatrices(
const LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields)
{
v_SetupInterpolationMatrices(pSession, pFields);
}
void Advection::Advect(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
......@@ -84,38 +63,5 @@ namespace Nektar
{
v_Advect(nConvectiveFields, fields, advVel, inarray, outarray);
}
void Advection::DivCFlux_1D(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, const NekDouble> &fluxX1,
const Array<OneD, const NekDouble> &numericalFlux,
Array<OneD, NekDouble> &divCFlux)
{
v_DivCFlux_1D(nConvectiveFields, fields, fluxX1, numericalFlux, divCFlux);
}
void Advection::DivCFlux_2D(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, const NekDouble> &fluxX1,
const Array<OneD, const NekDouble> &fluxX2,
const Array<OneD, const NekDouble> &numericalFlux,
Array<OneD, NekDouble> &divCFlux)
{
v_DivCFlux_2D(nConvectiveFields, fields, fluxX1, fluxX2, numericalFlux, divCFlux);
}
void Advection::DivCFlux_3D(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, const NekDouble> &fluxX1,
const Array<OneD, const NekDouble> &fluxX2,
const Array<OneD, const NekDouble> &fluxX3,
const Array<OneD, const NekDouble> &numericalFlux,
Array<OneD, NekDouble> &divCFlux)
{
v_DivCFlux_3D(nConvectiveFields, fields, fluxX1, fluxX2, fluxX3, numericalFlux, divCFlux);
}
}
}
......@@ -60,50 +60,14 @@ namespace Nektar
SOLVER_UTILS_EXPORT void InitObject(
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields);
SOLVER_UTILS_EXPORT void SetupMetrics(
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields);
SOLVER_UTILS_EXPORT void SetupCFunctions(
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields);
SOLVER_UTILS_EXPORT void SetupInterpolationMatrices(
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields);
SOLVER_UTILS_EXPORT void Advect(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &advVel,
const Array<OneD, Array<OneD, NekDouble> > &inarray,
Array<OneD, Array<OneD, NekDouble> > &outarray);
SOLVER_UTILS_EXPORT void DivCFlux_1D(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, const NekDouble> &fluxX1,
const Array<OneD, const NekDouble> &numericalFlux,
Array<OneD, NekDouble> &divCFlux);
SOLVER_UTILS_EXPORT void DivCFlux_2D(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, const NekDouble> &fluxX1,
const Array<OneD, const NekDouble> &fluxX2,
const Array<OneD, const NekDouble> &numericalFlux,
Array<OneD, NekDouble> &divCFlux);
SOLVER_UTILS_EXPORT void DivCFlux_3D(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, const NekDouble> &fluxX1,
const Array<OneD, const NekDouble> &fluxX2,
const Array<OneD, const NekDouble> &fluxX3,
const Array<OneD, const NekDouble> &numericalFlux,
Array<OneD, NekDouble> &divCFlux);
template<typename FuncPointerT, typename ObjectPointerT>
void SetFluxVector(FuncPointerT func, ObjectPointerT obj)
{
......@@ -122,28 +86,7 @@ namespace Nektar
{
};
virtual void v_SetupMetrics(
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields)
{
};
virtual void v_SetupCFunctions(
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields)
{
};
virtual void v_SetupInterpolationMatrices(
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields)
{
};
virtual void v_Advect(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
......@@ -151,39 +94,6 @@ namespace Nektar
const Array<OneD, Array<OneD, NekDouble> > &inarray,
Array<OneD, Array<OneD, NekDouble> > &outarray)=0;
virtual void v_DivCFlux_1D(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, const NekDouble> &fluxX1,
const Array<OneD, const NekDouble> &numericalFlux,
Array<OneD, NekDouble> &divCFlux)
{
};
virtual void v_DivCFlux_2D(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, const NekDouble> &fluxX1,
const Array<OneD, const NekDouble> &fluxX2,
const Array<OneD, const NekDouble> &numericalFlux,
Array<OneD, NekDouble> &divCFlux)
{
};
virtual void v_DivCFlux_3D(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, const NekDouble> &fluxX1,
const Array<OneD, const NekDouble> &fluxX2,
const Array<OneD, const NekDouble> &fluxX3,
const Array<OneD, const NekDouble> &numericalFlux,
Array<OneD, NekDouble> &divCFlux)
{
};
AdvectionFluxVecCB m_fluxVector;
RiemannSolverSharedPtr m_riemann;
};
......
......@@ -6,7 +6,6 @@ SET(SOLVER_UTILS_SOURCES
Diffusion/Diffusion.cpp
Diffusion/DiffusionLDG.cpp
Diffusion/DiffusionLDGNS.cpp
Diffusion/DiffusionLFR.cpp
Driver.cpp
DriverArnoldi.cpp
DriverModifiedArnoldi.cpp
......@@ -30,7 +29,6 @@ SET(SOLVER_UTILS_HEADERS
Diffusion/Diffusion.h
Diffusion/DiffusionLDG.h
Diffusion/DiffusionLDGNS.h
Diffusion/DiffusionLFR.h
Driver.h
DriverArnoldi.h
DriverModifiedArnoldi.h
......
......@@ -61,45 +61,5 @@ namespace Nektar
{
v_Diffuse(nConvectiveFields, fields, inarray, outarray);
}
void Diffusion::NumFluxforScalar(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &ufield,
Array<OneD, Array<OneD, Array<OneD, NekDouble> > >&uflux)
{
v_NumFluxforScalar(fields, ufield, uflux);
}
void Diffusion::WeakPenaltyforScalar(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const int var,
const Array<OneD, const NekDouble> &ufield,
Array<OneD, NekDouble> &penaltyflux,
NekDouble time)
{
v_WeakPenaltyforScalar(fields, var, ufield, penaltyflux, time);
}
void Diffusion::NumFluxforVector(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &ufield,
Array<OneD, Array<OneD, Array<OneD, NekDouble> > >&qfield,
Array<OneD, Array<OneD, NekDouble> > &qflux)
{
v_NumFluxforVector(fields, ufield, qfield, qflux);
}
void Diffusion::WeakPenaltyforVector(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const int var,
const int dir,
const Array<OneD, const NekDouble> &qfield,
Array<OneD, NekDouble> &penaltyflux,
NekDouble C11,
NekDouble time)
{
v_WeakPenaltyforVector(fields, var, dir, qfield,
penaltyflux, C11, time);
}
}
}
......@@ -76,33 +76,6 @@ namespace Nektar
const Array<OneD, Array<OneD, NekDouble> > &inarray,
Array<OneD, Array<OneD, NekDouble> > &outarray);
SOLVER_UTILS_EXPORT void NumFluxforScalar(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &ufield,
Array<OneD, Array<OneD, Array<OneD, NekDouble> > >&uflux);
SOLVER_UTILS_EXPORT void NumFluxforVector(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &ufield,
Array<OneD, Array<OneD, Array<OneD, NekDouble> > >&qfield,
Array<OneD, Array<OneD, NekDouble> > &qflux);
SOLVER_UTILS_EXPORT void WeakPenaltyforScalar(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const int var,
const Array<OneD, const NekDouble> &ufield,
Array<OneD, NekDouble> &penaltyflux,
NekDouble time);
SOLVER_UTILS_EXPORT void WeakPenaltyforVector(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const int var,
const int dir,
const Array<OneD, const NekDouble> &qfield,
Array<OneD, NekDouble> &penaltyflux,
NekDouble C11,
NekDouble time);
template<typename FuncPointerT, typename ObjectPointerT>
void SetFluxVector(FuncPointerT func, ObjectPointerT obj)
{
......@@ -132,45 +105,6 @@ namespace Nektar
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &inarray,
Array<OneD, Array<OneD, NekDouble> > &outarray)=0;
virtual void v_NumFluxforScalar(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &ufield,
Array<OneD, Array<OneD, Array<OneD, NekDouble> > >&uflux)
{
};
virtual void v_WeakPenaltyforScalar(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const int var,
const Array<OneD, const NekDouble> &ufield,
Array<OneD, NekDouble> &penaltyflux,
NekDouble time)
{
};
virtual void v_NumFluxforVector(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &ufield,
Array<OneD, Array<OneD, Array<OneD, NekDouble> > >&qfield,
Array<OneD, Array<OneD, NekDouble> > &qflux)
{
};
virtual void v_WeakPenaltyforVector(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const int var,
const int dir,
const Array<OneD, const NekDouble> &qfield,
Array<OneD, NekDouble> &penaltyflux,
NekDouble C11,
NekDouble time)
{
};
DiffusionFluxVecCB m_fluxVector;
DiffusionFluxVecCBNS m_fluxVectorNS;
......
......@@ -275,7 +275,9 @@ namespace Nektar
GetBoundaryConditionType() == SpatialDomains::eDirichlet)
{
Vmath::Vcopy(Nfps,
&BDphysics[id1], 1,
&(fields[var]->
GetBndCondExpansions()[i]->
UpdatePhys())[id1], 1,
&penaltyflux[id2], 1);
}
......@@ -375,7 +377,7 @@ namespace Nektar
// Imposing weak boundary condition with flux
if (fields[0]->GetBndCondExpansions().num_elements())
{
WeakPenaltyforVector(fields,
v_WeakPenaltyforVector(fields,
i, j,
qfield[j][i],
qfluxtemp,
......@@ -479,7 +481,9 @@ namespace Nektar
{
Vmath::Vmul(Nfps,
&m_traceNormals[dir][id2], 1,
&BDphysics[id1], 1,
&(fields[var]->
GetBndCondExpansions()[i]->
UpdatePhys())[id1], 1,
&penaltyflux[id2], 1);
}
}
......
This diff is collapsed.
///////////////////////////////////////////////////////////////////////////////
//
// File: DiffusionLFR.h
//
// For more information, please see: http://www.nektar.info
//
// The MIT License
//
// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
// Department of Aeronautics, Imperial College London (UK), and Scientific
// Computing and Imaging Institute, University of Utah (USA).
//
// License for the specific language governing rights and limitations under
// Permission is hereby granted, free of charge, to any person obtaining a
// copy of this software and associated documentation files (the "Software"),
// to deal in the Software without restriction, including without limitation
// the rights to use, copy, modify, merge, publish, distribute, sublicense,
// and/or sell copies of the Software, and to permit persons to whom the
// Software is furnished to do so, subject to the following conditions:
//
// The above copyright notice and this permission notice shall be included
// in all copies or substantial portions of the Software.
//
// THE SOFTWARE IS PROVIDED "AS IS", WITHOUzT WARRANTY OF ANY KIND, EXPRESS
// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
// DEALINGS IN THE SOFTWARE.
//
// Description: LDG diffusion class.
//
///////////////////////////////////////////////////////////////////////////////
#ifndef NEKTAR_SOLVERUTILS_DIFFUSIONLFR
#define NEKTAR_SOLVERUTILS_DIFFUSIONLFR
#include <SolverUtils/Diffusion/Diffusion.h>
namespace Nektar
{
namespace SolverUtils
{
class DiffusionLFR : public Diffusion
{
public:
static DiffusionSharedPtr create(std::string diffType)
{
return DiffusionSharedPtr(new DiffusionLFR());
}
static std::string type;
protected:
DiffusionLFR();
Array<OneD, Array<OneD, NekDouble> > m_traceNormals;
LibUtilities::SessionReaderSharedPtr m_session;
virtual void v_InitObject(
LibUtilities::SessionReaderSharedPtr pSession);
virtual void v_Diffuse(
const int nConvective,
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &inarray,
Array<OneD, Array<OneD, NekDouble> > &outarray);
virtual void v_NumFluxforScalar(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &ufield,
Array<OneD, Array<OneD, Array<OneD, NekDouble> > >&uflux);
virtual void v_WeakPenaltyforScalar(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const int var,
const Array<OneD, const NekDouble> &ufield,
Array<OneD, NekDouble> &penaltyflux,
NekDouble time);
virtual void v_NumFluxforVector(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &ufield,
Array<OneD, Array<OneD, Array<OneD, NekDouble> > >&qfield,
Array<OneD, Array<OneD, NekDouble> > &qflux);
virtual void v_WeakPenaltyforVector(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const int var,
const int dir,
const Array<OneD, const NekDouble> &qfield,
Array<OneD, NekDouble> &penaltyflux,
NekDouble C11,
NekDouble time);
};
}
}
#endif
......@@ -128,7 +128,7 @@ namespace Nektar
string advName;
string diffName;
string riemName;
m_session->LoadSolverInfo("AdvectionType", advName, "WeakDG");
m_session->LoadSolverInfo("DiffusionType", diffName, "LDGNS");
......@@ -148,6 +148,9 @@ namespace Nektar
m_riemannSolver = SolverUtils::GetRiemannSolverFactory()
.CreateInstance(riemName);
m_riemannSolverLDG = SolverUtils::GetRiemannSolverFactory()
.CreateInstance("Upwind");
m_riemannSolver->AddParam (
"gamma",
......@@ -159,9 +162,21 @@ namespace Nektar
"N",
&CompressibleFlowSystem::GetNormals, this);
m_advection->SetRiemannSolver(m_riemannSolver);
m_advection->InitObject (m_session, m_fields);
m_diffusion->InitObject (m_session);
m_riemannSolverLDG->AddParam (
"gamma",
&CompressibleFlowSystem::GetGamma, this);
m_riemannSolverLDG->AddScalar(
"velLoc",
&CompressibleFlowSystem::GetVelLoc, this);
m_riemannSolverLDG->AddVector(
"N",
&CompressibleFlowSystem::GetNormals, this);
m_advection->SetRiemannSolver (m_riemannSolver);
m_diffusion->SetRiemannSolver (m_riemannSolverLDG);
m_advection->InitObject (m_session, m_fields);
m_diffusion->InitObject (m_session);
break;
}
default:
......@@ -212,7 +227,7 @@ namespace Nektar
// user defined boundaries into account
int e, id1, id2, npts;
for(e = 0; e < m_fields[0]->
for (e = 0; e < m_fields[0]->
GetBndCondExpansions()[bcRegion]->GetExpSize(); ++e)
{
npts = m_fields[0]->GetBndCondExpansions()[bcRegion]->
......@@ -225,7 +240,7 @@ namespace Nektar
// For 2D/3D, define: v* = v - 2(v.n)n
Array<OneD,NekDouble> tmp(npts,0.0);
Array<OneD,NekDouble> tmp(npts, 0.0);
// Calculate (v.n)
for (i = 0; i < m_expdim; ++i)
......@@ -244,7 +259,6 @@ namespace Nektar
&tmp[0],1);
// Calculate v* = v - 2.0(v.n)n
for (i = 0; i < m_expdim; ++i)
{
Vmath::Vvtvp(npts,
......@@ -289,7 +303,7 @@ namespace Nektar
// user defined boundaries into account
int e, id1, id2, npts;
for(e = 0; e < m_fields[0]->
for (e = 0; e < m_fields[0]->
GetBndCondExpansions()[bcRegion]->GetExpSize(); ++e)
{
npts = m_fields[0]->GetBndCondExpansions()[bcRegion]->
......@@ -299,83 +313,10 @@ namespace Nektar
id2 = m_fields[0]->GetTrace()->GetPhys_Offset(
m_fields[0]->GetTraceMap()->
GetBndCondCoeffsToGlobalCoeffsMap(cnt+e));
/*
switch(m_expdim)
{
// Special case for 2D
case 2:
{
Array<OneD, NekDouble> tmp_n(npts);
Array<OneD, NekDouble> tmp_t(npts);
Vmath::Vmul(npts,
&Fwd[1][id2], 1,
&m_traceNormals[0][id2], 1,
&tmp_n[0], 1);
Vmath::Vvtvp(npts,
&Fwd[2][id2], 1,
&m_traceNormals[1][id2], 1,
&tmp_n[0], 1,
&tmp_n[0], 1);
Vmath::Vmul(npts,
&Fwd[1][id2], 1,
&m_traceNormals[1][id2], 1,
&tmp_t[0], 1);
Vmath::Vvtvm(npts,
&Fwd[2][id2], 1,
&m_traceNormals[0][id2], 1,
&tmp_t[0], 1,
&tmp_t[0], 1);
// negate fluxes
Vmath::Neg(npts, tmp_n, 1);
Vmath::Neg(npts, tmp_t, 1);
// rotate back to Cartesian
Vmath::Vmul(npts,
&tmp_t[0], 1,
&m_traceNormals[1][id2], 1,