Ignore:
Timestamp:
01/31/19 12:12:52 (5 years ago)
Author:
yushan
Message:

branch merged with trunk @1645. arch file (ep&mpi) added for ADA

File:
1 edited

Legend:

Unmodified
Added
Removed
  • XIOS/dev/dev_trunk_omp/src/transformation/domain_algorithm_interpolate.cpp

    r1630 r1646  
    2020#include "interpolate_domain.hpp" 
    2121#include "grid.hpp" 
     22#ifdef _usingEP 
    2223using namespace ep_lib; 
     24#endif 
    2325 
    2426namespace xios { 
     
    3234                                                                     std::map<int, int>& elementPositionInGridDst2AxisPosition, 
    3335                                                                     std::map<int, int>& elementPositionInGridDst2DomainPosition) 
     36TRY 
    3437{ 
    3538  std::vector<CDomain*> domainListDestP = gridDst->getDomains(); 
     
    4245  return (new CDomainAlgorithmInterpolate(domainListDestP[domainDstIndex], domainListSrcP[domainSrcIndex], interpolateDomain)); 
    4346} 
     47CATCH 
    4448 
    4549bool CDomainAlgorithmInterpolate::registerTrans() 
     50TRY 
    4651{ 
    4752  CGridTransformationFactory<CDomain>::registerTransformation(TRANS_INTERPOLATE_DOMAIN, create); 
    4853} 
     54CATCH 
    4955 
    5056CDomainAlgorithmInterpolate::CDomainAlgorithmInterpolate(CDomain* domainDestination, CDomain* domainSource, CInterpolateDomain* interpDomain) 
    5157: CDomainAlgorithmTransformation(domainDestination, domainSource), interpDomain_(interpDomain), writeToFile_(false), readFromFile_(false) 
     58TRY 
    5259{ 
    5360  CContext* context = CContext::getCurrent(); 
     
    94101     
    95102} 
     103CATCH 
    96104 
    97105/*! 
     
    99107*/ 
    100108void CDomainAlgorithmInterpolate::computeRemap() 
     109TRY 
    101110{ 
    102111  using namespace sphereRemap; 
     
    305314  CArray<double,2> boundsLonSrcUnmasked(nVertexSrc,nSrcLocalUnmasked); 
    306315  CArray<double,2> boundsLatSrcUnmasked(nVertexSrc,nSrcLocalUnmasked); 
     316  CArray<double,1> areaSrcUnmasked(nSrcLocalUnmasked); 
     317   
    307318  long int * globalSrcUnmasked = new long int [nSrcLocalUnmasked]; 
    308319 
    309320  nSrcLocalUnmasked=0 ; 
     321  bool hasSrcArea=domainSrc_->hasArea && !domainSrc_->radius.isEmpty() && !interpDomain_->use_area.isEmpty() && interpDomain_->use_area==true  ; 
     322  double srcAreaFactor ; 
     323  if (hasSrcArea) srcAreaFactor=1./(domainSrc_->radius*domainSrc_->radius) ; 
     324   
    310325  for (int idx=0 ; idx < nSrcLocal; idx++) 
    311326  { 
     
    317332        boundsLatSrcUnmasked(n,nSrcLocalUnmasked) = boundsLatSrc(n,idx) ; 
    318333      } 
     334      if (hasSrcArea) areaSrcUnmasked(nSrcLocalUnmasked) = domainSrc_->areavalue(idx)*srcAreaFactor ; 
    319335      globalSrcUnmasked[nSrcLocalUnmasked]=globalSrc[idx] ; 
    320336      ++nSrcLocalUnmasked ; 
    321337    } 
    322338  } 
    323  
     339  
    324340 
    325341  int nDstLocalUnmasked = 0 ; 
     
    328344  CArray<double,2> boundsLonDestUnmasked(nVertexDest,nDstLocalUnmasked); 
    329345  CArray<double,2> boundsLatDestUnmasked(nVertexDest,nDstLocalUnmasked); 
     346  CArray<double,1>   areaDstUnmasked(nDstLocalUnmasked); 
     347 
    330348  long int * globalDstUnmasked = new long int [nDstLocalUnmasked]; 
    331349 
    332350  nDstLocalUnmasked=0 ; 
     351  bool hasDstArea=domainDest_->hasArea && !domainDest_->radius.isEmpty() && !interpDomain_->use_area.isEmpty() && interpDomain_->use_area==true ; 
     352  double dstAreaFactor ; 
     353  if (hasDstArea) dstAreaFactor=1./(domainDest_->radius*domainDest_->radius) ; 
    333354  for (int idx=0 ; idx < nDstLocal; idx++) 
    334355  { 
     
    340361        boundsLatDestUnmasked(n,nDstLocalUnmasked) = boundsLatDest(n,idx) ; 
    341362      } 
     363      if (hasDstArea) areaDstUnmasked(nDstLocalUnmasked) = domainDest_->areavalue(idx)*dstAreaFactor ; 
    342364      globalDstUnmasked[nDstLocalUnmasked]=globalDst[idx] ; 
    343365      ++nDstLocalUnmasked ; 
     
    345367  } 
    346368 
    347   mapper.setSourceMesh(boundsLonSrcUnmasked.dataFirst(), boundsLatSrcUnmasked.dataFirst(), nVertexSrc, nSrcLocalUnmasked, &srcPole[0], globalSrcUnmasked); 
    348   mapper.setTargetMesh(boundsLonDestUnmasked.dataFirst(), boundsLatDestUnmasked.dataFirst(), nVertexDest, nDstLocalUnmasked, &dstPole[0], globalDstUnmasked); 
     369  double* ptAreaSrcUnmasked = NULL ; 
     370  if (hasSrcArea) ptAreaSrcUnmasked=areaSrcUnmasked.dataFirst() ; 
     371 
     372  double* ptAreaDstUnmasked = NULL ; 
     373  if (hasDstArea) ptAreaDstUnmasked=areaDstUnmasked.dataFirst() ; 
     374 
     375  mapper.setSourceMesh(boundsLonSrcUnmasked.dataFirst(), boundsLatSrcUnmasked.dataFirst(), ptAreaSrcUnmasked, nVertexSrc, nSrcLocalUnmasked, &srcPole[0], globalSrcUnmasked); 
     376  mapper.setTargetMesh(boundsLonDestUnmasked.dataFirst(), boundsLatDestUnmasked.dataFirst(), ptAreaDstUnmasked, nVertexDest, nDstLocalUnmasked, &dstPole[0], globalDstUnmasked); 
    349377 
    350378  std::vector<double> timings = mapper.computeWeights(orderInterp,renormalize,quantity); 
     
    400428 
    401429} 
     430CATCH 
    402431 
    403432void CDomainAlgorithmInterpolate::processPole(std::map<int,std::vector<std::pair<int,double> > >& interMapValuePole, 
    404433                                              int nbGlobalPointOnPole) 
     434TRY 
    405435{ 
    406436  CContext* context = CContext::getCurrent(); 
     
    468498 
    469499} 
     500CATCH 
    470501 
    471502/*! 
     
    473504*/ 
    474505void CDomainAlgorithmInterpolate::computeIndexSourceMapping_(const std::vector<CArray<double,1>* >& dataAuxInputs) 
     506TRY 
    475507{ 
    476508  if (readFromFile_)   
     
    481513  } 
    482514} 
     515CATCH 
    483516 
    484517void CDomainAlgorithmInterpolate::writeRemapInfo(std::map<int,std::vector<std::pair<int,double> > >& interpMapValue) 
     518TRY 
    485519 
    486520  writeInterpolationInfo(fileToReadWrite_, interpMapValue); 
    487521} 
     522CATCH 
    488523 
    489524void CDomainAlgorithmInterpolate::readRemapInfo() 
     525TRY 
    490526 
    491527  std::map<int,std::vector<std::pair<int,double> > > interpMapValue; 
     
    494530  exchangeRemapInfo(interpMapValue); 
    495531} 
     532CATCH 
    496533 
    497534void CDomainAlgorithmInterpolate::convertRemapInfo(std::map<int,std::vector<std::pair<int,double> > >& interpMapValue) 
     535TRY 
    498536{ 
    499537  CContext* context = CContext::getCurrent(); 
     
    520558  }       
    521559} 
     560CATCH 
    522561 
    523562/*! 
     
    525564*/ 
    526565void CDomainAlgorithmInterpolate::exchangeRemapInfo(std::map<int,std::vector<std::pair<int,double> > >& interpMapValue) 
     566TRY 
    527567{ 
    528568  CContext* context = CContext::getCurrent(); 
     
    706746  delete [] recvBuff; 
    707747} 
     748CATCH 
    708749  
    709750/*! Redefined some functions of CONetCDF4 to make use of them */ 
     
    716757int CDomainAlgorithmInterpolate::WriteNetCdf::addDimensionWrite(const StdString& name,  
    717758                                                                const StdSize size) 
     759TRY 
    718760{ 
    719761  return CONetCDF4::addDimension(name, size);   
    720762} 
     763CATCH 
    721764 
    722765int CDomainAlgorithmInterpolate::WriteNetCdf::addVariableWrite(const StdString& name, nc_type type, 
    723766                                                               const std::vector<StdString>& dim) 
     767TRY 
    724768{ 
    725769  return CONetCDF4::addVariable(name, type, dim); 
    726770} 
     771CATCH 
    727772 
    728773void CDomainAlgorithmInterpolate::WriteNetCdf::endDefinition() 
     774TRY 
    729775{ 
    730776  CONetCDF4::definition_end(); 
    731777} 
     778CATCH 
    732779 
    733780void CDomainAlgorithmInterpolate::WriteNetCdf::writeDataIndex(const CArray<int,1>& data, const StdString& name, 
     
    735782                                                              const std::vector<StdSize>* start, 
    736783                                                              const std::vector<StdSize>* count) 
     784TRY 
    737785{ 
    738786  CONetCDF4::writeData<int,1>(data, name, collective, record, start, count); 
    739787} 
     788CATCH 
    740789 
    741790void CDomainAlgorithmInterpolate::WriteNetCdf::writeDataIndex(const CArray<double,1>& data, const StdString& name, 
     
    743792                                                              const std::vector<StdSize>* start, 
    744793                                                              const std::vector<StdSize>* count) 
     794TRY 
    745795{ 
    746796  CONetCDF4::writeData<double,1>(data, name, collective, record, start, count); 
    747797} 
     798CATCH 
    748799 
    749800/* 
     
    754805void CDomainAlgorithmInterpolate::writeInterpolationInfo(std::string& filename, 
    755806                                                         std::map<int,std::vector<std::pair<int,double> > >& interpMapValue) 
     807TRY 
    756808{ 
    757809  CContext* context = CContext::getCurrent(); 
     
    806858  std::vector<StdSize> count(1, localNbWeight); 
    807859   
     860  int my_rank; 
     861  MPI_Comm_rank(client->intraComm, &my_rank); 
     862 
     863  #ifdef _usingEP 
    808864  int my_rank_loc = client->intraComm->ep_comm_ptr->size_rank_info[1].first; 
    809   int my_rank = client->intraComm->ep_comm_ptr->size_rank_info[0].first; 
    810    
     865  #elif _usingMPI 
     866  int my_rank_loc = 0; 
     867  #endif 
    811868   
    812869  
    813870  WriteNetCdf *netCdfWriter; 
    814871 
     872  #ifdef _usingEP 
    815873  MPI_Barrier_local(client->intraComm); 
    816    
     874  #endif 
     875 
    817876  if(my_rank_loc==0) 
    818877  { 
    819     info(100)<<"rank "<< my_rank <<" create weight info file"<< std::endl; 
     878    #pragma omp critical (_output) 
     879    { 
     880      info(100)<<"rank "<< my_rank <<" create weight info file"<< std::endl; 
     881    } 
    820882     
    821883    WriteNetCdf my_writer(filename, client->intraComm);   
    822     info(100)<<"rank "<< my_rank <<" file created"<< std::endl; 
     884    #pragma omp critical (_output) 
     885    { 
     886      info(100)<<"rank "<< my_rank <<" file created"<< std::endl; 
     887    } 
    823888    netCdfWriter = &my_writer;  
    824889   
     
    827892    netCdfWriter->addDimensionWrite("n_dst", n_dst); 
    828893    netCdfWriter->addDimensionWrite("n_weight", globalNbWeight); 
    829     info(100)<<"rank "<< my_rank <<" addDimensionWrite : n_src, n_dst, n_weight"<< std::endl; 
     894    #pragma omp critical (_output) 
     895    { 
     896      info(100)<<"rank "<< my_rank <<" addDimensionWrite : n_src, n_dst, n_weight"<< std::endl; 
     897    } 
    830898   
    831899    std::vector<StdString> dims(1,"n_weight"); 
     
    836904    netCdfWriter->addVariableWrite("weight", NC_DOUBLE, dims); 
    837905     
    838     info(100)<<"rank "<< my_rank <<" addVariableWrite : src_idx, dst_idx, weight"<< std::endl; 
     906    #pragma omp critical (_output) 
     907    { 
     908      info(100)<<"rank "<< my_rank <<" addVariableWrite : src_idx, dst_idx, weight"<< std::endl; 
     909    } 
    839910 
    840911    // End of definition 
    841912    netCdfWriter->endDefinition(); 
    842     info(100)<<"rank "<< my_rank <<" endDefinition"<< std::endl; 
     913    #pragma omp critical (_output) 
     914    { 
     915      info(100)<<"rank "<< my_rank <<" endDefinition"<< std::endl; 
     916    } 
    843917   
    844918    netCdfWriter->closeFile(); 
    845     info(100)<<"rank "<< my_rank <<" file closed"<< std::endl; 
    846   } 
    847    
     919    #pragma omp critical (_output) 
     920    { 
     921      info(100)<<"rank "<< my_rank <<" file closed"<< std::endl; 
     922    } 
     923  } 
     924   
     925  #ifdef _usingEP 
    848926  MPI_Barrier_local(client->intraComm); 
     927  #endif 
    849928   
    850929  #pragma omp critical (write_weight_data) 
    851930  { 
    852931    // open file 
    853     info(100)<<"rank "<< my_rank <<" writing in weight info file"<< std::endl; 
     932    #pragma omp critical (_output) 
     933    { 
     934      info(100)<<"rank "<< my_rank <<" writing in weight info file"<< std::endl; 
     935    } 
    854936     
    855937    WriteNetCdf my_writer(filename, true, client->intraComm);   
    856     info(100)<<"rank "<< my_rank <<" file opened"<< std::endl; 
     938    #pragma omp critical (_output) 
     939    { 
     940      info(100)<<"rank "<< my_rank <<" file opened"<< std::endl; 
     941    } 
    857942    netCdfWriter = &my_writer;  
    858943     
     
    864949      netCdfWriter->writeDataIndex(weights, "weight", false, 0, &start, &count); 
    865950       
    866       info(100)<<"rank "<< my_rank <<" WriteDataIndex : src_idx, dst_idx, weight"<< std::endl; 
     951      #pragma omp critical (_output) 
     952      { 
     953        info(100)<<"rank "<< my_rank <<" WriteDataIndex : src_idx, dst_idx, weight"<< std::endl; 
     954      } 
    867955    } 
    868956     
    869957    netCdfWriter->closeFile(); 
    870     info(100)<<"rank "<< my_rank <<" file closed"<< std::endl; 
     958    #pragma omp critical (_output) 
     959    { 
     960      info(100)<<"rank "<< my_rank <<" file closed"<< std::endl; 
     961    } 
    871962     
    872963  } 
    873964   
     965  #ifdef _usingEP 
    874966  MPI_Barrier_local(client->intraComm); 
    875    
    876  
    877 } 
     967  #endif 
     968   
     969 
     970} 
     971CATCH 
    878972 
    879973/*! 
     
    885979void CDomainAlgorithmInterpolate::readInterpolationInfo(std::string& filename, 
    886980                                                        std::map<int,std::vector<std::pair<int,double> > >& interpMapValue) 
     981TRY 
    887982{ 
    888983  int ncid ; 
     
    9461041  } 
    9471042} 
     1043CATCH 
    9481044 
    9491045void CDomainAlgorithmInterpolate::apply(const std::vector<std::pair<int,double> >& localIndex, 
     
    9521048                                            std::vector<bool>& flagInitial, 
    9531049                                            bool ignoreMissingValue, bool firstPass  ) 
     1050TRY 
    9541051{ 
    9551052  int nbLocalIndex = localIndex.size();    
     
    9931090  } 
    9941091} 
     1092CATCH 
    9951093 
    9961094void CDomainAlgorithmInterpolate::updateData(CArray<double,1>& dataOut) 
     1095TRY 
    9971096{ 
    9981097  if (detectMissingValue) 
     
    10011100    size_t nbIndex=dataOut.numElements() ;  
    10021101 
    1003     for (int idx = 0; idx < nbIndex; ++idx) 
    1004     { 
    1005       if (allMissing(idx)) dataOut(idx) = defaultValue; // If all data source are nan then data destination must be nan 
     1102    if (allMissing.numElements()>0) 
     1103    { 
     1104      for (int idx = 0; idx < nbIndex; ++idx) 
     1105      { 
     1106        if (allMissing(idx)) dataOut(idx) = defaultValue; // If all data source are nan then data destination must be nan 
     1107      } 
    10061108    } 
    10071109     
     
    10131115  } 
    10141116} 
    1015  
    1016 } 
     1117CATCH 
     1118 
     1119} 
Note: See TracChangeset for help on using the changeset viewer.