Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Static Public Member Functions | Static Public Attributes | Protected Member Functions | Static Protected Attributes | Friends | List of all members
Nektar::SolverUtils::DriverAdaptive Class Reference

Base class for the adaptive polynomial order driver. More...

#include <DriverAdaptive.h>

Inheritance diagram for Nektar::SolverUtils::DriverAdaptive:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SolverUtils::DriverAdaptive:
Collaboration graph
[legend]

Static Public Member Functions

static DriverSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession)
 Creates an instance of this class. More...
 

Static Public Attributes

static std::string className
 Name of the class. More...
 

Protected Member Functions

SOLVER_UTILS_EXPORT DriverAdaptive (const LibUtilities::SessionReaderSharedPtr pSession)
 Constructor. More...
 
virtual SOLVER_UTILS_EXPORT ~DriverAdaptive ()
 Destructor. More...
 
virtual SOLVER_UTILS_EXPORT void v_InitObject (std::ostream &out=std::cout)
 Second-stage initialisation. More...
 
virtual SOLVER_UTILS_EXPORT void v_Execute (std::ostream &out=std::cout)
 Virtual function for solve implementation. More...
 
SOLVER_UTILS_EXPORT void ReplaceExpansion (Array< OneD, MultiRegions::ExpListSharedPtr > &fields, std::map< int, int > deltaP)
 Update EXPANSIONS tag inside XML schema to reflect new polynomial order distribution. More...
 
- Protected Member Functions inherited from Nektar::SolverUtils::Driver
 Driver (const LibUtilities::SessionReaderSharedPtr pSession)
 Initialises EquationSystem class members. More...
 
virtual SOLVER_UTILS_EXPORT
Array< OneD, NekDouble
v_GetRealEvl (void)
 
virtual SOLVER_UTILS_EXPORT
Array< OneD, NekDouble
v_GetImagEvl (void)
 

Static Protected Attributes

static std::string driverLookupId
 
- Static Protected Attributes inherited from Nektar::SolverUtils::Driver
static std::string evolutionOperatorLookupIds []
 
static std::string evolutionOperatorDef
 
static std::string driverDefault
 

Friends

class MemoryManager< DriverAdaptive >
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::SolverUtils::Driver
virtual ~Driver ()
 Destructor. More...
 
SOLVER_UTILS_EXPORT void InitObject (std::ostream &out=std::cout)
 Initialise Object. More...
 
SOLVER_UTILS_EXPORT void Execute (std::ostream &out=std::cout)
 Execute driver. More...
 
SOLVER_UTILS_EXPORT Array
< OneD,
EquationSystemSharedPtr
GetEqu ()
 
SOLVER_UTILS_EXPORT Array
< OneD, NekDouble
GetRealEvl (void)
 
SOLVER_UTILS_EXPORT Array
< OneD, NekDouble
GetImagEvl (void)
 
- Protected Attributes inherited from Nektar::SolverUtils::Driver
LibUtilities::CommSharedPtr m_comm
 Communication object. More...
 
LibUtilities::SessionReaderSharedPtr m_session
 Session reader object. More...
 
LibUtilities::SessionReaderSharedPtr session_LinNS
 I the Coupling between SFD and arnoldi. More...
 
Array< OneD,
EquationSystemSharedPtr
m_equ
 Equation system to solve. More...
 
int m_nequ
 number of equations More...
 
enum EvolutionOperatorType m_EvolutionOperator
 Evolution Operator. More...
 

Detailed Description

Base class for the adaptive polynomial order driver.

Definition at line 48 of file DriverAdaptive.h.

Constructor & Destructor Documentation

Nektar::SolverUtils::DriverAdaptive::DriverAdaptive ( const LibUtilities::SessionReaderSharedPtr  pSession)
protected

Constructor.

Definition at line 62 of file DriverAdaptive.cpp.

64  : Driver(pSession)
65 {
66 }
Driver(const LibUtilities::SessionReaderSharedPtr pSession)
Initialises EquationSystem class members.
Definition: Driver.cpp:78
Nektar::SolverUtils::DriverAdaptive::~DriverAdaptive ( )
protectedvirtual

Destructor.

Definition at line 71 of file DriverAdaptive.cpp.

72 {
73 }

Member Function Documentation

static DriverSharedPtr Nektar::SolverUtils::DriverAdaptive::create ( const LibUtilities::SessionReaderSharedPtr pSession)
inlinestatic

Creates an instance of this class.

Definition at line 54 of file DriverAdaptive.h.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and CellMLToNektar.cellml_metadata::p.

56  {
59  p->InitObject();
60  return p;
61  }
boost::shared_ptr< Driver > DriverSharedPtr
A shared pointer to a Driver object.
Definition: Driver.h:52
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
void Nektar::SolverUtils::DriverAdaptive::ReplaceExpansion ( Array< OneD, MultiRegions::ExpListSharedPtr > &  fields,
std::map< int, int >  deltaP 
)
protected

Update EXPANSIONS tag inside XML schema to reflect new polynomial order distribution.

Parameters
fieldsInput fields.
deltaPMap of polynomial order expansions

Definition at line 445 of file DriverAdaptive.cpp.

References Nektar::LibUtilities::BasisTypeMap, Nektar::ParseUtils::GenerateSeqString(), Nektar::SolverUtils::Driver::m_equ, Nektar::SolverUtils::Driver::m_session, and Nektar::LibUtilities::ShapeTypeMap.

Referenced by v_Execute().

448 {
449  int nExp, nDim;
450  int expdim = m_equ[0]->UpdateFields()[0]->GetGraph()->GetMeshDimension();
451 
452  // Get field definitions
453  std::vector<LibUtilities::FieldDefinitionsSharedPtr> fielddefs =
454  fields[0]->GetFieldDefinitions();
455 
456  if (fielddefs[0]->m_numHomogeneousDir == 1)
457  {
458  nDim = expdim+1;
459  }
460  else
461  {
462  nDim = expdim;
463  }
464 
465  // Add variables to field definition
466  for (int i = 0; i < fielddefs.size(); ++i)
467  {
468  for (int j = 0; j < fields.num_elements(); ++j)
469  {
470  fielddefs[i]->m_fields.push_back(m_session->GetVariable(j));
471  }
472  }
473 
474  // Get tinyxml objects
475  TiXmlElement *exp_tag = m_session->GetElement("NEKTAR/EXPANSIONS");
476 
477  // Clear current expansions
478  exp_tag->Clear();
479 
480  // Write new expansion information
481  for (int f = 0; f < fielddefs.size(); ++f)
482  {
483  nExp = fielddefs[f]->m_elementIDs.size();
484 
485  //---------------------------------------------
486  // Write ELEMENTS
487  TiXmlElement *elemTag = new TiXmlElement("ELEMENTS");
488  exp_tag->LinkEndChild(elemTag);
489 
490  // Write FIELDS
491  std::string fieldsString;
492  {
493  std::stringstream fieldsStringStream;
494  bool first = true;
495  for (std::vector<int>::size_type i = 0;
496  i < fielddefs[f]->m_fields.size(); i++)
497  {
498  if (!first)
499  {
500  fieldsStringStream << ",";
501  }
502  fieldsStringStream << fielddefs[f]->m_fields[i];
503  first = false;
504  }
505  fieldsString = fieldsStringStream.str();
506  }
507  elemTag->SetAttribute("FIELDS", fieldsString);
508 
509  // Write SHAPE
510  std::string shapeString;
511  {
512  std::stringstream shapeStringStream;
513  shapeStringStream
514  << LibUtilities::ShapeTypeMap[fielddefs[f]->m_shapeType];
515  if (fielddefs[f]->m_numHomogeneousDir == 1)
516  {
517  shapeStringStream << "-HomogenousExp1D";
518  }
519  else if (fielddefs[f]->m_numHomogeneousDir == 2)
520  {
521  shapeStringStream << "-HomogenousExp2D";
522  }
523 
524  shapeString = shapeStringStream.str();
525  }
526  elemTag->SetAttribute("SHAPE", shapeString);
527 
528  // Write BASIS
529  std::string basisString;
530  {
531  std::stringstream basisStringStream;
532  bool first = true;
533  for (std::vector<LibUtilities::BasisType>::size_type i = 0;
534  i < fielddefs[f]->m_basis.size(); i++)
535  {
536  if (!first)
537  {
538  basisStringStream << ",";
539  }
540  basisStringStream
541  << LibUtilities::BasisTypeMap[fielddefs[f]->m_basis[i]];
542  first = false;
543  }
544  basisString = basisStringStream.str();
545  }
546  elemTag->SetAttribute("BASIS", basisString);
547 
548  // Write homogeneuous length details
549  if (fielddefs[f]->m_numHomogeneousDir)
550  {
551  std::string homoLenString;
552  {
553  std::stringstream homoLenStringStream;
554  bool first = true;
555  for (int i = 0; i < fielddefs[f]->m_numHomogeneousDir; ++i)
556  {
557  if (!first)
558  {
559  homoLenStringStream << ",";
560  }
561  homoLenStringStream
562  << fielddefs[f]->m_homogeneousLengths[i];
563  first = false;
564  }
565  homoLenString = homoLenStringStream.str();
566  }
567  elemTag->SetAttribute("HOMOGENEOUSLENGTHS", homoLenString);
568  }
569 
570  // Write homogeneuous planes/lines details
571  if (fielddefs[f]->m_numHomogeneousDir)
572  {
573  if (fielddefs[f]->m_homogeneousYIDs.size() > 0)
574  {
575  std::string homoYIDsString;
576  {
577  std::stringstream homoYIDsStringStream;
578  bool first = true;
579  for (int i = 0; i < fielddefs[f]->m_homogeneousYIDs.size();
580  i++)
581  {
582  if (!first)
583  {
584  homoYIDsStringStream << ",";
585  }
586  homoYIDsStringStream
587  << fielddefs[f]->m_homogeneousYIDs[i];
588  first = false;
589  }
590  homoYIDsString = homoYIDsStringStream.str();
591  }
592  elemTag->SetAttribute("HOMOGENEOUSYIDS", homoYIDsString);
593  }
594 
595  if (fielddefs[f]->m_homogeneousZIDs.size() > 0)
596  {
597  std::string homoZIDsString;
598  {
599  std::stringstream homoZIDsStringStream;
600  bool first = true;
601  for (int i = 0; i < fielddefs[f]->m_homogeneousZIDs.size();
602  i++)
603  {
604  if (!first)
605  {
606  homoZIDsStringStream << ",";
607  }
608  homoZIDsStringStream
609  << fielddefs[f]->m_homogeneousZIDs[i];
610  first = false;
611  }
612  homoZIDsString = homoZIDsStringStream.str();
613  }
614  elemTag->SetAttribute("HOMOGENEOUSZIDS", homoZIDsString);
615  }
616  }
617 
618  // Write NUMMODESPERDIR
619  std::string numModesString;
620  {
621  std::stringstream numModesStringStream;
622 
623  numModesStringStream << "MIXORDER:";
624  bool first = true;
625  int eId;
626  Array<OneD, int> order(nDim, 0);
627  for (int n = 0; n < nExp; n++)
628  {
629  eId = fielddefs[f]->m_elementIDs[n];
630 
631  for (int i = 0; i < expdim; i++)
632  {
633  order[i] = deltaP[eId];
634  }
635 
636  for (int i = 0; i < nDim; i++)
637  {
638  if (fielddefs[f]->m_uniOrder)
639  {
640  order[i] += fielddefs[f]->m_numModes[i];
641  }
642  else
643  {
644  order[i] += fielddefs[f]->m_numModes[n * nDim + i];
645  }
646 
647  if (!first)
648  {
649  numModesStringStream << ",";
650  }
651 
652  numModesStringStream << order[i];
653  first = false;
654  }
655  }
656  numModesString = numModesStringStream.str();
657  }
658  elemTag->SetAttribute("NUMMODESPERDIR", numModesString);
659 
660  // Write IDs. Should ideally look at ways of compressing this stream if
661  // just sequential.
662  elemTag->SetAttribute(
663  "ID", ParseUtils::GenerateSeqString(fielddefs[f]->m_elementIDs));
664  }
665 }
const char *const BasisTypeMap[]
Definition: Foundations.hpp:47
static std::string GenerateSeqString(const std::vector< unsigned int > &elmtids)
Definition: ParseUtils.hpp:159
const char *const ShapeTypeMap[]
Definition: ShapeType.hpp:66
Array< OneD, EquationSystemSharedPtr > m_equ
Equation system to solve.
Definition: Driver.h:94
LibUtilities::SessionReaderSharedPtr m_session
Session reader object.
Definition: Driver.h:88
void Nektar::SolverUtils::DriverAdaptive::v_Execute ( std::ostream &  out = std::cout)
protectedvirtual

