/[escript]/trunk/finley/src/CPPAdapter/MeshAdapter.cpp
ViewVC logotype

Diff of /trunk/finley/src/CPPAdapter/MeshAdapter.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1388 by trankine, Fri Jan 11 07:45:58 2008 UTC revision 1464 by gross, Tue Apr 1 23:27:09 2008 UTC
# Line 24  extern "C" { Line 24  extern "C" {
24  }  }
25  #include <vector>  #include <vector>
26    
27    #define IS_THERE_A_REASON_FOR_THIS_MAGIC_NAME_LENGTH  256
28    
29  using namespace std;  using namespace std;
30  using namespace escript;  using namespace escript;
31    
# Line 814  void MeshAdapter::addPDEToRHS( escript:: Line 816  void MeshAdapter::addPDEToRHS( escript::
816  //  //
817  void MeshAdapter::addPDEToTransportProblem(  void MeshAdapter::addPDEToTransportProblem(
818                       TransportProblemAdapter& tp, escript::Data& source, const escript::Data& M,                       TransportProblemAdapter& tp, escript::Data& source, const escript::Data& M,
819                       const escript::Data& A, const escript::Data& B, const escript::Data& C,const  escript::Data& D,const  escript::Data& X,const  escript::Data& Y,                       const escript::Data& A, const escript::Data& B, const escript::Data& C,
820                         const  escript::Data& D,const  escript::Data& X,const  escript::Data& Y,
821                       const escript::Data& d, const escript::Data& y,                       const escript::Data& d, const escript::Data& y,
822                       const escript::Data& d_contact,const escript::Data& y_contact) const                       const escript::Data& d_contact,const escript::Data& y_contact) const
823  {  {
824     DataArrayView::ShapeType shape;     DataArrayView::ShapeType shape;
825     source.expand();     source.expand();
    escript:: Data tmp(0.0,M.getDataPointShape(),tp.getFunctionSpace(),true);  
826     escriptDataC _source=source.getDataC();     escriptDataC _source=source.getDataC();
    escriptDataC _tmp=tmp.getDataC();  
827     escriptDataC _M=M.getDataC();     escriptDataC _M=M.getDataC();
828     escriptDataC _A=A.getDataC();     escriptDataC _A=A.getDataC();
829     escriptDataC _B=B.getDataC();     escriptDataC _B=B.getDataC();
# Line 838  void MeshAdapter::addPDEToTransportProbl Line 839  void MeshAdapter::addPDEToTransportProbl
839     Finley_Mesh* mesh=m_finleyMesh.get();     Finley_Mesh* mesh=m_finleyMesh.get();
840     Paso_FCTransportProblem* _tp = tp.getPaso_FCTransportProblem();     Paso_FCTransportProblem* _tp = tp.getPaso_FCTransportProblem();
841        
842       Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,_tp->mass_matrix, &_source, 0, 0, 0, &_M, 0, 0 );
    Finley_Assemble_LumpedSystem(mesh->Nodes,mesh->Elements,&_tmp, &_M);  
    checkFinleyError();  
    /* add mass matix to lumped mass matrix of transport problem */  
    double* tmp_prt=getSampleData(&_tmp,0);  
    int i;  
    int n=Paso_FCTransportProblem_getTotalNumRows(_tp);  
    #pragma omp parallel for private(i) schedule(static)  
    for (i=0;i<n ;++i) _tp->lumped_mass_matrix[i]+=tmp_prt[i];  
   
    Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,_tp->transport_matrix, &_source, &_A, 0, 0, &_D, &_X, &_Y );  
843     checkFinleyError();     checkFinleyError();
844    
845     Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,_tp->flux_matrix, &_source, 0, &_B, &_C, 0, 0, 0 );     Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,_tp->transport_matrix, &_source, &_A, &_B, &_C, &_D, &_X, &_Y );
846     checkFinleyError();     checkFinleyError();
847    
848     Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, _tp->transport_matrix, &_source, 0, 0, 0, &_d, 0, &_y );     Finley_Assemble_PDE(mesh->Nodes,mesh->FaceElements, _tp->transport_matrix, &_source, 0, 0, 0, &_d, 0, &_y );
# Line 1092  void MeshAdapter::interpolateOnDomain(es Line 1083  void MeshAdapter::interpolateOnDomain(es
1083            case(Elements):            case(Elements):
1084            case(ReducedElements):            case(ReducedElements):
1085                if (getMPISize()>1) {                if (getMPISize()>1) {
1086                   escript::Data temp=escript::Data( in,  continuousFunction(asAbstractContinuousDomain()) );                   escript::Data temp=escript::Data( in,  reducedContinuousFunction(asAbstractContinuousDomain()) );
1087                   escriptDataC _in2 = temp.getDataC();                   escriptDataC _in2 = temp.getDataC();
1088                   Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in2,&_target);                   Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in2,&_target);
1089                } else {                } else {
# Line 1102  void MeshAdapter::interpolateOnDomain(es Line 1093  void MeshAdapter::interpolateOnDomain(es
1093            case(FaceElements):            case(FaceElements):
1094            case(ReducedFaceElements):            case(ReducedFaceElements):
1095                if (getMPISize()>1) {                if (getMPISize()>1) {
1096                   escript::Data temp=escript::Data( in,  continuousFunction(asAbstractContinuousDomain()) );                   escript::Data temp=escript::Data( in,  reducedContinuousFunction(asAbstractContinuousDomain()) );
1097                   escriptDataC _in2 = temp.getDataC();                   escriptDataC _in2 = temp.getDataC();
1098                   Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in2,&_target);                   Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in2,&_target);
1099                } else {                } else {
# Line 1111  void MeshAdapter::interpolateOnDomain(es Line 1102  void MeshAdapter::interpolateOnDomain(es
1102               break;               break;
1103            case(Points):            case(Points):
1104                if (getMPISize()>1) {                if (getMPISize()>1) {
1105                   escript::Data temp=escript::Data( in,  continuousFunction(asAbstractContinuousDomain()) );                   escript::Data temp=escript::Data( in,  reducedContinuousFunction(asAbstractContinuousDomain()) );
1106                   escriptDataC _in2 = temp.getDataC();                   escriptDataC _in2 = temp.getDataC();
1107                   Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in2,&_target);                   Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in2,&_target);
1108                } else {                } else {
# Line 1123  void MeshAdapter::interpolateOnDomain(es Line 1114  void MeshAdapter::interpolateOnDomain(es
1114            case(ReducedContactElementsZero):            case(ReducedContactElementsZero):
1115            case(ReducedContactElementsOne):            case(ReducedContactElementsOne):
1116                if (getMPISize()>1) {                if (getMPISize()>1) {
1117                   escript::Data temp=escript::Data( in,  continuousFunction(asAbstractContinuousDomain()) );                   escript::Data temp=escript::Data( in,  reducedContinuousFunction(asAbstractContinuousDomain()) );
1118                   escriptDataC _in2 = temp.getDataC();                   escriptDataC _in2 = temp.getDataC();
1119                   Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in2,&_target);                   Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in2,&_target);
1120                } else {                } else {
# Line 1454  void MeshAdapter::setNewX(const escript: Line 1445  void MeshAdapter::setNewX(const escript:
1445  // saves a data array in openDX format:  // saves a data array in openDX format:
1446  void MeshAdapter::saveDX(const std::string& filename,const boost::python::dict& arg) const  void MeshAdapter::saveDX(const std::string& filename,const boost::python::dict& arg) const
1447  {  {
1448      int MAX_namelength=256;      unsigned int MAX_namelength=IS_THERE_A_REASON_FOR_THIS_MAGIC_NAME_LENGTH;
1449      const int num_data=boost::python::extract<int>(arg.attr("__len__")());      const int num_data=boost::python::extract<int>(arg.attr("__len__")());
1450    /* win32 refactor */    /* win32 refactor */
1451    char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;    char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
# Line 1502  void MeshAdapter::saveDX(const std::stri Line 1493  void MeshAdapter::saveDX(const std::stri
1493  // saves a data array in openVTK format:  // saves a data array in openVTK format:
1494  void MeshAdapter::saveVTK(const std::string& filename,const boost::python::dict& arg) const  void MeshAdapter::saveVTK(const std::string& filename,const boost::python::dict& arg) const
1495  {  {
1496      int MAX_namelength=256;      unsigned int MAX_namelength=IS_THERE_A_REASON_FOR_THIS_MAGIC_NAME_LENGTH;
1497      const int num_data=boost::python::extract<int>(arg.attr("__len__")());      const int num_data=boost::python::extract<int>(arg.attr("__len__")());
1498    /* win32 refactor */    /* win32 refactor */
1499    char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;    char* *names = (num_data>0) ? TMPMEMALLOC(num_data,char*) : (char**)NULL;
# Line 1601  SystemMatrixAdapter MeshAdapter::newSyst Line 1592  SystemMatrixAdapter MeshAdapter::newSyst
1592  // creates a TransportProblemAdapter  // creates a TransportProblemAdapter
1593  TransportProblemAdapter MeshAdapter::newTransportProblem(  TransportProblemAdapter MeshAdapter::newTransportProblem(
1594                        const double theta,                        const double theta,
                       const double dt_max,  
1595                        const int blocksize,                        const int blocksize,
1596                        const escript::FunctionSpace& functionspace,                        const escript::FunctionSpace& functionspace,
1597                        const int type) const                        const int type) const
# Line 1624  TransportProblemAdapter MeshAdapter::new Line 1614  TransportProblemAdapter MeshAdapter::new
1614      Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceOrder,reduceOrder);      Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceOrder,reduceOrder);
1615      checkFinleyError();      checkFinleyError();
1616      Paso_FCTransportProblem* transportProblem;      Paso_FCTransportProblem* transportProblem;
1617      transportProblem=Paso_FCTransportProblem_alloc(theta,dt_max,fsystemMatrixPattern,blocksize);      transportProblem=Paso_FCTransportProblem_alloc(theta,fsystemMatrixPattern,blocksize);
1618      checkPasoError();      checkPasoError();
1619      Paso_SystemMatrixPattern_free(fsystemMatrixPattern);      Paso_SystemMatrixPattern_free(fsystemMatrixPattern);
1620      return TransportProblemAdapter(transportProblem,theta,dt_max,blocksize,functionspace);      return TransportProblemAdapter(transportProblem,theta,blocksize,functionspace);
1621  }  }
1622    
1623  //  //
# Line 1863  escript::Data MeshAdapter::getSize() con Line 1853  escript::Data MeshAdapter::getSize() con
1853    
1854  int* MeshAdapter::borrowSampleReferenceIDs(int functionSpaceType) const  int* MeshAdapter::borrowSampleReferenceIDs(int functionSpaceType) const
1855  {  {
1856    int *out=0,i;    int *out = NULL;
1857    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
1858    switch (functionSpaceType) {    switch (functionSpaceType) {
1859      case(Nodes):      case(Nodes):

Legend:
Removed from v.1388  
changed lines
  Added in v.1464

  ViewVC Help
Powered by ViewVC 1.1.26