Changeset 2505 for XIOS3/branches

Timestamp:

05/23/23 16:20:10 (13 months ago)

Author:

jderouillat

Message:

Backport 2500, 2502 (convex cells) and 2503 (memory_report option) from XIOS2

Location:

XIOS3/branches/xios-3.0-beta

Files:

• extern/remap/src/meshutil.cpp (modified) (4 diffs)
• extern/remap/src/polyg.cpp (modified) (3 diffs)
• extern/remap/src/polyg.hpp (modified) (2 diffs)
• src/client.cpp (modified) (1 diff)
• src/cxios.cpp (modified) (2 diffs)
• src/cxios.hpp (modified) (1 diff)
• src/mem_checker.cpp (modified) (1 diff)
• src/node/variable.cpp (modified) (1 diff)
• src/server.cpp (modified) (1 diff)

XIOS3/branches/xios-3.0-beta/extern/remap/src/meshutil.cpp

@@ -11 +11 @@
 using namespace std;
 
-double computePolygoneArea(Elt& a, const Coord &pole)
+void computePolygonGeometry(Elt& a, const Coord &pole, double& area, Coord& bary)
 {
   using N = uint32_t;
@@ -49 +49 @@
   vector<N> indices_a_gno = mapbox::earcut<N>(polyline);
   
-  double area_a_gno=0 ;
+  area=0 ;
   for(int i=0;i<indices_a_gno.size()/3;++i)
     {
@@ -55 +55 @@
       Coord x1 = Ox * polyline[0][indices_a_gno[3*i+1]][0] + Oy* polyline[0][indices_a_gno[3*i+1]][1] + Oz ;
       Coord x2 = Ox * polyline[0][indices_a_gno[3*i+2]][0] + Oy* polyline[0][indices_a_gno[3*i+2]][1] + Oz ;
-      area_a_gno+=triarea(x0 * (1./norm(x0)),x1* (1./norm(x1)), x2* (1./norm(x2))) ;
-    }
-
+      area+=triarea(x0 * (1./norm(x0)),x1* (1./norm(x1)), x2* (1./norm(x2))) ;
+    }
+
+
+  bary = exact_barycentre(dstPolygon.data(),na) ;
+  
+  // check signed area of polygons on gnomonic plan => if <0 then switch barycenter and invert vertex numbering 
+  double signedArea = 0 ;
+  for(int n=0; n<na;n++) signedArea+= a_gno[n].x*a_gno[(n+1)%na].y-a_gno[(n+1)%na].x*a_gno[n].y ;
+  if (signedArea<0) 
+  {
+    bary = bary * (-1.) ;
+    switchOrientation(a.n, a.vertex,a.edge,a.d) ;
+  }
+  
   vect_points.clear();
   polyline.clear();
   indices_a_gno.clear();
-  delete [] a_gno ;
-  return area_a_gno ;
+
 }
 
@@ -68 +79 @@
 void cptEltGeom(Elt& elt, const Coord &pole)
 {
-  orient(elt.n, elt.vertex, elt.edge, elt.d, elt.x);
-  normals(elt, pole);
-  Coord gg;
-  elt.area = airbar(elt.n, elt.vertex, elt.edge, elt.d, pole, gg);
-  elt.x = gg;
+   orient(elt.n, elt.vertex, elt.edge, elt.d, elt.x);
+   normals(elt, pole);
+//  Coord gg;
+//  elt.area = airbar(elt.n, elt.vertex, elt.edge, elt.d, pole, gg);
+//  elt.x = gg;
 // overwrite area computation 
 
-  elt.area =  computePolygoneArea(elt, pole) ;
+   computePolygonGeometry(elt, pole, elt.area, elt.x) ;
+   normals(elt, pole);
+
 }
 

XIOS3/branches/xios-3.0-beta/extern/remap/src/polyg.cpp

@@ -18 +18 @@
   (because computing intersection area requires both polygons to have same orientation)
 */
+
 void orient(int N, Coord *vertex, Coord *edge, double *d, const Coord &g)
 {
@@ -24 +25 @@
   Coord vertical = crossprod(ga, gb);
   if (N > 2 && scalarprod(g, vertical) < 0)  // (GAxGB).G
-  {
-    for (int i = 0; i < N/2; i++)
-      swap(vertex[i], vertex[N-1-i]);
-
-    for (int i = 0; i < (N-1)/2; i++)
-    {
-      swap(edge[N-2-i], edge[i]);
-      swap(d[i], d[N-2-i]);
-    }
-  }
-}
-
+    switchOrientation(N, vertex, edge, d) ;
+  
+}
+
+void switchOrientation(int N, Coord *vertex, Coord *edge, double *d)
+{
+  for (int i = 0; i < N/2; i++) swap(vertex[i], vertex[N-1-i]);
+
+  for (int i = 0; i < (N-1)/2; i++)
+  {
+    swap(edge[N-2-i], edge[i]);
+    swap(d[i], d[N-2-i]);
+  }
+  
+}
 
 void normals(Coord *x, int n, Coord *a)
