Memosa-FVM  0.2
COMETDiscretizer< T > Class Template Reference

#include <COMETDiscretizer.h>

Collaboration diagram for COMETDiscretizer< T >:

Public Types

typedef NumTypeTraits< T >
::T_Scalar 
T_Scalar
 
typedef Kspace< T > Tkspace
 
typedef kvol< T > Tkvol
 
typedef pmode< T > Tmode
 
typedef Vector< T_Scalar, 3 > VectorT3
 
typedef Array< VectorT3VectorT3Array
 
typedef Array< T_ScalarTArray
 
typedef MatrixJML< T > TMatrix
 
typedef ArrowHeadMatrix< T > TArrow
 
typedef SquareMatrix< T > TSquare
 
typedef map< int, COMETBC< T > * > COMETBCMap
 
typedef COMETModelOptions< T > COpts
 
typedef Array< int > IntArray
 
typedef Array< bool > BoolArray
 
typedef Vector< int, 2 > VecInt2
 
typedef map< int, VecInt2FaceToFg
 
typedef Tmode::Refl_pair Refl_pair
 
typedef SquareTensor< T, 3 > T3Tensor
 
typedef KSConnectivity< T > TKConnectivity
 
typedef TKConnectivityTKCptr
 
typedef vector< TKCptrTKClist
 
typedef map< int, TKClistFgTKClistMap
 
typedef Gradient< T > GradType
 
typedef Array< GradTypeGradArray
 
typedef GradientModel< T > GradModelType
 
typedef
GradModelType::GradMatrixType 
GradMatrix
 

Public Member Functions

 COMETDiscretizer (const Mesh &mesh, const GeomFields &geomfields, PhononMacro &macro, Tkspace &kspace, COMETBCMap &bcMap, const IntArray &BCArray, const IntArray &BCfArray, COpts &options, const FgTKClistMap &FgToKsc)
 
void COMETSolveFine (const int dir, const int level)
 
void COMETSolveCoarse (const int dir, const int level)
 
void COMETSolveFull (const int dir, const int level)
 
void COMETConvectionFine (const int cell0, TArrow &Amat, TArray &BVec, const GradMatrix &gMat)
 
void COMETConvectionCoarse (const int cell0, TArrow &Amat, TArray &BVec)
 
void COMETConvection (const int cell, TSquare &Amat, TArray &BVec)
 
void COMETCollision (const int cell, TMatrix *Amat, TArray &BVec)
 
void COMETEquilibrium (const int cell, TMatrix *Amat, TArray &BVec)
 
void COMETSource (const int cell, TArray &BVec)
 
void COMETFullScatt (const int cell, TArray &s, TArray &BVec)
 
void ScatterPhonons (const int cell0)
 
void COMETShifted (const int cell, TMatrix *Amat, TArray &BVec)
 
void Distribute (const int cell, TArray &BVec, TArray &Rvec)
 
void findResidFine ()
 
void findResidCoarse (const bool plusFAS)
 
void findResidFull ()
 
TArray gatherResid (const int c)
 
getResidChange ()
 
getAveResid ()
 
void setfgFinder ()
 
int findFgId (const int faceIndex)
 
void ArrayAbs (TArray &o)
 
void makeValueArray (const int c, TArray &o)
 
void addFAS (const int c, TArray &bVec)
 
void updatee0 ()
 
void updatee0 (const int c)
 
void updateGhostFine (const int cell, const GradMatrix &gMat)
 
void updateGhostCoarse (const int cell)
 
void correctInterface (const int cell0, TArray &Bvec)
 
void updateeShifted ()
 
void outerProduct (const VectorT3 &v1, const VectorT3 &v2, T3Tensor &out)
 
scaledResid (const TArray &de, const int c)
 

Private Attributes

const Mesh_mesh
 
const GeomFields_geomFields
 
const StorageSite_cells
 
const StorageSite_faces
 
const CRConnectivity_cellFaces
 
const CRConnectivity_faceCells
 
const TArray_faceAreaMag
 
const VectorT3Array_faceArea
 
const TArray_cellVolume
 
const VectorT3Array_cellCoords
 
PhononMacro_macro
 
Tkspace_kspace
 
COMETBCMap_bcMap
 
const IntArray_BCArray
 
const IntArray_BCfArray
 
_aveResid
 
_residChange
 
FaceToFg _fgFinder
 
COpts _options
 
const FgTKClistMap_FaceToKSC
 
TArray_eArray
 
TArray_e0Array
 
TArray_resArray
 

Detailed Description

template<class T>
class COMETDiscretizer< T >

Definition at line 31 of file COMETDiscretizer.h.

Member Typedef Documentation

template<class T>
typedef Array<bool> COMETDiscretizer< T >::BoolArray

Definition at line 48 of file COMETDiscretizer.h.

template<class T>
typedef map<int,COMETBC<T>*> COMETDiscretizer< T >::COMETBCMap

Definition at line 45 of file COMETDiscretizer.h.

template<class T>
typedef COMETModelOptions<T> COMETDiscretizer< T >::COpts

Definition at line 46 of file COMETDiscretizer.h.

template<class T>
typedef map<int,VecInt2> COMETDiscretizer< T >::FaceToFg

Definition at line 50 of file COMETDiscretizer.h.

template<class T>
typedef map<int, TKClist> COMETDiscretizer< T >::FgTKClistMap

Definition at line 56 of file COMETDiscretizer.h.

template<class T>
typedef Array<GradType> COMETDiscretizer< T >::GradArray

Definition at line 58 of file COMETDiscretizer.h.

Definition at line 60 of file COMETDiscretizer.h.

template<class T>
typedef GradientModel<T> COMETDiscretizer< T >::GradModelType

Definition at line 59 of file COMETDiscretizer.h.

template<class T>
typedef Gradient<T> COMETDiscretizer< T >::GradType

Definition at line 57 of file COMETDiscretizer.h.

template<class T>
typedef Array<int> COMETDiscretizer< T >::IntArray

Definition at line 47 of file COMETDiscretizer.h.

template<class T>
typedef Tmode::Refl_pair COMETDiscretizer< T >::Refl_pair

Definition at line 51 of file COMETDiscretizer.h.

template<class T>
typedef SquareTensor<T,3> COMETDiscretizer< T >::T3Tensor

Definition at line 52 of file COMETDiscretizer.h.

template<class T>
typedef NumTypeTraits<T>::T_Scalar COMETDiscretizer< T >::T_Scalar

Definition at line 35 of file COMETDiscretizer.h.

template<class T>
typedef Array<T_Scalar> COMETDiscretizer< T >::TArray

Definition at line 41 of file COMETDiscretizer.h.

template<class T>
typedef ArrowHeadMatrix<T> COMETDiscretizer< T >::TArrow

Definition at line 43 of file COMETDiscretizer.h.

template<class T>
typedef vector<TKCptr> COMETDiscretizer< T >::TKClist

Definition at line 55 of file COMETDiscretizer.h.

template<class T>
typedef KSConnectivity<T> COMETDiscretizer< T >::TKConnectivity

Definition at line 53 of file COMETDiscretizer.h.

template<class T>
typedef TKConnectivity* COMETDiscretizer< T >::TKCptr

Definition at line 54 of file COMETDiscretizer.h.

template<class T>
typedef Kspace<T> COMETDiscretizer< T >::Tkspace

Definition at line 36 of file COMETDiscretizer.h.

template<class T>
typedef kvol<T> COMETDiscretizer< T >::Tkvol

Definition at line 37 of file COMETDiscretizer.h.

template<class T>
typedef MatrixJML<T> COMETDiscretizer< T >::TMatrix

Definition at line 42 of file COMETDiscretizer.h.

template<class T>
typedef pmode<T> COMETDiscretizer< T >::Tmode

Definition at line 38 of file COMETDiscretizer.h.

template<class T>
typedef SquareMatrix<T> COMETDiscretizer< T >::TSquare

Definition at line 44 of file COMETDiscretizer.h.

template<class T>
typedef Vector<int,2> COMETDiscretizer< T >::VecInt2

Definition at line 49 of file COMETDiscretizer.h.

template<class T>
typedef Vector<T_Scalar,3> COMETDiscretizer< T >::VectorT3

Definition at line 39 of file COMETDiscretizer.h.

template<class T>
typedef Array<VectorT3> COMETDiscretizer< T >::VectorT3Array

Definition at line 40 of file COMETDiscretizer.h.

Constructor & Destructor Documentation

template<class T>
COMETDiscretizer< T >::COMETDiscretizer ( const Mesh mesh,
const GeomFields geomfields,
PhononMacro macro,
Tkspace kspace,
COMETBCMap bcMap,
const IntArray BCArray,
const IntArray BCfArray,
COpts options,
const FgTKClistMap FgToKsc 
)
inline

Definition at line 62 of file COMETDiscretizer.h.

65  :
66  _mesh(mesh),
67  _geomFields(geomfields),
68  _cells(mesh.getCells()),
69  _faces(mesh.getFaces()),
70  _cellFaces(mesh.getCellFaces()),
72  _faceAreaMag((dynamic_cast<const TArray&>(_geomFields.areaMag[_faces]))),
73  _faceArea(dynamic_cast<const VectorT3Array&>(_geomFields.area[_faces])),
74  _cellVolume(dynamic_cast<const TArray&>(_geomFields.volume[_cells])),
75  _cellCoords(dynamic_cast<const VectorT3Array&>(_geomFields.coordinate[_cells])),
76  _macro(macro),
77  _kspace(kspace),
78  _bcMap(bcMap),
79  _BCArray(BCArray),
80  _BCfArray(BCfArray),
81  _aveResid(-1.),
82  _residChange(-1.),
83  _fgFinder(),
84  _options(options),
85  _FaceToKSC(FgToKsc),
86  _eArray(kspace.geteArray()),
87  _e0Array(kspace.gete0Array()),
88  _resArray(kspace.getResArray())
89  {}
const VectorT3Array & _cellCoords
PhononMacro & _macro
Field coordinate
Definition: GeomFields.h:19
const CRConnectivity & _faceCells
const IntArray & _BCArray
const StorageSite & _faces
const CRConnectivity & getAllFaceCells() const
Definition: Mesh.cpp:378
const CRConnectivity & getCellFaces() const
Definition: Mesh.cpp:454
const IntArray & _BCfArray
const TArray & _cellVolume
const StorageSite & getFaces() const
Definition: Mesh.h:108
const StorageSite & getCells() const
Definition: Mesh.h:109
Field volume
Definition: GeomFields.h:26
const CRConnectivity & _cellFaces
const TArray & _faceAreaMag
const GeomFields & _geomFields
const VectorT3Array & _faceArea
Field area
Definition: GeomFields.h:23
const StorageSite & _cells
Field areaMag
Definition: GeomFields.h:25
const FgTKClistMap & _FaceToKSC

Member Function Documentation

template<class T>
void COMETDiscretizer< T >::addFAS ( const int  c,
TArray bVec 
)
inline

Definition at line 1604 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, Kspace< T >::addFAS(), Kspace< T >::gettotmodes(), and PhononMacro::TlFASCorrection.

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), COMETDiscretizer< T >::COMETSolveFine(), COMETDiscretizer< T >::COMETSolveFull(), and COMETDiscretizer< T >::findResidCoarse().

1605  {
1606  const int totmodes=_kspace.gettotmodes();
1607  _kspace.addFAS(c,bVec);
1608  TArray& fasArray=dynamic_cast<TArray&>(_macro.TlFASCorrection[_cells]);
1609  bVec[totmodes]-=fasArray[c];
1610  }
PhononMacro & _macro
Array< T_Scalar > TArray
void addFAS(const int c, TArray &Bvec)
Definition: Kspace.h:1273
Field TlFASCorrection
Definition: PhononMacro.h:26
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
template<class T>
void COMETDiscretizer< T >::ArrayAbs ( TArray o)
inline

Definition at line 1576 of file COMETDiscretizer.h.

References fabs(), and Array< T >::getLength().

Referenced by COMETDiscretizer< T >::findResidCoarse(), COMETDiscretizer< T >::findResidFine(), and COMETDiscretizer< T >::findResidFull().

1577  {
1578  int length=o.getLength();
1579  for(int i=0;i<length;i++)
1580  o[i]=fabs(o[i]);
1581  }
Tangent fabs(const Tangent &a)
Definition: Tangent.h:312
template<class T>
void COMETDiscretizer< T >::COMETCollision ( const int  cell,
TMatrix Amat,
TArray BVec 
)
inline

Definition at line 996 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_cellVolume, COMETDiscretizer< T >::_e0Array, COMETDiscretizer< T >::_eArray, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, pmode< T >::calcde0dT(), MatrixJML< T >::getElement(), Kspace< T >::getGlobalIndex(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), Kspace< T >::getTau(), Kspace< T >::gettotmodes(), and PhononMacro::temperature.

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), COMETDiscretizer< T >::COMETSolveFine(), COMETDiscretizer< T >::COMETSolveFull(), COMETDiscretizer< T >::findResidCoarse(), COMETDiscretizer< T >::findResidFine(), and COMETDiscretizer< T >::findResidFull().

997  {
998  const int klen=_kspace.getlength();
999  const int totalmodes=_kspace.gettotmodes();
1000  const int order=totalmodes+1;
1001  TArray& Tlold=dynamic_cast<TArray&>(_macro.temperature[_cells]);
1002  T coeff;
1003  int cellIndex=_kspace.getGlobalIndex(cell,0);
1004 
1005  for(int k=0;k<klen;k++)
1006  {
1007  Tkvol& kvol=_kspace.getkvol(k);
1008  const int numModes=kvol.getmodenum();
1009  for(int m=0;m<numModes;m++)
1010  {
1011  Tmode& mode=kvol.getmode(m);
1012  const int count=mode.getIndex();
1013  const T tau=_kspace.getTau(cellIndex);
1014  T de0dT=mode.calcde0dT(Tlold[cell]);
1015  coeff=_cellVolume[cell]/tau;
1016  Amat->getElement(count,order)-=coeff*de0dT;
1017  Amat->getElement(count,count)+=coeff;
1018  BVec[count-1]-=coeff*_eArray[cellIndex];
1019  BVec[count-1]+=coeff*_e0Array[cellIndex];
1020  cellIndex++;
1021  }
1022  }
1023 
1024  }
PhononMacro & _macro
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
Field temperature
Definition: PhononMacro.h:21
int getmodenum()
Definition: kvol.h:43
Array< T_Scalar > TArray
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
const TArray & _cellVolume
const T getTau(const int index)
Definition: Kspace.h:1231
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
Definition: kvol.h:14
template<class T>
void COMETDiscretizer< T >::COMETConvection ( const int  cell,
TSquare Amat,
TArray BVec 
)
inline

