~vcs-imports/escript-finley/trunk

« back to all changes in this revision

Viewing changes to finley/src/CPPAdapter/MeshAdapter.cpp

  • Committer: jfenwick
  • Date: 2010-10-11 01:48:14 UTC
  • Revision ID: svn-v4:77569008-7704-0410-b7a0-a92fef0b09fd:trunk:3259
Merging dudley and scons updates from branches

Show diffs side-by-side

added added

removed removed

Lines of Context:
18
18
#ifdef USE_NETCDF
19
19
#include <netcdfcpp.h>
20
20
#endif
21
 
#ifdef PASO_MPI
 
21
#ifdef ESYS_MPI
22
22
#include <mpi.h>
23
 
#include "paso/Paso_MPI.h"
 
23
#include "esysUtils/Esys_MPI.h"
24
24
#endif
25
25
extern "C" {
26
26
#include "esysUtils/blocktimer.h"
85
85
}
86
86
void MeshAdapter::MPIBarrier() const
87
87
{
88
 
#ifdef PASO_MPI
 
88
#ifdef ESYS_MPI
89
89
   MPI_Barrier(m_finleyMesh.get()->MPIInfo->comm);
90
90
#endif
91
91
   return;
96
96
}
97
97
 
98
98
 
99
 
#ifdef PASO_MPI
 
99
#ifdef ESYS_MPI
100
100
  MPI_Comm
101
101
#else
102
102
  unsigned int