@@ -94 +98 @@
   {
     return  proj(gc_normalintegral(x, n));
+    //return new_barycentre(x,n) ;
   }
   else if (n == 0) return ORIGIN;
   else if (n == 2) return midpoint(x[0], x[1]);
   else if (n == 1) return x[0];
-
-  error_exit( "Missing return in : Coord exact_barycentre(const Coord *x, int n)" );
-  return ORIGIN;
-}
+}
+
+/* other methode to compute barycenter of spherical polygon
+   for a spherical polygon, the moment is half the sum of (a x b) / ||a x b|| * (angle between a and b) 
+   for each pair of consecutive vertices a,b.
+*/ 
+
+Coord new_barycentre(const Coord *x, int n)
+{
+  if (n >= 3)
+  {
+    Coord sum=ORIGIN ;
+    for (int i = 0; i < n; i++)
+    { 
+      sum=sum+crossprod(x[i],x[(i+1)%n])*(1./norm(crossprod(x[i],x[(i+1)%n]))*2.*atan2(norm(x[i]-x[(i+1)%n]),norm(x[i]+x[(i+1)%n]))) ; 
+    }
+    return proj(sum) ; 
+  }
+  else if (n == 0) return ORIGIN;
+  else if (n == 2) return midpoint(x[0], x[1]);
+  else if (n == 1) return x[0];
+}
+
+
 
 Coord sc_gc_moon_normalintegral(Coord a, Coord b, Coord pole)

XIOS3/branches/xios-3.0-beta/extern/remap/src/polyg.hpp

@@ -7 +7 @@
 
 void orient(int n, Coord *vertex, Coord *edge, double *d, const Coord &g);
-
+void switchOrientation(int N, Coord *vertex, Coord *edge, double *d) ;
 void normals(Coord *x, int n, Coord *a);
 
@@ -16 +16 @@
 double polygonarea(Coord *x, int n);
 Coord exact_barycentre(const Coord *x, int n) ;
+Coord new_barycentre(const Coord *x, int n) ;
 
 int packedPolygonSize(const Elt& e);

XIOS3/branches/xios-3.0-beta/src/client.cpp

@@ -520 +520 @@
       report(0)<< " Memory report : increasing it by a factor will increase performance, depending of the volume of data wrote in file at each time step of the file"<<endl ;
       report(100)<<CTimer::getAllCumulatedTime()<<endl ;
-      report(100)<<CMemChecker::getAllCumulatedMem()<<endl ;
+      if (CXios::reportMemory)
+      {
+        report(100)<<CMemChecker::getAllCumulatedMem()<<endl ;
+      }
       CWorkflowGraph::drawWorkFlowGraph_client();
 

XIOS3/branches/xios-3.0-beta/src/cxios.cpp

@@ -55 +55 @@
   bool CXios::checkSumSend=false ;
   bool CXios::logMemory=false ;
+  bool CXios::reportMemory=true ;
 
   CDaemonsManager*    CXios::daemonsManager_=nullptr ;
@@ -118 +119 @@
 
     logMemory = getin<bool>("log_memory", false);
+    reportMemory = getin<bool>("memory_report", true);
 
     globalComm=MPI_COMM_WORLD ;

XIOS3/branches/xios-3.0-beta/src/cxios.hpp

@@ -68 +68 @@
      static bool checkSumRecv; //!< For debugging, compute a checksum of fields received by the model through the XIOS client
 
-     static bool logMemory; //!< Activate memory monitoring for all XIOS process
+     static bool logMemory; //!< Activate memory monitoring for all XIOS process (generate CSV file for https://forge.ipsl.jussieu.fr/ioserver/chrome/site/XIOS_TOOLS/xios_memory.html)
+     static bool reportMemory; //!< Activate memory reporting for all XIOS process (report in log files)
 
      static const string defaultPoolId ;

XIOS3/branches/xios-3.0-beta/src/mem_checker.cpp

@@ -44 +44 @@
   {
     std::vector<double> memories(0);
+    if ( !CXios::reportMemory ) return memories ;
     memories.clear();
     

XIOS3/branches/xios-3.0-beta/src/node/variable.cpp

@@ -77 +77 @@
                                     "log_type",
                                     "max_buffer_size",
+                                    "memory_report",
                                     "memtrack_blocks",
                                     "memtrack_size",

XIOS3/branches/xios-3.0-beta/src/server.cpp

@@ -442 +442 @@
       report(0)<<"Performance report : Ratio : "<<CTimer::get("Process events").getCumulatedTime()/CTimer::get("XIOS server").getCumulatedTime()*100.<<"%"<<endl  ;
       report(100)<<CTimer::getAllCumulatedTime()<<endl ;
-      report(100)<<CMemChecker::getAllCumulatedMem()<<endl ;
+      if (CXios::reportMemory)
+      {
+        report(100)<<CMemChecker::getAllCumulatedMem()<<endl ;
+      }
       
       CWorkflowGraph::drawWorkFlowGraph_server();