Dendro  5.01
Dendro in Greek language means tree. The Dendro library is a large scale (262K cores on ORNL's Titan) distributed memory adaptive octree framework. The main goal of Dendro is to perform large scale multiphysics simulations efficeiently in mordern supercomputers. Dendro consists of efficient parallel data structures and algorithms to perform variational ( finite element) methods and finite difference mthods on 2:1 balanced arbitary adaptive octrees which enables the users to perform simulations raning from black holes (binary black hole mergers) to blood flow in human body, where applications ranging from relativity, astrophysics to biomedical engineering.
feVec.h
Go to the documentation of this file.
1 //
2 // Created by milinda on 10/30/18.
3 //
11 #ifndef DENDRO_5_0_FEVEC_H
12 #define DENDRO_5_0_FEVEC_H
13 
14 #include "oda.h"
15 
16 #ifdef BUILD_WITH_PETSC
17 #include "petscdmda.h"
18 #endif
19 
20 class feVec {
21 
22 protected:
25 
28 
31 
34 
35 #ifdef BUILD_WITH_PETSC
36 
37  DM m_uiPETSC_DA;
38 #endif
39 
40 public:
41 
46  {
47  m_uiOctDA=da;
48  }
49 
52  {
53 
54  }
59  //virtual void evalVec(VECType* out,double scale=1.0)=0;
60 
67  virtual void computeVec(const VECType* in,VECType* out,double scale=1.0)=0;
68 
69 
71  inline void setProblemDimensions(const Point& pt_min, const Point& pt_max)
72  {
73  m_uiPtMin=pt_min;
74  m_uiPtMax=pt_max;
75  }
76 
77  virtual void setPlaceholder(const double * v) { assert(false); }
78 
79 #ifdef BUILD_WITH_PETSC
80 // all PETSC function should go here.
81 
88  virtual void computeVec(const Vec& in,Vec& out,double scale=1.0)=0;
89 
90 #endif
91 
92 
93 };
94 
95 #endif //DENDRO_5_0_FEVEC_H
Point m_uiPtMax
problem domain max point
Definition: feVec.h:33
Point m_uiPtMin
problem domain min point
Definition: feVec.h:30
A point class.
Definition: point.h:36
Definition: oda.h:125
DAType
: DA type OCT: uses Dendro DA and PETSC uses petsc based DA
Definition: oda.h:105
~feVec()
deconstructor
Definition: feVec.h:51
virtual void computeVec(const VECType *in, VECType *out, double scale=1.0)=0
Evaluates the RHS of the PDE at specific points (for example evaluation at the quadrature points) ...
The class that manages the octree mesh that support FEM computations. Note that this file is a refact...
feVec(ot::DA *da)
: feVec constructor
Definition: feVec.h:45
void setProblemDimensions(const Point &pt_min, const Point &pt_max)
set the problem dimension
Definition: feVec.h:71
ot::DAType m_uiDaType
: type of the DA
Definition: feVec.h:27
Definition: feVec.h:20
ot::DA * m_uiOctDA
: pointer to OCT DA
Definition: feVec.h:24