103
103
#endif
104
104
MeshAdapter::getMPIComm() const
105
105
{
106
 
#ifdef PASO_MPI
 
106
#ifdef ESYS_MPI
107
107
        return m_finleyMesh->MPIInfo->comm;
108
108
#else
109
109
        return 0;
149
149
   int num_Elements_numNodes            = mesh->Elements->numNodes;
150
150
   int num_FaceElements_numNodes        = mesh->FaceElements->numNodes;
151
151
   int num_ContactElements_numNodes     = mesh->ContactElements->numNodes;
152
 
#ifdef PASO_MPI
 
152
#ifdef ESYS_MPI
153
153
   MPI_Status status;
154
154
#endif
155
155
 
156
156
/* Incoming token indicates it's my turn to write */
157
 
#ifdef PASO_MPI
 
157
#ifdef ESYS_MPI
158
158
   if (mpi_rank>0) MPI_Recv(&num_Tags, 0, MPI_INT, mpi_rank-1, 81800, mesh->MPIInfo->comm, &status);
159
159
#endif
160
160
 
161
 
   char *newFileName = Paso_MPI_appendRankToFileName(fileName.c_str(),
 
161
   char *newFileName = Esys_MPI_appendRankToFileName(fileName.c_str(),
162
162
                                                     mpi_size, mpi_rank);
163
163
 
164
164
   /* Figure out how much storage is required for tags */
529
529
   }
530
530
 
531
531
/* Send token to next MPI process so he can take his turn */
532
 
#ifdef PASO_MPI
 
532
#ifdef ESYS_MPI
533
533
   if (mpi_rank<mpi_size-1) MPI_Send(&num_Tags, 0, MPI_INT, mpi_rank+1, 81800, mesh->MPIInfo->comm);
534
534
#endif
535
535
 
1044
1044
      case(Elements):
1045
1045
      case(ReducedElements):
1046
1046
      if (getMPISize()>1) {
1047
 
         escript::Data temp=escript::Data( in,  continuousFunction(asAbstractContinuousDomain()) );
 
1047
         escript::Data temp=escript::Data( in,  continuousFunction(*this) );
1048
1048
         escriptDataC _in2 = temp.getDataC();
1049
1049
         Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in2,&_target);
1050
1050
      } else {
1054
1054
      case(FaceElements):
1055
1055
      case(ReducedFaceElements):
1056
1056
      if (getMPISize()>1) {
1057
 
         escript::Data temp=escript::Data( in,  continuousFunction(asAbstractContinuousDomain()) );
 
1057
         escript::Data temp=escript::Data( in,  continuousFunction(*this) );
1058
1058
         escriptDataC _in2 = temp.getDataC();
1059
1059
         Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in2,&_target);
1060
1060
   
1064
1064
      break;
1065
1065
      case(Points):
1066
1066
      if (getMPISize()>1) {
1067
 
         escript::Data temp=escript::Data( in,  continuousFunction(asAbstractContinuousDomain()) );
 
1067
         escript::Data temp=escript::Data( in,  continuousFunction(*this) );
1068
1068
         escriptDataC _in2 = temp.getDataC();
1069
1069
      } else {
1070
1070
         Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in,&_target);
1075
1075
      case(ReducedContactElementsZero):
1076
1076
      case(ReducedContactElementsOne):
1077
1077
      if (getMPISize()>1) {
1078
 
         escript::Data temp=escript::Data( in,  continuousFunction(asAbstractContinuousDomain()) );
 
1078
         escript::Data temp=escript::Data( in,  continuousFunction(*this) );
1079
1079
         escriptDataC _in2 = temp.getDataC();
1080
1080
         Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in2,&_target);
1081
1081
      } else {
1113
1113
      case(Elements):
1114
1114
      case(ReducedElements):
1115
1115
      if (getMPISize()>1) {
1116
 
         escript::Data temp=escript::Data( in,  reducedContinuousFunction(asAbstractContinuousDomain()) );
 
1116
         escript::Data temp=escript::Data( in,  reducedContinuousFunction(*this) );
1117
1117
         escriptDataC _in2 = temp.getDataC();
1118
1118
         Finley_Assemble_interpolate(mesh->Nodes,mesh->Elements,&_in2,&_target);
1119
1119
      } else {
1123
1123
      case(FaceElements):
1124
1124
      case(ReducedFaceElements):
1125
1125
      if (getMPISize()>1) {
1126
 
         escript::Data temp=escript::Data( in,  reducedContinuousFunction(asAbstractContinuousDomain()) );
 
1126
         escript::Data temp=escript::Data( in,  reducedContinuousFunction(*this) );
1127
1127
         escriptDataC _in2 = temp.getDataC();
1128
1128
         Finley_Assemble_interpolate(mesh->Nodes,mesh->FaceElements,&_in2,&_target);
1129
1129
      } else {
1132
1132
      break;
1133
1133
      case(Points):
1134
1134
      if (getMPISize()>1) {
1135
 
         escript::Data temp=escript::Data( in,  reducedContinuousFunction(asAbstractContinuousDomain()) );
 
1135
         escript::Data temp=escript::Data( in,  reducedContinuousFunction(*this) );
1136
1136
         escriptDataC _in2 = temp.getDataC();
1137
1137
         Finley_Assemble_interpolate(mesh->Nodes,mesh->Points,&_in2,&_target);
1138
1138
      } else {
1144
1144
      case(ReducedContactElementsZero):
1145
1145
      case(ReducedContactElementsOne):
1146
1146
      if (getMPISize()>1) {
1147
 
         escript::Data temp=escript::Data( in,  reducedContinuousFunction(asAbstractContinuousDomain()) );
 
1147
         escript::Data temp=escript::Data( in,  reducedContinuousFunction(*this) );
1148
1148
         escriptDataC _in2 = temp.getDataC();
1149
1149
         Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&_in2,&_target);
1150
1150
      } else {
1181
1181
      escriptDataC _arg=arg.getDataC();
1182
1182
      Finley_Assemble_NodeCoordinates(mesh->Nodes,&_arg);
1183
1183
   } else {
1184
 
      escript::Data tmp_data=Vector(0.0,continuousFunction(asAbstractContinuousDomain()),true);
 
1184
      escript::Data tmp_data=Vector(0.0,continuousFunction(*this),true);
1185
1185
      escriptDataC _tmp_data=tmp_data.getDataC();
1186
1186
      Finley_Assemble_NodeCoordinates(mesh->Nodes,&_tmp_data);
1187
1187
      // this is then interpolated onto arg:
1275
1275
   escriptDataC _arg=arg.getDataC();
1276
1276
   switch(arg.getFunctionSpace().getTypeCode()) {
1277
1277
   case(Nodes):
1278
 
   temp=escript::Data( arg, escript::function(asAbstractContinuousDomain()) );
 
1278
   temp=escript::Data( arg, escript::function(*this) );
1279
1279
   _temp=temp.getDataC();
1280
1280
   Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1281
1281
   break;
1282
1282
   case(ReducedNodes):
1283
 
   temp=escript::Data( arg, escript::function(asAbstractContinuousDomain()) );
 
1283
   temp=escript::Data( arg, escript::function(*this) );
1284
1284
   _temp=temp.getDataC();
1285
1285
   Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1286
1286
   break;
1312
1312
   Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&_arg,&integrals[0]);
1313
1313
   break;
1314
1314
   case(DegreesOfFreedom):
1315
 
   temp=escript::Data( arg, escript::function(asAbstractContinuousDomain()) );
 
1315
   temp=escript::Data( arg, escript::function(*this) );
1316
1316
   _temp=temp.getDataC();
1317
1317
   Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1318
1318
   break;
1319
1319
   case(ReducedDegreesOfFreedom):
1320
 
   temp=escript::Data( arg, escript::function(asAbstractContinuousDomain()) );
 
1320
   temp=escript::Data( arg, escript::function(*this) );
1321
1321
   _temp=temp.getDataC();
1322
1322
   Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&_temp,&integrals[0]);
1323
1323
   break;
1349
1349
   escript::Data temp;
1350
1350
   if (getMPISize()>1) {
1351
1351
      if( arg.getFunctionSpace().getTypeCode() == DegreesOfFreedom ) {
1352
 
         temp=escript::Data( arg,  continuousFunction(asAbstractContinuousDomain()) );
 
1352
         temp=escript::Data( arg,  continuousFunction(*this) );
1353
1353
         nodeDataC = temp.getDataC();
1354
1354
      } else if( arg.getFunctionSpace().getTypeCode() == ReducedDegreesOfFreedom ) {
1355
 
         temp=escript::Data( arg,  reducedContinuousFunction(asAbstractContinuousDomain()) );
 
1355
         temp=escript::Data( arg,  reducedContinuousFunction(*this) );
1356
1356
         nodeDataC = temp.getDataC();
1357
1357
      } else {
1358
1358
         nodeDataC = arg.getDataC();
1471
1471
   const MeshAdapter& newDomain=dynamic_cast<const MeshAdapter&>(*(new_x.getFunctionSpace().getDomain()));
1472
1472
   if (newDomain!=*this) 
1473
1473
      throw FinleyAdapterException("Error - Illegal domain of new point locations");
1474
 
   if ( new_x.getFunctionSpace() == continuousFunction(asAbstractContinuousDomain()) ) {
 
1474
   if ( new_x.getFunctionSpace() == continuousFunction(*this) ) {
1475
1475
       tmp = new_x.getDataC();
1476
1476
       Finley_Mesh_setCoordinates(mesh,&tmp);
1477
1477
   } else {
1478
 
       escript::Data new_x_inter=escript::Data( new_x,  continuousFunction(asAbstractContinuousDomain()) );
 
1478
       escript::Data new_x_inter=escript::Data( new_x,  continuousFunction(*this) );
1479
1479
       tmp = new_x_inter.getDataC();
1480
1480
       Finley_Mesh_setCoordinates(mesh,&tmp);
1481
1481
   }
1559
1559
 
1560
1560
bool MeshAdapter::ownSample(int fs_code, index_t id) const
1561
1561
{
1562
 
#ifdef PASO_MPI
 
1562
#ifdef ESYS_MPI
1563
1563
    index_t myFirstNode=0, myLastNode=0, k=0;
1564
1564
    index_t* globalNodeIndex=0;
1565
1565
    Finley_Mesh* mesh_p=m_finleyMesh.get();
1586
1586
//
1587
1587
// creates a SystemMatrixAdapter stiffness matrix an initializes it with zeros
1588
1588
//
1589
 
SystemMatrixAdapter MeshAdapter::newSystemMatrix(
 
1589
ASM_ptr MeshAdapter::newSystemMatrix(
1590
1590
                                                 const int row_blocksize,
1591
1591
                                                 const escript::FunctionSpace& row_functionspace,
1592
1592
                                                 const int column_blocksize,
1633
1633
   }
1634
1634
   checkPasoError();
1635
1635
   Paso_SystemMatrixPattern_free(fsystemMatrixPattern);
1636
 
   return SystemMatrixAdapter(fsystemMatrix,row_blocksize,row_functionspace,column_blocksize,column_functionspace);
 
1636
   SystemMatrixAdapter* sma=new SystemMatrixAdapter(fsystemMatrix, row_blocksize, row_functionspace, column_blocksize, column_functionspace);
 
1637
   return ASM_ptr(sma);
 
1638
//   return SystemMatrixAdapter(fsystemMatrix,row_blocksize,row_functionspace,column_blocksize,column_functionspace);
1637
1639
}
1638
1640
 
1639
1641
//
1640
1642
// creates a TransportProblemAdapter
1641
1643
//
1642
 
TransportProblemAdapter MeshAdapter::newTransportProblem(
 
1644
ATP_ptr MeshAdapter::newTransportProblem(
1643
1645
                                                         const bool useBackwardEuler,
1644
1646
                                                         const int blocksize,
1645
1647
                                                         const escript::FunctionSpace& functionspace,
1666
1668
   transportProblem=Paso_TransportProblem_alloc(useBackwardEuler,fsystemMatrixPattern,blocksize);
1667
1669
   checkPasoError();
1668
1670
   Paso_SystemMatrixPattern_free(fsystemMatrixPattern);
1669
 
   return TransportProblemAdapter(transportProblem,useBackwardEuler,blocksize,functionspace);
 
1671
   TransportProblemAdapter* tpa=new TransportProblemAdapter(transportProblem,useBackwardEuler,blocksize,functionspace);
 
1672
   return ATP_ptr(tpa);
 
1673
//   return TransportProblemAdapter(transportProblem,useBackwardEuler,blocksize,functionspace);
1670
1674
}
1671
1675
 
1672
1676
//
2039
2043
 
2040
2044
escript::Data MeshAdapter::getX() const
2041
2045
{
2042
 
   return continuousFunction(asAbstractContinuousDomain()).getX();
 
2046
   return continuousFunction(*this).getX();
2043
2047
}
2044
2048
 
2045
2049
escript::Data MeshAdapter::getNormal() const
2046
2050
{
2047
 
   return functionOnBoundary(asAbstractContinuousDomain()).getNormal();
 
2051
   return functionOnBoundary(*this).getNormal();
2048
2052
}
2049
2053
 
2050
2054
escript::Data MeshAdapter::getSize() const
2051
2055
{
2052
 
   return escript::function(asAbstractContinuousDomain()).getSize();
 
2056
   return escript::function(*this).getSize();
2053
2057
}
2054
2058
 
2055
2059
const int* MeshAdapter::borrowSampleReferenceIDs(int functionSpaceType) const
2397
2401
  return order;
2398
2402
}
2399
2403
 
2400
 
ReferenceElementSetWrapper::ReferenceElementSetWrapper(ElementTypeId id, index_t order, index_t reducedOrder)
 
2404
bool MeshAdapter::supportsContactElements() const
 
2405
{
 
2406
  return true;
 
2407
}
 
2408
 
 
2409
ReferenceElementSetWrapper::ReferenceElementSetWrapper(Finley_ElementTypeId id, index_t order, index_t reducedOrder)
2401
2410
{
2402
2411
  m_refSet = Finley_ReferenceElementSet_alloc(id, order, reducedOrder);
2403
2412
}