Changeset 844 for XIOS/trunk/extern/remap/src/mapper.cpp

Timestamp:

04/27/16 17:23:48 (8 years ago)

Author:

ymipsl

Message:

Management of masked cell for remapping algorithm.

YM

File:

• XIOS/trunk/extern/remap/src/mapper.cpp (modified) (8 diffs)

XIOS/trunk/extern/remap/src/mapper.cpp

@@ -109 +109 @@
 }
 
-vector<double> Mapper::computeWeights(int interpOrder)
+vector<double> Mapper::computeWeights(int interpOrder, bool renormalize)
 {
         vector<double> timings;
@@ -152 +152 @@
         if (mpiRank == 0 && verbose) cout << "Remapping..." << endl;
         tic = cputime();
-        nWeights = remap(&targetElements[0], targetElements.size(), interpOrder);
+        nWeights = remap(&targetElements[0], targetElements.size(), interpOrder, renormalize);
         timings.push_back(cputime() - tic);
 
@@ -166 +166 @@
    @param order is the order of interpolaton (must be 1 or 2).
 */
-int Mapper::remap(Elt *elements, int nbElements, int order)
+int Mapper::remap(Elt *elements, int nbElements, int order, bool renormalize)
 {
         int mpiSize, mpiRank;
@@ -281 +281 @@
                                 {
                                         sendGrad[rank][jj] = sstree.localElements[recvElement[rank][j]].grad;
+//          cout<<"grad  "<<jj<<"  "<<recvElement[rank][j]<<"  "<<sendGrad[rank][jj]<<" "<<sstree.localElements[recvElement[rank][j]].grad<<endl ;
                                         sendNeighIds[rank][jj] = sstree.localElements[recvElement[rank][j]].src_id;
                                         jj++;
@@ -286 +287 @@
                                         {
                                                 sendGrad[rank][jj] = sstree.localElements[recvElement[rank][j]].gradNeigh[i];
-                                                sendNeighIds[rank][jj] = sstree.localElements[recvElement[rank][j]].neighId[i];
+//            cout<<"grad  "<<jj<<"  "<<sendGrad[rank][jj]<<" "<<sstree.localElements[recvElement[rank][j]].grad<<endl ;
+            sendNeighIds[rank][jj] = sstree.localElements[recvElement[rank][j]].neighId[i];
                                                 jj++;
                                         }
@@ -368 +370 @@
                                         int kk = n1 * (NMAX + 1) + k;
                                         GloId neighID = recvNeighIds[rank][kk];
-                                        if (neighID.ind == -1) break;
-                                        wgt_map[neighID] +=
-                                                w * scalarprod(recvGrad[rank][kk], (*it)->x);
+                                        if (neighID.ind != -1)  wgt_map[neighID] += w * scalarprod(recvGrad[rank][kk], (*it)->x);
                                 }
 
                         }
                 }
-                for (map<GloId,double>::iterator it = wgt_map.begin(); it != wgt_map.end(); it++)
-                {
-                        this->remapMatrix[i] = it->second / e.area;
+
+    double renorm=0;
+    if (renormalize) 
+      for (map<GloId,double>::iterator it = wgt_map.begin(); it != wgt_map.end(); it++) renorm+=it->second / e.area;
+    else renorm=1. ;
+
+    for (map<GloId,double>::iterator it = wgt_map.begin(); it != wgt_map.end(); it++)
+                {
+                        this->remapMatrix[i] = (it->second / e.area) / renorm;
                         this->srcAddress[i] = it->first.ind;
                         this->srcRank[i] = it->first.rank;
@@ -669 +675 @@
                 }
 
+/*
                 for (int i = 0; i < elt->n; i++)
                 {
@@ -674 +681 @@
                                 error_exit("neighbour not found");
                 }
+*/
         }
 }