Virtual function for solve implementation.

Implements Nektar::SolverUtils::Driver.

Definition at line 83 of file DriverAdaptive.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL0, Nektar::LibUtilities::eHexahedron, Nektar::LibUtilities::eOrtho_A, Nektar::LibUtilities::eOrtho_B, Nektar::LibUtilities::eOrtho_C, Nektar::LibUtilities::ePrism, Nektar::LibUtilities::ePyramid, Nektar::LibUtilities::eQuadrilateral, Nektar::LibUtilities::eTetrahedron, Nektar::LibUtilities::eTriangle, Nektar::GlobalMapping::Mapping::Load(), Nektar::SolverUtils::Driver::m_comm, Nektar::SolverUtils::Driver::m_equ, Nektar::SolverUtils::Driver::m_session, Nektar::GlobalMapping::MappingSharedPtr, class_topology::P, Nektar::LibUtilities::ReduceMax, ReplaceExpansion(), Nektar::SolverUtils::Driver::v_InitObject(), Vmath::Vmul(), Vmath::Vsub(), and Vmath::Vsum().

84 {
85  time_t starttime, endtime;
86  NekDouble CPUtime;
87 
88  m_equ[0]->PrintSummary(out);
89 
90  // First run using original order
91  time(&starttime);
92  m_equ[0]->DoInitialise();
93 
94  // Obtain initial time in case a restart was used
95  NekDouble startTime = m_equ[0]->GetFinalTime();
96  m_equ[0]->DoSolve();
97 
98  // Load session parameters and solver info
99  bool isHomogeneous1D;
100  int nRuns, minP, maxP, sensorVar;
101  NekDouble lowerTol, upperTol;
102  m_session->LoadParameter ("NumRuns", nRuns, 1);
103  m_session->LoadParameter ("AdaptiveMinModes", minP, 4);
104  m_session->LoadParameter ("AdaptiveMaxModes", maxP, 12);
105  m_session->LoadParameter ("AdaptiveLowerTolerance", lowerTol, 1e-8);
106  m_session->LoadParameter ("AdaptiveUpperTolerance", upperTol, 1e-6);
107  m_session->LoadParameter ("AdaptiveSensorVariable", sensorVar, 0);
108  m_session->MatchSolverInfo("Homogeneous", "1D", isHomogeneous1D, false);
109 
110  // Get number of elements and planes
111  int nExp, nPlanes;
112  if (isHomogeneous1D)
113  {
114  nExp = m_equ[0]->UpdateFields()[0]->GetPlane(0)->GetExpSize();
115  nPlanes = m_equ[0]->UpdateFields()[0]->GetZIDs().num_elements();
116  }
117  else
118  {
119  nExp = m_equ[0]->UpdateFields()[0]->GetExpSize();
120  nPlanes = 1;
121  }
122  int expdim = m_equ[0]->UpdateFields()[0]->GetGraph()->GetMeshDimension();
123 
124  int nFields = m_equ[0]->UpdateFields().num_elements();
125  int numSteps = m_session->GetParameter("NumSteps");
126  NekDouble period = m_session->GetParameter("TimeStep") * numSteps;
127 
128  Array<OneD, NekDouble> coeffs, phys, physReduced, tmpArray;
129 
130  // Get mapping
133  m_equ[0]->UpdateFields());
134 
135  // Adaptive loop
136  Array<OneD, int> P(expdim);
137  Array<OneD, int> numPoints(expdim);
138  Array<OneD, LibUtilities::PointsKey> ptsKey(expdim);
140  for (int i = 1; i < nRuns; i++)
141  {
142  // Get field expansions
143  Array<OneD, MultiRegions::ExpListSharedPtr> fields =
144  m_equ[0]->UpdateFields();
145 
146  // Determine the change to be applied in the order
147  map<int, int> deltaP;
148  int offset = 0;
149  for (int n = 0; n < nExp; n++)
150  {
151  offset = fields[sensorVar]->GetPhys_Offset(n);
152  Exp = fields[sensorVar]->GetExp(n);
153 
154  for( int k = 0; k < expdim; ++k)
155  {
156  P[k] = Exp->GetBasis(k)->GetNumModes();
157  numPoints[k] = Exp->GetBasis(k)->GetNumPoints();
158  ptsKey[k] = LibUtilities::PointsKey (numPoints[k],
159  Exp->GetBasis(k)->GetPointsType());
160  }
161 
162  // Declare orthogonal basis.
164  switch (Exp->GetGeom()->GetShapeType())
165  {
167  {
168  LibUtilities::BasisKey Ba(LibUtilities::eOrtho_A, P[0] - 1,
169  ptsKey[0]);
170  LibUtilities::BasisKey Bb(LibUtilities::eOrtho_A, P[1] - 1,
171  ptsKey[1]);
172  OrthoExp = MemoryManager<
173  StdRegions::StdQuadExp>::AllocateSharedPtr(Ba, Bb);
174  break;
175  }
177  {
178  LibUtilities::BasisKey Ba(LibUtilities::eOrtho_A, P[0] - 1,
179  ptsKey[0]);
180  LibUtilities::BasisKey Bb(LibUtilities::eOrtho_B, P[1] - 1,
181  ptsKey[1]);
182  OrthoExp =
184  Ba, Bb);
185  break;
186  }
188  {
189  LibUtilities::BasisKey Ba(LibUtilities::eOrtho_A, P[0] - 1,
190  ptsKey[0]);
191  LibUtilities::BasisKey Bb(LibUtilities::eOrtho_B, P[1] - 1,
192  ptsKey[1]);
193  LibUtilities::BasisKey Bc(LibUtilities::eOrtho_C, P[2] - 1,
194  ptsKey[2]);
195  OrthoExp =
197  Ba, Bb, Bc);
198  break;
199  }
201  {
202  LibUtilities::BasisKey Ba(LibUtilities::eOrtho_A, P[0] - 1,
203  ptsKey[0]);
204  LibUtilities::BasisKey Bb(LibUtilities::eOrtho_A, P[1] - 1,
205  ptsKey[1]);
206  LibUtilities::BasisKey Bc(LibUtilities::eOrtho_C, P[2] - 1,
207  ptsKey[2]);
208  OrthoExp =
210  Ba, Bb, Bc);
211  break;
212  }
214  {
215  LibUtilities::BasisKey Ba(LibUtilities::eOrtho_A, P[0] - 1,
216  ptsKey[0]);
217  LibUtilities::BasisKey Bb(LibUtilities::eOrtho_A, P[1] - 1,
218  ptsKey[1]);
219  LibUtilities::BasisKey Bc(LibUtilities::eOrtho_B, P[2] - 1,
220  ptsKey[2]);
221  OrthoExp =
223  Ba, Bb, Bc);
224  break;
225  }
227  {
228  LibUtilities::BasisKey Ba(LibUtilities::eOrtho_A, P[0] - 1,
229  ptsKey[0]);
230  LibUtilities::BasisKey Bb(LibUtilities::eOrtho_A, P[1] - 1,
231  ptsKey[1]);
232  LibUtilities::BasisKey Bc(LibUtilities::eOrtho_A, P[2] - 1,
233  ptsKey[2]);
234  OrthoExp =
236  Ba, Bb, Bc);
237  break;
238  }
239  default:
240  ASSERTL0(false, "Shape not supported.");
241  break;
242  }
243 
244  int nq = OrthoExp->GetTotPoints();
245 
246  NekDouble error = 0;
247  NekDouble fac = 0;
248  NekDouble tmp = 0;
249 
250  coeffs = Array<OneD, NekDouble>(OrthoExp->GetNcoeffs());
251  physReduced = Array<OneD, NekDouble>(OrthoExp->GetTotPoints());
252  tmpArray = Array<OneD, NekDouble>(OrthoExp->GetTotPoints(), 0.0);
253 
254  // Refinement based only on one variable
255  for (int plane = 0; plane < nPlanes; plane++)
256  {
257  if (isHomogeneous1D)
258  {
259  phys =
260  fields[sensorVar]->GetPlane(plane)->GetPhys() + offset;
261  }
262  else
263  {
264  phys = fields[sensorVar]->GetPhys() + offset;
265  }
266 
267  // Project solution to lower order
268  OrthoExp->FwdTrans(phys, coeffs);
269  OrthoExp->BwdTrans(coeffs, physReduced);
270 
271  // Calculate error =||phys-physReduced||^2 / ||phys||^2
272  Vmath::Vsub(nq, phys, 1, physReduced, 1, tmpArray, 1);
273  Vmath::Vmul(nq, tmpArray, 1, tmpArray, 1, tmpArray, 1);
274  tmp = Exp->Integral(tmpArray);
275 
276  Vmath::Vmul(nq, phys, 1, phys, 1, tmpArray, 1);
277  fac = Exp->Integral(tmpArray);
278 
279  tmp = abs(tmp / fac);
280 
281  if (tmp != tmp)
282  {
283  ASSERTL0(false,
284  "Adaptive procedure encountered NaN value.");
285  }
286 
287  // Get maximum value along planes
288  error = (tmp > error) ? tmp : error;
289  }
290 
291  // Combine planes from different processes
292  m_session->GetComm()->GetColumnComm()->AllReduce(
293  error, LibUtilities::ReduceMax);
294 
295  // Override tolerances if function is defined
296  if (m_session->DefinesFunction("AdaptiveLowerTolerance"))
297  {
298  int nq = Exp->GetTotPoints();
299 
300  // Obtain points from the the element
301  Array<OneD, NekDouble> xc0, xc1, xc2;
302  xc0 = Array<OneD, NekDouble>(nq, 0.0);
303  xc1 = Array<OneD, NekDouble>(nq, 0.0);
304  xc2 = Array<OneD, NekDouble>(nq, 0.0);
305  Exp->GetCoords(xc0, xc1, xc2);
306 
307  // Evaluate function from session file
308  Array<OneD, NekDouble> tolerance(nq, 0.0);
310  m_session->GetFunction("AdaptiveLowerTolerance", 0);
311  ffunc->Evaluate(xc0, xc1, xc2, tolerance);
312  lowerTol = Vmath::Vsum(nq, tolerance, 1) / nq;
313  }
314 
315  if (m_session->DefinesFunction("AdaptiveUpperTolerance"))
316  {
317  int nq = Exp->GetTotPoints();
318 
319  // Obtain points from the the element
320  Array<OneD, NekDouble> xc0, xc1, xc2;
321  xc0 = Array<OneD, NekDouble>(nq, 0.0);
322  xc1 = Array<OneD, NekDouble>(nq, 0.0);
323  xc2 = Array<OneD, NekDouble>(nq, 0.0);
324  Exp->GetCoords(xc0, xc1, xc2);
325 
326  // Evaluate function from session file
327  Array<OneD, NekDouble> tolerance(nq, 0.0);
329  m_session->GetFunction("AdaptiveUpperTolerance", 0);
330  ffunc->Evaluate(xc0, xc1, xc2, tolerance);
331  upperTol = Vmath::Vsum(nq, tolerance, 1) / nq;
332  }
333 
334  // Determine what to do with the polynomial order
335  if ((error > upperTol) && (P[0] < maxP))
336  {
337  deltaP[Exp->GetGeom()->GetGlobalID()] = 1;
338  }
339  else if ((error < lowerTol) && P[0] > minP)
340  {
341  deltaP[Exp->GetGeom()->GetGlobalID()] = -1;
342  }
343  else
344  {
345  deltaP[Exp->GetGeom()->GetGlobalID()] = 0;
346  }
347  }
348 
349  // Write new expansion section to the session reader
350  ReplaceExpansion(fields, deltaP);
351 
352  // Reset GlobalLinSys Manager to avoid using too much memory
353  //
354  // @todo This could be made better by replacing individual matrices
355  // within the linear system.
356  if (LibUtilities::NekManager<MultiRegions::GlobalLinSysKey,
357  MultiRegions::GlobalLinSys>::
358  PoolCreated(std::string("GlobalLinSys")))
359  {
360  LibUtilities::NekManager<MultiRegions::GlobalLinSysKey,
361  MultiRegions::GlobalLinSys>::
362  ClearManager(std::string("GlobalLinSys"));
363  }
364 
365  int chkNumber = m_equ[0]->GetCheckpointNumber();
366  int chkSteps = m_equ[0]->GetCheckpointSteps();
367 
368  // Initialise driver again
370 
371  // Update mapping (must be before m_equ[0]->DoInitialise();)
372  mapping->ReplaceField(m_equ[0]->UpdateFields());
373 
374  // Set chkSteps to zero to avoid writing initial condition
375  m_equ[0]->SetCheckpointSteps(0);
376 
377  // Initialise equation
378  m_equ[0]->DoInitialise();
379  m_equ[0]->SetInitialStep(i * numSteps);
380  m_equ[0]->SetSteps(i * numSteps + numSteps);
381  m_equ[0]->SetTime(startTime + i * period);
382  m_equ[0]->SetBoundaryConditions(startTime + i * period);
383  m_equ[0]->SetCheckpointNumber(chkNumber);
384  m_equ[0]->SetCheckpointSteps(chkSteps);
385 
386  // Project solution to new expansion
387  for (int n = 0; n < nFields; n++)
388  {
389  m_equ[0]->UpdateFields()[n]->ExtractCoeffsToCoeffs(
390  fields[n], fields[n]->GetCoeffs(),
391  m_equ[0]->UpdateFields()[n]->UpdateCoeffs());
392  m_equ[0]->UpdateFields()[n]->BwdTrans_IterPerExp(
393  m_equ[0]->UpdateFields()[n]->GetCoeffs(),
394  m_equ[0]->UpdateFields()[n]->UpdatePhys());
395  }
396 
397  // Solve equation
398  m_equ[0]->DoSolve();
399  }
400 
401  time(&endtime);
402 
403  m_equ[0]->Output();
404 
405  if (m_comm->GetRank() == 0)
406  {
407  CPUtime = difftime(endtime, starttime);
408  cout << "-------------------------------------------" << endl;
409  cout << "Total Computation Time = " << CPUtime << "s" << endl;
410  cout << "-------------------------------------------" << endl;
411  }
412 
413  // Evaluate and output computation time and solution accuracy.
414  // The specific format of the error output is essential for the
415  // regression tests to work.
416 
417  // Evaluate L2 Error
418  for (int i = 0; i < m_equ[0]->GetNvariables(); ++i)
419  {
420  Array<OneD, NekDouble> exactsoln(m_equ[0]->GetTotPoints(), 0.0);
421 
422  // Evaluate "ExactSolution" function, or zero array
423  m_equ[0]->EvaluateExactSolution(i, exactsoln, m_equ[0]->GetFinalTime());
424 
425  NekDouble vL2Error = m_equ[0]->L2Error (i, exactsoln);
426  NekDouble vLinfError = m_equ[0]->LinfError(i, exactsoln);
427 
428  if (m_comm->GetRank() == 0)
429  {
430  out << "L 2 error (variable " << m_equ[0]->GetVariable(i)
431  << ") : " << vL2Error << endl;
432  out << "L inf error (variable " << m_equ[0]->GetVariable(i)
433  << ") : " << vLinfError << endl;
434  }
435  }
436 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
SOLVER_UTILS_EXPORT void ReplaceExpansion(Array< OneD, MultiRegions::ExpListSharedPtr > &fields, std::map< int, int > deltaP)
Update EXPANSIONS tag inside XML schema to reflect new polynomial order distribution.
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
Principle Orthogonal Functions .
Definition: BasisType.h:47
virtual SOLVER_UTILS_EXPORT void v_InitObject(std::ostream &out=std::cout)
Definition: Driver.cpp:92
Principle Orthogonal Functions .
Definition: BasisType.h:48
Principle Orthogonal Functions .
Definition: BasisType.h:46
double NekDouble
boost::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
GLOBAL_MAPPING_EXPORT typedef boost::shared_ptr< Mapping > MappingSharedPtr
A shared pointer to a Mapping object.
Definition: Mapping.h:51
boost::shared_ptr< Equation > EquationSharedPtr
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
Definition: Vmath.cpp:343
Array< OneD, EquationSystemSharedPtr > m_equ
Equation system to solve.
Definition: Driver.h:94
LibUtilities::CommSharedPtr m_comm
Communication object.
Definition: Driver.h:85
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
Definition: Vmath.cpp:737
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
static GLOBAL_MAPPING_EXPORT MappingSharedPtr Load(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
Return a pointer to the mapping, creating it on first call.
Definition: Mapping.cpp:265
LibUtilities::SessionReaderSharedPtr m_session
Session reader object.
Definition: Driver.h:88
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:183
void Nektar::SolverUtils::DriverAdaptive::v_InitObject ( std::ostream &  out = std::cout)
protectedvirtual

Second-stage initialisation.

Reimplemented from Nektar::SolverUtils::Driver.

Definition at line 78 of file DriverAdaptive.cpp.

References Nektar::SolverUtils::Driver::v_InitObject().

79 {
81 }
virtual SOLVER_UTILS_EXPORT void v_InitObject(std::ostream &out=std::cout)
Definition: Driver.cpp:92

Friends And Related Function Documentation

friend class MemoryManager< DriverAdaptive >
friend

Definition at line 51 of file DriverAdaptive.h.

Member Data Documentation

string Nektar::SolverUtils::DriverAdaptive::className
static
Initial value:

Name of the class.

Definition at line 64 of file DriverAdaptive.h.

string Nektar::SolverUtils::DriverAdaptive::driverLookupId
staticprotected
Initial value:

Definition at line 84 of file DriverAdaptive.h.