Commit 261d2b88 authored by Gianmarco Mengaldo's avatar Gianmarco Mengaldo

DiffusionLFR added. Some other minor modifications.

parent 707b0f9d
......@@ -2086,7 +2086,7 @@ namespace Nektar
}
}
/// Extract the data in fielddata into the coeffs
void ExpList::ExtractDataToCoeffs(
SpatialDomains::FieldDefinitionsSharedPtr &fielddef,
......@@ -2102,7 +2102,14 @@ namespace Nektar
v_ExtractCoeffsToCoeffs(fromExpList,fromCoeffs,toCoeffs);
}
void ExpList::v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field, Array<OneD, NekDouble> &coeffs)
/// NEW METHOD BY DAVE =================================================
void ExpList::v_ExtractDataToCoeffs(
SpatialDomains::FieldDefinitionsSharedPtr &fielddef,
std::vector<NekDouble> &fielddata, std::string &field,
Array<OneD, NekDouble> &coeffs)
{
int i, cnt;
int offset = 0;
......@@ -2144,8 +2151,19 @@ namespace Nektar
int n = vComm->GetSize();
int p = vComm->GetRank();
// Determine number of elements in fielddef inside this process.
for(cnt = i = 0; i < fielddef->m_elementIDs.size(); ++i)
{
if (ElmtID_to_ExpID.count(fielddef->m_elementIDs[i]) == 0)
{
continue;
}
++cnt;
}
Array<OneD, int> numEls(n, 0);
numEls[p] = m_exp->size();
numEls[p] = cnt;
vComm->AllReduce(numEls, LibUtilities::ReduceSum);
int totEls = Vmath::Vsum(n, numEls, 1);
......@@ -2156,8 +2174,8 @@ namespace Nektar
elOffsets[i] = elOffsets[i-1] + numEls[i-1];
}
Array<OneD, int> coeffsPerEl (totEls, 0.0);
Array<OneD, int> elmtGlobalIds(totEls, 0.0);
Array<OneD, int> coeffsPerEl (totEls, 0);
Array<OneD, int> elmtGlobalIds(totEls, 0);
// Determine number of coefficients in each local (to this
// partition) element.
......@@ -2171,7 +2189,7 @@ namespace Nektar
int eid = ElmtID_to_ExpID[fielddef->m_elementIDs[i]];
int datalen = (*m_exp)[eid]->CalcNumberOfCoefficients(
fielddef->m_numModes,modes_offset);
fielddef->m_numModes,modes_offset);
if(fielddef->m_uniOrder == true)
{
modes_offset = 0;
......@@ -2189,6 +2207,8 @@ namespace Nektar
for (i = 0; i < totEls; ++i)
{
ASSERTL0(coeffsElmtMap.count(elmtGlobalIds[i]) == 0,
"Error in communicating global ids!");
coeffsElmtMap[elmtGlobalIds[i]] = coeffsPerEl[i];
}
......@@ -2218,16 +2238,24 @@ namespace Nektar
else // unpack data to new order
{
(*m_exp)[eid]->ExtractDataToCoeffs(
&fielddata[offset], fielddef->m_numModes,
modes_offset, &coeffs[m_coeff_offset[eid]]);
&fielddata[offset], fielddef->m_numModes,
modes_offset, &coeffs[m_coeff_offset[eid]]);
}
offset += datalen;
}
}
}
/// END NEW METHOD BY DAVE =============================================
/*
void ExpList::v_ExtractDataToCoeffs(SpatialDomains::FieldDefinitionsSharedPtr &fielddef, std::vector<NekDouble> &fielddata, std::string &field, Array<OneD, NekDouble> &coeffs)
/// OLD METHOD NOT WORKING IN PARALLEL =================================
void ExpList::v_ExtractDataToCoeffs(
SpatialDomains::FieldDefinitionsSharedPtr &fielddef,
std::vector<NekDouble> &fielddata, std::string &field,
Array<OneD, NekDouble> &coeffs)
{
int i;
int offset = 0;
......@@ -2287,7 +2315,8 @@ namespace Nektar
}
}
}
*/
/// END OLD METHOD NOT WORKING IN PARALLEL =============================
*/
void ExpList::v_ExtractCoeffsToCoeffs(const boost::shared_ptr<ExpList> &fromExpList, const Array<OneD, const NekDouble> &fromCoeffs, Array<OneD, NekDouble> &toCoeffs)
{
......
......@@ -1201,7 +1201,7 @@ namespace Nektar
* direction X2.
* @param numericalFlux Riemann flux in the physical space.
* @param divCFlux Divergence of the corrective flux for 2D
* Problems.
* problems.
*
* \todo: Switch on shapes eventually here.
*/
......
......@@ -6,6 +6,7 @@ SET(SOLVER_UTILS_SOURCES
Diffusion/Diffusion.cpp
Diffusion/DiffusionLDG.cpp
Diffusion/DiffusionLDGNS.cpp
Diffusion/DiffusionLFR.cpp
Driver.cpp
DriverArnoldi.cpp
DriverModifiedArnoldi.cpp
......@@ -30,6 +31,7 @@ SET(SOLVER_UTILS_HEADERS
Diffusion/Diffusion.h
Diffusion/DiffusionLDG.h
Diffusion/DiffusionLDGNS.h
Diffusion/DiffusionLFR.h
Driver.h
DriverArnoldi.h
DriverModifiedArnoldi.h
......
......@@ -48,9 +48,10 @@ namespace Nektar
}
void Diffusion::InitObject(
const LibUtilities::SessionReaderSharedPtr pSession)
const LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields)
{
v_InitObject(pSession);
v_InitObject(pSession, pFields);
}
void Diffusion::Diffuse(
......
......@@ -67,7 +67,8 @@ namespace Nektar
{
public:
SOLVER_UTILS_EXPORT void InitObject(
LibUtilities::SessionReaderSharedPtr pSession);
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields);
SOLVER_UTILS_EXPORT void Diffuse(
const int nConvectiveFields,
......@@ -94,7 +95,8 @@ namespace Nektar
protected:
virtual void v_InitObject(
LibUtilities::SessionReaderSharedPtr pSession)
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields)
{
};
......
......@@ -60,7 +60,8 @@ namespace Nektar
LibUtilities::SessionReaderSharedPtr m_session;
virtual void v_InitObject(
LibUtilities::SessionReaderSharedPtr pSession);
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields);
virtual void v_Diffuse(
const int nConvective,
......@@ -77,8 +78,7 @@ namespace Nektar
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const int var,
const Array<OneD, const NekDouble> &ufield,
Array<OneD, NekDouble> &penaltyflux,
NekDouble time);
Array<OneD, NekDouble> &penaltyflux);
virtual void v_NumFluxforVector(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
......@@ -92,8 +92,7 @@ namespace Nektar
const int dir,
const Array<OneD, const NekDouble> &qfield,
Array<OneD, NekDouble> &penaltyflux,
NekDouble C11,
NekDouble time);
NekDouble C11);
};
}
}
......
......@@ -47,7 +47,8 @@ namespace Nektar
}
void DiffusionLDGNS::v_InitObject(
LibUtilities::SessionReaderSharedPtr pSession)
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields)
{
m_session = pSession;
}
......
......@@ -59,7 +59,8 @@ namespace Nektar
LibUtilities::SessionReaderSharedPtr m_session;
virtual void v_InitObject(
LibUtilities::SessionReaderSharedPtr pSession);
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields);
virtual void v_Diffuse(
const int nConvective,
......
This source diff could not be displayed because it is too large. You can view the blob instead.
///////////////////////////////////////////////////////////////////////////////
//
// 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", WITHOUT 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: LFR 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(diffType));
}
static std::string type[];
Array<OneD, Array<OneD, NekDouble> > m_Q2D_e0;
Array<OneD, Array<OneD, NekDouble> > m_Q2D_e1;
Array<OneD, Array<OneD, NekDouble> > m_Q2D_e2;
Array<OneD, Array<OneD, NekDouble> > m_Q2D_e3;
Array<OneD, Array<OneD, NekDouble> > m_dGL_xi1;
Array<OneD, Array<OneD, NekDouble> > m_dGR_xi1;
Array<OneD, Array<OneD, NekDouble> > m_dGL_xi2;
Array<OneD, Array<OneD, NekDouble> > m_dGR_xi2;
Array<OneD, Array<OneD, NekDouble> > m_dGL_xi3;
Array<OneD, Array<OneD, NekDouble> > m_dGR_xi3;
DNekMatSharedPtr m_Ixm;
DNekMatSharedPtr m_Ixp;
protected:
DiffusionLFR(std::string diffType);
Array<OneD, Array<OneD, NekDouble> > m_traceNormals;
Array<OneD, Array<OneD, Array<OneD,NekDouble> > > m_tanbasis;
LibUtilities::SessionReaderSharedPtr m_session;
std::string m_diffType;
virtual void v_InitObject(
LibUtilities::SessionReaderSharedPtr pSession,
Array<OneD, MultiRegions::ExpListSharedPtr> pFields);
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_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_DerCFlux_1D(
const int nConvectiveFields,
const Array<OneD, MultiRegions::ExpListSharedPtr>&fields,
const Array<OneD, const NekDouble> &flux,
const Array<OneD, const NekDouble> &iFlux,
Array<OneD, NekDouble> &derCFlux);
virtual void v_DerCFlux_2D(
const int nConvectiveFields,
const int direction,
const Array<OneD, MultiRegions::ExpListSharedPtr>&fields,
const Array<OneD, const NekDouble> &flux,
const Array<OneD, NekDouble> &iFlux,
Array<OneD, NekDouble> &derCFlux);
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);
/*
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);*/
};
}
}
#endif
......@@ -78,12 +78,14 @@ namespace Nektar
{
case MultiRegions::eDiscontinuous:
{
string diffName;
std::string diffName;
m_session->LoadSolverInfo("DiffusionType", diffName, "LDG");
m_diffusion = SolverUtils::GetDiffusionFactory().CreateInstance(diffName, diffName);
m_diffusion->SetFluxVector(&UnsteadyDiffusion::GetFluxVector, this);
m_diffusion->InitObject(m_session);
m_diffusion = SolverUtils::GetDiffusionFactory().
CreateInstance(diffName, diffName);
m_diffusion->SetFluxVector(&UnsteadyDiffusion::
GetFluxVector, this);
m_diffusion->InitObject(m_session, m_fields);
break;
}
......
......@@ -115,6 +115,9 @@ namespace Nektar
m_session->LoadParameter ("Twall", m_Twall, 300.15);
m_session->LoadSolverInfo("ViscosityType", m_ViscosityType, "Constant");
m_session->LoadParameter ("mu", m_mu, 1.78e-05);
m_session->LoadParameter ("thermalConductivity",
m_thermalConductivity, 0.0257);
// Type of advection class to be used
switch(m_projectionType)
......@@ -183,7 +186,7 @@ namespace Nektar
m_advection->SetRiemannSolver (m_riemannSolver);
m_diffusion->SetRiemannSolver (m_riemannSolverLDG);
m_advection->InitObject (m_session, m_fields);
m_diffusion->InitObject (m_session);
m_diffusion->InitObject (m_session, m_fields);
break;
}
default:
......@@ -506,15 +509,11 @@ namespace Nektar
for (i = 0; i < nvel; ++i)
{
Vmath::Vdiv(nTracePts, Fwd[i+1], 1, Fwd[0], 1, tmp1, 1);
Vmath::Vdiv(nTracePts, tmp1, 1, tmp1, 1, tmp1, 1);
Vmath::Vmul(nTracePts, tmp1, 1, tmp1, 1, tmp1, 1);
Vmath::Vadd(nTracePts, tmp1, 1, absVel, 1, absVel, 1);
}
}
Vmath::Vsqrt(nTracePts, absVel, 1, absVel, 1);
for (i = 0; i < nTracePts; ++i)
{
cout<<"i = "<< i << ";\t absVel = "<< absVel[i] << endl;
}
// Get speed of sound
Array<OneD, NekDouble > SoundSpeed(nTracePts);
Array<OneD, NekDouble > pressure (nTracePts);
......@@ -536,12 +535,7 @@ namespace Nektar
// The Mach number you should get here is the total Mach number
Vmath::Vdiv(nTracePts, absVel, 1, SoundSpeed, 1, Mach, 1);
for (i = 0; i < nTracePts; ++i)
{
cout<<"i = "<< i << ";\t Mach = "<< Mach[i] << endl;
}
// Auxiliary variables
int e, id1, id2, npts, pnt;
NekDouble cPlus, rPlus, cMinus, rMinus;
......@@ -944,7 +938,7 @@ namespace Nektar
int nPts = m_fields[0]->GetTotPoints();
// Note: the value below is referred to 20˚C
NekDouble thermalConductivity = /*0.0257*/0.0000257;
//NekDouble thermalConductivity = /*0.0257*/0.0000257;
// Stokes hypotesis
NekDouble lambda = -0.66666;
......@@ -1069,7 +1063,7 @@ namespace Nektar
&Sgg[0][0], 1, &STx[0], 1);
// k * dT/dx
Vmath::Smul(nPts, thermalConductivity,
Vmath::Smul(nPts, m_thermalConductivity,
&derivativesO1[0][2][0], 1,
&tmp1[0], 1);
......@@ -1092,7 +1086,7 @@ namespace Nektar
&Sxy[0], 1, &tmp1[0], 1);
// k * dT/dx
Vmath::Smul(nPts, thermalConductivity,
Vmath::Smul(nPts, m_thermalConductivity,
&derivativesO1[0][2][0], 1,
&tmp2[0], 1);
......@@ -1115,7 +1109,7 @@ namespace Nektar
&Sxy[0], 1, &tmp1[0], 1);
// k * dT/dy
Vmath::Smul(nPts, thermalConductivity,
Vmath::Smul(nPts, m_thermalConductivity,
&derivativesO1[1][2][0], 1,
&tmp2[0], 1);
......@@ -1144,7 +1138,7 @@ namespace Nektar
&Sxz[0], 1, &tmp2[0], 1);
// k * dT/dx
Vmath::Smul(nPts, thermalConductivity,
Vmath::Smul(nPts, m_thermalConductivity,
&derivativesO1[0][2][0], 1,
&tmp3[0], 1);
......@@ -1173,7 +1167,7 @@ namespace Nektar
&Syz[0], 1, &tmp2[0], 1);
// k * dT/dy
Vmath::Smul(nPts, thermalConductivity,
Vmath::Smul(nPts, m_thermalConductivity,
&derivativesO1[1][2][0], 1,
&tmp3[0], 1);
......@@ -1202,7 +1196,7 @@ namespace Nektar
&Syz[0], 1, &tmp2[0], 1);
// k * dT/dz
Vmath::Smul(nPts, thermalConductivity,
Vmath::Smul(nPts, m_thermalConductivity,
&derivativesO1[2][2][0], 1,
&tmp3[0], 1);
......
......@@ -102,6 +102,7 @@ namespace Nektar
NekDouble m_Twall;
std::string m_ViscosityType;
NekDouble m_mu;
NekDouble m_thermalConductivity;
CompressibleFlowSystem(
const LibUtilities::SessionReaderSharedPtr& pSession);
......
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment