Nektar++
CompressibleFlowSystemImplicit.h
Go to the documentation of this file.
1///////////////////////////////////////////////////////////////////////////////
2//
3// File: CompressibleFlowSystemImplicit.h
4//
5// For more information, please see: http://www.nektar.info
6//
7// The MIT License
8//
9// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10// Department of Aeronautics, Imperial College London (UK), and Scientific
11// Computing and Imaging Institute, University of Utah (USA).
12//
13// Permission is hereby granted, free of charge, to any person obtaining a
14// copy of this software and associated documentation files (the "Software"),
15// to deal in the Software without restriction, including without limitation
16// the rights to use, copy, modify, merge, publish, distribute, sublicense,
17// and/or sell copies of the Software, and to permit persons to whom the
18// Software is furnished to do so, subject to the following conditions:
19//
20// The above copyright notice and this permission notice shall be included
21// in all copies or substantial portions of the Software.
22//
23// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29// DEALINGS IN THE SOFTWARE.
30//
31// Description: Auxiliary functions for the incompressible
32// compressible flow system
33//
34///////////////////////////////////////////////////////////////////////////////
35
36#ifndef NEKTAR_SOLVERS_COMPRESSIBLEFLOWSOLVER_COMPRESSIBLEFLOWSYSTEMIMPLICIT_H
37#define NEKTAR_SOLVERS_COMPRESSIBLEFLOWSOLVER_COMPRESSIBLEFLOWSYSTEMIMPLICIT_H
38
43
44namespace Nektar
45{
46
47/**
48 *
49 */
50class CFSImplicit : virtual public CompressibleFlowSystem
51
52{
53public:
56
57 ~CFSImplicit() override = default;
58
59protected:
63
64 int m_nPadding = 1;
65
66 // Flag to control the update of preconditioning matrix.
68 int m_TotLinIts = 0;
70
71 /// Estimate the magnitude of each conserved varibles
73
75
76 /// coefff of spacial derivatives(rhs or m_F in GLM) in calculating the
77 /// residual of the whole equation(used in unsteady time integrations)
81
84
86
88
89 // flag to update shock capturing artificial viscosity. this flag is
90 // switched on/off in DoImplicitSolve() to ensure that the AV is only
91 // updated once every stage in a multi-stage time integration scheme
93
94 void v_InitObject(bool DeclareFields = true) override;
95
97
98 void v_DoSolve() override;
99
100 void v_PrintStatusInformation(const int step,
101 const NekDouble cpuTime) override;
102
103 void v_PrintSummaryStatistics(const NekDouble intTime) override;
104
105 void v_ALEInitObject(
106 int spaceDim,
108
111 const bool &flag);
112
114 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
115 Array<OneD, Array<OneD, NekDouble>> &out, const bool &flag);
116
117 void DoOdeImplicitRhs(
118 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
119 Array<OneD, Array<OneD, NekDouble>> &outarray, const NekDouble time);
120
121 void DoOdeRhsCoeff(const Array<OneD, const Array<OneD, NekDouble>> &inarray,
123 const NekDouble time);
124
125 void DoAdvectionCoeff(
126 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
127 Array<OneD, Array<OneD, NekDouble>> &outarray, const NekDouble time,
128 const Array<OneD, const Array<OneD, NekDouble>> &pFwd,
129 const Array<OneD, const Array<OneD, NekDouble>> &pBwd);
130
131 void DoDiffusionCoeff(
132 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
134 const Array<OneD, const Array<OneD, NekDouble>> &pFwd,
135 const Array<OneD, const Array<OneD, NekDouble>> &pBwd);
136
137 void DoImplicitSolve(
138 const Array<OneD, const Array<OneD, NekDouble>> &inpnts,
139 Array<OneD, Array<OneD, NekDouble>> &outpnt, const NekDouble time,
140 const NekDouble lambda);
141
143 const Array<OneD, const Array<OneD, NekDouble>> &inpnts,
144 const Array<OneD, const NekDouble> &inarray,
145 Array<OneD, NekDouble> &out, const NekDouble time,
146 const NekDouble lambda);
147
149 const Array<OneD, const NekDouble> &inarray,
150 Array<OneD, NekDouble> &out, const bool &centralDifferenceFlag);
151
152 void CalcRefValues(const Array<OneD, const NekDouble> &inarray);
153
154 void PreconCoeff(const Array<OneD, NekDouble> &inarray,
155 Array<OneD, NekDouble> &outarray, const bool &flag);
156
157 template <typename DataType, typename TypeNekBlkMatSharedPtr>
159 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
160 const Array<OneD, const TensorOfArray2D<NekDouble>> &qfield,
162 TensorOfArray4D<DataType> &StdMatDataDBB,
163 TensorOfArray5D<DataType> &StdMatDataDBDB);
164
165 template <typename DataType>
167 TensorOfArray5D<DataType> &StdMatDataDBDB);
168
169 template <typename DataType, typename TypeNekBlkMatSharedPtr>
171 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
176 Array<OneD, Array<OneD, DataType>> &TraceJacDerivSign,
177 TensorOfArray5D<DataType> &TraceIPSymJacArray);
178
179 template <typename DataType, typename TypeNekBlkMatSharedPtr>
180 void ElmtVarInvMtrx(
182 TypeNekBlkMatSharedPtr &gmtVar, const DataType &tmpDatatype);
183
184 template <typename DataType, typename TypeNekBlkMatSharedPtr>
185 void GetTraceJac(const Array<OneD, const Array<OneD, NekDouble>> &inarray,
189 Array<OneD, Array<OneD, DataType>> &TraceJacDerivSign,
190 TensorOfArray5D<DataType> &TraceIPSymJacArray);
191
192 template <typename DataType, typename TypeNekBlkMatSharedPtr>
194 const int nConvectiveFields,
196 const Array<OneD, const Array<OneD, NekDouble>> &AdvVel,
197 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
198 TensorOfArray3D<NekDouble> &qfield, const NekDouble &time,
199 const Array<OneD, const Array<OneD, NekDouble>> &Fwd,
200 const Array<OneD, const Array<OneD, NekDouble>> &Bwd,
201 TypeNekBlkMatSharedPtr &FJac, TypeNekBlkMatSharedPtr &BJac,
202 TensorOfArray5D<DataType> &TraceIPSymJacArray);
203
205 const Array<OneD, NekDouble> &normals,
206 DNekMatSharedPtr &FJac, const NekDouble efix,
207 const NekDouble fsw);
208
209 template <typename DataType, typename TypeNekBlkMatSharedPtr>
210 void TranSamesizeBlkDiagMatIntoArray(const TypeNekBlkMatSharedPtr &BlkMat,
211 TensorOfArray3D<DataType> &MatArray);
212
213 template <typename DataType, typename TypeNekBlkMatSharedPtr>
216 TensorOfArray4D<DataType> &TraceJacDerivArray);
217
218 template <typename DataType, typename TypeNekBlkMatSharedPtr>
221 const DataType valu);
222
223 template <typename DataType, typename TypeNekBlkMatSharedPtr>
225 Array<OneD, TypeNekBlkMatSharedPtr> &gmtxarray, const DataType valu);
226
228 const Array<OneD, unsigned int> nrow,
229 const Array<OneD, unsigned int> ncol)
230 {
231 mat =
233 SNekMatSharedPtr loc_matNvar;
234 for (int nelm = 0; nelm < nrow.size(); ++nelm)
235 {
236 int nrowsVars = nrow[nelm];
237 int ncolsVars = ncol[nelm];
238
240 nrowsVars, ncolsVars, 0.0);
241 mat->SetBlock(nelm, nelm, loc_matNvar);
242 }
243 }
244
246 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
250 Array<OneD, Array<OneD, NekSingle>> &TraceJacDerivSign,
251 TensorOfArray4D<NekSingle> &TraceJacArray,
252 TensorOfArray4D<NekSingle> &TraceJacDerivArray,
253 TensorOfArray5D<NekSingle> &TraceIPSymJacArray);
254
255 template <typename DataType, typename TypeNekBlkMatSharedPtr>
258 const NekDouble dtlamda);
259
261 const int nConvectiveFields, const int nElmtPnt,
262 const Array<OneD, const Array<OneD, NekDouble>> &locVars,
263 const Array<OneD, NekDouble> &normals, DNekMatSharedPtr &wspMat,
264 Array<OneD, Array<OneD, NekDouble>> &PntJacArray);
265
266 void GetFluxVectorJacPoint(const int nConvectiveFields,
267 const Array<OneD, NekDouble> &conservVar,
268 const Array<OneD, NekDouble> &normals,
269 DNekMatSharedPtr &fluxJac);
270
272 const int nConvectiveFields, const int nDim, const int nPts,
273 const int nTracePts, const NekDouble PenaltyFactor2,
275 const Array<OneD, const Array<OneD, NekDouble>> &AdvVel,
276 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
277 const NekDouble time, TensorOfArray3D<NekDouble> &qfield,
278 const Array<OneD, const Array<OneD, NekDouble>> &vFwd,
279 const Array<OneD, const Array<OneD, NekDouble>> &vBwd,
280 const Array<OneD, const TensorOfArray2D<NekDouble>> &qFwd,
281 const Array<OneD, const TensorOfArray2D<NekDouble>> &qBwd,
282 const Array<OneD, NekDouble> &MuVarTrace,
283 Array<OneD, int> &nonZeroIndex,
284 Array<OneD, Array<OneD, NekDouble>> &traceflux);
285
287 const int nConvectiveFields, const int nElmtPnt,
288 const Array<OneD, const Array<OneD, NekDouble>> &locVars,
289 const TensorOfArray3D<NekDouble> &locDerv,
290 const Array<OneD, NekDouble> &locmu,
291 const Array<OneD, NekDouble> &locDmuDT,
292 const Array<OneD, NekDouble> &normals, DNekMatSharedPtr &wspMat,
293 Array<OneD, Array<OneD, NekDouble>> &PntJacArray)
294 {
295 v_MinusDiffusionFluxJacPoint(nConvectiveFields, nElmtPnt, locVars,
296 locDerv, locmu, locDmuDT, normals, wspMat,
297 PntJacArray);
298 }
299
301 const MultiRegions::ExpListSharedPtr &explist,
302 const Array<OneD, const Array<OneD, NekDouble>> &normals,
303 const int nDervDir,
304 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
305 TensorOfArray5D<NekDouble> &ElmtJacArray, const int nFluxDir)
306 {
307 v_GetFluxDerivJacDirctn(explist, normals, nDervDir, inarray,
308 ElmtJacArray, nFluxDir);
309 }
310
312 const int nConvectiveFields, const int nElmtPnt, const int nDervDir,
313 const Array<OneD, const Array<OneD, NekDouble>> &locVars,
314 const Array<OneD, NekDouble> &locmu,
315 const Array<OneD, const Array<OneD, NekDouble>> &locnormal,
316 DNekMatSharedPtr &wspMat,
317 Array<OneD, Array<OneD, NekDouble>> &PntJacArray)
318 {
319 v_GetFluxDerivJacDirctnElmt(nConvectiveFields, nElmtPnt, nDervDir,
320 locVars, locmu, locnormal, wspMat,
321 PntJacArray);
322 }
323
325 const MultiRegions::ExpListSharedPtr &explist,
326 const Array<OneD, const Array<OneD, NekDouble>> &normals,
327 const int nDervDir,
328 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
330 {
331 v_GetFluxDerivJacDirctn(explist, normals, nDervDir, inarray, ElmtJac);
332 }
333
334 void CalcPhysDeriv(const Array<OneD, const Array<OneD, NekDouble>> &inarray,
336 {
337 v_CalcPhysDeriv(inarray, qfield);
338 }
339
340 void CalcMuDmuDT(const Array<OneD, const Array<OneD, NekDouble>> &inarray,
342 {
343 v_CalcMuDmuDT(inarray, mu, DmuDT);
344 }
345
346 virtual void v_DoDiffusionCoeff(
347 [[maybe_unused]] const Array<OneD, const Array<OneD, NekDouble>>
348 &inarray,
349 [[maybe_unused]] Array<OneD, Array<OneD, NekDouble>> &outarray,
350 [[maybe_unused]] const Array<OneD, const Array<OneD, NekDouble>> &pFwd,
351 [[maybe_unused]] const Array<OneD, const Array<OneD, NekDouble>> &pBwd)
352 {
353 }
354
355 virtual void v_CalcMuDmuDT(
356 [[maybe_unused]] const Array<OneD, const Array<OneD, NekDouble>>
357 &inarray,
358 [[maybe_unused]] Array<OneD, NekDouble> &mu,
359 [[maybe_unused]] Array<OneD, NekDouble> &DmuDT)
360 {
361 }
362
363 virtual void v_CalcPhysDeriv(
364 [[maybe_unused]] const Array<OneD, const Array<OneD, NekDouble>>
365 &inarray,
366 [[maybe_unused]] TensorOfArray3D<NekDouble> &qfield)
367 {
368 }
369
370 virtual void v_MinusDiffusionFluxJacPoint(
371 const int nConvectiveFields, const int nElmtPnt,
372 const Array<OneD, const Array<OneD, NekDouble>> &locVars,
373 const TensorOfArray3D<NekDouble> &locDerv,
374 const Array<OneD, NekDouble> &locmu,
375 const Array<OneD, NekDouble> &locDmuDT,
376 const Array<OneD, NekDouble> &normals, DNekMatSharedPtr &wspMat,
377 Array<OneD, Array<OneD, NekDouble>> &PntJacArray);
378
379 virtual void v_GetFluxDerivJacDirctn(
380 const MultiRegions::ExpListSharedPtr &explist,
381 const Array<OneD, const Array<OneD, NekDouble>> &normals,
382 const int nDervDir,
383 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
384 TensorOfArray5D<NekDouble> &ElmtJacArray, const int nFluxDir);
385
386 virtual void v_GetFluxDerivJacDirctnElmt(
387 const int nConvectiveFields, const int nElmtPnt, const int nDervDir,
388 const Array<OneD, const Array<OneD, NekDouble>> &locVars,
389 const Array<OneD, NekDouble> &locmu,
390 const Array<OneD, const Array<OneD, NekDouble>> &locnormal,
391 DNekMatSharedPtr &wspMat,
392 Array<OneD, Array<OneD, NekDouble>> &PntJacArray);
393
395 const MultiRegions::ExpListSharedPtr &explist,
396 const Array<OneD, const Array<OneD, NekDouble>> &normals,
397 const int nDervDir,
398 const Array<OneD, const Array<OneD, NekDouble>> &inarray,
400
401 bool v_UpdateTimeStepCheck() override;
402};
403
404} // namespace Nektar
405#endif
void CalcVolJacStdMat(TensorOfArray4D< DataType > &StdMatDataDBB, TensorOfArray5D< DataType > &StdMatDataDBDB)
void NonlinSysEvaluatorCoeff(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &out, const bool &flag)
void NumCalcRiemFluxJac(const int nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, const Array< OneD, NekDouble > > &AdvVel, const Array< OneD, const Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfield, const NekDouble &time, const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, TypeNekBlkMatSharedPtr &FJac, TypeNekBlkMatSharedPtr &BJac, TensorOfArray5D< DataType > &TraceIPSymJacArray)
void DoOdeImplicitRhs(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
void AddMatNSBlkDiagVol(const Array< OneD, const Array< OneD, NekDouble > > &inarray, const Array< OneD, const TensorOfArray2D< NekDouble > > &qfield, Array< OneD, Array< OneD, TypeNekBlkMatSharedPtr > > &gmtxarray, TensorOfArray4D< DataType > &StdMatDataDBB, TensorOfArray5D< DataType > &StdMatDataDBDB)
void GetFluxDerivJacDirctnElmt(const int nConvectiveFields, const int nElmtPnt, const int nDervDir, const Array< OneD, const Array< OneD, NekDouble > > &locVars, const Array< OneD, NekDouble > &locmu, const Array< OneD, const Array< OneD, NekDouble > > &locnormal, DNekMatSharedPtr &wspMat, Array< OneD, Array< OneD, NekDouble > > &PntJacArray)
void CalcTraceNumericalFlux(const int nConvectiveFields, const int nDim, const int nPts, const int nTracePts, const NekDouble PenaltyFactor2, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, const Array< OneD, NekDouble > > &AdvVel, const Array< OneD, const Array< OneD, NekDouble > > &inarray, const NekDouble time, TensorOfArray3D< NekDouble > &qfield, const Array< OneD, const Array< OneD, NekDouble > > &vFwd, const Array< OneD, const Array< OneD, NekDouble > > &vBwd, const Array< OneD, const TensorOfArray2D< NekDouble > > &qFwd, const Array< OneD, const TensorOfArray2D< NekDouble > > &qBwd, const Array< OneD, NekDouble > &MuVarTrace, Array< OneD, int > &nonZeroIndex, Array< OneD, Array< OneD, NekDouble > > &traceflux)
void DoImplicitSolve(const Array< OneD, const Array< OneD, NekDouble > > &inpnts, Array< OneD, Array< OneD, NekDouble > > &outpnt, const NekDouble time, const NekDouble lambda)
void DoImplicitSolveCoeff(const Array< OneD, const Array< OneD, NekDouble > > &inpnts, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out, const NekDouble time, const NekDouble lambda)
void GetFluxVectorJacDirElmt(const int nConvectiveFields, const int nElmtPnt, const Array< OneD, const Array< OneD, NekDouble > > &locVars, const Array< OneD, NekDouble > &normals, DNekMatSharedPtr &wspMat, Array< OneD, Array< OneD, NekDouble > > &PntJacArray)
void ElmtVarInvMtrx(Array< OneD, Array< OneD, TypeNekBlkMatSharedPtr > > &gmtxarray, TypeNekBlkMatSharedPtr &gmtVar, const DataType &tmpDatatype)
virtual void v_CalcMuDmuDT(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &mu, Array< OneD, NekDouble > &DmuDT)
void CalcPreconMatBRJCoeff(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, SNekBlkMatSharedPtr > > &gmtxarray, SNekBlkMatSharedPtr &gmtVar, Array< OneD, SNekBlkMatSharedPtr > &TraceJac, Array< OneD, SNekBlkMatSharedPtr > &TraceJacDeriv, Array< OneD, Array< OneD, NekSingle > > &TraceJacDerivSign, TensorOfArray4D< NekSingle > &TraceJacArray, TensorOfArray4D< NekSingle > &TraceJacDerivArray, TensorOfArray5D< NekSingle > &TraceIPSymJacArray)
virtual void v_MinusDiffusionFluxJacPoint(const int nConvectiveFields, const int nElmtPnt, const Array< OneD, const Array< OneD, NekDouble > > &locVars, const TensorOfArray3D< NekDouble > &locDerv, const Array< OneD, NekDouble > &locmu, const Array< OneD, NekDouble > &locDmuDT, const Array< OneD, NekDouble > &normals, DNekMatSharedPtr &wspMat, Array< OneD, Array< OneD, NekDouble > > &PntJacArray)
void PointFluxJacobianPoint(const Array< OneD, NekDouble > &Fwd, const Array< OneD, NekDouble > &normals, DNekMatSharedPtr &FJac, const NekDouble efix, const NekDouble fsw)
void v_DoSolve() override
Solves an unsteady problem.
~CFSImplicit() override=default
CFSImplicit(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
void GetTraceJac(const Array< OneD, const Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfield, Array< OneD, TypeNekBlkMatSharedPtr > &TraceJac, Array< OneD, TypeNekBlkMatSharedPtr > &TraceJacDeriv, Array< OneD, Array< OneD, DataType > > &TraceJacDerivSign, TensorOfArray5D< DataType > &TraceIPSymJacArray)
void PreconCoeff(const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const bool &flag)
void MultiplyElmtInvMassPlusSource(Array< OneD, Array< OneD, TypeNekBlkMatSharedPtr > > &gmtxarray, const NekDouble dtlamda)
void GetFluxDerivJacDirctn(const MultiRegions::ExpListSharedPtr &explist, const Array< OneD, const Array< OneD, NekDouble > > &normals, const int nDervDir, const Array< OneD, const Array< OneD, NekDouble > > &inarray, TensorOfArray5D< NekDouble > &ElmtJacArray, const int nFluxDir)
void NonlinSysEvaluatorCoeff1D(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out, const bool &flag)
void CalcRefValues(const Array< OneD, const NekDouble > &inarray)
void TransTraceJacMatToArray(const Array< OneD, TypeNekBlkMatSharedPtr > &TraceJac, TensorOfArray4D< DataType > &TraceJacDerivArray)
void v_PrintStatusInformation(const int step, const NekDouble cpuTime) override
Print Status Information.
void CalcPhysDeriv(const Array< OneD, const Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfield)
void MinusDiffusionFluxJacPoint(const int nConvectiveFields, const int nElmtPnt, const Array< OneD, const Array< OneD, NekDouble > > &locVars, const TensorOfArray3D< NekDouble > &locDerv, const Array< OneD, NekDouble > &locmu, const Array< OneD, NekDouble > &locDmuDT, const Array< OneD, NekDouble > &normals, DNekMatSharedPtr &wspMat, Array< OneD, Array< OneD, NekDouble > > &PntJacArray)
void CalcMuDmuDT(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &mu, Array< OneD, NekDouble > &DmuDT)
Array< OneD, Array< OneD, NekDouble > > m_solutionPhys
void DoOdeRhsCoeff(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
Compute the right-hand side.
void DoDiffusionCoeff(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, const Array< OneD, NekDouble > > &pFwd, const Array< OneD, const Array< OneD, NekDouble > > &pBwd)
Add the diffusions terms to the right-hand side Similar to DoDiffusion() but with outarray in coeffic...
void AllocateNekBlkMatDig(SNekBlkMatSharedPtr &mat, const Array< OneD, unsigned int > nrow, const Array< OneD, unsigned int > ncol)
virtual void v_GetFluxDerivJacDirctnElmt(const int nConvectiveFields, const int nElmtPnt, const int nDervDir, const Array< OneD, const Array< OneD, NekDouble > > &locVars, const Array< OneD, NekDouble > &locmu, const Array< OneD, const Array< OneD, NekDouble > > &locnormal, DNekMatSharedPtr &wspMat, Array< OneD, Array< OneD, NekDouble > > &PntJacArray)
void v_ALEInitObject(int spaceDim, Array< OneD, MultiRegions::ExpListSharedPtr > &fields) override
TensorOfArray4D< NekSingle > m_stdSMatDataDBB
void TranSamesizeBlkDiagMatIntoArray(const TypeNekBlkMatSharedPtr &BlkMat, TensorOfArray3D< DataType > &MatArray)
virtual void v_GetFluxDerivJacDirctn(const MultiRegions::ExpListSharedPtr &explist, const Array< OneD, const Array< OneD, NekDouble > > &normals, const int nDervDir, const Array< OneD, const Array< OneD, NekDouble > > &inarray, TensorOfArray5D< NekDouble > &ElmtJacArray, const int nFluxDir)
void v_PrintSummaryStatistics(const NekDouble intTime) override
Print Summary Statistics.
virtual void v_DoDiffusionCoeff(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, const Array< OneD, NekDouble > > &pFwd, const Array< OneD, const Array< OneD, NekDouble > > &pBwd)
void GetFluxDerivJacDirctn(const MultiRegions::ExpListSharedPtr &explist, const Array< OneD, const Array< OneD, NekDouble > > &normals, const int nDervDir, const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, DNekMatSharedPtr > > &ElmtJac)
void GetFluxVectorJacPoint(const int nConvectiveFields, const Array< OneD, NekDouble > &conservVar, const Array< OneD, NekDouble > &normals, DNekMatSharedPtr &fluxJac)
TensorOfArray5D< NekSingle > m_stdSMatDataDBDB
NekDouble m_TimeIntegLambda
coefff of spacial derivatives(rhs or m_F in GLM) in calculating the residual of the whole equation(us...
void Fill2DArrayOfBlkDiagonalMat(Array< OneD, Array< OneD, TypeNekBlkMatSharedPtr > > &gmtxarray, const DataType valu)
Array< OneD, NekDouble > m_magnitdEstimat
Estimate the magnitude of each conserved varibles.
void AddMatNSBlkDiagBnd(const Array< OneD, const Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfield, TensorOfArray2D< TypeNekBlkMatSharedPtr > &gmtxarray, Array< OneD, TypeNekBlkMatSharedPtr > &TraceJac, Array< OneD, TypeNekBlkMatSharedPtr > &TraceJacDeriv, Array< OneD, Array< OneD, DataType > > &TraceJacDerivSign, TensorOfArray5D< DataType > &TraceIPSymJacArray)
void MatrixMultiplyMatrixFreeCoeff(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out, const bool &centralDifferenceFlag)
LibUtilities::NekNonlinSysIterSharedPtr m_nonlinsol
virtual void v_CalcPhysDeriv(const Array< OneD, const Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfield)
void Fill1DArrayOfBlkDiagonalMat(Array< OneD, TypeNekBlkMatSharedPtr > &gmtxarray, const DataType valu)
virtual void v_GetFluxDerivJacDirctn(const MultiRegions::ExpListSharedPtr &explist, const Array< OneD, const Array< OneD, NekDouble > > &normals, const int nDervDir, const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, DNekMatSharedPtr > > &ElmtJac)
void v_InitObject(bool DeclareFields=true) override
Initialization object for CFSImplicit class.
void DoAdvectionCoeff(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time, const Array< OneD, const Array< OneD, NekDouble > > &pFwd, const Array< OneD, const Array< OneD, NekDouble > > &pBwd)
Compute the advection terms for the right-hand side.
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
std::shared_ptr< NekNonlinSysIter > NekNonlinSysIterSharedPtr
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
std::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:174
std::shared_ptr< SNekBlkMat > SNekBlkMatSharedPtr
std::shared_ptr< SNekMat > SNekMatSharedPtr
std::shared_ptr< PreconCfs > PreconCfsSharedPtr
Definition: PreconCfs.h:126
std::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:75
double NekDouble