/[escript]/branches/trilinos_from_5897/dudley/src/Assemble_CopyNodalData.cpp
ViewVC logotype

Diff of /branches/trilinos_from_5897/dudley/src/Assemble_CopyNodalData.cpp

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

revision 6078 by caltinay, Wed Mar 2 04:13:26 2016 UTC revision 6079 by caltinay, Mon Mar 21 12:22:38 2016 UTC
# Line 14  Line 14 
14  *  *
15  *****************************************************************************/  *****************************************************************************/
16    
 /****************************************************************************  
   
   Assemblage routines: copies data between different types nodal  
   representations  
   
 *****************************************************************************/  
   
17  #include "Assemble.h"  #include "Assemble.h"
18  #include "Util.h"  #include "Util.h"
19    
20  namespace dudley {  namespace dudley {
21    
22  void Assemble_CopyNodalData(Dudley_NodeFile* nodes, escript::Data* out, const escript::Data* in)  void Assemble_CopyNodalData(const NodeFile* nodes, escript::Data& out,
23                                const escript::Data& in)
24  {  {
25      if (!nodes)      if (!nodes)
26          return;          return;
27    
28      const int mpiSize = nodes->MPIInfo->size;      const int mpiSize = nodes->MPIInfo->size;
29      const int numComps = out->getDataPointSize();      const int numComps = out.getDataPointSize();
30      const int in_data_type = in->getFunctionSpace().getTypeCode();      const int in_data_type = in.getFunctionSpace().getTypeCode();
31      const int out_data_type = out->getFunctionSpace().getTypeCode();      const int out_data_type = out.getFunctionSpace().getTypeCode();
32    
33      // check out and in      // check out and in
34      if (numComps != in->getDataPointSize()) {      if (numComps != in.getDataPointSize()) {
35          throw DudleyException("Assemble_CopyNodalData: number of components of input and output Data do not match.");          throw DudleyException("Assemble_CopyNodalData: number of components of input and output Data do not match.");
36      } else if (!out->actsExpanded()) {      } else if (!out.actsExpanded()) {
37          throw DudleyException("Assemble_CopyNodalData: expanded Data object is expected for output data.");          throw DudleyException("Assemble_CopyNodalData: expanded Data object is expected for output data.");
38      }      }
39    
40      // more sophisticated test needed for overlapping node/DOF counts      // more sophisticated test needed for overlapping node/DOF counts
41      if (in_data_type == DUDLEY_NODES) {      if (in_data_type == DUDLEY_NODES) {
42          if (!in->numSamplesEqual(1, Dudley_NodeFile_getNumNodes(nodes))) {          if (!in.numSamplesEqual(1, nodes->getNumNodes())) {
             throw DudleyException("Assemble_CopyNodalData: illegal number of samples of input Data object");  
         }  
     } else if (in_data_type == DUDLEY_REDUCED_NODES) {  
         if (!in->numSamplesEqual(1, Dudley_NodeFile_getNumReducedNodes(nodes))) {  
43              throw DudleyException("Assemble_CopyNodalData: illegal number of samples of input Data object");              throw DudleyException("Assemble_CopyNodalData: illegal number of samples of input Data object");
44          }          }
45      } else if (in_data_type == DUDLEY_DEGREES_OF_FREEDOM) {      } else if (in_data_type == DUDLEY_DEGREES_OF_FREEDOM) {
46          if (!in->numSamplesEqual(1, Dudley_NodeFile_getNumDegreesOfFreedom(nodes))) {          if (!in.numSamplesEqual(1, nodes->getNumDegreesOfFreedom())) {
47              throw DudleyException("Assemble_CopyNodalData: illegal number of samples of input Data object");              throw DudleyException("Assemble_CopyNodalData: illegal number of samples of input Data object");
48          }          }
49          if ((((out_data_type == DUDLEY_NODES) || (out_data_type == DUDLEY_DEGREES_OF_FREEDOM)) && !in->actsExpanded() && (mpiSize > 1))) {          if ((((out_data_type == DUDLEY_NODES) || (out_data_type == DUDLEY_DEGREES_OF_FREEDOM)) && !in.actsExpanded() && (mpiSize > 1))) {
50    
51              throw DudleyException("Assemble_CopyNodalData: DUDLEY_DEGREES_OF_FREEDOM to DUDLEY_NODES or DUDLEY_DEGREES_OF_FREEDOM requires expanded input data on more than one processor.");              throw DudleyException("Assemble_CopyNodalData: DUDLEY_DEGREES_OF_FREEDOM to DUDLEY_NODES or DUDLEY_DEGREES_OF_FREEDOM requires expanded input data on more than one processor.");
52          }          }
     } else if (in_data_type == DUDLEY_REDUCED_DEGREES_OF_FREEDOM) {  
         if (!in->numSamplesEqual(1, Dudley_NodeFile_getNumReducedDegreesOfFreedom(nodes))) {  
             throw DudleyException("Assemble_CopyNodalData: illegal number of samples of input Data object");  
         }  
         if ((out_data_type == DUDLEY_DEGREES_OF_FREEDOM) && !in->actsExpanded() && (mpiSize > 1)) {  
             throw DudleyException("Assemble_CopyNodalData: DUDLEY_REDUCED_DEGREES_OF_FREEDOM to DUDLEY_DEGREES_OF_FREEDOM requires expanded input data on more than one processor.");  
         }  
53      } else {      } else {
54          throw DudleyException("Assemble_CopyNodalData: illegal function space type for target object");          throw DudleyException("Assemble_CopyNodalData: illegal function space type for target object");
55      }      }
56    
57      if (out_data_type == DUDLEY_NODES) {      dim_t numOut = 0;
58          if (!out->numSamplesEqual(1, Dudley_NodeFile_getNumNodes(nodes))) {      switch (out_data_type) {
59              throw DudleyException("Assemble_CopyNodalData: illegal number of samples of output Data object");          case DUDLEY_NODES:
60          }              numOut = nodes->getNumNodes();
61      } else if (out_data_type == DUDLEY_REDUCED_NODES) {              break;
62          if (!out->numSamplesEqual(1, Dudley_NodeFile_getNumReducedNodes(nodes))) {  
63              throw DudleyException("Assemble_CopyNodalData: illegal number of samples of output Data object");          case DUDLEY_DEGREES_OF_FREEDOM:
64          }              numOut = nodes->getNumDegreesOfFreedom();
65      } else if (out_data_type == DUDLEY_DEGREES_OF_FREEDOM) {              break;
66          if (!out->numSamplesEqual(1, Dudley_NodeFile_getNumDegreesOfFreedom(nodes))) {  
67              throw DudleyException("Assemble_CopyNodalData: illegal number of samples of output Data object");          default:
68          }              throw escript::ValueError("Assemble_CopyNodalData: illegal function space type for source object");
     } else if (out_data_type == DUDLEY_REDUCED_DEGREES_OF_FREEDOM) {  
         if (!out->numSamplesEqual(1, Dudley_NodeFile_getNumReducedDegreesOfFreedom(nodes))) {  
             throw DudleyException("Assemble_CopyNodalData: illegal number of samples of output Data object");  
         }  
     } else {  
         throw DudleyException("Assemble_CopyNodalData: illegal function space type for source object");  
69      }      }
70    
71      size_t numComps_size = numComps * sizeof(double);      if (!out.numSamplesEqual(1, numOut)) {
72            throw escript::ValueError("Assemble_CopyNodalData: illegal number of samples of output Data object");
73        }
74    
75        const size_t numComps_size = numComps * sizeof(double);
76    
77      /**************************** DUDLEY_NODES ******************************/      /**************************** DUDLEY_NODES ******************************/
78      if (in_data_type == DUDLEY_NODES) {      if (in_data_type == DUDLEY_NODES) {
79          out->requireWrite();          out.requireWrite();
         if (out_data_type == DUDLEY_NODES) {  
 #pragma omp parallel for  
             for (index_t n = 0; n < nodes->nodesMapping->numNodes; n++) {  
                 memcpy(out->getSampleDataRW(n), in->getSampleDataRO(n), numComps_size);  
             }  
         } else if (out_data_type == DUDLEY_REDUCED_NODES) {  
 #pragma omp parallel for  
             for (index_t n = 0; n < nodes->reducedNodesMapping->numTargets; n++) {  
                 memcpy(out->getSampleDataRW(n),  
                        in->getSampleDataRO(nodes->reducedNodesMapping->map[n]),  
                        numComps_size);  
             }  
         } else if (out_data_type == DUDLEY_DEGREES_OF_FREEDOM) {  
             const dim_t nComps = nodes->degreesOfFreedomDistribution->getMyNumComponents();  
 #pragma omp parallel for  
             for (index_t n = 0; n < nComps; n++) {  
                 memcpy(out->getSampleDataRW(n),  
                        in->getSampleDataRO(nodes->degreesOfFreedomMapping->map[n]), numComps_size);  
             }  
         } else if (out_data_type == DUDLEY_REDUCED_DEGREES_OF_FREEDOM) {  
             const dim_t nComps = nodes->reducedDegreesOfFreedomDistribution->getMyNumComponents();  
 #pragma omp parallel for  
             for (index_t n = 0; n < nComps; n++) {  
                 memcpy(out->getSampleDataRW(n),  
                        in->getSampleDataRO(nodes->reducedDegreesOfFreedomMapping->map[n]), numComps_size);  
             }  
         }  
     /************************ DUDLEY_REDUCED_NODES **************************/  
     } else if (in_data_type == DUDLEY_REDUCED_NODES) {  
         out->requireWrite();  
80          if (out_data_type == DUDLEY_NODES) {          if (out_data_type == DUDLEY_NODES) {
             throw DudleyException("Assemble_CopyNodalData: cannot copy from reduced nodes to nodes.");  
         } else if (out_data_type == DUDLEY_REDUCED_NODES) {  
81  #pragma omp parallel for  #pragma omp parallel for
82              for (index_t n = 0; n < nodes->reducedNodesMapping->numNodes; n++) {              for (index_t n = 0; n < numOut; n++) {
83                  memcpy(out->getSampleDataRW(n), in->getSampleDataRO(n), numComps_size);                  memcpy(out.getSampleDataRW(n), in.getSampleDataRO(n), numComps_size);
84              }              }
85          } else if (out_data_type == DUDLEY_DEGREES_OF_FREEDOM) {          } else if (out_data_type == DUDLEY_DEGREES_OF_FREEDOM) {
86              throw DudleyException("Assemble_CopyNodalData: cannot copy from reduced nodes to degrees of freedom.");              const index_t* map = nodes->borrowDegreesOfFreedomTarget();
         } else if (out_data_type == DUDLEY_REDUCED_DEGREES_OF_FREEDOM) {  
             const dim_t nComps = nodes->reducedDegreesOfFreedomDistribution->getMyNumComponents();  
87  #pragma omp parallel for  #pragma omp parallel for
88              for (index_t n = 0; n < nComps; n++) {              for (index_t n = 0; n < numOut; n++) {
89                  const dim_t k = nodes->reducedDegreesOfFreedomMapping->map[n];                  memcpy(out.getSampleDataRW(n), in.getSampleDataRO(map[n]),
90                  memcpy(out->getSampleDataRW(n),                         numComps_size);
                        in->getSampleDataRO(nodes->reducedNodesMapping->target[k]), numComps_size);  
91              }              }
92          }          }
93      /********************** DUDLEY_DEGREES_OF_FREEDOM ***********************/      /********************** DUDLEY_DEGREES_OF_FREEDOM ***********************/
94      } else if (in_data_type == DUDLEY_DEGREES_OF_FREEDOM) {      } else if (in_data_type == DUDLEY_DEGREES_OF_FREEDOM) {
95          out->requireWrite();          out.requireWrite();
96          if (out_data_type == DUDLEY_NODES) {          if (out_data_type == DUDLEY_NODES) {
97              paso::Coupler_ptr coupler(new paso::Coupler(nodes->degreesOfFreedomConnector, numComps));              paso::Coupler_ptr coupler(new paso::Coupler(nodes->degreesOfFreedomConnector, numComps));
98              // safe provided coupler->copyAll is called before the pointer              // safe provided coupler->copyAll is called before the pointer
99              // in "in" is invalidated              // in "in" is invalidated
100              const_cast<escript::Data*>(in)->resolve();              const_cast<escript::Data*>(&in)->resolve();
101              coupler->startCollect(in->getDataRO());                coupler->startCollect(in.getDataRO());  
102              const double* recv_buffer = coupler->finishCollect();              const double* recv_buffer = coupler->finishCollect();
103              const index_t upperBound = nodes->degreesOfFreedomDistribution->getMyNumComponents();              const index_t upperBound = nodes->getNumDegreesOfFreedom();
104                const index_t* target = nodes->borrowTargetDegreesOfFreedom();
105  #pragma omp parallel for  #pragma omp parallel for
106              for (index_t n = 0; n < nodes->numNodes; n++) {              for (index_t n = 0; n < numOut; n++) {
107                  const index_t k = nodes->degreesOfFreedomMapping->target[n];                  const index_t k = target[n];
108                  if (k < upperBound) {                  if (k < upperBound) {
109                      memcpy(out->getSampleDataRW(n), in->getSampleDataRO(k), numComps_size);                      memcpy(out.getSampleDataRW(n), in.getSampleDataRO(k),
110                               numComps_size);
111                  } else {                  } else {
112                      memcpy(out->getSampleDataRW(n),                      memcpy(out.getSampleDataRW(n),
113                             &recv_buffer[(k - upperBound) * numComps], numComps_size);                             &recv_buffer[(k - upperBound) * numComps],
114                  }                             numComps_size);
             }  
         } else if (out_data_type == DUDLEY_REDUCED_NODES) {  
             paso::Coupler_ptr coupler(new paso::Coupler(nodes->degreesOfFreedomConnector, numComps));  
             // safe provided coupler->copyAll is called before the pointer  
             // in "in" is invalidated  
             const_cast<escript::Data*>(in)->resolve();  
             coupler->startCollect(in->getDataRO());    
             const double* recv_buffer = coupler->finishCollect();  
             const index_t upperBound = nodes->degreesOfFreedomDistribution->getMyNumComponents();  
   
 #pragma omp parallel for  
             for (index_t n = 0; n < nodes->reducedNodesMapping->numTargets; n++) {  
                 const index_t l = nodes->reducedNodesMapping->map[n];  
                 const index_t k = nodes->degreesOfFreedomMapping->target[l];  
                 if (k < upperBound) {  
                     memcpy(out->getSampleDataRW(n), in->getSampleDataRO(k), numComps_size);  
                 } else {  
                     memcpy(out->getSampleDataRW(n),  
                            &recv_buffer[(k - upperBound)*numComps], numComps_size);  
115                  }                  }
116              }              }
117          } else if (out_data_type == DUDLEY_DEGREES_OF_FREEDOM) {          } else if (out_data_type == DUDLEY_DEGREES_OF_FREEDOM) {
             const dim_t nComps = nodes->degreesOfFreedomDistribution->getMyNumComponents();  
118  #pragma omp parallel for  #pragma omp parallel for
119              for (index_t n = 0; n < nComps; n++) {              for (index_t n = 0; n < numOut; n++) {
120                  memcpy(out->getSampleDataRW(n), in->getSampleDataRO(n), numComps_size);                  memcpy(out.getSampleDataRW(n), in.getSampleDataRO(n),
121              }                         numComps_size);
         } else if (out_data_type == DUDLEY_REDUCED_DEGREES_OF_FREEDOM) {  
             const dim_t nComps = nodes->reducedDegreesOfFreedomDistribution->getMyNumComponents();  
 #pragma omp parallel for  
             for (index_t n = 0; n < nComps; n++) {  
                 const dim_t k = nodes->reducedDegreesOfFreedomMapping->map[n];  
                 memcpy(out->getSampleDataRW(n),  
                        in->getSampleDataRO(nodes->degreesOfFreedomMapping->target[k]), numComps_size);  
             }  
         }  
   
     /****************** DUDLEY_REDUCED_DEGREES_OF_FREEDOM *******************/  
     } else if (in_data_type == DUDLEY_REDUCED_DEGREES_OF_FREEDOM) {  
         if (out_data_type == DUDLEY_NODES) {  
             throw DudleyException("Assemble_CopyNodalData: cannot copy from reduced degrees of freedom to nodes.");  
         } else if (out_data_type == DUDLEY_REDUCED_NODES) {  
             paso::Coupler_ptr coupler(new paso::Coupler(nodes->reducedDegreesOfFreedomConnector, numComps));  
             const dim_t upperBound = nodes->reducedDegreesOfFreedomDistribution->getMyNumComponents();  
             // safe provided coupler->copyAll is called before the pointer  
             // in "in" is invalidated  
             const_cast<escript::Data*>(in)->resolve();  
             coupler->startCollect(in->getDataRO());    
             out->requireWrite();  
             const double* recv_buffer = coupler->finishCollect();  
 #pragma omp parallel for  
             for (index_t n = 0; n < nodes->reducedNodesMapping->numTargets; n++) {  
                 const index_t l = nodes->reducedNodesMapping->map[n];  
                 const index_t k = nodes->reducedDegreesOfFreedomMapping->target[l];  
                 if (k < upperBound) {  
                     memcpy(out->getSampleDataRW(n), in->getSampleDataRO(k), numComps_size);  
                 } else {  
                     memcpy(out->getSampleDataRW(n),  
                            &recv_buffer[(k - upperBound)*numComps], numComps_size);  
                 }  
             }  
         } else if (out_data_type == DUDLEY_REDUCED_DEGREES_OF_FREEDOM) {  
             const dim_t nComps = nodes->reducedDegreesOfFreedomDistribution->getMyNumComponents();  
             out->requireWrite();  
 #pragma omp parallel for  
             for (index_t n = 0; n < nComps; n++) {  
                 memcpy(out->getSampleDataRW(n), in->getSampleDataRO(n), numComps_size);  
122              }              }
         } else if (out_data_type == DUDLEY_DEGREES_OF_FREEDOM) {  
             throw DudleyException("Assemble_CopyNodalData: cannot copy from reduced degrees of freedom to degrees of freedom.");  
123          }          }
124      }      } // in_data_type
125  }  }
126    
127  } // namespace dudley  } // namespace dudley

Legend:
Removed from v.6078  
changed lines
  Added in v.6079

  ViewVC Help
Powered by ViewVC 1.1.26