Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
DiffusionLDGNS.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File: DiffusionLDGNS.cpp
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 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: LDGNS diffusion class.
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #include <SolverUtils/Diffusion/DiffusionLDGNS.h>
37 #include <iostream>
38 #include <iomanip>
39 
40 #include <LocalRegions/Expansion2D.h>
41 
42 namespace Nektar
43 {
44  namespace SolverUtils
45  {
46  std::string DiffusionLDGNS::type = GetDiffusionFactory().
47  RegisterCreatorFunction("LDGNS", DiffusionLDGNS::create);
48 
49  DiffusionLDGNS::DiffusionLDGNS()
50  {
51  }
52 
53  void DiffusionLDGNS::v_InitObject(
54  LibUtilities::SessionReaderSharedPtr pSession,
55  Array<OneD, MultiRegions::ExpListSharedPtr> pFields)
56  {
57  m_session = pSession;
58  m_session->LoadParameter ("Gamma", m_gamma, 1.4);
59  m_session->LoadParameter ("GasConstant", m_gasConstant, 287.058);
60  m_session->LoadParameter ("Twall", m_Twall, 300.15);
61  m_session->LoadSolverInfo("ViscosityType", m_ViscosityType,
62  "Constant");
63  m_session->LoadParameter ("mu", m_mu, 1.78e-05);
64  m_session->LoadParameter ("thermalConductivity",
65  m_thermalConductivity, 0.0257);
66  m_session->LoadParameter ("rhoInf", m_rhoInf, 1.225);
67  m_session->LoadParameter ("pInf", m_pInf, 101325);
68 
69  // Setting up the normals
70  int i;
71  int nDim = pFields[0]->GetCoordim(0);
72  int nTracePts = pFields[0]->GetTrace()->GetTotPoints();
73 
74  m_spaceDim = nDim;
75  if (pSession->DefinesSolverInfo("HOMOGENEOUS"))
76  {
77  m_spaceDim = 3;
78  }
79 
80  m_diffDim = m_spaceDim - nDim;
81 
82  m_traceVel = Array<OneD, Array<OneD, NekDouble> >(m_spaceDim);
83  m_traceNormals = Array<OneD, Array<OneD, NekDouble> >(m_spaceDim);
84  for(i = 0; i < m_spaceDim; ++i)
85  {
86  m_traceVel[i] = Array<OneD, NekDouble> (nTracePts, 0.0);
87  m_traceNormals[i] = Array<OneD, NekDouble> (nTracePts);
88  }
89  pFields[0]->GetTrace()->GetNormals(m_traceNormals);
90  }
91 
92  /**
93  * @brief Calculate weak DG Diffusion in the LDG form for the
94  * Navier-Stokes (NS) equations:
95  *
96  * \f$ \langle\psi, \hat{u}\cdot n\rangle
97  * - \langle\nabla\psi \cdot u\rangle
98  * \langle\phi, \hat{q}\cdot n\rangle -
99  * (\nabla \phi \cdot q) \rangle \f$
100  *
101  * The equations that need a diffusion operator are those related
102  * with the velocities and with the energy.
103  *
104  */
105  void DiffusionLDGNS::v_Diffuse(
106  const int nConvectiveFields,
107  const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
108  const Array<OneD, Array<OneD, NekDouble> > &inarray,
109  Array<OneD, Array<OneD, NekDouble> > &outarray)
110  {
111  int i, j;
112  int nDim = fields[0]->GetCoordim(0);
113  int nScalars = inarray.num_elements();
114  int nPts = fields[0]->GetTotPoints();
115  int nCoeffs = fields[0]->GetNcoeffs();
116  int nTracePts = fields[0]->GetTrace()->GetTotPoints();
117 
118  Array<OneD, NekDouble> tmp1(nCoeffs);
119  Array<OneD, Array<OneD, NekDouble> > tmp2(nConvectiveFields);
120 
121  Array<OneD, Array<OneD, Array<OneD, NekDouble> > >
122  numericalFluxO1(m_spaceDim);
123  Array<OneD, Array<OneD, Array<OneD, NekDouble> > >
124  derivativesO1(m_spaceDim);
125 
126  Array<OneD, Array<OneD, NekDouble> > fluxvector(m_spaceDim);
127 
128  for (j = 0; j < m_spaceDim; ++j)
129  {
130  numericalFluxO1[j] = Array<OneD, Array<OneD, NekDouble> >(
131  nScalars);
132  derivativesO1[j] = Array<OneD, Array<OneD, NekDouble> >(
133  nScalars);
134 
135  for (i = 0; i < nScalars; ++i)
136  {
137  numericalFluxO1[j][i] = Array<OneD, NekDouble>(
138  nTracePts, 0.0);
139  derivativesO1[j][i] = Array<OneD, NekDouble>(nPts, 0.0);
140  }
141  }
142 
143  // Compute the numerical fluxes for the first order derivatives
144  v_NumericalFluxO1(fields, inarray, numericalFluxO1);
145 
146  for (j = 0; j < nDim; ++j)
147  {
148  for (i = 0; i < nScalars; ++i)
149  {
150  fields[i]->IProductWRTDerivBase (j, inarray[i], tmp1);
151  Vmath::Neg (nCoeffs, tmp1, 1);
152  fields[i]->AddTraceIntegral (numericalFluxO1[j][i],
153  tmp1);
154  fields[i]->SetPhysState (false);
155  fields[i]->MultiplyByElmtInvMass(tmp1, tmp1);
156  fields[i]->BwdTrans (tmp1, derivativesO1[j][i]);
157  }
158  }
159 
160  // For 3D Homogeneous 1D only take derivatives in 3rd direction
161  if (m_diffDim == 1)
162  {
163  for (i = 0; i < nScalars; ++i)
164  {
165  derivativesO1[2][i] = m_homoDerivs[i];
166  }
167  }
168 
169  // Initialisation viscous tensor
170  m_viscTensor = Array<OneD, Array<OneD, Array<OneD, NekDouble> > >
171  (m_spaceDim);
172  Array<OneD, Array<OneD, NekDouble> > viscousFlux(nConvectiveFields);
173 
174  for (j = 0; j < m_spaceDim; ++j)
175  {
176  m_viscTensor[j] = Array<OneD, Array<OneD, NekDouble> >(
177  nScalars+1);
178  for (i = 0; i < nScalars+1; ++i)
179  {
180  m_viscTensor[j][i] = Array<OneD, NekDouble>(nPts, 0.0);
181  }
182  }
183 
184  for (i = 0; i < nConvectiveFields; ++i)
185  {
186  viscousFlux[i] = Array<OneD, NekDouble>(nTracePts, 0.0);
187  }
188 
189  m_fluxVectorNS(inarray, derivativesO1, m_viscTensor);
190 
191  // Compute u from q_{\eta} and q_{\xi}
192  // Obtain numerical fluxes
193  v_NumericalFluxO2(fields, inarray, m_viscTensor, viscousFlux);
194 
195  for (i = 0; i < nConvectiveFields; ++i)
196  {
197  tmp2[i] = Array<OneD, NekDouble>(nCoeffs, 0.0);
198 
199  for (j = 0; j < nDim; ++j)
200  {
201  fields[i]->IProductWRTDerivBase(j, m_viscTensor[j][i], tmp1);
202  Vmath::Vadd(nCoeffs, tmp1, 1, tmp2[i], 1, tmp2[i], 1);
203  }
204 
205  // Evaulate <\phi, \hat{F}\cdot n> - outarray[i]
206  Vmath::Neg (nCoeffs, tmp2[i], 1);
207  fields[i]->AddTraceIntegral (viscousFlux[i], tmp2[i]);
208  fields[i]->SetPhysState (false);
209  fields[i]->MultiplyByElmtInvMass(tmp2[i], tmp2[i]);
210  fields[i]->BwdTrans (tmp2[i], outarray[i]);
211  }
212  }
213 
214  /**
215  * @brief Builds the numerical flux for the 1st order derivatives
216  *
217  */
218  void DiffusionLDGNS::v_NumericalFluxO1(
219  const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
220  const Array<OneD, Array<OneD, NekDouble> > &inarray,
221  Array<OneD, Array<OneD, Array<OneD, NekDouble> > >
222  &numericalFluxO1)
223  {
224  int i, j;
225  int nTracePts = fields[0]->GetTrace()->GetTotPoints();
226  int nScalars = inarray.num_elements();
227  int nDim = fields[0]->GetCoordim(0);
228 
229  Array<OneD, NekDouble > Vn (nTracePts, 0.0);
230  Array<OneD, NekDouble > fluxtemp(nTracePts, 0.0);
231 
232  // Get the normal velocity Vn
233  for(i = 0; i < nDim; ++i)
234  {
235  Vmath::Svtvp(nTracePts, 1.0, m_traceNormals[i], 1,
236  Vn, 1, Vn, 1);
237  }
238 
239  // Store forwards/backwards space along trace space
240  Array<OneD, Array<OneD, NekDouble> > Fwd (nScalars);
241  Array<OneD, Array<OneD, NekDouble> > Bwd (nScalars);
242  Array<OneD, Array<OneD, NekDouble> > numflux(nScalars);
243 
244  for (i = 0; i < nScalars; ++i)
245  {
246  Fwd[i] = Array<OneD, NekDouble>(nTracePts);
247  Bwd[i] = Array<OneD, NekDouble>(nTracePts);
248  numflux[i] = Array<OneD, NekDouble>(nTracePts);
249  fields[i]->GetFwdBwdTracePhys(inarray[i], Fwd[i], Bwd[i]);
250  fields[0]->GetTrace()->Upwind(Vn, Fwd[i], Bwd[i], numflux[i]);
251  }
252 
253  // Modify the values in case of boundary interfaces
254  if (fields[0]->GetBndCondExpansions().num_elements())
255  {
256  v_WeakPenaltyO1(fields, inarray, numflux);
257  }
258 
259  // Splitting the numerical flux into the dimensions
260  for (j = 0; j < m_spaceDim; ++j)
261  {
262  for (i = 0; i < nScalars; ++i)
263  {
264  Vmath::Vmul(nTracePts, m_traceNormals[j], 1,
265  numflux[i], 1, numericalFluxO1[j][i], 1);
266  }
267  }
268  }
269 
270  /**
271  * @brief Imposes appropriate bcs for the 1st order derivatives
272  *
273  */
274  void DiffusionLDGNS::v_WeakPenaltyO1(
275  const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
276  const Array<OneD, Array<OneD, NekDouble> > &inarray,
277  Array<OneD, Array<OneD, NekDouble> > &penaltyfluxO1)
278  {
279  int cnt;
280  int i, j, e;
281  int id1, id2;
282 
283  int nBndEdgePts, nBndEdges, nBndRegions;
284 
285  int nTracePts = fields[0]->GetTrace()->GetTotPoints();
286  int nScalars = inarray.num_elements();
287 
288  Array<OneD, NekDouble> tmp1(nTracePts, 0.0);
289  Array<OneD, NekDouble> tmp2(nTracePts, 0.0);
290  Array<OneD, NekDouble> Tw(nTracePts, m_Twall);
291 
292  Array< OneD, Array<OneD, NekDouble > > scalarVariables(nScalars);
293  Array< OneD, Array<OneD, NekDouble > > uplus(nScalars);
294 
295  // Extract internal values of the scalar variables for Neumann bcs
296  for (i = 0; i < nScalars; ++i)
297  {
298  scalarVariables[i] = Array<OneD, NekDouble>(nTracePts, 0.0);
299 
300  uplus[i] = Array<OneD, NekDouble>(nTracePts, 0.0);
301  fields[i]->ExtractTracePhys(inarray[i], uplus[i]);
302  }
303 
304  // Compute boundary conditions for velocities
305  for (i = 0; i < nScalars-1; ++i)
306  {
307  // Note that cnt has to loop on nBndRegions and nBndEdges
308  // and has to be reset to zero per each equation
309  cnt = 0;
310  nBndRegions = fields[i+1]->
311  GetBndCondExpansions().num_elements();
312  for (j = 0; j < nBndRegions; ++j)
313  {
314  nBndEdges = fields[i+1]->
315  GetBndCondExpansions()[j]->GetExpSize();
316  for (e = 0; e < nBndEdges; ++e)
317  {
318  nBndEdgePts = fields[i+1]->
319  GetBndCondExpansions()[j]->GetExp(e)->GetTotPoints();
320 
321  id1 = fields[i+1]->
322  GetBndCondExpansions()[j]->GetPhys_Offset(e);
323 
324  id2 = fields[0]->GetTrace()->
325  GetPhys_Offset(fields[0]->GetTraceMap()->
326  GetBndCondTraceToGlobalTraceMap(cnt++));
327 
328  // Reinforcing bcs for velocity in case of Wall bcs
329  if (fields[i]->GetBndConditions()[j]->
330  GetUserDefined() ==
331  SpatialDomains::eWallViscous)
332  {
333  Vmath::Zero(nBndEdgePts,
334  &scalarVariables[i][id2], 1);
335  }
336 
337  // Imposing velocity bcs if not Wall
338  else if (fields[i]->GetBndConditions()[j]->
339  GetBoundaryConditionType() ==
340  SpatialDomains::eDirichlet)
341  {
342  Vmath::Vdiv(nBndEdgePts,
343  &(fields[i+1]->GetBndCondExpansions()[j]->
344  UpdatePhys())[id1], 1,
345  &(fields[0]->GetBndCondExpansions()[j]->
346  UpdatePhys())[id1], 1,
347  &scalarVariables[i][id2], 1);
348  }
349 
350  // For Dirichlet boundary condition: uflux = u_bcs
351  if (fields[i]->GetBndConditions()[j]->
352  GetBoundaryConditionType() ==
353  SpatialDomains::eDirichlet)
354  {
355  Vmath::Vcopy(nBndEdgePts,
356  &scalarVariables[i][id2], 1,
357  &penaltyfluxO1[i][id2], 1);
358  }
359 
360  // For Neumann boundary condition: uflux = u_+
361  else if ((fields[i]->GetBndConditions()[j])->
362  GetBoundaryConditionType() ==
363  SpatialDomains::eNeumann)
364  {
365  Vmath::Vcopy(nBndEdgePts,
366  &uplus[i][id2], 1,
367  &penaltyfluxO1[i][id2], 1);
368  }
369 
370  // Building kinetic energy to be used for T bcs
371  Vmath::Vmul(nBndEdgePts,
372  &scalarVariables[i][id2], 1,
373  &scalarVariables[i][id2], 1,
374  &tmp1[id2], 1);
375 
376  Vmath::Smul(nBndEdgePts, 0.5,
377  &tmp1[id2], 1,
378  &tmp1[id2], 1);
379 
380  Vmath::Vadd(nBndEdgePts,
381  &tmp2[id2], 1,
382  &tmp1[id2], 1,
383  &tmp2[id2], 1);
384  }
385  }
386  }
387 
388  // Compute boundary conditions for temperature
389  cnt = 0;
390  nBndRegions = fields[nScalars]->
391  GetBndCondExpansions().num_elements();
392  for (j = 0; j < nBndRegions; ++j)
393  {
394  nBndEdges = fields[nScalars]->
395  GetBndCondExpansions()[j]->GetExpSize();
396  for (e = 0; e < nBndEdges; ++e)
397  {
398  nBndEdgePts = fields[nScalars]->
399  GetBndCondExpansions()[j]->GetExp(e)->GetTotPoints();
400 
401  id1 = fields[nScalars]->
402  GetBndCondExpansions()[j]->GetPhys_Offset(e);
403 
404  id2 = fields[0]->GetTrace()->
405  GetPhys_Offset(fields[0]->GetTraceMap()->
406  GetBndCondTraceToGlobalTraceMap(cnt++));
407 
408  // Imposing Temperature Twall at the wall
409  if (fields[i]->GetBndConditions()[j]->
410  GetUserDefined() ==
411  SpatialDomains::eWallViscous)
412  {
413  Vmath::Vcopy(nBndEdgePts,
414  &Tw[0], 1,
415  &scalarVariables[nScalars-1][id2], 1);
416  }
417  // Imposing Temperature through condition on the Energy
418  // for no wall boundaries (e.g. farfield)
419  else if (fields[i]->GetBndConditions()[j]->
420  GetBoundaryConditionType() ==
421  SpatialDomains::eDirichlet)
422  {
423  // Divide E by rho
424  Vmath::Vdiv(nBndEdgePts,
425  &(fields[nScalars]->
426  GetBndCondExpansions()[j]->
427  GetPhys())[id1], 1,
428  &(fields[0]->
429  GetBndCondExpansions()[j]->
430  GetPhys())[id1], 1,
431  &scalarVariables[nScalars-1][id2], 1);
432 
433  // Subtract kinetic energy to E/rho
434  Vmath::Vsub(nBndEdgePts,
435  &scalarVariables[nScalars-1][id2], 1,
436  &tmp2[id2], 1,
437  &scalarVariables[nScalars-1][id2], 1);
438 
439  // Multiply by constant factor (gamma-1)/R
440  Vmath::Smul(nBndEdgePts, (m_gamma - 1)/m_gasConstant,
441  &scalarVariables[nScalars-1][id2], 1,
442  &scalarVariables[nScalars-1][id2], 1);
443  }
444 
445  // For Dirichlet boundary condition: uflux = u_bcs
446  if (fields[nScalars]->GetBndConditions()[j]->
447  GetBoundaryConditionType() ==
448  SpatialDomains::eDirichlet)
449  {
450  Vmath::Vcopy(nBndEdgePts,
451  &scalarVariables[nScalars-1][id2], 1,
452  &penaltyfluxO1[nScalars-1][id2], 1);
453  }
454 
455  // For Neumann boundary condition: uflux = u_+
456  else if ((fields[nScalars]->GetBndConditions()[j])->
457  GetBoundaryConditionType() ==
458  SpatialDomains::eNeumann)
459  {
460  Vmath::Vcopy(nBndEdgePts,
461  &uplus[nScalars-1][id2], 1,
462  &penaltyfluxO1[nScalars-1][id2], 1);
463  }
464  }
465  }
466  }
467 
468  /**
469  * @brief Build the numerical flux for the 2nd order derivatives
470  *
471  */
472  void DiffusionLDGNS::v_NumericalFluxO2(
473  const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
474  const Array<OneD, Array<OneD, NekDouble> > &ufield,
475  Array<OneD, Array<OneD, Array<OneD, NekDouble> > > &qfield,
476  Array<OneD, Array<OneD, NekDouble> > &qflux)
477  {
478  int i, j;
479  int nTracePts = fields[0]->GetTrace()->GetTotPoints();
480  int nVariables = fields.num_elements();
481  int nDim = fields[0]->GetCoordim(0);
482 
483  Array<OneD, NekDouble > Fwd(nTracePts);
484  Array<OneD, NekDouble > Bwd(nTracePts);
485  Array<OneD, NekDouble > Vn (nTracePts, 0.0);
486 
487  Array<OneD, NekDouble > qFwd (nTracePts);
488  Array<OneD, NekDouble > qBwd (nTracePts);
489  Array<OneD, NekDouble > qfluxtemp(nTracePts, 0.0);
490 
491  // Get the normal velocity Vn
492  for(i = 0; i < nDim; ++i)
493  {
494  Vmath::Svtvp(nTracePts, 1.0, m_traceNormals[i], 1,
495  Vn, 1, Vn, 1);
496  }
497 
498  // Evaulate Riemann flux
499  // qflux = \hat{q} \cdot u = q \cdot n
500  // Notice: i = 1 (first row of the viscous tensor is zero)
501  for (i = 1; i < nVariables; ++i)
502  {
503  qflux[i] = Array<OneD, NekDouble> (nTracePts, 0.0);
504  for (j = 0; j < nDim; ++j)
505  {
506  // Compute qFwd and qBwd value of qfield in position 'ji'
507  fields[i]->GetFwdBwdTracePhys(qfield[j][i], qFwd, qBwd);
508 
509  // Get Riemann flux of qflux --> LDG implies upwind
510  fields[i]->GetTrace()->Upwind(Vn, qBwd, qFwd, qfluxtemp);
511 
512  // Multiply the Riemann flux by the trace normals
513  Vmath::Vmul(nTracePts, m_traceNormals[j], 1, qfluxtemp, 1,
514  qfluxtemp, 1);
515 
516  // Impose weak boundary condition with flux
517  if (fields[0]->GetBndCondExpansions().num_elements())
518  {
519  v_WeakPenaltyO2(fields, i, j, qfield[j][i], qfluxtemp);
520  }
521 
522  // Store the final flux into qflux
523  Vmath::Vadd(nTracePts, qfluxtemp, 1, qflux[i], 1,
524  qflux[i], 1);
525  }
526  }
527  }
528 
529 
530  /**
531  * @brief Imposes appropriate bcs for the 2nd order derivatives
532  *
533  */
534  void DiffusionLDGNS::v_WeakPenaltyO2(
535  const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
536  const int var,
537  const int dir,
538  const Array<OneD, const NekDouble> &qfield,
539  Array<OneD, NekDouble> &penaltyflux)
540  {
541  int cnt = 0;
542  int nBndEdges, nBndEdgePts;
543  int i, e;
544  int id2;
545 
546  int nTracePts = fields[0]->GetTrace()->GetTotPoints();
547  int nBndRegions = fields[var]->GetBndCondExpansions().num_elements();
548 
549  Array<OneD, NekDouble > uterm(nTracePts);
550  Array<OneD, NekDouble > qtemp(nTracePts);
551 
552  // Extract the physical values of the solution at the boundaries
553  fields[var]->ExtractTracePhys(qfield, qtemp);
554 
555  // Loop on the boundary regions to apply appropriate bcs
556  for (i = 0; i < nBndRegions; ++i)
557  {
558  // Number of boundary regions related to region 'i'
559  nBndEdges = fields[var]->
560  GetBndCondExpansions()[i]->GetExpSize();
561 
562  // Weakly impose bcs by modifying flux values
563  for (e = 0; e < nBndEdges; ++e)
564  {
565  nBndEdgePts = fields[var]->
566  GetBndCondExpansions()[i]->GetExp(e)->GetTotPoints();
567 
568  id2 = fields[0]->GetTrace()->
569  GetPhys_Offset(fields[0]->GetTraceMap()->
570  GetBndCondTraceToGlobalTraceMap(cnt++));
571 
572  // In case of Dirichlet bcs:
573  // uflux = gD
574  if(fields[var]->GetBndConditions()[i]->
575  GetBoundaryConditionType() == SpatialDomains::eDirichlet)
576  {
577  Vmath::Vmul(nBndEdgePts,
578  &m_traceNormals[dir][id2], 1,
579  &qtemp[id2], 1,
580  &penaltyflux[id2], 1);
581  }
582  // 3.4) In case of Neumann bcs:
583  // uflux = u+
584  else if((fields[var]->GetBndConditions()[i])->
585  GetBoundaryConditionType() == SpatialDomains::eNeumann)
586  {
587  ASSERTL0(false,
588  "Neumann bcs not implemented for LDGNS");
589 
590  /*
591  Vmath::Vmul(nBndEdgePts,
592  &m_traceNormals[dir][id2], 1,
593  &(fields[var]->
594  GetBndCondExpansions()[i]->
595  UpdatePhys())[id1], 1,
596  &penaltyflux[id2], 1);
597  */
598  }
599  }
600  }
601  }
602  }//end of namespace SolverUtils
603 }//end of namespace Nektar
Generated on Sun Mar 15 2015 00:11:30 for Nektar++ by   doxygen 1.8.1.2