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.
heatMat.h
1 //
2 // Created by milinda on 11/21/18.
3 //
4 
5 #ifndef DENDRO_5_0_HEATMAT_H
6 #define DENDRO_5_0_HEATMAT_H
7 
8 #include "oda.h"
9 #include "feMatrix.h"
10 
11 namespace HeatEq
12 {
13  class HeatMat : public feMatrix<HeatMat>{
14 
15  private:
16  // some additional work space variables to perform elemental MatVec
17  double* imV1;
18  double* imV2;
19  double* Qx;
20  double* Qy;
21  double* Qz;
22 
23 
24  public:
26  HeatMat(ot::DA* da,unsigned int dof=1);
27 
29  ~HeatMat();
30 
32  virtual void elementalMatVec(const VECType* in,VECType* out, double*coords=NULL,double scale=1.0);
33 
35  bool preMatVec(const VECType* in,VECType* out,double scale=1.0);
36 
38  bool postMatVec(const VECType* in,VECType* out,double scale=1.0);
39 
41  double gridX_to_X(double x);
43  double gridY_to_Y(double y);
45  double gridZ_to_Z(double z);
46 
47  int cgSolve(double * x ,double * b,int max_iter, double& tol,unsigned int var=0);
48 
49 
50 
51  };
52 
53 
54 }
55 
56 
57 
58 
59 #endif //DENDRO_5_0_HEATMAT_H
Definition: heatMat.h:11
double gridZ_to_Z(double z)
octree grid z to domin z
Definition: heatMat.cpp:158
double gridX_to_X(double x)
octree grid x to domin x
Definition: heatMat.cpp:145
HeatMat(ot::DA *da, unsigned int dof=1)
: constructor
Definition: heatMat.cpp:7
bool postMatVec(const VECType *in, VECType *out, double scale=1.0)
things need to be performed after matvec (i.e. coords transform)
Definition: heatMat.cpp:129
Definition: oda.h:125
The class that manages the octree mesh that support FEM computations. Note that this file is a refact...
virtual void elementalMatVec(const VECType *in, VECType *out, double *coords=NULL, double scale=1.0)
Definition: heatMat.cpp:39
bool preMatVec(const VECType *in, VECType *out, double scale=1.0)
things need to be performed before matvec (i.e. coords transform)
Definition: heatMat.cpp:112
~HeatMat()
default destructor
Definition: heatMat.cpp:19
Definition: heatMat.h:13
double gridY_to_Y(double y)
octree grid y to domin y
Definition: heatMat.cpp:151
Definition: feMatrix.h:85