Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members | Related Pages
examples / src / backend

omgNeumann.C File Reference

#include <mpi.h>
#include <cstdio>
#include "oda.h"
#include "omg.h"
#include "Point.h"
#include "parUtils.h"
#include "octUtils.h"
#include "TreeNode.h"
#include "handleStencils.h"
#include <cstdlib>
#include "sys.h"
#include <sstream>
#include "omgNeumann.h"
#include "dendro.h"

Go to the source code of this file.

Classes

struct  Jac2MFreeData

Defines

#define ASSIGN_MAT_PROP_FINE_TO_COARSE_ELEM_BLOCK
#define BUILD_FULL_JAC_TYPE2_BLOCK(B)
#define CHK_AND_SCATTER_FINE_TO_COARSE_BLOCK
#define FINE_TO_COARSE_BLOCK
#define JAC_TYPE2_DIAG_BLOCK
#define JAC_TYPE2_ELEM_DIAG_BLOCK
#define JAC_TYPE2_ELEM_MULT_BLOCK
#define JAC_TYPE2_MULT_BLOCK

Functions

void CalculateCenters (ot::DA *da, std::vector< double > &centers)
PetscErrorCode ComputeJacobian2 (ot::DAMG damg, Mat J, Mat B)
PetscErrorCode ComputeRHS_omgNeumann (ot::DAMG damg, Vec in)
PetscErrorCode CreateJacobian2 (ot::DAMG damg, Mat *jac)
void DestroyUserContexts (ot::DAMG *damg)
void getActiveStateAndActiveCommForKSP_Shell_Jac2or3 (Mat mat, bool &activeState, MPI_Comm &activeComm)
void getPrivateMatricesForKSP_Shell_Jac2 (Mat mat, Mat *AmatPrivate, Mat *PmatPrivate, MatStructure *pFlag)
PetscErrorCode Jacobian2MatDestroy (Mat J)
PetscErrorCode Jacobian2MatGetDiagonal (Mat J, Vec diag)
PetscErrorCode Jacobian2MatMult (Mat J, Vec in, Vec out)
PetscErrorCode Jacobian2ShellMatMult (Mat J, Vec in, Vec out)
void SetPDECoefFromPts (ot::DAMG *damg, const std::vector< double > &centers, void(*CalcCoef)(const std::vector< double > &pts, std::vector< double > &values))
void solve_neumann (std::vector< double > &pts, void(*CalcCoef)(const std::vector< double > &pts, std::vector< double > &values), void(*CalcRHS)(const std::vector< double > &pts, std::vector< double > &values), int numMultigridLevels, Vec &sol, ot::DAMG *&damg)
 This function calculates the approximate solution to the scalar equation -div(alpha*grad u) + beta * u = f in the unit cube subject to homogeneous Neumann boundary conditions.
void solve_neumann_oct (std::vector< ot::TreeNode > &octs, void(*CalcCoef)(const std::vector< double > &pts, std::vector< double > &values), void(*CalcRHS)(const std::vector< double > &pts, std::vector< double > &values), int numMultigridLevels, Vec &sol, ot::DAMG *&damg)
 This function calculates the approximate solution to the scalar equation -div(alpha*grad u) + beta * u = f in the unit cube subject to homogeneous Neumann boundary conditions.

Variables

std::vector< double > force_values
double **** LaplacianType2Stencil
double **** MassType2Stencil
double *** RHSType2Stencil

Detailed Description

Author:
Ilya Lashuk, ilya.lashuk@gmail.com

Definition in file omgNeumann.C.

Define Documentation

#define ASSIGN_MAT_PROP_FINE_TO_COARSE_ELEM_BLOCK
 

Definition at line 73 of file omgNeumann.C.

#define BUILD_FULL_JAC_TYPE2_BLOCK (  B   ) 
 

Definition at line 560 of file omgNeumann.C.

#define CHK_AND_SCATTER_FINE_TO_COARSE_BLOCK
 

Definition at line 40 of file omgNeumann.C.

#define FINE_TO_COARSE_BLOCK
 

Definition at line 115 of file omgNeumann.C.

#define JAC_TYPE2_DIAG_BLOCK
 

Value:

{\
  ot::DA* da = damg->da;\
  iC(VecZeroEntries(diag));\
  double *matPropArr;\
  /*Elem,Ghost,Read-only,1 dof*/\
  da->vecGetBuffer<double>(*(data->matProp),matPropArr,true,true,true,2);\
  PetscScalar *diagArr;\
  /*Nodal,Non-Ghosted,Write,1 dof*/\
  da->vecGetBuffer(diag,diagArr,false,false,false,1);\
  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>()) {\
      JAC_TYPE2_ELEM_DIAG_BLOCK\
    } /*end i*/\
  } /*end if active*/\
  da->vecRestoreBuffer<double>(*(data->matProp),matPropArr,true,true,true,2);\
  da->vecRestoreBuffer(diag,diagArr,false,false,false,1);\
  PetscLogFlops(44*(da->getGhostedElementSize()));\
}

Definition at line 416 of file omgNeumann.C.

#define JAC_TYPE2_ELEM_DIAG_BLOCK
 

Value:

{\
  unsigned int idx = da->curr();\
  unsigned int lev = da->getLevel(idx);\
  double h = hFac*(1u << (maxD - lev));\
  double matP1 = matPropArr[2*idx];\
  double matP2 = matPropArr[2*idx+1];\
  double fac1 = matP1*h/2.0;\
  double fac2 = matP2*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++) {\
    diagArr[indices[k]] +=  ((fac1*(LaplacianType2Stencil[childNum][elemType][k][k])) +\
        (fac2*(MassType2Stencil[childNum][elemType][k][k])));\
  } /*end k*/\
}

Definition at line 396 of file omgNeumann.C.

#define JAC_TYPE2_ELEM_MULT_BLOCK
 

Value:

{\
  unsigned int idx = da->curr();\
  unsigned int lev = da->getLevel(idx);\
  double h = hFac*(1u << (maxD - lev));\
  double matP1 = matPropArr[2*idx];\
  double matP2 = matPropArr[2*idx+1];\
  double fac1 = matP1*h/2.0;\
  double fac2 = matP2*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++) {\
    for(int j=0;j<8;j++) {\
      outArr[indices[k]] +=  (((fac1*(LaplacianType2Stencil[childNum][elemType][k][j])) +\
            (fac2*(MassType2Stencil[childNum][elemType][k][j])))*inArr[indices[j]]);\
    }/*end for j*/\
  }/*end for k*/\
}

Definition at line 468 of file omgNeumann.C.

#define JAC_TYPE2_MULT_BLOCK
 

Definition at line 491 of file omgNeumann.C.

Function Documentation

void CalculateCenters (  ot::DA *  da,
std::vector< double > &  centers
)  [static]
 

Definition at line 788 of file omgNeumann.C.

References ot::DA::createVector(), ot::DA::curr(), ot::DA::end(), ot::DA::getCurrentOffset(), ot::DA::getLevel(), ot::DA::getMaxDepth(), ot::DA::iAmActive(), ot::DA::init(), ot::DA::next(), Point::xint(), Point::yint(), and Point::zint().

Referenced by solve_neumann(), and solve_neumann_oct().

PetscErrorCode ComputeJacobian2 (  ot::DAMG  damg,
Mat  J,
Mat  B
) 
 

Definition at line 616 of file omgNeumann.C.

References BUILD_FULL_JAC_TYPE2_BLOCK, ot::DAMG, Jac2MFreeData::Jmat_private, and ot::_p_DAMG::user.

Referenced by solve_neumann(), and solve_neumann_oct().

PetscErrorCode ComputeRHS_omgNeumann (  ot::DAMG  damg,
Vec  in
) 
 

Definition at line 749 of file omgNeumann.C.

References ot::DA::curr(), ot::_p_DAMG::da, ot::DAMG, ot::DA::end(), force_values, GET_ETYPE_BLOCK, ot::DA::getChildNumber(), ot::DA::getHangingNodeIndex(), ot::DA::getLevel(), ot::DA::getNodeIndices(), ot::DA::iAmActive(), ot::DA::init(), ot::DA::next(), RHSType2Stencil, ot::DA::vecGetBuffer(), and ot::DA::vecRestoreBuffer().

Referenced by solve_neumann(), and solve_neumann_oct().

PetscErrorCode CreateJacobian2 (  ot::DAMG  damg,
Mat *  jac
) 
 

Definition at line 643 of file omgNeumann.C.

References ot::_p_DAMG::comm, ot::DA::computedLocalToGlobal(), ot::DA::computeLocalToGlobalMappings(), ot::DA::createActiveMatrix(), ot::DA::createMatrix(), ot::_p_DAMG::da, ot::DAMG, ot::DA::getComm(), ot::DA::getCommActive(), ot::DA::getNodeSize(), ot::DA::getNpesActive(), ot::DA::getNpesAll(), ot::DA::iAmActive(), Jac2MFreeData::inTmp, Jacobian2MatDestroy(), Jacobian2MatGetDiagonal(), Jacobian2MatMult(), Jacobian2ShellMatMult(), Jac2MFreeData::Jmat_private, ot::_p_DAMG::nlevels, Jac2MFreeData::outTmp, ot::_p_DAMG::totalLevels, and ot::_p_DAMG::user.

Referenced by solve_neumann(), and solve_neumann_oct().

void DestroyUserContexts (  ot::DAMG *  damg  ) 
 

Definition at line 300 of file omgNeumann.C.

References ot::DAMG, Jac2MFreeData::inTmp, Jac2MFreeData::Jmat_private, Jac2MFreeData::matProp, ot::_p_DAMG::nlevels, Jac2MFreeData::outTmp, and ot::_p_DAMG::user.

Referenced by main(), solve_neumann(), and solve_neumann_oct().

void getActiveStateAndActiveCommForKSP_Shell_Jac2or3 (  Mat  mat,
bool &  activeState,
MPI_Comm &  activeComm
) 
 

Definition at line 357 of file omgNeumann.C.

References ot::_p_DAMG::da, ot::DAMG, ot::DA::getCommActive(), and ot::DA::iAmActive().

void getPrivateMatricesForKSP_Shell_Jac2 (  Mat  mat,
Mat *  AmatPrivate,
Mat *  PmatPrivate,
MatStructure *  pFlag
) 
 

Definition at line 370 of file omgNeumann.C.

References ot::DAMG, and Jac2MFreeData::Jmat_private.

PetscErrorCode Jacobian2MatDestroy (  Mat  J  ) 
 

Definition at line 387 of file omgNeumann.C.

Referenced by CreateJacobian2().

PetscErrorCode Jacobian2MatGetDiagonal (  Mat  J,
Vec  diag
) 
 

Definition at line 442 of file omgNeumann.C.

References ot::DAMG, iC, Jac2MFreeData::isFinestLevel, and JAC_TYPE2_DIAG_BLOCK.

Referenced by CreateJacobian2().

PetscErrorCode Jacobian2MatMult (  Mat  J,
Vec  in,
Vec  out
) 
 

Definition at line 533 of file omgNeumann.C.

References ot::DAMG, iC, Jac2MFreeData::isFinestLevel, and JAC_TYPE2_MULT_BLOCK.

Referenced by CreateJacobian2().

PetscErrorCode Jacobian2ShellMatMult (  Mat  J,
Vec  in,
Vec  out
) 
 

Definition at line 324 of file omgNeumann.C.

References ot::DAMG, Jac2MFreeData::inTmp, Jac2MFreeData::Jmat_private, and Jac2MFreeData::outTmp.

Referenced by CreateJacobian2().

void SetPDECoefFromPts (  ot::DAMG *  damg,
const std::vector< double > &  centers,
void(*)(const std::vector< double > &pts, std::vector< double > &values)  CalcCoef
) 
 

Definition at line 188 of file omgNeumann.C.

References CHK_AND_SCATTER_FINE_TO_COARSE_BLOCK, ot::DA::createVector(), ot::DA::curr(), ot::_p_DAMG::da, ot::DAMG, ot::DA::end(), FINE_TO_COARSE_BLOCK, ot::DA::getCommActive(), ot::DA::iAmActive(), ot::DA::init(), ot::DA::next(), ot::_p_DAMG::nlevels, ot::_p_DAMG::user, ot::DA::vecGetBuffer(), and ot::DA::vecRestoreBuffer().

Referenced by solve_neumann(), and solve_neumann_oct().

void solve_neumann (  std::vector< double > &  pts,
void(*)(const std::vector< double > &pts, std::vector< double > &values)  CalcCoef,
void(*)(const std::vector< double > &pts, std::vector< double > &values)  CalcRHS,
int  numMultigridLevels,
Vec &  sol,
ot::DAMG *&  damg
) 
 

This function calculates the approximate solution to the scalar equation -div(alpha*grad u) + beta * u = f in the unit cube subject to homogeneous Neumann boundary conditions.

Author:
Ilya Lashuk
Parameters:
pts coordinates of points. Each point must be in [0,1)^3. Octree will be built based on these points. No octant of the octree will contain more than 1 point. Coordinates are stored in a sequence X1, Y1, Z1, X2, Y2, Z2, ...
CalcCoef function pointer to function which will calculate the values of "alpha" and "beta" (see the equation above) for the center of each octant. Coordinates of all centers are provided in the first parameter (layout is X1,Y1,Z1,X2,Y2,Z2,...). Values should be returned in second parameter. The following layout should be used: alpha1,beta1,alpha2,beta2,.... Caller will set the array size for the second parameter (to be equal to 2/3 size of the first parameter).
CalcRHS function pointer to function which will calculate the value of "f" (see the equation above) for the center of each octant. Coordinates of all centers are provided in the first parameter (layout is the same as above). Values of "f" should be returned in second parameter. Caller will set the array size for the second parameter (to be equal to 1/3 size of the first parameter).
numMultigridLevels desired number of multigrid levels.
sol in this parameter the solution is returned. This object is automatically destroyed when the "damg" object (see next parameter) is destroyed.
damg in this parameter the multigrid context (used to calculate solution) is returned. User must destroy this object when no longer needed by calling DAMGDestroy

Definition at line 832 of file omgNeumann.C.

References ot::balanceOctree(), CalculateCenters(), ComputeJacobian2(), ComputeRHS_omgNeumann(), CreateJacobian2(), createLmatType2(), createMmatType2(), createRHSType2(), ot::DAMG, ot::DAMGCreateAndSetDA(), DAMGGetx, ot::DAMGSetKSP(), ot::DAMGSolve(), DendroIntL, destroyLmatType2(), destroyMmatType2(), destroyRHSType2(), DestroyUserContexts(), force_values, LaplacianType2Stencil, MassType2Stencil, ot::points2Octree(), RHSType2Stencil, and SetPDECoefFromPts().

Referenced by main().

void solve_neumann_oct (  std::vector< ot::TreeNode > &  octs,
void(*)(const std::vector< double > &pts, std::vector< double > &values)  CalcCoef,
void(*)(const std::vector< double > &pts, std::vector< double > &values)  CalcRHS,
int  numMultigridLevels,
Vec &  sol,
ot::DAMG *&  damg
) 
 

This function calculates the approximate solution to the scalar equation -div(alpha*grad u) + beta * u = f in the unit cube subject to homogeneous Neumann boundary conditions.

Author:
Ilya Lashuk
Parameters:
octs sequence of octants. Octree will be built by completing this sequence of octants.
CalcCoef function pointer to function which will calculate the values of "alpha" and "beta" (see the equation above) for the center of each octant. Coordinates of all centers are provided in the first parameter (layout is X1,Y1,Z1,X2,Y2,Z2,...). Values should be returned in second parameter. The following layout should be used: alpha1,beta1,alpha2,beta2,.... Caller will set the array size for the second parameter (to be equal to 2/3 size of the first parameter).
CalcRHS function pointer to function which will calculate the value of "f" (see the equation above) for the center of each octant. Coordinates of all centers are provided in the first parameter (layout is the same as above). Values of "f" should be returned in second parameter. Caller will set the array size for the second parameter (to be equal to 1/3 size of the first parameter).
numMultigridLevels desired number of multigrid levels.
sol in this parameter the solution is returned. This object is automatically destroyed when the "damg" object (see next parameter) is destroyed.
damg in this parameter the multigrid context (used to calculate solution) is returned. User must destroy this object when no longer needed by calling DAMGDestroy

Definition at line 949 of file omgNeumann.C.

References ot::balanceOctree(), CalculateCenters(), ot::completeOctree(), ComputeJacobian2(), ComputeRHS_omgNeumann(), CreateJacobian2(), createLmatType2(), createMmatType2(), createRHSType2(), ot::DAMG, ot::DAMGCreateAndSetDA(), DAMGGetx, ot::DAMGSetKSP(), ot::DAMGSolve(), DendroIntL, destroyLmatType2(), destroyMmatType2(), destroyRHSType2(), DestroyUserContexts(), force_values, ot::DA::getMaxDepth(), LaplacianType2Stencil, MassType2Stencil, RHSType2Stencil, and SetPDECoefFromPts().

Referenced by main().

Variable Documentation

std::vector<double> force_values [static]
 

Definition at line 30 of file omgNeumann.C.

Referenced by ComputeRHS_omgNeumann(), solve_neumann(), and solve_neumann_oct().

double**** LaplacianType2Stencil [static]
 

Definition at line 27 of file omgNeumann.C.

double**** MassType2Stencil [static]
 

Definition at line 28 of file omgNeumann.C.

double*** RHSType2Stencil [static]
 

Definition at line 29 of file omgNeumann.C.

Referenced by ComputeRHS_omgNeumann(), solve_neumann(), and solve_neumann_oct().

Generated on Tue Mar 24 16:14:09 2009 for DENDRO by  doxygen 1.3.9.1