Definition at line 867 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_BCfArray, COMETDiscretizer< T >::_bcMap, COMETDiscretizer< T >::_cellFaces, COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_faceArea, COMETDiscretizer< T >::_faceCells, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::findFgId(), CRConnectivity::getCount(), kvol< T >::getdk3(), pmode< T >::getfield(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), pmode< T >::getReflpair(), and pmode< T >::getv().

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), COMETDiscretizer< T >::COMETSolveFine(), COMETDiscretizer< T >::COMETSolveFull(), COMETDiscretizer< T >::findResidCoarse(), COMETDiscretizer< T >::findResidFine(), and COMETDiscretizer< T >::findResidFull().

868  {
869  const int neibcount=_cellFaces.getCount(cell);
870 
871  for(int j=0;j<neibcount;j++)
872  {
873  const int f=_cellFaces(cell,j);
874  int cell2=_faceCells(f,1);
875  const int klen=_kspace.getlength();
876  VectorT3 Af=_faceArea[f];
877 
878  T flux;
879 
880  if(cell2==cell)
881  {
882  Af=Af*(-1.);
883  cell2=_faceCells(f,0);
884  }
885 
886  if(_BCfArray[f]==2) //If the face in question is an implicit reflecting face
887  {
888  int Fgid=findFgId(f);
889  T refl=(*(_bcMap[Fgid]))["specifiedReflection"];
890  T oneMinusRefl=1.-refl;
891 
892  //first sweep - have to make sumVdotA
893  T sumVdotA=0.;
894  for(int k=0;k<klen;k++)
895  {
897  T dk3=kvol.getdk3();
898  const int numModes=kvol.getmodenum();
899  for(int m=0;m<numModes;m++)
900  {
901  Tmode& mode=kvol.getmode(m);
902  const int count=mode.getIndex();
903  VectorT3 vg=mode.getv();
904  Field& efield=mode.getfield();
905  TArray& eArray=dynamic_cast<TArray&>(efield[_cells]);
906  flux=vg[0]*Af[0]+vg[1]*Af[1]+vg[2]*Af[2];
907  if(flux>T_Scalar(0))
908  {
909  Amat(count,count)-=flux;
910  BVec[count-1]-=flux*eArray[cell];
911  }
912  else
913  {
914  sumVdotA-=flux*dk3;
915  Refl_pair& rpairs=mode.getReflpair(Fgid);
916  const int kk=rpairs.second.second;
917  Tkvol& kkvol=_kspace.getkvol(kk);
918  const int ccount=kkvol.getmode(m).getIndex();
919  Field& eefield=kkvol.getmode(m).getfield();
920  TArray& eeArray=dynamic_cast<TArray&>(eefield[_cells]);
921  Amat(count,ccount)-=flux*refl;
922  BVec[count-1]-=eeArray[cell]*refl*flux;
923  }
924  }
925  }
926 
927  T inverseSumVdotA=1./sumVdotA;
928 
929  //Second sweep
930  for(int k=0;k<klen;k++)
931  {
932  Tkvol& kvol=_kspace.getkvol(k);
933  const int numModes=kvol.getmodenum();
934  for(int m=0;m<numModes;m++)
935  {
936  const int count=kvol.getmode(m).getIndex();
937  VectorT3 vg=kvol.getmode(m).getv();
938  flux=vg[0]*Af[0]+vg[1]*Af[1]+vg[2]*Af[2];
939  if(flux<T_Scalar(0))
940  {
941  for(int kk=0;kk<klen;kk++)
942  {
943  Tkvol& kkvol=_kspace.getkvol(kk);
944  const int numMODES=kkvol.getmodenum();
945  T ddk3=kkvol.getdk3();
946  for(int mm=0;mm<numMODES;mm++)
947  {
948  VectorT3 vvg=kkvol.getmode(mm).getv();
949  T VdotA=vvg[0]*Af[0]+vvg[1]*Af[1]+vvg[2]*Af[2];
950  if(VdotA>T_Scalar(0))
951  {
952  const int ccount=kkvol.getmode(mm).getIndex();
953  Field& eefield=kkvol.getmode(mm).getfield();
954  TArray& eeArray=dynamic_cast<TArray&>(eefield[_cells]);
955  Amat(count,ccount)-=flux*VdotA*ddk3*oneMinusRefl*inverseSumVdotA;
956  BVec[count-1]-=flux*VdotA*ddk3
957  *eeArray[cell]*inverseSumVdotA*oneMinusRefl;
958  }
959  }
960  }
961  }
962  }
963  }
964 
965 
966  }
967  else //if the face in question is not implicit
968  {
969  for(int k=0;k<klen;k++)
970  {
971  Tkvol& kvol=_kspace.getkvol(k);
972  const int numModes=kvol.getmodenum();
973  for(int m=0;m<numModes;m++)
974  {
975  Tmode& mode=kvol.getmode(m);
976  const int count=mode.getIndex();
977  VectorT3 vg=mode.getv();
978  Field& efield=mode.getfield();
979  TArray& eArray=dynamic_cast<TArray&>(efield[_cells]);
980  flux=vg[0]*Af[0]+vg[1]*Af[1]+vg[2]*Af[2];
981  if(flux>T_Scalar(0))
982  {
983  Amat(count,count)-=flux;
984  BVec[count-1]-=flux*eArray[cell];
985  }
986  else
987  BVec[count-1]-=flux*eArray[cell2];
988  }
989  }
990  }
991 
992  }
993 
994  }
int getCount(const int i) const
const CRConnectivity & _faceCells
Definition: Field.h:14
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
Tmode & getmode(int n) const
Definition: kvol.h:44
int getlength() const
Definition: Kspace.h:391
NumTypeTraits< T >::T_Scalar T_Scalar
int getmodenum()
Definition: kvol.h:43
Tmode::Refl_pair Refl_pair
Array< T_Scalar > TArray
int getIndex()
Definition: pmode.h:73
const IntArray & _BCfArray
const CRConnectivity & _cellFaces
const VectorT3Array & _faceArea
const StorageSite & _cells
Definition: kvol.h:14
T getdk3()
Definition: kvol.h:42
int findFgId(const int faceIndex)
template<class T>
void COMETDiscretizer< T >::COMETConvectionCoarse ( const int  cell0,
TArrow Amat,
TArray BVec 
)
inline

Definition at line 819 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cellFaces, COMETDiscretizer< T >::_eArray, COMETDiscretizer< T >::_faceArea, COMETDiscretizer< T >::_faceCells, COMETDiscretizer< T >::_kspace, CRConnectivity::getCount(), ArrowHeadMatrix< X, K >::getElement(), Kspace< T >::getGlobalIndex(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), and pmode< T >::getv().

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), COMETDiscretizer< T >::COMETSolveFull(), COMETDiscretizer< T >::findResidCoarse(), and COMETDiscretizer< T >::findResidFull().

820  {
821  /* This is the COMET discretization for cells that do not
822  do not have a face which is a reflecting boundary. When
823  there is a face with a reflecting boundary, we can no longer
824  use an arrowhead matrix as the structure becomes unknown a priori.
825  */
826  const int neibcount=_cellFaces.getCount(cell0);
827 
828  for(int j=0;j<neibcount;j++)
829  {
830  const int f=_cellFaces(cell0,j);
831  int cell1=_faceCells(f,1);
832  VectorT3 Af=_faceArea[f];
833  if(cell1==cell0)
834  {
835  cell1=_faceCells(f,0);
836  Af*=-1.;
837  }
838  const int klen=_kspace.getlength();
839 
840  T flux;
841 
842  for(int k=0;k<klen;k++)
843  {
844 
846  const int numModes=kvol.getmodenum();
847  for(int m=0;m<numModes;m++)
848  {
849  Tmode& mode=kvol.getmode(m);
850  const int count=mode.getIndex();
851  const VectorT3 vg=mode.getv();
852  flux=vg[0]*Af[0]+vg[1]*Af[1]+vg[2]*Af[2];
853 
854  if(flux>T_Scalar(0))
855  {
856  Amat.getElement(count,count)+=flux;
857  BVec[count-1]-=flux*_eArray[_kspace.getGlobalIndex(cell0,count-1)];
858  }
859  else
860  BVec[count-1]-=flux*_eArray[_kspace.getGlobalIndex(cell1,count-1)];
861  }
862  }
863 
864  }
865  }
int getCount(const int i) const
const CRConnectivity & _faceCells
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
NumTypeTraits< T >::T_Scalar T_Scalar
int getmodenum()
Definition: kvol.h:43
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
const CRConnectivity & _cellFaces
const VectorT3Array & _faceArea
Definition: kvol.h:14
template<class T>
void COMETDiscretizer< T >::COMETConvectionFine ( const int  cell0,
TArrow Amat,
TArray BVec,
const GradMatrix gMat 
)
inline

Definition at line 579 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_BCfArray, COMETDiscretizer< T >::_cellCoords, COMETDiscretizer< T >::_cellFaces, COMETDiscretizer< T >::_eArray, COMETDiscretizer< T >::_faceArea, COMETDiscretizer< T >::_faceCells, COMETDiscretizer< T >::_faces, COMETDiscretizer< T >::_geomFields, COMETDiscretizer< T >::_kspace, computeLimitCoeff2(), GeomFields::coordinate, GradientMatrix< T_Scalar >::getCoeff(), CRConnectivity::getCount(), ArrowHeadMatrix< X, K >::getElement(), Kspace< T >::getGlobalIndex(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), Kspace< T >::gettotmodes(), pmode< T >::getv(), and Array< T >::zero().

Referenced by COMETDiscretizer< T >::COMETSolveFine(), and COMETDiscretizer< T >::findResidFine().

581  {
582 
583  const int neibcount=_cellFaces.getCount(cell0);
584  GradArray Grads(_kspace.gettotmodes());
585  TArray pointMin(_kspace.gettotmodes());
586  pointMin=-1;
587  TArray pointMax(_kspace.gettotmodes());
588  pointMax.zero();
589  TArray neibMin(_kspace.gettotmodes());
590  neibMin=-1;
591  TArray neibMax(_kspace.gettotmodes());
592  neibMax.zero();
593  TArray pointLim(_kspace.gettotmodes());
594  pointLim=100.;
595  TArray neibLim(_kspace.gettotmodes());
596  neibLim=100.;
597 
598  for(int k=0;k<_kspace.gettotmodes();k++)
599  Grads[k].zero();
600 
601  const VectorT3Array& faceCoords=
602  dynamic_cast<const VectorT3Array&>(_geomFields.coordinate[_faces]);
603 
604  for(int j=0;j<neibcount;j++) //first loop to get grad and max/min vals
605  {
606  const int f=_cellFaces(cell0,j);
607  int cell1=_faceCells(f,1);
608  if(cell1==cell0)
609  cell1=_faceCells(f,0);
610  const VectorT3 Gcoeff=gMat.getCoeff(cell0, cell1);
611 
612  int c0ind=_kspace.getGlobalIndex(cell0,0);
613  int c1ind=_kspace.getGlobalIndex(cell1,0);
614 
615  for(int k=0;k<_kspace.gettotmodes();k++)
616  {
617  const T e1=_eArray[c1ind];
618  const T e0=_eArray[c0ind];
619  T& curMax=pointMax[k];
620  T& curMin=pointMin[k];
621  Grads[k].accumulate(Gcoeff, e1-e0);
622 
623  if(e1>curMax)
624  curMax=e1;
625 
626  if(curMin==-1)
627  curMin=e1;
628  else
629  {
630  if(e1<curMin)
631  curMin=e1;
632  }
633 
634  c0ind++;
635  c1ind++;
636  }
637  }
638 
639  for(int j=0;j<neibcount;j++) //second loop to calculate limiting coeff
640  {
641  const int f=_cellFaces(cell0,j);
642  int cell1=_faceCells(f,1);
643  if(cell1==cell0)
644  cell1=_faceCells(f,0);
645 
646  const VectorT3 dr0(faceCoords[f]-_cellCoords[cell0]);
647  int c0ind=_kspace.getGlobalIndex(cell0,0);
648  vanLeer lf;
649 
650  for(int k=0;k<_kspace.gettotmodes();k++)
651  {
652  const T minVal=pointMin[k];
653  const T maxVal=pointMax[k];
654  const T de0(Grads[k]*dr0);
655  T& cl=pointLim[k];
656  computeLimitCoeff2(cl, _eArray[c0ind], de0, minVal, maxVal, lf);
657  c0ind++;
658  }
659  }
660 
661  for(int j=0;j<neibcount;j++) //loop to construct point matrix
662  {
663  const int f=_cellFaces(cell0,j);
664  int cell1=_faceCells(f,1);
665  VectorT3 Af=_faceArea[f];
666  if(cell1==cell0)
667  {
668  cell1=_faceCells(f,0);
669  Af*=-1.;
670  }
671 
672  if(_BCfArray[f]==0) //interior face
673  {
674  const int klen=_kspace.getlength();
675 
676  GradArray NeibGrads(_kspace.gettotmodes());
677 
678  for(int k=0;k<_kspace.gettotmodes();k++)
679  NeibGrads[k].zero();
680 
681  neibMax.zero();
682  neibMin=-1.;
683  const int neibcount1=_cellFaces.getCount(cell1);
684  for(int nj=0;nj<neibcount1;nj++) //loop to make gradients and min/max
685  {
686  const int f1=_cellFaces(cell1,nj);
687  int cell11=_faceCells(f1,1);
688  if(cell1==cell11)
689  cell11=_faceCells(f1,0);
690 
691  const VectorT3 Gcoeff=gMat.getCoeff(cell1, cell11);
692 
693  int c1ind=_kspace.getGlobalIndex(cell1,0);
694  int c11ind=_kspace.getGlobalIndex(cell11,0);
695 
696  for(int k=0;k<_kspace.gettotmodes();k++)
697  {
698  const T e1=_eArray[c1ind];
699  const T e11=_eArray[c11ind];
700  T& curMax=neibMax[k];
701  T& curMin=neibMin[k];
702  NeibGrads[k].accumulate(Gcoeff,e11-e1);
703 
704  if(e11>curMax)
705  curMax=e11;
706 
707  if(curMin==-1)
708  curMin=e11;
709  else
710  {
711  if(e11<curMin)
712  curMin=e11;
713  }
714 
715  c11ind++;
716  c1ind++;
717  }
718  }
719 
720  neibLim=100.;
721  for(int nj=0;nj<neibcount1;nj++) //second loop to calculate limiting coeff
722  {
723  const int f1=_cellFaces(cell1,nj);
724  int cell11=_faceCells(f1,1);
725  if(cell1==cell11)
726  cell11=_faceCells(f1,0);
727 
728  const VectorT3 dr1(faceCoords[f1]-_cellCoords[cell1]);
729  int c1ind=_kspace.getGlobalIndex(cell1,0);
730  vanLeer lf;
731 
732  for(int k=0;k<_kspace.gettotmodes();k++)
733  {
734  const T minVal=neibMin[k];
735  const T maxVal=neibMax[k];
736  const T de1(NeibGrads[k]*dr1);
737  T& cl=neibLim[k];
738  computeLimitCoeff2(cl, _eArray[c1ind], de1, minVal, maxVal, lf);
739  if(_BCfArray[f]!=0)
740  NeibGrads[k].zero();
741  c1ind++;
742  }
743  }
744 
745 
746  T flux;
747 
748  for(int k=0;k<klen;k++)
749  {
751  const int numModes=kvol.getmodenum();
752  for(int m=0;m<numModes;m++)
753  {
754  Tmode& mode=kvol.getmode(m);
755  const int count=mode.getIndex();
756  VectorT3 vg=mode.getv();
757  //T tau=mode.gettau();
758  flux=vg[0]*Af[0]+vg[1]*Af[1]+vg[2]*Af[2];
759 
760  const int c0ind=_kspace.getGlobalIndex(cell0,count-1);
761  const int c1ind=_kspace.getGlobalIndex(cell1,count-1);
762 
763  const GradType& grad=Grads[count-1];
764  const GradType& neibGrad=NeibGrads[count-1];
765  //vanLeer vl;
766 
767  if(flux>T_Scalar(0))
768  {
769  const VectorT3 rVec=_cellCoords[cell1]-_cellCoords[cell0];
770  const VectorT3 fVec=faceCoords[f]-_cellCoords[cell0];
771  //const T r=gMat.computeR(grad,_eArray,rVec,c0ind,c1ind);
772 
773  T SOU=(fVec[0]*grad[0]+fVec[1]*grad[1]+
774  fVec[2]*grad[2])*pointLim[count-1];
775  Amat.getElement(count,count)-=flux;
776  BVec[count-1]-=flux*_eArray[c0ind]-flux*SOU;
777  }
778  else
779  {
780  const VectorT3 rVec=_cellCoords[cell0]-_cellCoords[cell1];
781  const VectorT3 fVec=faceCoords[f]-_cellCoords[cell1];
782  //const T r=gMat.computeR(grad,_eArray,rVec,c1ind,c0ind);
783 
784  T SOU=(fVec[0]*neibGrad[0]+fVec[1]*neibGrad[1]+
785  fVec[2]*neibGrad[2])*neibLim[count-1];
786  BVec[count-1]-=flux*_eArray[c1ind]-flux*SOU;
787  }
788 
789  }
790  }
791  }
792  else
793  {
794  const int klen=_kspace.getlength();
795  T flux(0);
796  for(int k=0;k<klen;k++)
797  {
798 
799  Tkvol& kvol=_kspace.getkvol(k);
800  const int numModes=kvol.getmodenum();
801  for(int m=0;m<numModes;m++)
802  {
803  Tmode& mode=kvol.getmode(m);
804  const int count=mode.getIndex();
805  const VectorT3 vg=mode.getv();
806  flux=vg[0]*Af[0]+vg[1]*Af[1]+vg[2]*Af[2];
807 
808  if(flux>T_Scalar(0))
809  Amat.getElement(count,count)-=flux;
810 
811  BVec[count-1]-=flux*_eArray[_kspace.getGlobalIndex(cell1,count-1)];
812  }
813  }
814  }
815 
816  }
817  }
int getCount(const int i) const
const VectorT3Array & _cellCoords
Field coordinate
Definition: GeomFields.h:19
const CRConnectivity & _faceCells
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
const StorageSite & _faces
NumTypeTraits< T >::T_Scalar T_Scalar
int getmodenum()
Definition: kvol.h:43
Array< GradType > GradArray
Array< T_Scalar > TArray
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
const IntArray & _BCfArray
const CRConnectivity & _cellFaces
void computeLimitCoeff2(T &lc, const T x, const T &dx, const T &min, const T &max, const LimitFunc &f)
Definition: FluxLimiters.h:85
const GeomFields & _geomFields
int gettotmodes()
Definition: Kspace.h:393
const VectorT3Array & _faceArea
Gradient< T > GradType
Definition: kvol.h:14
template<class T>
void COMETDiscretizer< T >::COMETEquilibrium ( const int  cell,
TMatrix Amat,
TArray BVec 
)
inline

Definition at line 1026 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_e0Array, COMETDiscretizer< T >::_eArray, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, Kspace< T >::getde0taudT(), kvol< T >::getdk3(), Kspace< T >::getDK3(), MatrixJML< T >::getElement(), Kspace< T >::getGlobalIndex(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), Kspace< T >::getTau(), Kspace< T >::gettotmodes(), and PhononMacro::temperature.

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), COMETDiscretizer< T >::COMETSolveFine(), COMETDiscretizer< T >::COMETSolveFull(), COMETDiscretizer< T >::findResidCoarse(), COMETDiscretizer< T >::findResidFine(), and COMETDiscretizer< T >::findResidFull().

1027  {
1028  const int klen=_kspace.getlength();
1029  const int totalmodes=_kspace.gettotmodes();
1030  const int order=totalmodes+1;
1031  TArray& Tlold=dynamic_cast<TArray&>(_macro.temperature[_cells]);
1032  const T tauTot=_kspace.getde0taudT(cell,Tlold[cell]);
1033  T coeff;
1034  const T DK3=_kspace.getDK3();
1035  int cellIndex=_kspace.getGlobalIndex(cell,0);
1036  const T hbar=6.582119e-16; // (eV s)
1037 
1038  for(int k=0;k<klen;k++)
1039  {
1040  Tkvol& kvol=_kspace.getkvol(k);
1041  const T dk3=kvol.getdk3();
1042  const int numModes=kvol.getmodenum();
1043  for(int m=0;m<numModes;m++)
1044  {
1045  Tmode& mode=kvol.getmode(m);
1046  const int count=mode.getIndex();
1047  const T tau=_kspace.getTau(cellIndex);
1048  //const T energy=mode.getomega()*hbar;
1049  coeff=(dk3/DK3)/tau;
1050  Amat->getElement(order,count)-=coeff;
1051  BVec[totalmodes]+=coeff*_eArray[cellIndex];
1052  BVec[totalmodes]-=coeff*_e0Array[cellIndex];
1053  cellIndex++;
1054  }
1055  }
1056  Amat->getElement(order,order)=tauTot;
1057  }
PhononMacro & _macro
T getde0taudT(const int c, T Tl)
Definition: Kspace.h:646
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
Field temperature
Definition: PhononMacro.h:21
T getDK3() const
Definition: Kspace.h:408
Array< T_Scalar > TArray
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
const T getTau(const int index)
Definition: Kspace.h:1231
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
Definition: kvol.h:14
T getdk3()
Definition: kvol.h:42
template<class T>
void COMETDiscretizer< T >::COMETFullScatt ( const int  cell,
TArray s,
TArray BVec 
)
inline

Definition at line 1062 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cellVolume, COMETDiscretizer< T >::_kspace, Kspace< T >::getGlobalIndex(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), Kspace< T >::getTau(), and Kspace< T >::gettotmodes().

Referenced by COMETDiscretizer< T >::findResidFull().

1063  {
1064  const int klen=_kspace.getlength();
1065  const int totalmodes=_kspace.gettotmodes();
1066  //const int order=totalmodes+1;
1067  //TArray& Tlold=dynamic_cast<TArray&>(_macro.temperature[_cells]);
1068  //const T DK3=_kspace.getDK3();
1069  int cellIndex=_kspace.getGlobalIndex(cell,0);
1070  TArray ds(totalmodes);
1071  //_kspace.getSourceTerm(cell,s,ds);
1072 
1073  const T relFac(1);
1074  //const T impl(1e0);
1075 
1076  for(int k=0;k<klen;k++)
1077  {
1078  Tkvol& kvol=_kspace.getkvol(k);
1079  //const T dk3=kvol.getdk3();
1080  const int numModes=kvol.getmodenum();
1081  for(int m=0;m<numModes;m++)
1082  {
1083  Tmode& mode=kvol.getmode(m);
1084  const T tau=_kspace.getTau(cellIndex);
1085  const int count=mode.getIndex();
1086  //const T coeff=_cellVolume[cell]/tau;
1087 
1088  BVec[count-1]+=(s[count-1])*_cellVolume[cell];
1089  //-(_e0Array[cellIndex]-_eArray[cellIndex])/tau)*_cellVolume[cell];
1090 
1091  BVec[count-1]*=relFac;
1092  cellIndex++;
1093  }
1094  }
1095  }
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
int getmodenum()
Definition: kvol.h:43
Array< T_Scalar > TArray
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
const TArray & _cellVolume
const T getTau(const int index)
Definition: Kspace.h:1231
int gettotmodes()
Definition: Kspace.h:393
Definition: kvol.h:14
template<class T>
void COMETDiscretizer< T >::COMETShifted ( const int  cell,
TMatrix Amat,
TArray BVec 
)
inline

Definition at line 1153 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_cellVolume, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, Kspace< T >::getde0taudT(), MatrixJML< T >::getElement(), pmode< T >::geteShifted(), pmode< T >::getfield(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), pmode< T >::gettauN(), Kspace< T >::gettotmodes(), and PhononMacro::temperature.

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), COMETDiscretizer< T >::COMETSolveFine(), and COMETDiscretizer< T >::COMETSolveFull().

1154  { //adds to collision and equilibrium
1155  const int klen=_kspace.getlength();
1156  const int totalmodes=_kspace.gettotmodes();
1157  const int order=totalmodes+1;
1158  TArray& Tlold=dynamic_cast<TArray&>(_macro.temperature[_cells]);
1159  const T tauTot=_kspace.getde0taudT(cell,Tlold[cell]);
1160  T coeff;
1161 
1162  for(int k=0;k<klen;k++)
1163  {
1164  Tkvol& kvol=_kspace.getkvol(k);
1165  const int numModes=kvol.getmodenum();
1166  for(int m=0;m<numModes;m++)
1167  {
1168  Tmode& mode=kvol.getmode(m);
1169  const int count=mode.getIndex();
1170  T tau=mode.gettauN();
1171  Field& efield=mode.getfield();
1172  Field& eShiftedfield=mode.geteShifted();
1173  TArray& eArray=dynamic_cast<TArray&>(efield[_cells]);
1174  TArray& eShifted=dynamic_cast<TArray&>(eShiftedfield[_cells]);
1175  coeff=_cellVolume[cell]/tau;
1176  Amat->getElement(count,count)-=coeff;
1177  Amat->getElement(order,count)+=coeff/tauTot;
1178  BVec[count-1]-=coeff*eArray[cell];
1179  BVec[count-1]+=coeff*eShifted[cell];
1180  BVec[totalmodes]+=coeff*eArray[cell]/tauTot;
1181  BVec[totalmodes]-=coeff*eShifted[cell]/tauTot;
1182  }
1183  }
1184 
1185  }
PhononMacro & _macro
T getde0taudT(const int c, T Tl)
Definition: Kspace.h:646
Definition: Field.h:14
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
Field temperature
Definition: PhononMacro.h:21
int getmodenum()
Definition: kvol.h:43
Array< T_Scalar > TArray
const TArray & _cellVolume
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
Definition: kvol.h:14
template<class T>
void COMETDiscretizer< T >::COMETSolveCoarse ( const int  dir,
const int  level 
)
inline

Definition at line 244 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_BCArray, COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_FaceToKSC, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, COMETDiscretizer< T >::_options, COMETDiscretizer< T >::addFAS(), COMETDiscretizer< T >::COMETCollision(), COMETDiscretizer< T >::COMETConvection(), COMETDiscretizer< T >::COMETConvectionCoarse(), COMETDiscretizer< T >::COMETEquilibrium(), COMETDiscretizer< T >::COMETShifted(), COMETDiscretizer< T >::COMETSource(), COMETDiscretizer< T >::correctInterface(), COMETDiscretizer< T >::Distribute(), fabs(), StorageSite::getSelfCount(), Kspace< T >::gettotmodes(), COMETModelOptions< T >::maxNewton, COMETModelOptions< T >::minNewton, COMETModelOptions< T >::NewtonTol, COMETDiscretizer< T >::scaledResid(), SquareMatrix< T, N >::Solve(), ArrowHeadMatrix< X, K >::Solve(), PhononMacro::temperature, COMETDiscretizer< T >::updatee0(), COMETDiscretizer< T >::updateGhostCoarse(), COMETModelOptions< T >::withNormal, SquareMatrix< T, N >::zero(), ArrowHeadMatrix< X, K >::zero(), and Array< T >::zero().

Referenced by COMETModel< T >::smooth().

245  {
246  const int cellcount=_cells.getSelfCount();
247  int start;
248 
249  if(dir==1)
250  start=0;
251  if(dir==-1)
252  start=cellcount-1;
253  const int totalmodes=_kspace.gettotmodes();
254  TArray Bvec(totalmodes+1);
255  TArray Resid(totalmodes+1);
256  TArrow AMat(totalmodes+1);
257  const T newTol=_options.NewtonTol;
258  const int maxNew=_options.maxNewton;
259  const int minNew=_options.minNewton;
260  TArray& Tl=dynamic_cast<TArray&>(_macro.temperature[_cells]);
261 
262  for(int c=start;((c<cellcount)&&(c>-1));c+=dir)
263  {
264 
265  if(_BCArray[c]==0) //no reflections at all--interior cell or temperature boundary
266  {
267  T dt=1;
268  T rscaled=10;
269  int NewtonIters=0;
270 
271  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
272  {
273 
274  Bvec.zero();
275  Resid.zero();
276  AMat.zero();
277 
278  COMETConvectionCoarse(c,AMat,Bvec);
279  COMETCollision(c,&AMat,Bvec);
280  COMETEquilibrium(c,&AMat,Bvec);
281  COMETSource(c,Bvec);
282 
283  if(_options.withNormal)
284  COMETShifted(c,&AMat,Bvec);
285 
286  if(level>0)
287  addFAS(c,Bvec);
288 
289  Resid=Bvec;
290  AMat.Solve(Bvec);
291  rscaled=scaledResid(Bvec,c);
292 
293  Distribute(c,Bvec,Resid);
294  updatee0(c);
295  if(!_FaceToKSC.empty())
296  correctInterface(c,Bvec);
297  NewtonIters++;
298  dt=0.5*fabs(Bvec[totalmodes]/Tl[c])+rscaled*0.5;
299  }
300 
301  }
302  else if(_BCArray[c]==1) //Implicit reflecting boundary only
303  {
304  T dt=1;
305  int NewtonIters=0;
306  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
307  {
308  TSquare AMatS(totalmodes+1);
309  Bvec.zero();
310  Resid.zero();
311  AMatS.zero();
312 
313  COMETConvection(c,AMatS,Bvec);
314  COMETCollision(c,&AMatS,Bvec);
315  COMETEquilibrium(c,&AMatS,Bvec);
316 
317  if(level>0)
318  addFAS(c,Bvec);
319 
320  Resid=Bvec;
321  AMatS.Solve(Bvec);
322 
323  Distribute(c,Bvec,Resid);
324  updatee0(c);
325  NewtonIters++;
326  dt=fabs(Bvec[totalmodes]);
327  }
328  }
329  else if(_BCArray[c]==2) //Explicit boundary only
330  {
331  T dt=1;
332  T rscaled=10;
333  int NewtonIters=0;
335 
336  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
337  {
338 
339  Bvec.zero();
340  Resid.zero();
341  AMat.zero();
342 
343  COMETConvectionCoarse(c,AMat,Bvec);
344  COMETCollision(c,&AMat,Bvec);
345  COMETEquilibrium(c,&AMat,Bvec);
346  COMETSource(c,Bvec);
347 
348  if(_options.withNormal)
349  COMETShifted(c,&AMat,Bvec);
350 
351  if(level>0)
352  addFAS(c,Bvec);
353 
354  Resid=Bvec;
355  AMat.Solve(Bvec);
356  rscaled=scaledResid(Bvec,c);
357 
358  Distribute(c,Bvec,Resid);
359  dt=0.5*fabs(Bvec[totalmodes]/Tl[c])+0.5*rscaled;
360  updatee0(c);
362  if(!_FaceToKSC.empty())
363  correctInterface(c,Bvec);
364  NewtonIters++;
365  }
366  }
367  else if(_BCArray[c]==3) //Mix Implicit/Explicit reflecting boundary
368  {
369  T dt=1;
370  int NewtonIters=0;
372  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
373  {
374  TSquare AMatS(totalmodes+1);
375  Bvec.zero();
376  Resid.zero();
377  AMatS.zero();
378 
379  COMETConvection(c,AMatS,Bvec);
380  COMETCollision(c,&AMatS,Bvec);
381  COMETEquilibrium(c,&AMatS,Bvec);
382 
383  if(level>0)
384  addFAS(c,Bvec);
385 
386  Resid=Bvec;
387  AMatS.Solve(Bvec);
388 
389  Distribute(c,Bvec,Resid);
390  updatee0(c);
392  NewtonIters++;
393  dt=fabs(Bvec[totalmodes]);
394  }
395  }
396  else
397  throw CException("Unexpected value for boundary cell map.");
398 
399  }
400 
401  }
int getSelfCount() const
Definition: StorageSite.h:40
PhononMacro & _macro
void COMETConvection(const int cell, TSquare &Amat, TArray &BVec)
ArrowHeadMatrix< T > TArrow
const IntArray & _BCArray
Field temperature
Definition: PhononMacro.h:21
void Distribute(const int cell, TArray &BVec, TArray &Rvec)
void COMETSource(const int cell, TArray &BVec)
T scaledResid(const TArray &de, const int c)
void COMETCollision(const int cell, TMatrix *Amat, TArray &BVec)
Array< T_Scalar > TArray
void addFAS(const int c, TArray &bVec)
SquareMatrix< T > TSquare
void correctInterface(const int cell0, TArray &Bvec)
void COMETEquilibrium(const int cell, TMatrix *Amat, TArray &BVec)
Tangent fabs(const Tangent &a)
Definition: Tangent.h:312
void COMETConvectionCoarse(const int cell0, TArrow &Amat, TArray &BVec)
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
void updateGhostCoarse(const int cell)
void COMETShifted(const int cell, TMatrix *Amat, TArray &BVec)
const FgTKClistMap & _FaceToKSC
template<class T>
void COMETDiscretizer< T >::COMETSolveFine ( const int  dir,
const int  level 
)
inline

Definition at line 91 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_BCArray, COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_FaceToKSC, COMETDiscretizer< T >::_geomFields, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_mesh, COMETDiscretizer< T >::_options, COMETDiscretizer< T >::addFAS(), COMETDiscretizer< T >::COMETCollision(), COMETDiscretizer< T >::COMETConvection(), COMETDiscretizer< T >::COMETConvectionFine(), COMETDiscretizer< T >::COMETEquilibrium(), COMETDiscretizer< T >::COMETShifted(), COMETDiscretizer< T >::correctInterface(), COMETDiscretizer< T >::Distribute(), fabs(), GradientModel< T >::getGradientMatrix(), StorageSite::getSelfCount(), Kspace< T >::gettotmodes(), COMETModelOptions< T >::maxNewton, COMETModelOptions< T >::minNewton, COMETModelOptions< T >::NewtonTol, SquareMatrix< T, N >::Solve(), ArrowHeadMatrix< X, K >::Solve(), COMETDiscretizer< T >::updatee0(), COMETDiscretizer< T >::updateGhostFine(), COMETModelOptions< T >::withNormal, SquareMatrix< T, N >::zero(), ArrowHeadMatrix< X, K >::zero(), and Array< T >::zero().

Referenced by COMETModel< T >::smooth().

92  {
93  const int cellcount=_cells.getSelfCount();
94  int start;
95 
96  if(dir==1)
97  start=0;
98  if(dir==-1)
99  start=cellcount-1;
100  const int totalmodes=_kspace.gettotmodes();
101  TArray Bvec(totalmodes+1);
102  TArray Resid(totalmodes+1);
103  TArrow AMat(totalmodes+1);
104  const T newTol=_options.NewtonTol;
105  const int maxNew=_options.maxNewton;
106  const int minNew=_options.minNewton;
107 
109 
110  for(int c=start;((c<cellcount)&&(c>-1));c+=dir)
111  {
112 
113  if(_BCArray[c]==0) //no reflections at all--interior cell or temperature boundary
114  {
115  T dt=1;
116  int NewtonIters=0;
117 
118  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
119  {
120 
121  Bvec.zero();
122  Resid.zero();
123  AMat.zero();
124 
125  updateGhostFine(c, gradMatrix);
126  COMETConvectionFine(c,AMat,Bvec,gradMatrix);
127  COMETCollision(c,&AMat,Bvec);
128  COMETEquilibrium(c,&AMat,Bvec);
129 
130  if(_options.withNormal)
131  COMETShifted(c,&AMat,Bvec);
132 
133  if(level>0)
134  addFAS(c,Bvec);
135 
136  Resid=Bvec;
137  AMat.Solve(Bvec);
138 
139  Distribute(c,Bvec,Resid);
140  updatee0(c);
141  NewtonIters++;
142  dt=fabs(Bvec[totalmodes]);
143  }
144  }
145  else if(_BCArray[c]==1) //Implicit reflecting boundary only
146  {
147  T dt=1;
148  int NewtonIters=0;
149  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
150  {
151  TSquare AMatS(totalmodes+1);
152  Bvec.zero();
153  Resid.zero();
154  AMatS.zero();
155 
156  COMETConvection(c,AMatS,Bvec);
157  COMETCollision(c,&AMatS,Bvec);
158  COMETEquilibrium(c,&AMatS,Bvec);
159 
160  if(level>0)
161  addFAS(c,Bvec);
162 
163  Resid=Bvec;
164  AMatS.Solve(Bvec);
165 
166  Distribute(c,Bvec,Resid);
167  updatee0(c);
168  NewtonIters++;
169  dt=fabs(Bvec[totalmodes]);
170  }
171  }
172  else if(_BCArray[c]==2) //Explicit boundary only
173  {
174  T dt=1;
175  int NewtonIters=0;
176  updateGhostFine(c, gradMatrix);
177  //updateGhostCoarse(c);
178 
179  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
180  {
181 
182  Bvec.zero();
183  Resid.zero();
184  AMat.zero();
185 
186  COMETConvectionFine(c,AMat,Bvec,gradMatrix);
187  COMETCollision(c,&AMat,Bvec);
188  COMETEquilibrium(c,&AMat,Bvec);
189 
190  if(_options.withNormal)
191  COMETShifted(c,&AMat,Bvec);
192 
193  if(level>0)
194  addFAS(c,Bvec);
195 
196  Resid=Bvec;
197  AMat.Solve(Bvec);
198 
199  Distribute(c,Bvec,Resid);
200  dt=fabs(Bvec[totalmodes]);
201  updatee0(c);
202  updateGhostFine(c, gradMatrix);
203  if(!_FaceToKSC.empty())
204  correctInterface(c,Bvec);
205  NewtonIters++;
206  }
207  }
208  else if(_BCArray[c]==3) //Mix Implicit/Explicit reflecting boundary
209  {
210  T dt=1;
211  int NewtonIters=0;
212  updateGhostFine(c,gradMatrix);
213  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
214  {
215  TSquare AMatS(totalmodes+1);
216  Bvec.zero();
217  Resid.zero();
218  AMatS.zero();
219 
220  COMETConvection(c,AMatS,Bvec);
221  COMETCollision(c,&AMatS,Bvec);
222  COMETEquilibrium(c,&AMatS,Bvec);
223 
224  if(level>0)
225  addFAS(c,Bvec);
226 
227  Resid=Bvec;
228  AMatS.Solve(Bvec);
229 
230  Distribute(c,Bvec,Resid);
231  updatee0(c);
232  updateGhostFine(c, gradMatrix);
233  NewtonIters++;
234  dt=fabs(Bvec[totalmodes]);
235  }
236  }
237  else
238  throw CException("Unexpected value for boundary cell map.");
239 
240  }
241 
242  }
int getSelfCount() const
Definition: StorageSite.h:40
void COMETConvection(const int cell, TSquare &Amat, TArray &BVec)
ArrowHeadMatrix< T > TArrow
const IntArray & _BCArray
void COMETConvectionFine(const int cell0, TArrow &Amat, TArray &BVec, const GradMatrix &gMat)
void Distribute(const int cell, TArray &BVec, TArray &Rvec)
void COMETCollision(const int cell, TMatrix *Amat, TArray &BVec)
Array< T_Scalar > TArray
void addFAS(const int c, TArray &bVec)
SquareMatrix< T > TSquare
void correctInterface(const int cell0, TArray &Bvec)
static GradMatrixType & getGradientMatrix(const Mesh &mesh, const GeomFields &geomFields)
void COMETEquilibrium(const int cell, TMatrix *Amat, TArray &BVec)
Tangent fabs(const Tangent &a)
Definition: Tangent.h:312
const GeomFields & _geomFields
int gettotmodes()
Definition: Kspace.h:393
void updateGhostFine(const int cell, const GradMatrix &gMat)
const StorageSite & _cells
GradModelType::GradMatrixType GradMatrix
void COMETShifted(const int cell, TMatrix *Amat, TArray &BVec)
const FgTKClistMap & _FaceToKSC
template<class T>
void COMETDiscretizer< T >::COMETSolveFull ( const int  dir,
const int  level 
)
inline

Definition at line 403 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_BCArray, COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_FaceToKSC, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_options, COMETDiscretizer< T >::addFAS(), COMETDiscretizer< T >::COMETCollision(), COMETDiscretizer< T >::COMETConvection(), COMETDiscretizer< T >::COMETConvectionCoarse(), COMETDiscretizer< T >::COMETEquilibrium(), COMETDiscretizer< T >::COMETShifted(), COMETDiscretizer< T >::correctInterface(), COMETDiscretizer< T >::Distribute(), fabs(), StorageSite::getSelfCount(), Kspace< T >::gettotmodes(), COMETModelOptions< T >::maxNewton, COMETModelOptions< T >::minNewton, COMETModelOptions< T >::NewtonTol, COMETDiscretizer< T >::ScatterPhonons(), SquareMatrix< T, N >::Solve(), ArrowHeadMatrix< X, K >::Solve(), COMETDiscretizer< T >::updatee0(), COMETDiscretizer< T >::updateGhostCoarse(), COMETModelOptions< T >::withNormal, SquareMatrix< T, N >::zero(), ArrowHeadMatrix< X, K >::zero(), and Array< T >::zero().

Referenced by COMETModel< T >::smooth().

404  {
405  const int cellcount=_cells.getSelfCount();
406  int start;
407 
408  if(dir==1)
409  start=0;
410  if(dir==-1)
411  start=cellcount-1;
412  const int totalmodes=_kspace.gettotmodes();
413  TArray Bvec(totalmodes+1);
414  TArray Resid(totalmodes+1);
415  TArray s(totalmodes);
416  TArray ds(totalmodes);
417  TArrow AMat(totalmodes+1);
418  const T newTol=_options.NewtonTol;
419  const int maxNew=_options.maxNewton;
420  const int minNew=_options.minNewton;
421 
422  for(int c=start;((c<cellcount)&&(c>-1));c+=dir)
423  {
424 
425  if(_BCArray[c]==0) //no reflections at all--interior cell or temperature boundary
426  {
427  T dt=1;
428  int NewtonIters=0;
429 
430  for(int srcIts=0;srcIts<1;srcIts++)
431  {
432 
433  s.zero();
434  //_kspace.getSourceTerm(c,s,ds);
435 
436  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
437  {
438 
439  Bvec.zero();
440  Resid.zero();
441  AMat.zero();
442 
443  COMETConvectionCoarse(c,AMat,Bvec);
444  COMETCollision(c,&AMat,Bvec);
445  COMETEquilibrium(c,&AMat,Bvec);
446  //COMETFullScatt(c,s,Bvec);
447 
448  if(_options.withNormal)
449  COMETShifted(c,&AMat,Bvec);
450 
451  if(level>0)
452  addFAS(c,Bvec);
453 
454  Resid=Bvec;
455  AMat.Solve(Bvec);
456  //AMat.SolveDiag(Bvec);
457  //AMat.smoothGS(Bvec);
458 
459  Distribute(c,Bvec,Resid);
460  updatee0(c);
461  NewtonIters++;
462  dt=fabs(Bvec[totalmodes]);
463  }
464 
465  ScatterPhonons(c);
466 
467  }
468  }
469  else if(_BCArray[c]==1) //Implicit reflecting boundary only
470  {
471  T dt=1;
472  int NewtonIters=0;
473  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
474  {
475  TSquare AMatS(totalmodes+1);
476  Bvec.zero();
477  Resid.zero();
478  AMatS.zero();
479 
480  COMETConvection(c,AMatS,Bvec);
481  COMETCollision(c,&AMatS,Bvec);
482  COMETEquilibrium(c,&AMatS,Bvec);
483 
484  if(level>0)
485  addFAS(c,Bvec);
486 
487  Resid=Bvec;
488  AMatS.Solve(Bvec);
489 
490  Distribute(c,Bvec,Resid);
491  updatee0(c);
492  NewtonIters++;
493  dt=fabs(Bvec[totalmodes]);
494  }
495  }
496  else if(_BCArray[c]==2) //Explicit boundary only
497  {
498  T dt=1;
499  int NewtonIters=0;
501 
502  for(int srcIts=0;srcIts<1;srcIts++)
503  {
504  s.zero();
505  //_kspace.getSourceTerm(c,s,ds);
506 
507  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
508  {
509 
510  Bvec.zero();
511  Resid.zero();
512  AMat.zero();
513 
514  COMETConvectionCoarse(c,AMat,Bvec);
515  COMETCollision(c,&AMat,Bvec);
516  COMETEquilibrium(c,&AMat,Bvec);
517  //COMETFullScatt(c,s,Bvec);
518 
519  if(_options.withNormal)
520  COMETShifted(c,&AMat,Bvec);
521 
522  if(level>0)
523  addFAS(c,Bvec);
524 
525  Resid=Bvec;
526  AMat.Solve(Bvec);
527  //AMat.SolveDiag(Bvec);
528  //AMat.smoothGS(Bvec);
529 
530  Distribute(c,Bvec,Resid);
531  dt=fabs(Bvec[totalmodes]);
532  updatee0(c);
534  if(!_FaceToKSC.empty())
535  correctInterface(c,Bvec);
536  NewtonIters++;
537  }
538 
539  ScatterPhonons(c);
540 
541  }
542  }
543  else if(_BCArray[c]==3) //Mix Implicit/Explicit reflecting boundary
544  {
545  T dt=1;
546  int NewtonIters=0;
548  while((dt>newTol && NewtonIters<maxNew) || NewtonIters<minNew)
549  {
550  TSquare AMatS(totalmodes+1);
551  Bvec.zero();
552  Resid.zero();
553  AMatS.zero();
554 
555  COMETConvection(c,AMatS,Bvec);
556  COMETCollision(c,&AMatS,Bvec);
557  COMETEquilibrium(c,&AMatS,Bvec);
558 
559  if(level>0)
560  addFAS(c,Bvec);
561 
562  Resid=Bvec;
563  AMatS.Solve(Bvec);
564 
565  Distribute(c,Bvec,Resid);
566  updatee0(c);
568  NewtonIters++;
569  dt=fabs(Bvec[totalmodes]);
570  }
571  }
572  else
573  throw CException("Unexpected value for boundary cell map.");
574 
575  }
576 
577  }
int getSelfCount() const
Definition: StorageSite.h:40
void COMETConvection(const int cell, TSquare &Amat, TArray &BVec)
ArrowHeadMatrix< T > TArrow
const IntArray & _BCArray
void ScatterPhonons(const int cell0)
void Distribute(const int cell, TArray &BVec, TArray &Rvec)
void COMETCollision(const int cell, TMatrix *Amat, TArray &BVec)
Array< T_Scalar > TArray
void addFAS(const int c, TArray &bVec)
SquareMatrix< T > TSquare
void correctInterface(const int cell0, TArray &Bvec)
void COMETEquilibrium(const int cell, TMatrix *Amat, TArray &BVec)
Tangent fabs(const Tangent &a)
Definition: Tangent.h:312
void COMETConvectionCoarse(const int cell0, TArrow &Amat, TArray &BVec)
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
void updateGhostCoarse(const int cell)
void COMETShifted(const int cell, TMatrix *Amat, TArray &BVec)
const FgTKClistMap & _FaceToKSC
template<class T>
void COMETDiscretizer< T >::COMETSource ( const int  cell,
TArray BVec 
)
inline

Definition at line 1059 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cellVolume, COMETDiscretizer< T >::_kspace, and Kspace< T >::addSource().

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), and COMETDiscretizer< T >::findResidCoarse().

1060  {_kspace.addSource(cell, BVec,_cellVolume[cell]);}
const TArray & _cellVolume
void addSource(const int c, TArray &BVec, const T cv)
Definition: Kspace.h:1693
template<class T>
void COMETDiscretizer< T >::correctInterface ( const int  cell0,
TArray Bvec 
)
inline

Definition at line 2032 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_BCfArray, COMETDiscretizer< T >::_cellFaces, COMETDiscretizer< T >::_faceCells, COMETDiscretizer< T >::_FaceToKSC, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_mesh, COMETDiscretizer< T >::Distribute(), COMETDiscretizer< T >::findFgId(), CRConnectivity::getCount(), Kspace< T >::getDK3(), Mesh::getFaceGroup(), StorageSite::getOffset(), Kspace< T >::gettotmodes(), KSConnectivity< T >::multiplySelf(), FaceGroup::site, and Array< T >::zero().

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), COMETDiscretizer< T >::COMETSolveFine(), and COMETDiscretizer< T >::COMETSolveFull().

2033  {
2034  const int klen=_kspace.gettotmodes();
2035  const T DK3=_kspace.getDK3();
2036  TArray Correction(klen+1);
2037 
2038  const int neibcount=_cellFaces.getCount(cell0);
2039  for(int j=0;j<neibcount;j++)
2040  {
2041  const int f=_cellFaces(cell0,j);
2042  if(_BCfArray[f]==4)
2043  {
2044  int Fgid=findFgId(f);
2045  const FaceGroup& fg=_mesh.getFaceGroup(Fgid);
2046  const StorageSite& faces0=fg.site;
2047  const int offset=faces0.getOffset();
2048  int cell1=_faceCells(f,1);
2049  if(cell1==cell0)
2050  cell1=_faceCells(f,0);
2051  const TKConnectivity& TKC=*(_FaceToKSC.find(Fgid)->second)[f-offset];
2052  TKC.multiplySelf(Bvec,Correction,DK3);
2053  Bvec.zero();
2054  Distribute(cell1,Correction, Bvec);
2055  }
2056  }
2057 
2058  }
int getCount(const int i) const
const CRConnectivity & _faceCells
Definition: Mesh.h:28
KSConnectivity< T > TKConnectivity
void Distribute(const int cell, TArray &BVec, TArray &Rvec)
const FaceGroup & getFaceGroup(const int fgId) const
Definition: Mesh.cpp:1570
T getDK3() const
Definition: Kspace.h:408
Array< T_Scalar > TArray
const IntArray & _BCfArray
const CRConnectivity & _cellFaces
int getOffset() const
Definition: StorageSite.h:87
int gettotmodes()
Definition: Kspace.h:393
const FgTKClistMap & _FaceToKSC
StorageSite site
Definition: Mesh.h:40
int findFgId(const int faceIndex)
template<class T>
void COMETDiscretizer< T >::Distribute ( const int  cell,
TArray BVec,
TArray Rvec 
)
inline

Definition at line 1187 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_eArray, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, COMETDiscretizer< T >::_resArray, PhononMacro::deltaT, Kspace< T >::getGlobalIndex(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), Kspace< T >::gettotmodes(), PhononMacro::temperature, and PhononMacro::TlResidual.

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), COMETDiscretizer< T >::COMETSolveFine(), COMETDiscretizer< T >::COMETSolveFull(), COMETDiscretizer< T >::correctInterface(), COMETDiscretizer< T >::findResidCoarse(), COMETDiscretizer< T >::findResidFine(), and COMETDiscretizer< T >::findResidFull().

1188  {
1189  const int klen=_kspace.getlength();
1190  const int totalmodes=_kspace.gettotmodes();
1191  int cellIndex=_kspace.getGlobalIndex(cell,0);
1192 
1193  for(int k=0;k<klen;k++)
1194  {
1195  Tkvol& kvol=_kspace.getkvol(k);
1196  const int numModes=kvol.getmodenum();
1197  for(int m=0;m<numModes;m++)
1198  {
1199  Tmode& mode=kvol.getmode(m);
1200  const int count=mode.getIndex()-1;
1201  _eArray[cellIndex]+=BVec[count];
1202  _resArray[cellIndex]=-Rvec[count];
1203  cellIndex++;
1204  }
1205  }
1206 
1207  TArray& TlArray=dynamic_cast<TArray&>(_macro.temperature[_cells]);
1208  TlArray[cell]+=BVec[totalmodes];
1209  TArray& deltaTArray=dynamic_cast<TArray&>(_macro.deltaT[_cells]);
1210  deltaTArray[cell]=BVec[totalmodes];
1211  TArray& TlResArray=dynamic_cast<TArray&>(_macro.TlResidual[_cells]);
1212  TlResArray[cell]=-Rvec[totalmodes];
1213  }
PhononMacro & _macro
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
Field deltaT
Definition: PhononMacro.h:22
Field temperature
Definition: PhononMacro.h:21
int getmodenum()
Definition: kvol.h:43
Array< T_Scalar > TArray
Field TlResidual
Definition: PhononMacro.h:24
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
Definition: kvol.h:14
template<class T>
int COMETDiscretizer< T >::findFgId ( const int  faceIndex)
inline

Definition at line 1563 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_fgFinder.

Referenced by COMETDiscretizer< T >::COMETConvection(), COMETDiscretizer< T >::correctInterface(), COMETDiscretizer< T >::updateGhostCoarse(), and COMETDiscretizer< T >::updateGhostFine().

1564  {
1565  FaceToFg::iterator id;
1566  for(id=_fgFinder.begin();id!=_fgFinder.end();id++)
1567  {
1568  if(id->second[0]<=faceIndex && id->second[1]>faceIndex)
1569  return id->first;
1570  }
1571  cout<<"Face index: "<<faceIndex<<endl;
1572  throw CException("Didn't find matching FaceGroup!");
1573  return -1;
1574  }
template<class T>
void COMETDiscretizer< T >::findResidCoarse ( const bool  plusFAS)
inline

Definition at line 1314 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_aveResid, COMETDiscretizer< T >::_BCArray, COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_residChange, COMETDiscretizer< T >::addFAS(), COMETDiscretizer< T >::ArrayAbs(), COMETDiscretizer< T >::COMETCollision(), COMETDiscretizer< T >::COMETConvection(), COMETDiscretizer< T >::COMETConvectionCoarse(), COMETDiscretizer< T >::COMETEquilibrium(), COMETDiscretizer< T >::COMETSource(), COMETDiscretizer< T >::Distribute(), fabs(), StorageSite::getSelfCount(), Kspace< T >::gettotmodes(), ArrowHeadMatrix< X, K >::getTraceAbs(), SquareMatrix< T, N >::getTraceAbs(), COMETDiscretizer< T >::makeValueArray(), ArrowHeadMatrix< X, K >::multiply(), SquareMatrix< T, N >::multiply(), COMETDiscretizer< T >::updateGhostCoarse(), SquareMatrix< T, N >::zero(), ArrowHeadMatrix< X, K >::zero(), and Array< T >::zero().

Referenced by COMETModel< T >::updateResid().

1315  {
1316  const int cellcount=_cells.getSelfCount();
1317  const int totalmodes=_kspace.gettotmodes();
1318  TArray ResidSum(totalmodes+1);
1319  TArray Bsum(totalmodes+1);
1320  T ResidScalar=0.;
1321  T traceSum=0.;
1322  TArray Bvec(totalmodes+1);
1323  TArray Resid(totalmodes+1);
1324  TArray Dummy(totalmodes+1);
1325  TArrow AMat(totalmodes+1);
1326  ResidSum.zero();
1327  Bsum.zero();
1328  Dummy.zero();
1329 
1330  for(int c=0;c<cellcount;c++)
1331  {
1332  Bvec.zero();
1333  Resid.zero();
1334  Dummy.zero();
1335 
1336  if(_BCArray[c]==0 || _BCArray[c]==2) //Arrowhead
1337  {
1338  if(_BCArray[c]==2)
1339  updateGhostCoarse(c);
1340 
1341  AMat.zero();
1342 
1343  COMETConvectionCoarse(c,AMat,Bvec);
1344  COMETCollision(c,&AMat,Bvec);
1345  COMETEquilibrium(c,&AMat,Bvec);
1346  COMETSource(c,Bvec);
1347 
1348  if(plusFAS)
1349  addFAS(c,Bvec);
1350 
1351  traceSum+=AMat.getTraceAbs();
1352 
1353  Resid=Bvec;
1354  Bvec.zero();
1355  Distribute(c,Bvec,Resid);
1356  ArrayAbs(Resid);
1357  ResidSum+=Resid;
1358 
1359  AMat.multiply(Resid,Bvec);
1360  Resid=Bvec;
1361 
1362  makeValueArray(c,Bvec);
1363  AMat.multiply(Bvec,Dummy);
1364  Bvec=Dummy;
1365  ArrayAbs(Bvec);
1366  ArrayAbs(Resid);
1367  Bsum+=Bvec;
1368  }
1369  else if(_BCArray[c]==1 || _BCArray[c]==3) //General Dense
1370  {
1371  TSquare AMatS(totalmodes+1);
1372  AMatS.zero();
1373 
1374  COMETConvection(c,AMatS,Bvec);
1375  COMETCollision(c,&AMatS,Bvec);
1376  COMETEquilibrium(c,&AMatS,Bvec);
1377 
1378  if(plusFAS)
1379  addFAS(c,Bvec);
1380 
1381  traceSum+=AMatS.getTraceAbs();
1382  Resid=Bvec;
1383  Bvec.zero();
1384  Distribute(c,Bvec,Resid);
1385 
1386  AMatS.multiply(Resid,Bvec);
1387  Resid=Bvec;
1388 
1389  makeValueArray(c,Bvec);
1390  ArrayAbs(Resid);
1391  ArrayAbs(Bvec);
1392  Bsum+=Bvec;
1393  ResidSum+=Resid;
1394  }
1395  else
1396  throw CException("Unexpected value for boundary cell map.");
1397  }
1398 
1399  //traceSum=0; //added
1400  for(int o=0;o<totalmodes+1;o++)
1401  {
1402  ResidScalar+=ResidSum[o];
1403  //traceSum+=Bsum[o]; //added
1404  }
1405 
1406  ResidScalar/=traceSum;
1407 
1408  if(_aveResid==-1)
1409  {_aveResid=ResidScalar;}
1410  else
1411  {
1412  _residChange=fabs(_aveResid-ResidScalar)/_aveResid;
1413  _aveResid=ResidScalar;
1414  }
1415  }
void makeValueArray(const int c, TArray &o)
int getSelfCount() const
Definition: StorageSite.h:40
void COMETConvection(const int cell, TSquare &Amat, TArray &BVec)
ArrowHeadMatrix< T > TArrow
const IntArray & _BCArray
void Distribute(const int cell, TArray &BVec, TArray &Rvec)
void COMETSource(const int cell, TArray &BVec)
void COMETCollision(const int cell, TMatrix *Amat, TArray &BVec)
void ArrayAbs(TArray &o)
Array< T_Scalar > TArray
void addFAS(const int c, TArray &bVec)
SquareMatrix< T > TSquare
void COMETEquilibrium(const int cell, TMatrix *Amat, TArray &BVec)
Tangent fabs(const Tangent &a)
Definition: Tangent.h:312
void COMETConvectionCoarse(const int cell0, TArrow &Amat, TArray &BVec)
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
void updateGhostCoarse(const int cell)
template<class T>
void COMETDiscretizer< T >::findResidFine ( )
inline

Definition at line 1215 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_aveResid, COMETDiscretizer< T >::_BCArray, COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_geomFields, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_mesh, COMETDiscretizer< T >::_residChange, COMETDiscretizer< T >::ArrayAbs(), COMETDiscretizer< T >::COMETCollision(), COMETDiscretizer< T >::COMETConvection(), COMETDiscretizer< T >::COMETConvectionFine(), COMETDiscretizer< T >::COMETEquilibrium(), COMETDiscretizer< T >::Distribute(), fabs(), GradientModel< T >::getGradientMatrix(), StorageSite::getSelfCount(), Kspace< T >::gettotmodes(), ArrowHeadMatrix< X, K >::getTraceAbs(), SquareMatrix< T, N >::getTraceAbs(), COMETDiscretizer< T >::makeValueArray(), ArrowHeadMatrix< X, K >::multiply(), SquareMatrix< T, N >::multiply(), SquareMatrix< T, N >::zero(), ArrowHeadMatrix< X, K >::zero(), and Array< T >::zero().

Referenced by COMETModel< T >::updateResid().

1216  {
1217  const int cellcount=_cells.getSelfCount();
1218  const int totalmodes=_kspace.gettotmodes();
1219  TArray ResidSum(totalmodes+1);
1220  TArray Bsum(totalmodes+1);
1221  T ResidScalar=0.;
1222  T traceSum=0.;
1223  TArray Bvec(totalmodes+1);
1224  TArray Resid(totalmodes+1);
1225  TArray Dummy(totalmodes+1);
1226  TArrow AMat(totalmodes+1);
1227  ResidSum.zero();
1228  Bsum.zero();
1229  Dummy.zero();
1230 
1232 
1233  for(int c=0;c<cellcount;c++)
1234  {
1235  Bvec.zero();
1236  Resid.zero();
1237  Dummy.zero();
1238 
1239  if(_BCArray[c]==0 || _BCArray[c]==2) //Arrowhead
1240  {
1241 
1242  //updateGhostFine(c, gradMatrix);
1243  //updateGhostCoarse(c);
1244 
1245  AMat.zero();
1246 
1247  COMETConvectionFine(c,AMat,Bvec,gradMatrix);
1248  COMETCollision(c,&AMat,Bvec);
1249  COMETEquilibrium(c,&AMat,Bvec);
1250 
1251  traceSum+=AMat.getTraceAbs();
1252 
1253  Resid=Bvec;
1254  Bvec.zero();
1255  Distribute(c,Bvec,Resid);
1256  ArrayAbs(Resid);
1257  ResidSum+=Resid;
1258 
1259  AMat.multiply(Resid,Bvec);
1260  Resid=Bvec;
1261 
1262  makeValueArray(c,Bvec);
1263  AMat.multiply(Bvec,Dummy);
1264  Bvec=Dummy;
1265  ArrayAbs(Bvec);
1266  ArrayAbs(Resid);
1267  Bsum+=Bvec;
1268  }
1269  else if(_BCArray[c]==1 || _BCArray[c]==3) //General Dense
1270  {
1271  TSquare AMatS(totalmodes+1);
1272  AMatS.zero();
1273 
1274  COMETConvection(c,AMatS,Bvec);
1275  COMETCollision(c,&AMatS,Bvec);
1276  COMETEquilibrium(c,&AMatS,Bvec);
1277 
1278  traceSum+=AMatS.getTraceAbs();
1279  Resid=Bvec;
1280  Bvec.zero();
1281  Distribute(c,Bvec,Resid);
1282 
1283  AMatS.multiply(Resid,Bvec);
1284  Resid=Bvec;
1285 
1286  makeValueArray(c,Bvec);
1287  ArrayAbs(Resid);
1288  ArrayAbs(Bvec);
1289  Bsum+=Bvec;
1290  ResidSum+=Resid;
1291  }
1292  else
1293  throw CException("Unexpected value for boundary cell map.");
1294  }
1295 
1296  //traceSum=0; //added
1297  for(int o=0;o<totalmodes+1;o++)
1298  {
1299  ResidScalar+=ResidSum[o];
1300  //traceSum+=Bsum[o]; //added
1301  }
1302 
1303  ResidScalar/=traceSum;
1304 
1305  if(_aveResid==-1)
1306  {_aveResid=ResidScalar;}
1307  else
1308  {
1309  _residChange=fabs(_aveResid-ResidScalar)/_aveResid;
1310  _aveResid=ResidScalar;
1311  }
1312  }
void makeValueArray(const int c, TArray &o)
int getSelfCount() const
Definition: StorageSite.h:40
void COMETConvection(const int cell, TSquare &Amat, TArray &BVec)
ArrowHeadMatrix< T > TArrow
const IntArray & _BCArray
void COMETConvectionFine(const int cell0, TArrow &Amat, TArray &BVec, const GradMatrix &gMat)
void Distribute(const int cell, TArray &BVec, TArray &Rvec)
void COMETCollision(const int cell, TMatrix *Amat, TArray &BVec)
void ArrayAbs(TArray &o)
Array< T_Scalar > TArray
SquareMatrix< T > TSquare
static GradMatrixType & getGradientMatrix(const Mesh &mesh, const GeomFields &geomFields)
void COMETEquilibrium(const int cell, TMatrix *Amat, TArray &BVec)
Tangent fabs(const Tangent &a)
Definition: Tangent.h:312
const GeomFields & _geomFields
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
GradModelType::GradMatrixType GradMatrix
template<class T>
void COMETDiscretizer< T >::findResidFull ( )
inline

Definition at line 1417 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_aveResid, COMETDiscretizer< T >::_BCArray, COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_residChange, COMETDiscretizer< T >::ArrayAbs(), COMETDiscretizer< T >::COMETCollision(), COMETDiscretizer< T >::COMETConvection(), COMETDiscretizer< T >::COMETConvectionCoarse(), COMETDiscretizer< T >::COMETEquilibrium(), COMETDiscretizer< T >::COMETFullScatt(), COMETDiscretizer< T >::Distribute(), fabs(), StorageSite::getSelfCount(), Kspace< T >::getSourceTerm(), Kspace< T >::gettotmodes(), ArrowHeadMatrix< X, K >::getTraceAbs(), SquareMatrix< T, N >::getTraceAbs(), COMETDiscretizer< T >::makeValueArray(), ArrowHeadMatrix< X, K >::multiply(), SquareMatrix< T, N >::multiply(), COMETDiscretizer< T >::updateGhostCoarse(), SquareMatrix< T, N >::zero(), ArrowHeadMatrix< X, K >::zero(), and Array< T >::zero().

Referenced by COMETModel< T >::updateResid().

1418  {
1419  const int cellcount=_cells.getSelfCount();
1420  const int totalmodes=_kspace.gettotmodes();
1421  TArray ResidSum(totalmodes+1);
1422  TArray Bsum(totalmodes+1);
1423  T ResidScalar=0.;
1424  T traceSum=0.;
1425  TArray Bvec(totalmodes+1);
1426  TArray Resid(totalmodes+1);
1427  TArray Dummy(totalmodes+1);
1428  TArray s(totalmodes);
1429  TArray ds(totalmodes);
1430  TArrow AMat(totalmodes+1);
1431  ResidSum.zero();
1432  Bsum.zero();
1433  Dummy.zero();
1434 
1435  for(int c=0;c<cellcount;c++)
1436  {
1437  Bvec.zero();
1438  Resid.zero();
1439  Dummy.zero();
1440  s.zero();
1441  ds.zero();
1442 
1443  _kspace.getSourceTerm(c,s,ds);
1444 
1445  if(_BCArray[c]==0 || _BCArray[c]==2) //Arrowhead
1446  {
1447  if(_BCArray[c]==2)
1448  updateGhostCoarse(c);
1449 
1450  AMat.zero();
1451 
1452  COMETConvectionCoarse(c,AMat,Bvec);
1453  COMETCollision(c,&AMat,Bvec);
1454  // COMETEquilibrium(c,&AMat,Bvec);
1455  COMETFullScatt(c,s,Bvec);
1456 
1457  traceSum+=AMat.getTraceAbs();
1458 
1459  Resid=Bvec;
1460  Bvec.zero();
1461  Distribute(c,Bvec,Resid);
1462  ArrayAbs(Resid);
1463  ResidSum+=Resid;
1464 
1465  AMat.multiply(Resid,Bvec);
1466  Resid=Bvec;
1467 
1468  makeValueArray(c,Bvec);
1469  AMat.multiply(Bvec,Dummy);
1470  Bvec=Dummy;
1471  ArrayAbs(Bvec);
1472  ArrayAbs(Resid);
1473  Bsum+=Bvec;
1474  }
1475  else if(_BCArray[c]==1 || _BCArray[c]==3) //General Dense
1476  {
1477  TSquare AMatS(totalmodes+1);
1478  AMatS.zero();
1479 
1480  COMETConvection(c,AMatS,Bvec);
1481  COMETCollision(c,&AMatS,Bvec);
1482  COMETEquilibrium(c,&AMatS,Bvec);
1483 
1484  traceSum+=AMatS.getTraceAbs();
1485  Resid=Bvec;
1486  Bvec.zero();
1487  Distribute(c,Bvec,Resid);
1488 
1489  AMatS.multiply(Resid,Bvec);
1490  Resid=Bvec;
1491 
1492  makeValueArray(c,Bvec);
1493  ArrayAbs(Resid);
1494  ArrayAbs(Bvec);
1495  Bsum+=Bvec;
1496  ResidSum+=Resid;
1497  }
1498  else
1499  throw CException("Unexpected value for boundary cell map.");
1500  }
1501 
1502  //traceSum=0; //added
1503  for(int o=0;o<totalmodes+1;o++)
1504  {
1505  ResidScalar+=ResidSum[o];
1506  //traceSum+=Bsum[o]; //added
1507  }
1508 
1509  ResidScalar/=traceSum;
1510 
1511  if(_aveResid==-1)
1512  {_aveResid=ResidScalar;}
1513  else
1514  {
1515  _residChange=fabs(_aveResid-ResidScalar)/_aveResid;
1516  _aveResid=ResidScalar;
1517  }
1518  }
void makeValueArray(const int c, TArray &o)
int getSelfCount() const
Definition: StorageSite.h:40
void COMETConvection(const int cell, TSquare &Amat, TArray &BVec)
ArrowHeadMatrix< T > TArrow
const IntArray & _BCArray
void Distribute(const int cell, TArray &BVec, TArray &Rvec)
void COMETFullScatt(const int cell, TArray &s, TArray &BVec)
void COMETCollision(const int cell, TMatrix *Amat, TArray &BVec)
void ArrayAbs(TArray &o)
Array< T_Scalar > TArray
SquareMatrix< T > TSquare
void COMETEquilibrium(const int cell, TMatrix *Amat, TArray &BVec)
Tangent fabs(const Tangent &a)
Definition: Tangent.h:312
void COMETConvectionCoarse(const int cell0, TArrow &Amat, TArray &BVec)
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
void getSourceTerm(const int c, TArray &s, TArray &ds)
Definition: Kspace.h:1468
void updateGhostCoarse(const int cell)
template<class T>
TArray COMETDiscretizer< T >::gatherResid ( const int  c)
inline

Definition at line 1520 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), pmode< T >::getresid(), Kspace< T >::gettotmodes(), and PhononMacro::TlResidual.

1521  {
1522  const int klen=_kspace.getlength();
1523  const int totalmodes=_kspace.gettotmodes();
1524  const int order=totalmodes+1;
1525 
1526  TArray rArray(order);
1527 
1528  for(int k=0;k<klen;k++)
1529  {
1530  Tkvol& kvol=_kspace.getkvol(k);
1531  const int numModes=kvol.getmodenum();
1532  for(int m=0;m<numModes;m++)
1533  {
1534  Tmode& mode=kvol.getmode(m);
1535  const int count=mode.getIndex();
1536  Field& resfield=mode.getresid();
1537  TArray& resArray=dynamic_cast<TArray&>(resfield[_cells]);
1538  rArray[count]=resArray[c];
1539  }
1540  }
1541  rArray[order]=_macro.TlResidual[c];
1542  }
PhononMacro & _macro
Definition: Field.h:14
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
int getmodenum()
Definition: kvol.h:43
Array< T_Scalar > TArray
Field TlResidual
Definition: PhononMacro.h:24
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
Definition: kvol.h:14
template<class T>
T COMETDiscretizer< T >::getAveResid ( )
inline

Definition at line 1545 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_aveResid.

Referenced by COMETModel< T >::updateResid().

1545 {return _aveResid;}
template<class T>
T COMETDiscretizer< T >::getResidChange ( )
inline

Definition at line 1544 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_residChange.

1544 {return _residChange;}
template<class T>
void COMETDiscretizer< T >::makeValueArray ( const int  c,
TArray o 
)
inline

Definition at line 1583 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_eArray, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, Kspace< T >::getGlobalIndex(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), Kspace< T >::gettotmodes(), and PhononMacro::temperature.

Referenced by COMETDiscretizer< T >::findResidCoarse(), COMETDiscretizer< T >::findResidFine(), and COMETDiscretizer< T >::findResidFull().

1584  {
1585  int klen=_kspace.getlength();
1586  const int totmodes=_kspace.gettotmodes();
1587  int cellIndex=_kspace.getGlobalIndex(c,0);
1588  for(int k=0;k<klen;k++)
1589  {
1590  Tkvol& kvol=_kspace.getkvol(k);
1591  const int numModes=kvol.getmodenum();
1592  for(int m=0;m<numModes;m++)
1593  {
1594  Tmode& mode=kvol.getmode(m);
1595  const int count=mode.getIndex()-1;
1596  o[count]=_eArray[cellIndex];
1597  cellIndex++;
1598  }
1599  }
1600  TArray& Tl=dynamic_cast<TArray&>(_macro.temperature[_cells]);
1601  o[totmodes]=Tl[c];
1602  }
PhononMacro & _macro
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
Field temperature
Definition: PhononMacro.h:21
int getmodenum()
Definition: kvol.h:43
Array< T_Scalar > TArray
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
int gettotmodes()
Definition: Kspace.h:393
const StorageSite & _cells
Definition: kvol.h:14
template<class T>
void COMETDiscretizer< T >::outerProduct ( const VectorT3 v1,
const VectorT3 v2,
T3Tensor out 
)
inline

Definition at line 2124 of file COMETDiscretizer.h.

Referenced by COMETDiscretizer< T >::updateeShifted().

2125  {
2126  for(int i=0;i<3;i++)
2127  for(int j=0;j<3;j++)
2128  out(i,j)=v1[i]*v2[j];
2129  }
template<class T>
T COMETDiscretizer< T >::scaledResid ( const TArray de,
const int  c 
)
inline

Definition at line 2131 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_eArray, COMETDiscretizer< T >::_kspace, fabs(), Kspace< T >::getGlobalIndex(), and Array< T >::getLength().

Referenced by COMETDiscretizer< T >::COMETSolveCoarse().

2132  {
2133  T scaled(0);
2134  int cellIndex(_kspace.getGlobalIndex(c,0));
2135  for(int i=0;i<de.getLength()-1;i++)
2136  {
2137  scaled+=fabs(de[i])/_eArray[cellIndex];
2138  cellIndex++;
2139  }
2140 
2141  return scaled/de.getLength();
2142  }
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
Tangent fabs(const Tangent &a)
Definition: Tangent.h:312
int getLength() const
Definition: Array.h:87
template<class T>
void COMETDiscretizer< T >::ScatterPhonons ( const int  cell0)
inline

