#include "petsc.h"#include "omg.h"#include "oda.h"#include "vecMass.h"#include "elasticityJac.h"Go to the source code of this file.
Defines | |
| #define | BUILD_FULL_VEC_MASS_BLOCK |
| #define | VEC_MASS_DIAG_BLOCK |
| #define | VEC_MASS_ELEM_DIAG_BLOCK |
| #define | VEC_MASS_ELEM_MULT_BLOCK |
| #define | VEC_MASS_MULT_BLOCK |
Functions | |
| PetscErrorCode | ComputeConstVecMass (ot::DAMG damg, Mat J, Mat B) |
| PetscErrorCode | ConstVecMassMatDestroy (Mat J) |
| PetscErrorCode | ConstVecMassMatGetDiagonal (Mat J, Vec diag) |
| PetscErrorCode | ConstVecMassMatMult (Mat J, Vec in, Vec out) |
| PetscErrorCode | CreateConstVecMass (ot::DAMG damg, Mat *jac) |
Variables | |
| double **** | MassType2Stencil |
Definition in file vecMass.C.
|
|
Definition at line 247 of file vecMass.C. Referenced by ComputeConstVecMass(). |
|
|
Value: {\
ot::DA* da = damg->da;\
iC(VecZeroEntries(diag));\
PetscScalar *diagArr;\
unsigned char* bdyArr = data->bdyArr;\
/*Nodal,Non-Ghosted,Write,3 dof*/\
da->vecGetBuffer(diag, diagArr, false, false, false, 3);\
unsigned int maxD;\
double hFac;\
if(da->iAmActive()) {\
maxD = (da->getMaxDepth());\
hFac = 1.0/((double)(1u << (maxD-1)));\
/*Loop through All Elements including ghosted*/\
for(da->init<ot::DA_FLAGS::ALL>(); da->curr() < da->end<ot::DA_FLAGS::ALL>(); da->next<ot::DA_FLAGS::ALL>()) {\
VEC_MASS_ELEM_DIAG_BLOCK\
} /*end i*/\
} /*end if active*/\
da->vecRestoreBuffer(diag, diagArr, false, false, false, 3);\
/*2 IOP = 1 FLOP. Loop counters are included too.*/\
PetscLogFlops(93*(da->getGhostedElementSize()));\
}
Definition at line 94 of file vecMass.C. Referenced by ConstVecMassMatGetDiagonal(). |
|
|
Value: {\
unsigned int idx = da->curr();\
unsigned int lev = da->getLevel(idx);\
double h = hFac*(1u << (maxD - lev));\
double fac = h*h*h/8.0;\
unsigned int indices[8];\
da->getNodeIndices(indices);\
unsigned char childNum = da->getChildNumber();\
unsigned char hnMask = da->getHangingNodeIndex(idx);\
unsigned char elemType = 0;\
GET_ETYPE_BLOCK(elemType,hnMask,childNum)\
for(int k = 0; k < 8; k++) {\
if(bdyArr[indices[k]]) {\
/*Dirichlet Node*/\
for(int dof = 0; dof < 3; dof++) {\
diagArr[(3*indices[k])+dof] = 1.0;\
} /*end dof*/\
} else { \
for(int dof = 0; dof < 3; dof++) {\
diagArr[(3*indices[k])+dof] += (fac*(MassType2Stencil[childNum][elemType][k][k]));\
} /*end dof*/\
}\
} /*end k*/\
}
|
|
|
Value: {\
unsigned int idx = da->curr();\
unsigned int lev = da->getLevel(idx);\
double h = hFac*(1u << (maxD - lev));\
double fac = h*h*h/8.0;\
unsigned int indices[8];\
da->getNodeIndices(indices);\
unsigned char childNum = da->getChildNumber();\
unsigned char hnMask = da->getHangingNodeIndex(idx);\
unsigned char elemType = 0;\
GET_ETYPE_BLOCK(elemType, hnMask, childNum)\
for(int k = 0;k < 8;k++) {\
if(bdyArr[indices[k]]) {\
/*Dirichlet Node Row*/\
for(int dof = 0; dof < 3; dof++) {\
outArr[(3*indices[k]) + dof] = inArr[(3*indices[k]) + dof];\
}/*end for dof*/\
} else {\
for(int j=0;j<8;j++) {\
/*Avoid Dirichlet Node Columns*/\
if(!(bdyArr[indices[j]])) {\
for(int dof = 0; dof < 3; dof++) {\
outArr[(3*indices[k]) + dof] += (fac*(MassType2Stencil[childNum][elemType][k][j])\
*inArr[(3*indices[j]) + dof]);\
}/*end for dof*/\
}\
}/*end for j*/\
}\
}/*end for k*/\
}
|
|
|
Definition at line 174 of file vecMass.C. Referenced by ConstVecMassMatMult(). |
|
||||||||||||||||
|
Definition at line 337 of file vecMass.C. References BUILD_FULL_VEC_MASS_BLOCK, ot::DAMG, and ot::_p_DAMG::user. Referenced by ComputeElasticityRHS(). |
|
|
Definition at line 63 of file vecMass.C. Referenced by CreateConstVecMass(). |
|
||||||||||||
|
Definition at line 116 of file vecMass.C. References ot::DAMG, iC, and VEC_MASS_DIAG_BLOCK. Referenced by CreateConstVecMass(). |
|
||||||||||||||||
|
Definition at line 219 of file vecMass.C. References ot::DAMG, iC, and VEC_MASS_MULT_BLOCK. Referenced by CreateConstVecMass(). |
|
||||||||||||
|
Definition at line 22 of file vecMass.C. References ot::_p_DAMG::comm, ot::DA::computedLocalToGlobal(), ot::DA::computeLocalToGlobalMappings(), ConstVecMassMatDestroy(), ConstVecMassMatGetDiagonal(), ConstVecMassMatMult(), ot::DA::createMatrix(), ot::_p_DAMG::da, ot::DAMG, ot::DA::getNodeSize(), and ot::_p_DAMG::nlevels. Referenced by ComputeElasticityRHS(). |
|
|
Definition at line 45 of file checkError.C. |
1.3.9.1