Changeset 1220


Ignore:
Timestamp:
07/20/17 09:18:34 (7 years ago)
Author:
yushan
Message:

test_remap_omp tested on ADA except two fields

Location:
XIOS/dev/branch_openmp
Files:
18 edited

Legend:

Unmodified
Added
Removed
  • XIOS/dev/branch_openmp/extern/remap/src/gridRemap.hpp

    r1205 r1220  
    1414Coord readPole(std::istream&); 
    1515 
    16 extern CRemapGrid srcGrid; 
    17 extern CRemapGrid tgtGrid; 
     16 
    1817 
    1918} 
  • XIOS/dev/branch_openmp/extern/remap/src/intersect.cpp

    r1205 r1220  
    1515namespace sphereRemap { 
    1616 
     17extern CRemapGrid srcGrid; 
     18#pragma omp threadprivate(srcGrid) 
     19 
     20extern CRemapGrid tgtGrid; 
     21#pragma omp threadprivate(tgtGrid) 
     22 
    1723using namespace std; 
    1824 
     
    2127int neighbour_idx(const Elt& a, const Elt& b) 
    2228{ 
    23         for (int i = 0; i < a.n; i++) 
    24         { 
    25                 for (int j = 0; j < b.n; j++) 
    26                 { 
    27                         assert(squaredist(a.vertex[ i       ], b.vertex[ j       ]) > EPS*EPS || 
    28                                squaredist(a.vertex[(i+1)%a.n], b.vertex[(j+1)%b.n]) > EPS*EPS); 
    29                         if (   squaredist(a.vertex[ i       ], b.vertex[ j           ]) < 1e-13*1e-13 && 
    30                                squaredist(a.vertex[(i+1)%a.n], b.vertex[(j+b.n-1)%b.n]) < 1e-13*1e-13) 
    31                         { 
    32                                 return i; 
    33                         } 
    34                 } 
    35         } 
    36         return NOT_FOUND; 
     29  for (int i = 0; i < a.n; i++) 
     30  { 
     31    for (int j = 0; j < b.n; j++) 
     32    { 
     33      assert(squaredist(a.vertex[ i       ], b.vertex[ j       ]) > EPS*EPS || 
     34             squaredist(a.vertex[(i+1)%a.n], b.vertex[(j+1)%b.n]) > EPS*EPS); 
     35      if (   squaredist(a.vertex[ i       ], b.vertex[ j           ]) < 1e-13*1e-13 && 
     36             squaredist(a.vertex[(i+1)%a.n], b.vertex[(j+b.n-1)%b.n]) < 1e-13*1e-13) 
     37      { 
     38        return i; 
     39      } 
     40    } 
     41  } 
     42  return NOT_FOUND; 
    3743} 
    3844 
     
    205211bool isNeighbour(Elt& a, const Elt& b) 
    206212{ 
    207         // return neighbour_idx(a, b) != NOT_FOUND; 
     213  // return neighbour_idx(a, b) != NOT_FOUND; 
    208214  return insertNeighbour(a,b,false) ; 
    209215} 
     
    213219void intersect(Elt *a, Elt *b) 
    214220{ 
    215         int na = a->n; /* vertices of a */ 
    216         int nb = b->n; /* vertices of b */ 
    217         Coord *c   = new Coord[na+nb]; 
    218         Coord *c2  = new Coord[na+nb]; 
    219         Coord *xc  = new Coord[na+nb]; 
    220         Coord *xc2 = new Coord[na+nb]; 
    221         Coord gc, gc2; 
    222         double *d = new double[na+nb]; 
    223         double *d2 = new double[na+nb]; 
    224         double are, are2; 
    225         Ipt ipt[NMAX*NMAX]; 
    226         Ipt ipt2[NMAX*NMAX]; 
    227         ptsec(a, b, ipt); 
    228         /* make ipt2 transpose of ipt */ 
    229         for (int ii = 0; ii < na; ii++) 
    230                 for (int jj = 0; jj < nb; jj++) 
    231                         ipt2[jj*na+ii] = ipt[ii*nb+jj]; 
    232         list<Sgm> iscot; 
    233         recense(a, b, ipt, iscot, 0); 
    234         recense(b, a, ipt2, iscot, 1); 
    235  
    236         int nseg = iscot.size(); 
    237         int nc = 0; 
    238         int nc2 = 0; 
    239         while (iscot.size() && nc < 2) 
    240                 nc = assemble(iscot, c, d, xc); 
    241         while (iscot.size() && nc2 < 2) 
    242                 nc2 = assemble(iscot, c2, d2, xc2); 
    243 //      assert(nseg == nc + nc2 || nseg == 1); // unused segment 
     221  int na = a->n; /* vertices of a */ 
     222  int nb = b->n; /* vertices of b */ 
     223  Coord *c   = new Coord[na+nb]; 
     224  Coord *c2  = new Coord[na+nb]; 
     225  Coord *xc  = new Coord[na+nb]; 
     226  Coord *xc2 = new Coord[na+nb]; 
     227  Coord gc, gc2; 
     228  double *d = new double[na+nb]; 
     229  double *d2 = new double[na+nb]; 
     230  double are, are2; 
     231  Ipt ipt[NMAX*NMAX]; 
     232  Ipt ipt2[NMAX*NMAX]; 
     233  ptsec(a, b, ipt); 
     234  /* make ipt2 transpose of ipt */ 
     235  for (int ii = 0; ii < na; ii++) 
     236    for (int jj = 0; jj < nb; jj++) 
     237      ipt2[jj*na+ii] = ipt[ii*nb+jj]; 
     238  list<Sgm> iscot; 
     239  recense(a, b, ipt, iscot, 0); 
     240  recense(b, a, ipt2, iscot, 1); 
     241 
     242  int nseg = iscot.size(); 
     243  int nc = 0; 
     244  int nc2 = 0; 
     245  while (iscot.size() && nc < 2) 
     246    nc = assemble(iscot, c, d, xc); 
     247  while (iscot.size() && nc2 < 2) 
     248    nc2 = assemble(iscot, c2, d2, xc2); 
     249//  assert(nseg == nc + nc2 || nseg == 1); // unused segment 
    244250 
    245251        if (!(nseg == nc + nc2 || nseg == 1)) 
     
    260266 
    261267//    intersect_ym(a,b) ; 
    262         if (nc == 1) nc = 0; 
    263         if (nc2 == 1) nc2 = 0; 
    264         gc = barycentre(xc, nc); 
    265         gc2 = barycentre(xc2, nc2); 
    266         orient(nc, xc, c, d, gc); 
    267  
    268         Coord pole = srcGrid.pole; 
    269         if (pole == ORIGIN) pole = tgtGrid.pole; 
    270         const double MINBASE = 1e-11; 
    271         if (nc == 2) /* nc is the number of vertices of super mesh element */ 
    272         { 
    273                 double base = arcdist(xc[0], xc[1]); 
     268  if (nc == 1) nc = 0; 
     269  if (nc2 == 1) nc2 = 0; 
     270  gc = barycentre(xc, nc); 
     271  gc2 = barycentre(xc2, nc2); 
     272  orient(nc, xc, c, d, gc); 
     273 
     274  Coord pole = srcGrid.pole; 
     275  if (pole == ORIGIN) pole = tgtGrid.pole; 
     276  const double MINBASE = 1e-11; 
     277  if (nc == 2) /* nc is the number of vertices of super mesh element */ 
     278  { 
     279    double base = arcdist(xc[0], xc[1]); 
    274280cerr << "DID ARRIVE " << base << xc[0] << xc[1] << endl; 
    275                 gc = midpoint(gc, midpointSC(xc[0], xc[1])); 
    276                 /* intersection area `are` must be zero here unless we have one great and one small circle */ 
    277                 are = alun(base, fabs(scalarprod(xc[0], pole))); 
    278         } 
    279         else 
    280         { 
    281                 are = airbar(nc, xc, c, d, pole, gc); 
    282         } 
    283         if (nc2 == 2) 
    284         { 
    285                 double base = arcdist(xc2[0], xc2[1]); 
     281    gc = midpoint(gc, midpointSC(xc[0], xc[1])); 
     282    /* intersection area `are` must be zero here unless we have one great and one small circle */ 
     283    are = alun(base, fabs(scalarprod(xc[0], pole))); 
     284  } 
     285  else 
     286  { 
     287    are = airbar(nc, xc, c, d, pole, gc); 
     288  } 
     289  if (nc2 == 2) 
     290  { 
     291    double base = arcdist(xc2[0], xc2[1]); 
    286292cerr << "DID ARRIVE " << base << xc2[0] << xc2[1] << endl; 
    287                 assert(base > MINBASE); 
    288                 gc2 = midpoint(gc2, midpointSC(xc2[0], xc2[1])); 
    289                 are2 = alun(base, fabs(scalarprod(xc2[0], pole))); // 0 
    290         } 
    291         else 
    292         { 
    293                 are2 = airbar(nc2, xc2, c2, d2, pole, gc2); 
    294         } 
     293    assert(base > MINBASE); 
     294    gc2 = midpoint(gc2, midpointSC(xc2[0], xc2[1])); 
     295    are2 = alun(base, fabs(scalarprod(xc2[0], pole))); // 0 
     296  } 
     297  else 
     298  { 
     299    are2 = airbar(nc2, xc2, c2, d2, pole, gc2); 
     300  } 
    295301 
    296302//  double ym_area=intersect_ym(a,b) ; 
    297303 
    298304  if (nc > 1) 
    299         { 
    300                 /* create one super mesh polygon that src and dest point to */ 
    301                 Polyg *is = new Polyg; 
    302                 is->x = gc; 
    303                 is->area = are; 
    304                 is->id = b->id; 
    305                 is->src_id = b->src_id; 
    306                 is->n = nc; 
    307                 (a->is).push_back(is); 
    308                 (b->is).push_back(is); 
     305  { 
     306    /* create one super mesh polygon that src and dest point to */ 
     307    Polyg *is = new Polyg; 
     308    is->x = gc; 
     309    is->area = are; 
     310    is->id = b->id; 
     311    is->src_id = b->src_id; 
     312    is->n = nc; 
     313    (a->is).push_back(is); 
     314    (b->is).push_back(is); 
    309315/* 
    310           if (  2*fabs(are-ym_area)/(are+ym_area) > 1.1 && ym_area>1e-8) 
     316    if (  2*fabs(are-ym_area)/(are+ym_area) > 1.1 && ym_area>1e-8) 
    311317    { 
    312318      cout<<"Big area difference : "<<are<<"  "<<ym_area<<endl ; 
     
    314320    } 
    315321*/ 
    316 //              cout<<"intersection : "<<are<<" "<< ym_area<<"  diff : "<<fabs(are-ym_area)<<"  ratio : "<<fabs(are-ym_area)/(0.5*(are+ym_area))<<endl ; 
    317         } 
    318         if (nc2 > 1) 
    319         { 
    320                 Polyg *is = new Polyg; 
    321                 is->x = gc2; 
    322                 is->area = are2; 
    323                 is->id = b->id; /* intersection holds id of corresponding source element (see Elt class definition for details about id) */ 
    324                 is->src_id = b->src_id; 
    325                 is->n = nc2; 
    326                 (a->is).push_back(is); 
    327                 (b->is).push_back(is); 
     322//    cout<<"intersection : "<<are<<" "<< ym_area<<"  diff : "<<fabs(are-ym_area)<<"  ratio : "<<fabs(are-ym_area)/(0.5*(are+ym_area))<<endl ; 
     323  } 
     324  if (nc2 > 1) 
     325  { 
     326    Polyg *is = new Polyg; 
     327    is->x = gc2; 
     328    is->area = are2; 
     329    is->id = b->id; /* intersection holds id of corresponding source element (see Elt class definition for details about id) */ 
     330    is->src_id = b->src_id; 
     331    is->n = nc2; 
     332    (a->is).push_back(is); 
     333    (b->is).push_back(is); 
    328334/* 
    329     if (        2*fabs(are-ym_area)/(are+ym_area) > 1.1 && ym_area>1e-8 ) 
     335    if (  2*fabs(are-ym_area)/(are+ym_area) > 1.1 && ym_area>1e-8 ) 
    330336    { 
    331337      cout<<"Big area difference : "<<are<<"  "<<ym_area<<endl ; 
     
    333339    } 
    334340*/ 
    335 //              cout<<"intersection : "<<are2<<" "<< ym_area<<"  diff : "<<fabs(are-ym_area)<<"  ratio : "<<fabs(are-ym_area)/(0.5*(are+ym_area))<<endl ; 
    336         } 
     341//    cout<<"intersection : "<<are2<<" "<< ym_area<<"  diff : "<<fabs(are-ym_area)<<"  ratio : "<<fabs(are-ym_area)/(0.5*(are+ym_area))<<endl ; 
     342  } 
    337343/* 
    338344  if (nc<=1 && nc2<=1) 
     
    344350  } 
    345351*/ 
    346         delete [] c; 
    347         delete [] c2; 
    348         delete [] xc; 
    349         delete [] xc2; 
    350         delete [] d; 
    351         delete [] d2; 
    352 } 
    353  
    354 } 
     352  delete [] c; 
     353  delete [] c2; 
     354  delete [] xc; 
     355  delete [] xc2; 
     356  delete [] d; 
     357  delete [] d2; 
     358} 
     359 
     360} 
  • XIOS/dev/branch_openmp/extern/remap/src/intersection_ym.cpp

    r1205 r1220  
    1313 
    1414namespace sphereRemap { 
     15 
     16extern CRemapGrid srcGrid; 
     17#pragma omp threadprivate(srcGrid) 
     18 
     19extern CRemapGrid tgtGrid; 
     20#pragma omp threadprivate(tgtGrid) 
     21 
    1522 
    1623using namespace std; 
  • XIOS/dev/branch_openmp/extern/remap/src/mapper.cpp

    r1205 r1220  
    1414 
    1515namespace sphereRemap { 
     16 
     17extern CRemapGrid srcGrid; 
     18#pragma omp threadprivate(srcGrid) 
     19 
     20extern CRemapGrid tgtGrid; 
     21#pragma omp threadprivate(tgtGrid) 
     22 
    1623 
    1724/* A subdivition of an array into N sub-arays 
  • XIOS/dev/branch_openmp/extern/remap/src/parallel_tree.cpp

    r1205 r1220  
    1414 
    1515namespace sphereRemap { 
     16 
     17extern CRemapGrid srcGrid; 
     18#pragma omp threadprivate(srcGrid) 
     19 
     20extern CRemapGrid tgtGrid; 
     21#pragma omp threadprivate(tgtGrid) 
    1622 
    1723static const int assignLevel = 2; 
  • XIOS/dev/branch_openmp/extern/remap/src/polyg.cpp

    r1205 r1220  
    77 
    88#include "polyg.hpp" 
     9#include <stdio.h> 
    910 
    1011namespace sphereRemap { 
     
    4546Coord barycentre(const Coord *x, int n) 
    4647{ 
     48 
    4749        if (n == 0) return ORIGIN; 
    4850        Coord bc = ORIGIN; 
    4951        for (int i = 0; i < n; i++) 
     52        { 
    5053                bc = bc + x[i]; 
     54        }        
    5155        /* both distances can be equal down to roundoff when norm(bc) < mashineepsilon  
    5256           which can occur when weighted with tiny area */ 
    53  
    54         assert(squaredist(bc, proj(bc)) <= squaredist(bc, proj(bc * (-1.0)))); 
    55         //if (squaredist(bc, proj(bc)) > squaredist(bc, proj(bc * (-1.0)))) return proj(bc * (-1.0)); 
     57   
     58         
     59  //assert(squaredist(bc, proj(bc)) <= squaredist(bc, proj(bc * (-1.0))));       
    5660 
    5761        return proj(bc); 
     
    173177                t[1] = x[i]; 
    174178                t[2] = x[ii]; 
    175                 double sc=scalarprod(crossprod(t[1] - t[0], t[2] - t[0]), t[0]) ; 
     179    double sc=scalarprod(crossprod(t[1] - t[0], t[2] - t[0]), t[0]) ; 
    176180                assert(sc >= -1e-10); // Error: tri a l'env (wrong orientation) 
    177181                double area_gc = triarea(t[0], t[1], t[2]); 
     182                if(area_gc<=0) printf("area_gc = %e\n", area_gc); 
    178183                double area_sc_gc_moon = 0; 
    179184                if (d[i]) /* handle small circle case */ 
     
    183188                        char sgl = (mext > 0) ? -1 : 1; 
    184189                        area_sc_gc_moon = sgl * alun(arcdist(t[1], t[2]), fabs(scalarprod(t[1], pole))); 
     190                        //if(area_sc_gc_moon<=0) printf("area_sc_gc_moon = %e\n", area_sc_gc_moon); 
    185191                        gg_exact = gg_exact + sc_gc_moon_normalintegral(t[1], t[2], pole); 
    186192                } 
    187193                area += area_gc + area_sc_gc_moon; /* for "spherical circle segment" sum triangular part (at) and "small moon" and => account for small circle */ 
    188194                g[i] = barycentre(t, 3) * (area_gc + area_sc_gc_moon); 
     195                //printf("g[%d] = (%e,%e,%e) * (%e+%e) = (%e,%e,%e) norm = %e\n", i, barycentre(t, 3).x, barycentre(t, 3).y, barycentre(t, 3).z, area_gc,  area_sc_gc_moon, g[i].x, g[i].y, g[i].z, norm(g[i])); 
    189196        } 
    190197        gg = barycentre(g, N); 
  • XIOS/dev/branch_openmp/extern/src_ep_dev/ep_lib.cpp

    r1196 r1220  
    77using namespace std; 
    88 
    9 std::list< ep_lib::MPI_Request* > ** EP_PendingRequests = 0; 
     9std::list< ep_lib::MPI_Request* > * EP_PendingRequests = 0; 
    1010#pragma omp threadprivate(EP_PendingRequests) 
    1111 
  • XIOS/dev/branch_openmp/extern/src_ep_dev/ep_message.cpp

    r1196 r1220  
    5353      } 
    5454      #elif _intelmpi 
    55       //#pragma omp critical (_mpi_call) 
    56       //::MPI_Improbe(MPI_ANY_SOURCE, MPI_ANY_TAG, mpi_comm, &flag, &message, &status);  
    57       #pragma omp critical (_mpi_call) 
    58       { 
    59         ::MPI_Iprobe(MPI_ANY_SOURCE, MPI_ANY_TAG, mpi_comm, &flag, &status); 
    60         if(flag) 
    61         { 
    62           Debug("find message in mpi comm \n"); 
    63           mpi_source = status.MPI_SOURCE; 
    64           int tag = status.MPI_TAG; 
    65           ::MPI_Mprobe(mpi_source, tag, mpi_comm, &message, &status); 
    66  
    67         } 
    68       } 
     55      #pragma omp critical (_mpi_call) 
     56      ::MPI_Improbe(MPI_ANY_SOURCE, MPI_ANY_TAG, mpi_comm, &flag, &message, &status);  
    6957      #endif 
    7058       
     
    7260      { 
    7361 
    74         MPI_Message *msg_block = new MPI_Message; //printf("myRank = %d, new ok\n", myRank); 
    75         msg_block->mpi_message = message;  //printf("myRank = %d, msg_block->mpi_message = %d\n", myRank, msg_block->mpi_message); 
    76         msg_block->ep_tag = bitset<15>(status.MPI_TAG >> 16).to_ulong();  //printf("myRank = %d, msg_block->ep_tag = %d\n", myRank, msg_block->ep_tag); 
    77         int src_loc       = bitset<8> (status.MPI_TAG >> 8) .to_ulong();  //printf("myRank = %d, src_loc = %d\n",           myRank, src_loc); 
    78         int dest_loc      = bitset<8> (status.MPI_TAG)      .to_ulong();  //printf("myRank = %d, dest_loc = %d\n",          myRank, dest_loc); 
    79         int src_mpi       = status.MPI_SOURCE;                            //printf("myRank = %d, src_mpi = %d\n",           myRank, src_mpi); 
     62        MPI_Message *msg_block = new MPI_Message;  
     63        msg_block->mpi_message = message;   
     64        msg_block->ep_tag = bitset<15>(status.MPI_TAG >> 16).to_ulong();  
     65        int src_loc       = bitset<8> (status.MPI_TAG >> 8) .to_ulong();  
     66        int dest_loc      = bitset<8> (status.MPI_TAG)      .to_ulong(); 
     67        int src_mpi       = status.MPI_SOURCE; 
    8068              
    81         msg_block->ep_src  = get_ep_rank(comm, src_loc,  src_mpi);        //printf("myRank = %d, msg_block->ep_src = %d\n",  myRank, msg_block->ep_src); 
     69        msg_block->ep_src  = get_ep_rank(comm, src_loc,  src_mpi);        
    8270        int dest_mpi = comm.ep_comm_ptr->size_rank_info[2].first; 
    8371        int ep_dest = get_ep_rank(comm, dest_loc, dest_mpi); 
    84         //printf("myRank = %d, probed one message, ep_src = %d, ep_dest = %d, tag = %d, message = %d\n", myRank, msg_block->ep_src, ep_dest, msg_block->ep_tag, msg_block->mpi_message); 
    8572        msg_block->mpi_status = new ::MPI_Status(status); 
    8673 
     
    9380          #pragma omp flush 
    9481          ptr_comm_target->ep_comm_ptr->message_queue->push_back(*msg_block);   
    95           //printf("myRank = %d, push_back OK, ep_src = %d, ep_tag = %d, dest = %d(%d)\n", myRank,  
    96           //                                                                                   ptr_comm_target->ep_comm_ptr->message_queue->back().ep_src, 
    97           //                                                                                   ptr_comm_target->ep_comm_ptr->message_queue->back().ep_tag, 
    98           //                                                                                   ep_dest, dest_loc); 
    9982     
    10083          #pragma omp flush 
    10184        } 
    10285         
    103         delete msg_block; 
    104         //printf("myRank = %d, delete msg_block, queue size = %lu\n", myRank, ptr_comm_target->ep_comm_ptr->message_queue->size()); 
    105          
     86        delete msg_block;         
    10687      } 
    10788 
     
    146127      #elif _intelmpi 
    147128      #pragma omp critical (_mpi_call) 
    148       { 
    149         ::MPI_Iprobe(MPI_ANY_SOURCE, MPI_ANY_TAG, mpi_comm, &flag, &status); 
    150         if(flag) 
    151         { 
    152           Debug("find message in mpi comm \n"); 
    153           mpi_source = status.MPI_SOURCE; 
    154           int tag = status.MPI_TAG; 
    155           ::MPI_Mprobe(mpi_source, tag, mpi_comm, &message, &status); 
    156  
    157         } 
    158       } 
    159       //::MPI_Improbe(MPI_ANY_SOURCE, MPI_ANY_TAG, mpi_comm, &flag, &message, &status);        
     129      ::MPI_Improbe(MPI_ANY_SOURCE, MPI_ANY_TAG, mpi_comm, &flag, &message, &status);        
    160130      #endif 
    161131 
     
    213183      #elif _intelmpi 
    214184      #pragma omp critical (_mpi_call) 
    215       { 
    216         ::MPI_Iprobe(MPI_ANY_SOURCE, MPI_ANY_TAG, mpi_comm, &flag, &status); 
    217         if(flag) 
    218         { 
    219           Debug("find message in mpi comm \n"); 
    220           mpi_source = status.MPI_SOURCE; 
    221           int tag = status.MPI_TAG; 
    222           ::MPI_Mprobe(mpi_source, tag, mpi_comm, &message, &status); 
    223  
    224         } 
    225       } 
    226       //::MPI_Improbe(MPI_ANY_SOURCE, MPI_ANY_TAG, mpi_comm, &flag, &message, &status);        
     185      ::MPI_Improbe(MPI_ANY_SOURCE, MPI_ANY_TAG, mpi_comm, &flag, &message, &status);        
    227186      #endif 
    228187 
     
    251210          #pragma omp flush 
    252211          ptr_comm_target->ep_comm_ptr->message_queue->push_back(*msg_block); 
    253           //printf("probed one message, ep_src = %d, tag = %d, mpi_status = %p (%p), message = %d\n", msg_block->ep_src, msg_block->ep_tag, msg_block->mpi_status, &status, msg_block->mpi_message); 
    254212          #pragma omp flush 
    255213        } 
     
    286244        MPI_Imrecv((*it)->buf, recv_count, (*it)->ep_datatype, &message, *it); 
    287245        (*it)->type = 3; 
    288         //printf("request add = %p, mpi_request=%d\n", *it, (*it)->mpi_request); 
    289246        EP_PendingRequests->erase(it); 
    290247        it = EP_PendingRequests->begin(); 
  • XIOS/dev/branch_openmp/extern/src_ep_dev/ep_probe.cpp

    r1185 r1220  
    3030    #pragma omp flush 
    3131 
    32      
    33     //printf("iprobe, message queue size = %lu, queue = %p\n", comm.ep_comm_ptr->message_queue->size(), comm.ep_comm_ptr->message_queue); 
    34  
    35     #pragma omp flush 
    36  
    3732    #pragma omp critical (_query) 
    3833    if(!comm.ep_comm_ptr->message_queue->empty()) 
     
    4641        { 
    4742          Debug("find message\n"); 
    48           // printf("iprobe, find in local message queue %p, src = %d, tag = %d\n", comm.ep_comm_ptr->message_queue, it->ep_src, it->ep_tag); 
    49  
    5043          *flag = true; 
    5144 
     
    131124          #pragma omp critical (_query2) 
    132125          {               
    133             //printf("local message erased. src = %d, dest = %d, tag = %d\n", it->ep_src, it->ep_dest, it->ep_tag);      
    134126            delete it->mpi_status; 
    135127            comm.ep_comm_ptr->message_queue->erase(it); 
  • XIOS/dev/branch_openmp/extern/src_ep_dev/ep_recv.cpp

    r1196 r1220  
    8888    { 
    8989      EP_PendingRequests = new std::list< MPI_Request* >; 
    90       //printf("proc %d(%d) : EP_PendingRequests allocated, add = %p\n", dest_rank, world_rank, EP_PendingRequests);   
    9190    } 
    9291 
     
    9594 
    9695    Request_Check(); 
    97     //printf("proc %d(%d) : EP_PendingRequests insert one request, src = %d(%d), tag = %d(%d), size = %d; request add = %p\n",  
    98     //        dest_rank, world_rank, EP_PendingRequests->back()->ep_src, request->ep_src,  
    99     //        EP_PendingRequests->back()->ep_tag, request->ep_tag,  
    100     //        EP_PendingRequests->size(), request); 
    101      
    102     // check all EP_PendingRequests       
    103     //for(std::list<MPI_Request* >::iterator it = EP_PendingRequests->begin(); it!=EP_PendingRequests->end(); ) 
    104     //{ 
    105     //if((*it)->type == 3)  
    106     //{ 
    107     //    EP_PendingRequests->erase(it); 
    108    //     it = EP_PendingRequests->begin(); 
    109     //    continue; 
    110      // } 
    111          
    112       //int probed = false; 
    113       //MPI_Message pending_message; 
    114       //MPI_Status pending_status; 
    115      
    116       //MPI_Improbe((*it)->ep_src, (*it)->ep_tag, (*it)->comm, &probed, &pending_message, &pending_status); 
    117      
    118       //if(probed)  
    119       //{  
    120         //int count; 
    121         //MPI_Get_count(&pending_status, (*it)->ep_datatype, &count); 
    122         //MPI_Imrecv((*it)->buf, count, (*it)->ep_datatype, &pending_message, *it); 
    12396 
    124         //EP_PendingRequests->erase(it); 
    125         //if(EP_PendingRequests->empty()) return 0; 
    126          
    127         //it = EP_PendingRequests->begin(); 
    128         //continue; 
    129      // } 
    130  
    131       //it++; 
    132    // } 
    133      
    13497    return 0; 
    13598  } 
  • XIOS/dev/branch_openmp/extern/src_ep_dev/ep_test.cpp

    r1185 r1220  
    4141    if(request->type == 2)   // irecv message not probed 
    4242    { 
    43       Message_Check(request->comm); 
     43      Request_Check(); 
     44       
    4445      #pragma omp flush 
    45       MPI_Message message; 
    46       MPI_Improbe(request->ep_src, request->ep_tag, request->comm, flag, &message, status); 
    47       if(*flag) 
    48       { 
    49          
    50         int count; 
    51         MPI_Get_count(status, request->ep_datatype, &count); 
    52         MPI_Imrecv(request->buf, count, request->ep_datatype, &message, request); 
    53         printf("in ep_test, found message src = %d, tag = %d, type = %d\n", request->ep_src, request->ep_tag, request->type); 
    54         MPI_Test(request, flag, status); 
    55       } 
    56       return 0; 
     46       
    5747    } 
    5848 
     
    6252      ::MPI_Status mpi_status; 
    6353       
    64       ::MPI_Errhandler_set(MPI_COMM_WORLD_STD, MPI_ERRORS_RETURN); 
    65       int error_code = ::MPI_Test(mpi_request, flag, &mpi_status); 
    66       if (error_code != MPI_SUCCESS) { 
     54      ::MPI_Test(mpi_request, flag, &mpi_status); 
    6755       
    68          char error_string[BUFSIZ]; 
    69          int length_of_error_string, error_class; 
    70        
    71          ::MPI_Error_class(error_code, &error_class); 
    72          ::MPI_Error_string(error_class, error_string, &length_of_error_string); 
    73          printf("%s\n", error_string); 
    74       } 
    7556       
    7657      if(*flag) 
     
    8364        //MPI_Get_count(status, request->ep_datatype, &count); 
    8465        //check_sum_recv(request->buf, count, request->ep_datatype, request->ep_src, request->ep_tag, request->comm, 2); 
    85       } 
    86  
    87       status->ep_src = request->ep_src; 
    88       status->ep_tag = request->ep_tag; 
    89       status->ep_datatype = request->ep_datatype; 
    90  
    91        
     66      }   
    9267 
    9368      return 0; 
  • XIOS/dev/branch_openmp/extern/src_ep_dev/ep_type.hpp

    r1196 r1220  
    174174    } 
    175175 
    176     // ep_intercomm(ep_intercomm &ref) 
    177     // { 
    178     //   printf("calling copy Constructor of ep_intercomm\n"); 
    179     //   ep_intercomm return_intercomm; 
    180  
    181     //   return_intercomm.mpi_inter_comm = ref.mpi_inter_comm; 
    182     //   return_intercomm.intercomm_rank_map = ref.intercomm_rank_map; 
    183     //   return_intercomm.local_rank_map = ref.local_rank_map; 
    184     //   return_intercomm.remote_rank_map = ref.remote_rank_map; 
    185     //   return_intercomm.size_rank_info[0] = ref.size_rank_info[0]; 
    186     //   return_intercomm.size_rank_info[1] = ref.size_rank_info[1]; 
    187     //   return_intercomm.size_rank_info[2] = ref.size_rank_info[2]; 
    188     //   return_intercomm.local_comm = ref.local_comm; 
    189     //   return_intercomm.intercomm_tag = ref.intercomm_tag; 
    190     // } 
    191176 
    192177    bool operator == (ep_intercomm right) 
     
    212197    } 
    213198 
    214     // ep_intercomm operator = (ep_intercomm ref) 
    215     // { 
    216     //   printf("calling = operator of ep_intercomm\n"); 
    217     //   ep_intercomm return_intercomm; 
    218  
    219     //   return_intercomm.mpi_inter_comm = ref.mpi_inter_comm; 
    220     //   return_intercomm.intercomm_rank_map = ref.intercomm_rank_map; 
    221     //   return_intercomm.local_rank_map = ref.local_rank_map; 
    222     //   return_intercomm.remote_rank_map = ref.remote_rank_map; 
    223     //   return_intercomm.size_rank_info[0] = ref.size_rank_info[0]; 
    224     //   return_intercomm.size_rank_info[1] = ref.size_rank_info[1]; 
    225     //   return_intercomm.size_rank_info[2] = ref.size_rank_info[2]; 
    226     //   return_intercomm.local_comm = ref.local_comm; 
    227     //   return_intercomm.intercomm_tag = ref.intercomm_tag; 
    228     // } 
     199 
    229200  }; 
    230201 
     
    271242    } 
    272243 
    273     // ep_communicator operator = (ep_communicator ref) 
    274     // { 
    275     //   printf("calling = operator of ep_communicator\n"); 
    276     //   ep_communicator return_ep; 
    277  
    278     //   return_ep.intercomm = ref.intercomm; 
    279     //   return_ep.comm_label = ref.comm_label; 
    280     //   return_ep.message_queue = ref.message_queue; 
    281     //   return_ep.comm_list = ref.comm_list; 
    282     //   return_ep.size_rank_info[0] = ref.size_rank_info[0]; 
    283     //   return_ep.size_rank_info[1] = ref.size_rank_info[1]; 
    284     //   return_ep.size_rank_info[2] = ref.size_rank_info[2]; 
    285     // } 
    286244  }; 
    287245 
     
    387345    } 
    388346 
    389     // MPI_Comm operator = (MPI_Comm ref) 
    390     // { 
    391     //   printf("calling = operator of MPI_Comm\n"); 
    392     //   MPI_Comm return_comm; 
    393  
    394     //   return_comm.mpi_comm = ref.mpi_comm; 
    395     //   return_comm.is_ep = ref.is_ep; 
    396     //   return_comm.is_intercomm = ref.is_intercomm; 
    397     //   return_comm.my_buffer = ref.my_buffer; 
    398     //   return_comm.ep_barrier = ref.ep_barrier; 
    399     //   return_comm.rank_map = ref.rank_map; 
    400     //   return_comm.ep_comm_ptr = ref.ep_comm_ptr; 
    401     // } 
    402347  }; 
    403348 
     
    502447            //    <MPI_Fint,thread_num>   EP_Comm 
    503448 
    504   //static std::list<MPI_Request * > *EP_PendingRequests = 0;  
    505449} 
    506450 
  • XIOS/dev/branch_openmp/extern/src_ep_dev/ep_wait.cpp

    r1196 r1220  
    5959    ::MPI_Status* mpi_status = new ::MPI_Status[count]; 
    6060 
    61     //if(EP_PendingRequests != 0) printf("pending size = %d, add = %p\n", EP_PendingRequests->size(), EP_PendingRequests); 
    62  
    6361    while(std::accumulate(finished.begin(), finished.end(), 0) < count) 
    6462    { 
  • XIOS/dev/branch_openmp/inputs/REMAP/iodef.xml

    r1203 r1220  
    3636       
    3737      <file_group id="read_then_write_files" enabled=".TRUE."> 
    38        <file id="output_regular_pole" name="output_dst_regular" enabled=".FALSE."  > 
    39           <field field_ref="tmp_field_0" name="field_regular_0" enabled=".FALSE."/> 
    40           <field field_ref="dst_field_regular" name="field_regular" enabled=".FALSE."/> 
     38       <file id="output_regular_pole" name="output_dst_regular" enabled=".TRUE."  > 
     39          <field field_ref="tmp_field_0" name="field_regular_0" enabled=".TRUE."/> 
     40          <field field_ref="dst_field_regular" name="field_regular" enabled=".TRUE."/> 
    4141          <field field_ref="dst_field_regular_pole_0" name="field_regular_pole_0" enabled=".FALSE." /> 
    4242          <field field_ref="dst_field_regular_pole_1" name="field_regular_pole_1" enabled=".FALSE." /> 
    4343       </file> 
    44        <file id="output_dst_curvilinear" name="output_dst_curvilinear" enabled=".FALSE." > 
     44       <file id="output_dst_curvilinear" name="output_dst_curvilinear" enabled=".TRUE." > 
    4545          <field field_ref="tmp_field_1" operation="instant"/> 
    4646       </file> 
    47        <file id="output_dst_unstructured" name="output_dst_unstructured" enabled=".FALSE." > 
     47       <file id="output_dst_unstructured" name="output_dst_unstructured" enabled=".TRUE." > 
    4848          <field field_ref="tmp_field_2" operation="instant"/> 
    4949       </file> 
     
    5252      <file_group id="write_files" > 
    5353        <file id="output_2D" name="output_2D" enabled=".TRUE."> 
    54           <field field_ref="src_field_2D" name="field_src" enabled=".FALSE."/> 
    55           <field field_ref="src_field_2D_clone" name="field_src_clone" default_value="100000" enabled=".FALSE."/> 
    56           <field field_ref="src_field_2D" name="field_dst_regular_0"  domain_ref="dst_domain_regular_pole" enabled=".FALSE."/> 
     54          <field field_ref="src_field_2D" name="field_src" enabled=".TRUE."/> 
     55          <field field_ref="src_field_2D_clone" name="field_src_clone" default_value="100000" enabled=".TRUE."/> 
     56          <field field_ref="src_field_2D" name="field_dst_regular_0"  domain_ref="dst_domain_regular_pole" enabled=".TRUE."/> 
    5757          <field field_ref="dst_field_2D" name="field_dst_regular_1" enabled=".TRUE." /> 
    58           <field field_ref="dst_field_2D_regular_pole" name="field_dst_regular_2" enabled=".FALSE."/> 
    59           <field field_ref="dst_field_2D_clone" name="field_dst_regular_3" detect_missing_value=".false." default_value="100000" enabled=".FALSE."/> 
    60           <field field_ref="dst_field_2D_extract" name="field_dst_regular_4" enabled=".FALSE."/> 
     58          <field field_ref="dst_field_2D_regular_pole" name="field_dst_regular_2" enabled=".TRUE."/> 
     59          <field field_ref="dst_field_2D_clone" name="field_dst_regular_3" detect_missing_value=".false." default_value="100000" enabled=".TRUE."/> 
     60          <field field_ref="dst_field_2D_extract" name="field_dst_regular_4" enabled=".TRUE."/> 
    6161        </file>  
    62        <file id="output_3D" name="output_3D" enabled=".FALSE."> 
     62       <file id="output_3D" name="output_3D" enabled=".TRUE."> 
    6363          <field field_ref="src_field_3D" name="field_src" /> 
    6464          <field field_ref="src_field_3D_pression" name="field" /> 
    6565          <field field_ref="dst_field_3D_interp" name="field_dst_interp_domain" /> 
    66           <field field_ref="dst_field_3D_interp" name="field_dst_interp_domain_axis" domain_ref="dst_domain_regular_pole" enabled=".FALSE."/>   
    67        </file> 
    68        <file id="output_4D" name="output_4D" enabled=".FALSE."> 
     66          <field field_ref="dst_field_3D_interp" name="field_dst_interp_domain_axis" domain_ref="dst_domain_regular_pole" enabled=".TRUE."/>   
     67       </file> 
     68       <file id="output_4D" name="output_4D" enabled=".TRUE."> 
    6969          <field field_ref="src_field_4D" name="field_4D" /> 
    70           <field field_ref="dst_field_4D_extract" name="field_4D_extract" enabled=".FALSE."/> 
     70          <field field_ref="dst_field_4D_extract" name="field_4D_extract" enabled=".TRUE."/> 
    7171        </file> 
    7272     </file_group> 
    73      <file_group id="read_files" enabled=".FALSE." > 
     73     <file_group id="read_files" enabled=".TRUE." > 
    7474       <file id="output_src_regular" name="output_src_regular" mode="read" > 
    7575          <field id="src_field_regular" name="field" grid_ref="src_grid_regular_read" operation="instant"/> 
     
    212212 
    213213</simulation> 
    214  
  • XIOS/dev/branch_openmp/inputs/Unstruct/iodef.xml

    r1209 r1220  
    5353       <generate_rectilinear_domain id="domain_regular_pole3"/> 
    5454       <interpolate_domain write_weight="false" order="1" renormalize="true"/>        
    55        <zoom_domain ibegin="0" ni="70" jbegin="0" nj="70" />  
     55       <zoom_domain ibegin="0" ni="40" jbegin="0" nj="40" />  
    5656     </domain> 
    5757 
  • XIOS/dev/branch_openmp/src/client_server_mapping.cpp

    r843 r1220  
    6464  MPI_Allgather(&nbConnectedServer,1,MPI_INT,recvCount,1,MPI_INT,clientIntraComm) ; 
    6565 
     66   
     67  for(int i=0; i<nbClient; i++) 
     68    printf("MPI_Allgather : recvCount[%d] = %d\n", i, recvCount[i]); 
     69 
    6670  displ[0]=0 ; 
    6771  for(int n=1;n<nbClient;n++) displ[n]=displ[n-1]+recvCount[n-1] ; 
     
    7175 
    7276  MPI_Allgatherv(sendBuff,nbConnectedServer,MPI_INT,recvBuff,recvCount,displ,MPI_INT,clientIntraComm) ; 
     77 
     78  for(int i=0; i<recvSize; i++) 
     79    printf("MPI_Allgatherv : recvBuff[%d] = %d\n", i, recvBuff[i]); 
     80 
     81 
    7382  for(int n=0;n<recvSize;n++) clientRes[recvBuff[n]]++ ; 
    7483 
  • XIOS/dev/branch_openmp/src/filter/spatial_transform_filter.cpp

    r1203 r1220  
    5353    CSpatialTransformFilterEngine* spaceFilter = static_cast<CSpatialTransformFilterEngine*>(engine); 
    5454    CDataPacketPtr outputPacket = spaceFilter->applyFilter(data, outputDefaultValue); 
     55    printf("spalceFilter applied\n"); 
    5556    if (outputPacket) 
    5657      onOutputReady(outputPacket); 
     
    122123    double defaultValue = std::numeric_limits<double>::quiet_NaN(); 
    123124    if (0 != dataDest.numElements()) ignoreMissingValue = NumTraits<double>::isnan(dataDest(0)); 
    124  
     125     
     126    const std::list<CGridTransformation::SendingIndexGridSourceMap> *listLocalIndexSend_ptr = & (gridTransformation->getLocalIndexToSendFromGridSource()); 
     127     
    125128    const std::list<CGridTransformation::SendingIndexGridSourceMap>& listLocalIndexSend = gridTransformation->getLocalIndexToSendFromGridSource(); 
    126129    const std::list<CGridTransformation::RecvIndexGridDestinationMap>& listLocalIndexToReceive = gridTransformation->getLocalIndexToReceiveOnGridDest(); 
     
    131134    CArray<double,1> dataCurrentDest(dataSrc.copy()); 
    132135 
    133     std::list<CGridTransformation::SendingIndexGridSourceMap>::const_iterator itListSend  = listLocalIndexSend.begin(), 
    134                                                                               iteListSend = listLocalIndexSend.end(); 
     136    std::list<CGridTransformation::SendingIndexGridSourceMap>::const_iterator itListSend  = listLocalIndexSend_ptr->begin(), 
     137                                                                              iteListSend = listLocalIndexSend_ptr->end(); 
    135138    std::list<CGridTransformation::RecvIndexGridDestinationMap>::const_iterator itListRecv = listLocalIndexToReceive.begin(); 
    136139    std::list<size_t>::const_iterator itNbListRecv = listNbLocalIndexToReceive.begin(); 
     
    186189        int srcRank = itRecv->first; 
    187190        int countSize = itRecv->second.size(); 
    188          
    189191        MPI_Irecv(recvBuff + currentBuff, countSize, MPI_DOUBLE, srcRank, 12, client->intraComm, &sendRecvRequest[position]); 
    190192        position++; 
  • XIOS/dev/branch_openmp/src/test/test_remap_omp.f90

    r1203 r1220  
    222222  CALL xios_close_context_definition() 
    223223 
    224 !  CALL xios_get_domain_attr("src_domain_regular_read", ni=src_tmp_ni, nj=src_tmp_nj) 
    225 !  ALLOCATE(tmp_field_0(src_tmp_ni*src_tmp_nj)) 
    226  
    227 !  CALL xios_get_axis_attr("src_axis_curvilinear_read", n=src_tmp_n) 
    228 !  CALL xios_get_domain_attr("src_domain_curvilinear_read", ni=src_tmp_ni, nj=src_tmp_nj) 
    229 !  ALLOCATE(tmp_field_1(src_tmp_ni*src_tmp_nj*src_tmp_n)) 
    230  
    231 !  CALL xios_get_domain_attr("src_domain_unstructured_read", ni=src_tmp_ni, nj=src_tmp_nj) 
    232 !  ALLOCATE(tmp_field_2(src_tmp_ni*src_tmp_nj)) 
    233    
    234 !  CALL xios_recv_field("src_field_regular", tmp_field_0) 
    235 !  CALL xios_recv_field("src_field_curvilinear", tmp_field_1) 
    236 !  CALL xios_recv_field("src_field_unstructured", tmp_field_2) 
     224  CALL xios_get_domain_attr("src_domain_regular_read", ni=src_tmp_ni, nj=src_tmp_nj) 
     225  ALLOCATE(tmp_field_0(src_tmp_ni*src_tmp_nj)) 
     226 
     227  CALL xios_get_axis_attr("src_axis_curvilinear_read", n=src_tmp_n) 
     228  CALL xios_get_domain_attr("src_domain_curvilinear_read", ni=src_tmp_ni, nj=src_tmp_nj) 
     229  ALLOCATE(tmp_field_1(src_tmp_ni*src_tmp_nj*src_tmp_n)) 
     230 
     231  CALL xios_get_domain_attr("src_domain_unstructured_read", ni=src_tmp_ni, nj=src_tmp_nj) 
     232  ALLOCATE(tmp_field_2(src_tmp_ni*src_tmp_nj)) 
     233   
     234  CALL xios_recv_field("src_field_regular", tmp_field_0) 
     235  CALL xios_recv_field("src_field_curvilinear", tmp_field_1) 
     236  CALL xios_recv_field("src_field_unstructured", tmp_field_2) 
    237237 
    238238  DO ts=1,10 
     
    252252    CALL xios_send_field("src_field_4D",src_field_4D) 
    253253    CALL xios_send_field("src_field_3D_pression",src_field_pression) 
    254  !   CALL xios_send_field("tmp_field_0",tmp_field_0) 
    255  !   CALL xios_send_field("tmp_field_1",tmp_field_1) 
    256  !   CALL xios_send_field("tmp_field_2",tmp_field_2) 
     254    CALL xios_send_field("tmp_field_0",tmp_field_0) 
     255    CALL xios_send_field("tmp_field_1",tmp_field_1) 
     256    CALL xios_send_field("tmp_field_2",tmp_field_2) 
    257257    CALL wait_us(5000) ; 
    258258   ENDDO 
     
    262262  DEALLOCATE(src_lon, src_lat, src_boundslon,src_boundslat, src_field_2D) 
    263263  DEALLOCATE(dst_lon, dst_lat, dst_boundslon,dst_boundslat) 
    264   !DEALLOCATE(tmp_field_0, tmp_field_1, tmp_field_2) 
     264  DEALLOCATE(tmp_field_0, tmp_field_1, tmp_field_2) 
    265265   
    266266  CALL xios_finalize() 
     
    293293 
    294294 
    295  
Note: See TracChangeset for help on using the changeset viewer.