Definition at line 1097 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cellFaces, COMETDiscretizer< T >::_cellVolume, COMETDiscretizer< T >::_e0Array, COMETDiscretizer< T >::_eArray, COMETDiscretizer< T >::_faceArea, COMETDiscretizer< T >::_faceCells, COMETDiscretizer< T >::_kspace, CRConnectivity::getCount(), Kspace< T >::getGlobalIndex(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), Kspace< T >::gettotmodes(), pmode< T >::getv(), Kspace< T >::ScatterPhonons(), and Array< T >::zero().

Referenced by COMETDiscretizer< T >::COMETSolveFull().

1098  {
1099 
1100  TArray C(_kspace.gettotmodes());
1101  TArray B(_kspace.gettotmodes());
1102  TArray V(_kspace.gettotmodes());
1103  TArray newE(_kspace.gettotmodes());
1104  C.zero();
1105  B.zero();
1106  V.zero();
1107  newE.zero();
1108 
1109  const int neibcount=_cellFaces.getCount(cell0);
1110 
1111  for(int j=0;j<neibcount;j++)
1112  {
1113  const int f=_cellFaces(cell0,j);
1114  int cell1=_faceCells(f,1);
1115  VectorT3 Af=_faceArea[f];
1116  if(cell1==cell0)
1117  {
1118  cell1=_faceCells(f,0);
1119  Af*=-1.;
1120  }
1121  const int klen=_kspace.getlength();
1122 
1123  T flux;
1124  for(int k=0;k<klen;k++)
1125  {
1126 
1127  Tkvol& kvol=_kspace.getkvol(k);
1128  const int numModes=kvol.getmodenum();
1129  for(int m=0;m<numModes;m++)
1130  {
1131  Tmode& mode=kvol.getmode(m);
1132  const int count=mode.getIndex();
1133  const VectorT3 vg=mode.getv();
1134  flux=vg[0]*Af[0]+vg[1]*Af[1]+vg[2]*Af[2];
1135 
1136  if(flux>T_Scalar(0)) //outgoing
1137  {
1138  V[count-1]+=flux/_cellVolume[cell0];
1139  B[count-1]+=flux*_e0Array[_kspace.getGlobalIndex(cell0,count-1)]/_cellVolume[cell0];
1140  C[count-1]+=flux*_eArray[_kspace.getGlobalIndex(cell0,count-1)]/_cellVolume[cell0];
1141  }
1142  else //incoming
1143  C[count-1]+=flux*_eArray[_kspace.getGlobalIndex(cell1,count-1)]/_cellVolume[cell0];
1144 
1145  }
1146  }
1147  }
1148 
1149  _kspace.ScatterPhonons(cell0, 50, C, B, V, newE, _cellVolume[cell0]);
1150 
1151  }
int getCount(const int i) const
const CRConnectivity & _faceCells
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
NumTypeTraits< T >::T_Scalar T_Scalar
int getmodenum()
Definition: kvol.h:43
Array< T_Scalar > TArray
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
const TArray & _cellVolume
const CRConnectivity & _cellFaces
void ScatterPhonons(const int c, const int totIts, TArray &C, TArray &B, const TArray &V, TArray &newE, const T cv)
Definition: Kspace.h:1476
int gettotmodes()
Definition: Kspace.h:393
const VectorT3Array & _faceArea
Definition: kvol.h:14
template<class T>
void COMETDiscretizer< T >::setfgFinder ( )
inline

Definition at line 1547 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_fgFinder, COMETDiscretizer< T >::_mesh, Mesh::getBoundaryFaceGroups(), StorageSite::getCount(), StorageSite::getOffset(), FaceGroup::id, and FaceGroup::site.

Referenced by COMETModel< T >::smooth(), and COMETModel< T >::updateResid().

1548  {
1549  foreach(const FaceGroupPtr fgPtr, _mesh.getBoundaryFaceGroups())
1550  {
1551  const FaceGroup& fg=*fgPtr;
1552  const int off=fg.site.getOffset();
1553  const int cnt=fg.site.getCount();
1554  const int id=fg.id;
1555  VecInt2 BegEnd;
1556 
1557  BegEnd[0]=off;
1558  BegEnd[1]=off+cnt;
1559  _fgFinder[id]=BegEnd;
1560  }
1561  }
const FaceGroupList & getBoundaryFaceGroups() const
Definition: Mesh.h:187
shared_ptr< FaceGroup > FaceGroupPtr
Definition: Mesh.h:46
Vector< int, 2 > VecInt2
Definition: Mesh.h:28
const int id
Definition: Mesh.h:41
int getOffset() const
Definition: StorageSite.h:87
int getCount() const
Definition: StorageSite.h:39
StorageSite site
Definition: Mesh.h:40
template<class T>
void COMETDiscretizer< T >::updatee0 ( )
inline

Definition at line 1612 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_e0Array, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, pmode< T >::calce0(), Kspace< T >::getGlobalIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), StorageSite::getSelfCount(), and PhononMacro::temperature.

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), COMETDiscretizer< T >::COMETSolveFine(), and COMETDiscretizer< T >::COMETSolveFull().

1613  {
1614  TArray& Tl=dynamic_cast<TArray&>(_macro.temperature[_cells]);
1615 
1616  const T cellCount=_cells.getSelfCount();
1617  int klen=_kspace.getlength();
1618  for(int c=0;c<cellCount;c++)
1619  {
1620  int cellIndex=_kspace.getGlobalIndex(c,0);
1621  for(int k=0;k<klen;k++)
1622  {
1623  Tkvol& kvol=_kspace.getkvol(k);
1624  const int numModes=kvol.getmodenum();
1625  for(int m=0;m<numModes;m++)
1626  {
1627  Tmode& mode=kvol.getmode(m);
1628  _e0Array[cellIndex]=mode.calce0(Tl[c]);
1629  cellIndex++;
1630  }
1631  }
1632  }
1633  }
int getSelfCount() const
Definition: StorageSite.h:40
PhononMacro & _macro
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
Field temperature
Definition: PhononMacro.h:21
int getmodenum()
Definition: kvol.h:43
Array< T_Scalar > TArray
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
const StorageSite & _cells
Definition: kvol.h:14
template<class T>
void COMETDiscretizer< T >::updatee0 ( const int  c)
inline

Definition at line 1635 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_e0Array, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, pmode< T >::calce0(), Kspace< T >::getGlobalIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), PhononMacro::temperature, and Kspace< T >::updateTau().

1636  {
1637  TArray& Tl=dynamic_cast<TArray&>(_macro.temperature[_cells]);
1638 
1639  int klen=_kspace.getlength();
1640  int cellIndex=_kspace.getGlobalIndex(c,0);
1641  for(int k=0;k<klen;k++)
1642  {
1643  Tkvol& kvol=_kspace.getkvol(k);
1644  const int numModes=kvol.getmodenum();
1645  for(int m=0;m<numModes;m++)
1646  {
1647  Tmode& mode=kvol.getmode(m);
1648  _e0Array[cellIndex]=mode.calce0(Tl[c]);
1649  _kspace.updateTau(c,Tl[c]);
1650  cellIndex++;
1651  }
1652  }
1653  }
PhononMacro & _macro
void updateTau(const int c, const T Tl)
Definition: Kspace.h:1491
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
Field temperature
Definition: PhononMacro.h:21
int getmodenum()
Definition: kvol.h:43
Array< T_Scalar > TArray
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
const StorageSite & _cells
Definition: kvol.h:14
template<class T>
void COMETDiscretizer< T >::updateeShifted ( )
inline

Definition at line 2060 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, kvol< T >::getdk3(), pmode< T >::geteShifted(), pmode< T >::getfield(), kvol< T >::getkvec(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), StorageSite::getSelfCount(), inverse(), PhononMacro::lam, COMETDiscretizer< T >::outerProduct(), PhononMacro::temperature, Vector< T, N >::zero(), and SquareTensor< T, N >::zero().

2061  {
2062  const int cellcount=_cells.getSelfCount();
2063  TArray& Tl=dynamic_cast<TArray&>(_macro.temperature[_cells]);
2064  VectorT3Array& lam=dynamic_cast<VectorT3Array&>(_macro.lam[_cells]);
2065 
2066  for(int c=0;c<cellcount;c++)
2067  {
2068  T3Tensor TensorSum;
2069  VectorT3 VectorSum;
2070  VectorT3 lambda;
2071  T TpreFactor=0.;
2072  T VpreFactor=0.;
2073 
2074  int klen=_kspace.getlength();
2075  for(int k=0;k<klen;k++)
2076  {
2077  Tkvol& kvol=_kspace.getkvol(k);
2078  const int numModes=kvol.getmodenum();
2079  for(int m=0;m<numModes;m++)
2080  {
2081  Tmode& mode=kvol.getmode(m);
2082  Field& eField=mode.getfield();
2083  TArray& eArray=dynamic_cast<TArray&>(eField[_cells]);
2084  TpreFactor+=mode.calcTensorPrefactor(Tl[c]);
2085  VpreFactor+=mode.calcVectorPrefactor(Tl[c],eArray[c]);
2086  }
2087  }
2088 
2089  TensorSum.zero();
2090  VectorSum.zero();
2091 
2092  for(int k=0;k<klen;k++)
2093  {
2094  Tkvol& kvol=_kspace.getkvol(k);
2095  T dk3=kvol.getdk3();
2096  VectorT3 Kvec=kvol.getkvec();
2097  T3Tensor TempTensor;
2098  outerProduct(Kvec,Kvec,TempTensor);
2099  TensorSum+=TempTensor*dk3*TpreFactor;
2100  VectorSum+=Kvec*dk3*VpreFactor;
2101  }
2102 
2103  lambda=inverse(TensorSum)*VectorSum;
2104  lam[c]=lambda;
2105 
2106  for(int k=0;k<klen;k++)
2107  {
2108  Tkvol& kvol=_kspace.getkvol(k);
2109  VectorT3 Kvec=kvol.getkvec();
2110  T shift=Kvec[0]*lambda[0]+Kvec[1]*lambda[1]+Kvec[2]*lambda[2];
2111  const int numModes=kvol.getmodenum();
2112  for(int m=0;m<numModes;m++)
2113  {
2114  Tmode& mode=kvol.getmode(m);
2115  Field& eShiftedField=mode.geteShifted();
2116  TArray& eShifted=dynamic_cast<TArray&>(eShiftedField[_cells]);
2117  eShifted[c]=mode.calcShifted(Tl[c],shift);
2118  }
2119  }
2120 
2121  }
2122  }
int getSelfCount() const
Definition: StorageSite.h:40
PhononMacro & _macro
Definition: Field.h:14
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
Field temperature
Definition: PhononMacro.h:21
int getmodenum()
Definition: kvol.h:43
SquareTensor< T, 3 > T3Tensor
Array< T_Scalar > TArray
void outerProduct(const VectorT3 &v1, const VectorT3 &v2, T3Tensor &out)
SquareMatrix< T, 2 > inverse(SquareMatrix< T, 2 > &a)
void zero()
Definition: Vector.h:156
const StorageSite & _cells
Definition: kvol.h:14
T getdk3()
Definition: kvol.h:42
template<class T>
void COMETDiscretizer< T >::updateGhostCoarse ( const int  cell)
inline

Definition at line 1932 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_BCfArray, COMETDiscretizer< T >::_bcMap, COMETDiscretizer< T >::_cellFaces, COMETDiscretizer< T >::_cells, COMETDiscretizer< T >::_eArray, COMETDiscretizer< T >::_faceArea, COMETDiscretizer< T >::_faceCells, COMETDiscretizer< T >::_kspace, COMETDiscretizer< T >::_macro, pmode< T >::calce0(), Kspace< T >::calcTemp(), COMETDiscretizer< T >::findFgId(), CRConnectivity::getCount(), kvol< T >::getdk3(), Kspace< T >::getDK3(), Kspace< T >::getGlobalIndex(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), pmode< T >::getReflpair(), pmode< T >::getv(), and PhononMacro::temperature.

Referenced by COMETDiscretizer< T >::COMETSolveCoarse(), COMETDiscretizer< T >::COMETSolveFull(), COMETDiscretizer< T >::findResidCoarse(), and COMETDiscretizer< T >::findResidFull().

1933  {
1934  const int neibcount=_cellFaces.getCount(cell);
1935  const int klen=_kspace.getlength();
1936  TArray& Tl=dynamic_cast<TArray&>(_macro.temperature[_cells]);
1937  T DK3=_kspace.getDK3();
1938 
1939  for(int j=0;j<neibcount;j++)
1940  {
1941 
1942  const int f=_cellFaces(cell,j);
1943  if(_BCfArray[f]==3)
1944  {
1945  T SumEVdotA(0);
1946  int Fgid=findFgId(f);
1947  const T refl=(*(_bcMap[Fgid]))["specifiedReflection"];
1948  int cell2=_faceCells(f,1);
1949  VectorT3 Af=_faceArea[f];
1950 
1951  if(cell2==cell)
1952  {
1953  Af=Af*-1.;
1954  cell2=_faceCells(f,0);
1955  }
1956 
1957  const int cellstart=_kspace.getGlobalIndex(cell,0);
1958  int cellIndex=_kspace.getGlobalIndex(cell,0);
1959  int cell2Index=_kspace.getGlobalIndex(cell2,0);
1960 
1961  //first loop to sum up the diffuse reflection
1962  for(int k=0;k<klen;k++)
1963  {
1964  Tkvol& kvol=_kspace.getkvol(k);
1965  const int numModes=kvol.getmodenum();
1966  T dk3=kvol.getdk3();
1967  for(int m=0;m<numModes;m++)
1968  {
1969  Tmode& mode=kvol.getmode(m);
1970  VectorT3 vg=mode.getv();
1971  const T VdotA=Af[0]*vg[0]+Af[1]*vg[1]+Af[2]*vg[2];
1972 
1973  if(VdotA>0)
1974  {
1975  _eArray[cell2Index]=_eArray[cellIndex]*refl;
1976  SumEVdotA+=_eArray[cellIndex]*VdotA*(dk3/DK3);
1977  }
1978  else
1979  {
1980  if(refl==1.)
1981  {
1982  Refl_pair& rpairs=mode.getReflpair(Fgid);
1983  const int kk=rpairs.second.second;
1984  Tkvol& kkvol=_kspace.getkvol(kk);
1985  Tmode& refMode=kkvol.getmode(m);
1986  const int rIndex=cellstart+refMode.getIndex()-1;
1987  _eArray[cell2Index]=_eArray[rIndex]*refl;
1988  }
1989  else
1990  _eArray[cell2Index]=0;
1991  }
1992  cellIndex++;
1993  cell2Index++;
1994  }
1995  }
1996 
1997  //second loop to add in the diffuse reflection
1998  if(refl==0.)
1999  {
2000 
2001  T wallTemp=Tl[cell2];
2002  T esum=SumEVdotA*DK3;
2003  _kspace.calcTemp(wallTemp,esum,Af);
2004  SumEVdotA=wallTemp;
2005  Tl[cell2]=wallTemp;
2006 
2007  for(int k=0;k<klen;k++)
2008  {
2009  Tkvol& kvol=_kspace.getkvol(k);
2010  const int numModes=kvol.getmodenum();
2011  for(int m=0;m<numModes;m++)
2012  {
2013  Tmode& mode=kvol.getmode(m);
2014  VectorT3 vg=mode.getv();
2015  const int Index=mode.getIndex()-1;
2016  const int cell2Index=_kspace.getGlobalIndex(cell2,Index);
2017 
2018  const T oneMinusRefl=1.-refl;
2019  //const T VdotA=Af[0]*vg[0]+Af[1]*vg[1]+Af[2]*vg[2];
2020 
2021  _eArray[cell2Index]+=mode.calce0(Tl[cell2])*oneMinusRefl;
2022 
2023  }
2024  }
2025  }
2026 
2027  }
2028  }
2029 
2030  }
int getCount(const int i) const
PhononMacro & _macro
const CRConnectivity & _faceCells
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
Field temperature
Definition: PhononMacro.h:21
int getmodenum()
Definition: kvol.h:43
T getDK3() const
Definition: Kspace.h:408
void calcTemp(T &guess, const T e_sum, const Tvec An)
Definition: Kspace.h:443
Tmode::Refl_pair Refl_pair
Array< T_Scalar > TArray
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
const IntArray & _BCfArray
const CRConnectivity & _cellFaces
const VectorT3Array & _faceArea
const StorageSite & _cells
Definition: kvol.h:14
int findFgId(const int faceIndex)
template<class T>
void COMETDiscretizer< T >::updateGhostFine ( const int  cell,
const GradMatrix gMat 
)
inline

Definition at line 1655 of file COMETDiscretizer.h.

References COMETDiscretizer< T >::_BCfArray, COMETDiscretizer< T >::_bcMap, COMETDiscretizer< T >::_cellCoords, COMETDiscretizer< T >::_cellFaces, COMETDiscretizer< T >::_eArray, COMETDiscretizer< T >::_faceArea, COMETDiscretizer< T >::_faceCells, COMETDiscretizer< T >::_faces, COMETDiscretizer< T >::_geomFields, COMETDiscretizer< T >::_kspace, computeLimitCoeff2(), GeomFields::coordinate, COMETDiscretizer< T >::findFgId(), GradientMatrix< T_Scalar >::getCoeff(), CRConnectivity::getCount(), Kspace< T >::getGlobalIndex(), pmode< T >::getIndex(), Kspace< T >::getkvol(), Kspace< T >::getlength(), kvol< T >::getmode(), kvol< T >::getmodenum(), pmode< T >::getReflpair(), Kspace< T >::gettotmodes(), pmode< T >::getv(), and Array< T >::zero().

Referenced by COMETDiscretizer< T >::COMETSolveFine().

1656  {
1657  const int neibcount=_cellFaces.getCount(cell);
1658  const int klen=_kspace.getlength();
1659  //TArray SumEVdotA(neibcount);
1660  //TArray& Tl=dynamic_cast<TArray&>(_macro.temperature[_cells]);
1661  //T DK3=_kspace.getDK3();
1662 
1663  //SumEVdotA.zero();
1664 
1665  const VectorT3Array& faceCoords=
1666  dynamic_cast<const VectorT3Array&>(_geomFields.coordinate[_faces]);
1667 
1668  TArray pointMin(_kspace.gettotmodes());
1669  pointMin=-1;
1670  TArray pointMax(_kspace.gettotmodes());
1671  pointMax.zero();
1672  TArray pointLim(_kspace.gettotmodes());
1673  pointLim=100.;
1674 
1675  GradArray Grads(_kspace.gettotmodes());
1676  Grads.zero();
1677 
1678  VectorT3 Gcoeff;
1679 
1680  for(int j=0;j<neibcount;j++) //first loop to get grad and max/min vals
1681  {
1682  const int f=_cellFaces(cell,j);
1683  int cell1=_faceCells(f,1);
1684  if(cell1==cell)
1685  cell1=_faceCells(f,0);
1686  const VectorT3 Gcoeff=gMat.getCoeff(cell, cell1);
1687 
1688  int c0ind=_kspace.getGlobalIndex(cell,0);
1689  int c1ind=_kspace.getGlobalIndex(cell1,0);
1690 
1691  for(int k=0;k<_kspace.gettotmodes();k++)
1692  {
1693  const T e1=_eArray[c1ind];
1694  const T e0=_eArray[c0ind];
1695  T& curMax=pointMax[k];
1696  T& curMin=pointMin[k];
1697  Grads[k].accumulate(Gcoeff, e1-e0);
1698 
1699  if(e1>curMax)
1700  curMax=e1;
1701 
1702  if(curMin==-1)
1703  curMin=e1;
1704  else
1705  {
1706  if(e1<curMin)
1707  curMin=e1;
1708  }
1709 
1710  c0ind++;
1711  c1ind++;
1712  }
1713  }
1714 
1715  for(int j=0;j<neibcount;j++) //second loop to calculate limiting coeff
1716  {
1717  const int f=_cellFaces(cell,j);
1718  int cell1=_faceCells(f,1);
1719  if(cell1==cell)
1720  cell1=_faceCells(f,0);
1721 
1722  const VectorT3 dr0(faceCoords[f]-_cellCoords[cell]);
1723  int c0ind=_kspace.getGlobalIndex(cell,0);
1724  vanLeer lf;
1725 
1726  for(int k=0;k<_kspace.gettotmodes();k++)
1727  {
1728  const T minVal=pointMin[k];
1729  const T maxVal=pointMax[k];
1730  const T de0(Grads[k]*dr0);
1731  T& cl=pointLim[k];
1732  computeLimitCoeff2(cl, _eArray[c0ind], de0, minVal, maxVal, lf);
1733  c0ind++;
1734  }
1735  }
1736 
1737 
1738  for(int j=0;j<neibcount;j++)
1739  {
1740 
1741  const int f=_cellFaces(cell,j);
1742 
1743  if(_BCfArray[f]==1)
1744  {//temperature
1745  //int Fgid=findFgId(f);
1746  int cell2=_faceCells(f,1);
1747  VectorT3 Af=_faceArea[f];
1748 
1749  if(cell2==cell)
1750  {
1751  Af=Af*-1.;
1752  cell2=_faceCells(f,0);
1753  }
1754 
1755  int cellIndex=_kspace.getGlobalIndex(cell,0);
1756  int cell2Index=_kspace.getGlobalIndex(cell2,0);
1757  //vanLeer vl;
1758 
1759  for(int k=0;k<klen;k++)
1760  {
1761  Tkvol& kvol=_kspace.getkvol(k);
1762  const int numModes=kvol.getmodenum();
1763  //T dk3=kvol.getdk3();
1764  for(int m=0;m<numModes;m++)
1765  {
1766  Tmode& mode=kvol.getmode(m);
1767  VectorT3 vg=mode.getv();
1768  const T VdotA=Af[0]*vg[0]+Af[1]*vg[1]+Af[2]*vg[2];
1769  const int count=mode.getIndex();
1770 
1771  GradType& grad=Grads[count-1];
1772 
1773  if(VdotA>0)
1774  {
1775  const VectorT3 fVec=faceCoords[f]-_cellCoords[cell];
1776  const VectorT3 rVec=_cellCoords[cell2]-_cellCoords[cell];
1777  //const T r=gMat.computeR(grad,_eArray,rVec,cellIndex,cell2Index);
1778 
1779  T SOU=(fVec[0]*grad[0]+fVec[1]*grad[1]+
1780  fVec[2]*grad[2])*pointLim[count-1];
1781  _eArray[cell2Index]=_eArray[cellIndex]+SOU;
1782  }
1783  cellIndex++;
1784  cell2Index++;
1785  }
1786  }
1787 
1788  }
1789  else if(_BCfArray[f]==3)
1790  {//reflecting
1791  int Fgid=findFgId(f);
1792  const T refl=(*(_bcMap[Fgid]))["specifiedReflection"];
1793  int cell2=_faceCells(f,1);
1794  VectorT3 Af=_faceArea[f];
1795 
1796  if(cell2==cell)
1797  {
1798  Af=Af*-1.;
1799  cell2=_faceCells(f,0);
1800  }
1801 
1802  //const int cellstart=_kspace.getGlobalIndex(cell,0);
1803  const int cell2start=_kspace.getGlobalIndex(cell2,0);
1804  int cellIndex=_kspace.getGlobalIndex(cell,0);
1805  int cell2Index=_kspace.getGlobalIndex(cell2,0);
1806  vanLeer vl;
1807 
1808  //first loop to sum up the diffuse reflection
1809  for(int k=0;k<klen;k++)
1810  {
1811  Tkvol& kvol=_kspace.getkvol(k);
1812  const int numModes=kvol.getmodenum();
1813  //T dk3=kvol.getdk3();
1814  for(int m=0;m<numModes;m++)
1815  {
1816  Tmode& mode=kvol.getmode(m);
1817  VectorT3 vg=mode.getv();
1818  const T VdotA=Af[0]*vg[0]+Af[1]*vg[1]+Af[2]*vg[2];
1819  const int count=mode.getIndex();
1820 
1821  GradType& grad=Grads[count-1];
1822 
1823  if(VdotA>0)
1824  {
1825  VectorT3 rVec=_cellCoords[cell2]-_cellCoords[cell];
1826  VectorT3 fVec=faceCoords[f]-_cellCoords[cell];
1827  //T r=gMat.computeR(grad,_eArray,rVec,cellIndex,cell2Index);
1828  T SOU=(grad[0]*fVec[0]+grad[1]*fVec[1]+grad[2]*fVec[2])*pointLim[count-1];
1829  _eArray[cell2Index]=(_eArray[cellIndex]+SOU)*refl;
1830  }
1831  cellIndex++;
1832  cell2Index++;
1833  }
1834  }
1835 
1836  cellIndex=_kspace.getGlobalIndex(cell,0);
1837  cell2Index=_kspace.getGlobalIndex(cell2,0);
1838 
1839  //first loop to sum up the diffuse reflection
1840  for(int k=0;k<klen;k++)
1841  {
1842  Tkvol& kvol=_kspace.getkvol(k);
1843  const int numModes=kvol.getmodenum();
1844  //T dk3=kvol.getdk3();
1845  for(int m=0;m<numModes;m++)
1846  {
1847  Tmode& mode=kvol.getmode(m);
1848  VectorT3 vg=mode.getv();
1849  const T VdotA=Af[0]*vg[0]+Af[1]*vg[1]+Af[2]*vg[2];
1850  //const int count=mode.getIndex();
1851 
1852  if(VdotA<0)
1853  {
1854  Refl_pair& rpairs=mode.getReflpair(Fgid);
1855  const int kk=rpairs.second.second;
1856  Tkvol& kkvol=_kspace.getkvol(kk);
1857  Tmode& refMode=kkvol.getmode(m);
1858  const int rIndex=cell2start+refMode.getIndex()-1;
1859  _eArray[cell2Index]=_eArray[rIndex]*refl;
1860  }
1861  cellIndex++;
1862  cell2Index++;
1863  }
1864  }
1865 
1866  }
1867  }
1868 
1869  /*
1870  T wallTemp(Tl[cell]);
1871  for(int j=0;j<neibcount;j++)
1872  {
1873  const int f=_cellFaces(cell,j);
1874  if(_BCfArray[f]==3)
1875  {
1876  int cell2=_faceCells(f,1);
1877  VectorT3 Af=_faceArea[f];
1878  if(cell2==cell)
1879  {
1880  Af=Af*-1.;
1881  cell2=_faceCells(f,0);
1882  }
1883  T esum=SumEVdotA[j]*DK3;
1884  _kspace.calcTemp(wallTemp,esum,Af);
1885  SumEVdotA[j]=wallTemp;
1886  Tl[cell2]=wallTemp;
1887  }
1888  }
1889 
1890  //second loop to add in the diffuse reflection
1891  for(int k=0;k<klen;k++)
1892  {
1893  Tkvol& kvol=_kspace.getkvol(k);
1894  const int numModes=kvol.getmodenum();
1895  for(int m=0;m<numModes;m++)
1896  {
1897  Tmode& mode=kvol.getmode(m);
1898  VectorT3 vg=mode.getv();
1899  Field& eField=mode.getfield();
1900  TArray& eArray=dynamic_cast<TArray&>(eField[_cells]);
1901  for(int j=0;j<neibcount;j++)
1902  {
1903  const int f=_cellFaces(cell,j);
1904  if(_BCfArray[f]==3)
1905  {
1906  int Fgid=findFgId(f);
1907  const T refl=(*(_bcMap[Fgid]))["specifiedReflection"];
1908  const T oneMinusRefl=1.-refl;
1909  int cell2=_faceCells(f,1);
1910  VectorT3 Af=_faceArea[f];
1911 
1912  if(cell2==cell)
1913  {
1914  Af=Af*-1.;
1915  cell2=_faceCells(f,0);
1916  }
1917 
1918  const T VdotA=Af[0]*vg[0]+Af[1]*vg[1]+Af[2]*vg[2];
1919 
1920  if(VdotA<0)
1921  eArray[cell2]+=mode.calce0(SumEVdotA[j])*oneMinusRefl;
1922  else
1923  eArray[cell2]+=mode.calce0(SumEVdotA[j])*oneMinusRefl;
1924 
1925  }
1926  }
1927  }
1928  }
1929  */
1930  }
int getCount(const int i) const
const VectorT3Array & _cellCoords
Field coordinate
Definition: GeomFields.h:19
const CRConnectivity & _faceCells
Tkvol & getkvol(int n) const
Definition: Kspace.h:390
int getlength() const
Definition: Kspace.h:391
const StorageSite & _faces
int getmodenum()
Definition: kvol.h:43
Array< GradType > GradArray
Tmode::Refl_pair Refl_pair
Array< T_Scalar > TArray
int getGlobalIndex(const int cell, const int count)
Definition: Kspace.h:1216
const IntArray & _BCfArray
const CRConnectivity & _cellFaces
void computeLimitCoeff2(T &lc, const T x, const T &dx, const T &min, const T &max, const LimitFunc &f)
Definition: FluxLimiters.h:85
const GeomFields & _geomFields
int gettotmodes()
Definition: Kspace.h:393
const VectorT3Array & _faceArea
Gradient< T > GradType
Definition: kvol.h:14
int findFgId(const int faceIndex)

Member Data Documentation

template<class T>
const VectorT3Array& COMETDiscretizer< T >::_cellCoords
private
template<class T>
const TArray& COMETDiscretizer< T >::_faceAreaMag
private

Definition at line 2152 of file COMETDiscretizer.h.

template<class T>
const StorageSite& COMETDiscretizer< T >::_faces
private
template<class T>
FaceToFg COMETDiscretizer< T >::_fgFinder
private
template<class T>
TArray& COMETDiscretizer< T >::_resArray
private

Definition at line 2168 of file COMETDiscretizer.h.

Referenced by COMETDiscretizer< T >::Distribute().


The documentation for this class was generated from the following file: