source: XIOS/dev/dev_ym/XIOS_COUPLING/src/node/context.cpp @ 1869

Last change on this file since 1869 was 1869, checked in by ymipsl, 19 months ago

Some update...

YM

  • Property copyright set to
    Software name : XIOS (Xml I/O Server)
    http://forge.ipsl.jussieu.fr/ioserver
    Creation date : January 2009
    Licence : CeCCIL version2
    see license file in root directory : Licence_CeCILL_V2-en.txt
    or http://www.cecill.info/licences/Licence_CeCILL_V2-en.html
    Holder : CEA/LSCE (Laboratoire des Sciences du CLimat et de l'Environnement)
    CNRS/IPSL (Institut Pierre Simon Laplace)
    Project Manager : Yann Meurdesoif
    yann.meurdesoif@cea.fr
File size: 93.3 KB
Line 
1#include "context.hpp"
2#include "attribute_template.hpp"
3#include "object_template.hpp"
4#include "group_template.hpp"
5
6#include "calendar_type.hpp"
7#include "duration.hpp"
8
9#include "context_client.hpp"
10#include "context_server.hpp"
11#include "nc4_data_output.hpp"
12#include "node_type.hpp"
13#include "message.hpp"
14#include "type.hpp"
15#include "xios_spl.hpp"
16#include "timer.hpp"
17#include "memtrack.hpp"
18#include <limits>
19#include <fstream>
20#include "server.hpp"
21#include "distribute_file_server2.hpp"
22#include "services_manager.hpp"
23#include "contexts_manager.hpp"
24#include "cxios.hpp"
25#include "client.hpp"
26#include "coupler_in.hpp"
27#include "coupler_out.hpp"
28
29namespace xios {
30
31  std::shared_ptr<CContextGroup> CContext::root;
32
33   /// ////////////////////// Définitions ////////////////////// ///
34
35   CContext::CContext(void)
36      : CObjectTemplate<CContext>(), CContextAttributes()
37      , calendar(), hasClient(false), hasServer(false)
38      , isPostProcessed(false), finalized(false)
39      , client(nullptr), server(nullptr)
40      , allProcessed(false), countChildContextFinalized_(0), isProcessingEvent_(false)
41
42   { /* Ne rien faire de plus */ }
43
44   CContext::CContext(const StdString & id)
45      : CObjectTemplate<CContext>(id), CContextAttributes()
46      , calendar(), hasClient(false), hasServer(false)
47      , isPostProcessed(false), finalized(false)
48      , client(nullptr), server(nullptr)
49      , allProcessed(false), countChildContextFinalized_(0), isProcessingEvent_(false)
50   { /* Ne rien faire de plus */ }
51
52   CContext::~CContext(void)
53   {
54     delete client;
55     delete server;
56     for (std::vector<CContextClient*>::iterator it = clientPrimServer.begin(); it != clientPrimServer.end(); it++)  delete *it;
57     for (std::vector<CContextServer*>::iterator it = serverPrimServer.begin(); it != serverPrimServer.end(); it++)  delete *it;
58
59   }
60
61   //----------------------------------------------------------------
62   //! Get name of context
63   StdString CContext::GetName(void)   { return (StdString("context")); }
64   StdString CContext::GetDefName(void){ return (CContext::GetName()); }
65   ENodeType CContext::GetType(void)   { return (eContext); }
66
67   //----------------------------------------------------------------
68
69   /*!
70   \brief Get context group (context root)
71   \return Context root
72   */
73   CContextGroup* CContext::getRoot(void)
74   TRY
75   {
76      if (root.get()==NULL) root=std::shared_ptr<CContextGroup>(new CContextGroup(xml::CXMLNode::GetRootName()));
77      return root.get();
78   }
79   CATCH
80
81   //----------------------------------------------------------------
82
83   /*!
84   \brief Get calendar of a context
85   \return Calendar
86   */
87   std::shared_ptr<CCalendar> CContext::getCalendar(void) const
88   TRY
89   {
90      return (this->calendar);
91   }
92   CATCH
93
94   //----------------------------------------------------------------
95
96   /*!
97   \brief Set a context with a calendar
98   \param[in] newCalendar new calendar
99   */
100   void CContext::setCalendar(std::shared_ptr<CCalendar> newCalendar)
101   TRY
102   {
103      this->calendar = newCalendar;
104   }
105   CATCH_DUMP_ATTR
106
107   //----------------------------------------------------------------
108   /*!
109   \brief Parse xml file and write information into context object
110   \param [in] node xmld node corresponding in xml file
111   */
112   void CContext::parse(xml::CXMLNode & node)
113   TRY
114   {
115      CContext::SuperClass::parse(node);
116
117      // PARSING POUR GESTION DES ENFANTS
118      xml::THashAttributes attributes = node.getAttributes();
119
120      if (attributes.end() != attributes.find("src"))
121      {
122         StdIFStream ifs ( attributes["src"].c_str() , StdIFStream::in );
123         if ( (ifs.rdstate() & std::ifstream::failbit ) != 0 )
124            ERROR("void CContext::parse(xml::CXMLNode & node)",
125                  <<endl<< "Can not open <"<<attributes["src"].c_str()<<"> file" );
126         if (!ifs.good())
127            ERROR("CContext::parse(xml::CXMLNode & node)",
128                  << "[ filename = " << attributes["src"] << " ] Bad xml stream !");
129         xml::CXMLParser::ParseInclude(ifs, attributes["src"], *this);
130      }
131
132      if (node.getElementName().compare(CContext::GetName()))
133         DEBUG("Le noeud is wrong defined but will be considered as a context !");
134
135      if (!(node.goToChildElement()))
136      {
137         DEBUG("Le context ne contient pas d'enfant !");
138      }
139      else
140      {
141         do { // Parcours des contextes pour traitement.
142
143            StdString name = node.getElementName();
144            attributes.clear();
145            attributes = node.getAttributes();
146
147            if (attributes.end() != attributes.find("id"))
148            {
149              DEBUG(<< "Definition node has an id,"
150                    << "it will not be taking account !");
151            }
152
153#define DECLARE_NODE(Name_, name_)    \
154   if (name.compare(C##Name_##Definition::GetDefName()) == 0) \
155   { C##Name_##Definition::create(C##Name_##Definition::GetDefName()) -> parse(node); continue; }
156#define DECLARE_NODE_PAR(Name_, name_)
157#include "node_type.conf"
158
159            DEBUG(<< "The element \'"     << name
160                  << "\' in the context \'" << CContext::getCurrent()->getId()
161                  << "\' is not a definition !");
162
163         } while (node.goToNextElement());
164
165         node.goToParentElement(); // Retour au parent
166      }
167   }
168   CATCH_DUMP_ATTR
169
170   //----------------------------------------------------------------
171   //! Show tree structure of context
172   void CContext::ShowTree(StdOStream & out)
173   TRY
174   {
175      StdString currentContextId = CContext::getCurrent() -> getId();
176      std::vector<CContext*> def_vector =
177         CContext::getRoot()->getChildList();
178      std::vector<CContext*>::iterator
179         it = def_vector.begin(), end = def_vector.end();
180
181      out << "<? xml version=\"1.0\" ?>" << std::endl;
182      out << "<"  << xml::CXMLNode::GetRootName() << " >" << std::endl;
183
184      for (; it != end; it++)
185      {
186         CContext* context = *it;
187         CContext::setCurrent(context->getId());
188         out << *context << std::endl;
189      }
190
191      out << "</" << xml::CXMLNode::GetRootName() << " >" << std::endl;
192      CContext::setCurrent(currentContextId);
193   }
194   CATCH
195
196   //----------------------------------------------------------------
197
198   //! Convert context object into string (to print)
199   StdString CContext::toString(void) const
200   TRY
201   {
202      StdOStringStream oss;
203      oss << "<" << CContext::GetName()
204          << " id=\"" << this->getId() << "\" "
205          << SuperClassAttribute::toString() << ">" << std::endl;
206      if (!this->hasChild())
207      {
208         //oss << "<!-- No definition -->" << std::endl; // fait planter l'incrémentation
209      }
210      else
211      {
212
213#define DECLARE_NODE(Name_, name_)    \
214   if (C##Name_##Definition::has(C##Name_##Definition::GetDefName())) \
215   oss << * C##Name_##Definition::get(C##Name_##Definition::GetDefName()) << std::endl;
216#define DECLARE_NODE_PAR(Name_, name_)
217#include "node_type.conf"
218
219      }
220      oss << "</" << CContext::GetName() << " >";
221      return (oss.str());
222   }
223   CATCH
224
225   //----------------------------------------------------------------
226
227   /*!
228   \brief Find all inheritace among objects in a context.
229   \param [in] apply (true) write attributes of parent into ones of child if they are empty
230                     (false) write attributes of parent into a new container of child
231   \param [in] parent unused
232   */
233   void CContext::solveDescInheritance(bool apply, const CAttributeMap * const UNUSED(parent))
234   TRY
235   {
236#define DECLARE_NODE(Name_, name_)    \
237   if (C##Name_##Definition::has(C##Name_##Definition::GetDefName())) \
238     C##Name_##Definition::get(C##Name_##Definition::GetDefName())->solveDescInheritance(apply);
239#define DECLARE_NODE_PAR(Name_, name_)
240#include "node_type.conf"
241   }
242   CATCH_DUMP_ATTR
243
244   //----------------------------------------------------------------
245
246   //! Verify if all root definition in the context have child.
247   bool CContext::hasChild(void) const
248   TRY
249   {
250      return (
251#define DECLARE_NODE(Name_, name_)    \
252   C##Name_##Definition::has(C##Name_##Definition::GetDefName())   ||
253#define DECLARE_NODE_PAR(Name_, name_)
254#include "node_type.conf"
255      false);
256}
257   CATCH
258
259   //----------------------------------------------------------------
260
261   void CContext::CleanTree(void)
262   TRY
263   {
264#define DECLARE_NODE(Name_, name_) C##Name_##Definition::ClearAllAttributes();
265#define DECLARE_NODE_PAR(Name_, name_)
266#include "node_type.conf"
267   }
268   CATCH
269
270   ///---------------------------------------------------------------
271
272
273    /*!
274    Sets client buffers.
275    \param [in] contextClient
276    \param [in] bufferForWriting True if buffers are used for sending data for writing
277    This flag is only true for client and server-1 for communication with server-2
278  */
279   void CContext::setClientServerBuffer(CContextClient* contextClient, bool bufferForWriting)
280   TRY
281   {
282      // Estimated minimum event size for small events (20 is an arbitrary constant just for safety)
283     const size_t minEventSize = CEventClient::headerSize + 20 * sizeof(int);
284
285      // Ensure there is at least some room for 20 of such events in the buffers
286      size_t minBufferSize = std::max(CXios::minBufferSize, 20 * minEventSize);
287
288#define DECLARE_NODE(Name_, name_)    \
289     if (minBufferSize < sizeof(C##Name_##Definition)) minBufferSize = sizeof(C##Name_##Definition);
290#define DECLARE_NODE_PAR(Name_, name_)
291#include "node_type.conf"
292#undef DECLARE_NODE
293#undef DECLARE_NODE_PAR
294
295     // Compute the buffer sizes needed to send the attributes and data corresponding to fields
296     std::map<int, StdSize> maxEventSize;
297     std::map<int, StdSize> bufferSize = getAttributesBufferSize(maxEventSize, contextClient, bufferForWriting);
298     std::map<int, StdSize> dataBufferSize = getDataBufferSize(maxEventSize, contextClient, bufferForWriting);
299
300     std::map<int, StdSize>::iterator it, ite = dataBufferSize.end();
301     for (it = dataBufferSize.begin(); it != ite; ++it)
302       if (it->second > bufferSize[it->first]) bufferSize[it->first] = it->second;
303
304     // Apply the buffer size factor, check that we are above the minimum buffer size and below the maximum size
305     ite = bufferSize.end();
306     for (it = bufferSize.begin(); it != ite; ++it)
307     {
308       it->second *= CXios::bufferSizeFactor;
309       if (it->second < minBufferSize) it->second = minBufferSize;
310       if (it->second > CXios::maxBufferSize) it->second = CXios::maxBufferSize;
311     }
312
313     // Leaders will have to send some control events so ensure there is some room for those in the buffers
314     if (contextClient->isServerLeader())
315     {
316       const std::list<int>& ranks = contextClient->getRanksServerLeader();
317       for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
318       {
319         if (!bufferSize.count(*itRank))
320         {
321           bufferSize[*itRank] = minBufferSize;
322           maxEventSize[*itRank] = minEventSize;
323         }
324       }
325     }
326     contextClient->setBufferSize(bufferSize, maxEventSize);
327   }
328   CATCH_DUMP_ATTR
329
330   //! Verify whether a context is initialized
331   bool CContext::isInitialized(void)
332   TRY
333   {
334     return hasClient;
335   }
336   CATCH_DUMP_ATTR
337
338
339   void CContext::init(CServerContext* parentServerContext, MPI_Comm intraComm, int serviceType)
340   TRY
341   {
342     parentServerContext_ = parentServerContext ;
343     if (serviceType==CServicesManager::CLIENT) 
344       initClient(intraComm, serviceType) ;
345     else
346       initServer(intraComm, serviceType) ;
347    }
348    CATCH_DUMP_ATTR
349
350
351
352//! Initialize client side
353   void CContext::initClient(MPI_Comm intraComm, int serviceType)
354   TRY
355   {
356      intraComm_=intraComm ;
357      MPI_Comm_rank(intraComm_, &intraCommRank_) ;
358      MPI_Comm_size(intraComm_, &intraCommSize_) ;
359
360      serviceType_ = CServicesManager::CLIENT ;
361      if (serviceType_==CServicesManager::CLIENT)
362      {
363        hasClient=true ;
364        hasServer=false ;
365      }
366      contextId_ = getId() ;
367     
368      attached_mode=true ;
369      if (!CXios::isUsingServer()) attached_mode=false ;
370
371
372      string contextRegistryId=getId() ;
373      registryIn=new CRegistry(intraComm);
374      registryIn->setPath(contextRegistryId) ;
375     
376      int commRank ;
377      MPI_Comm_rank(intraComm_,&commRank) ;
378      if (commRank==0) registryIn->fromFile("xios_registry.bin") ;
379      registryIn->bcastRegistry() ;
380      registryOut=new CRegistry(intraComm_) ;
381      registryOut->setPath(contextRegistryId) ;
382     
383   }
384   CATCH_DUMP_ATTR
385
386   
387   void CContext::initServer(MPI_Comm intraComm, int serviceType)
388   TRY
389   {
390     hasServer=true;
391     intraComm_=intraComm ;
392     MPI_Comm_rank(intraComm_, &intraCommRank_) ;
393     MPI_Comm_size(intraComm_, &intraCommSize_) ;
394
395     serviceType_=serviceType ;
396
397     if (serviceType_==CServicesManager::GATHERER)
398     {
399       hasClient=true ;
400       hasServer=true ;
401     }
402     else if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::OUT_SERVER)
403     {
404       hasClient=false ;
405       hasServer=true ;
406     }
407
408     CXios::getContextsManager()->getContextId(getId(), contextId_, intraComm) ;
409     
410     registryIn=new CRegistry(intraComm);
411     registryIn->setPath(contextId_) ;
412     
413     int commRank ;
414     MPI_Comm_rank(intraComm_,&commRank) ;
415     if (commRank==0) registryIn->fromFile("xios_registry.bin") ;
416   
417     registryIn->bcastRegistry() ;
418     registryOut=new CRegistry(intraComm) ;
419     registryOut->setPath(contextId_) ;
420
421   }
422   CATCH_DUMP_ATTR
423
424
425  void CContext::createClientInterComm(MPI_Comm interCommClient, MPI_Comm interCommServer) // for servers
426  TRY
427  {
428    MPI_Comm intraCommClient ;
429    MPI_Comm_dup(intraComm_, &intraCommClient);
430    comms.push_back(intraCommClient);
431    // attached_mode=parentServerContext_->isAttachedMode() ; //ym probably inherited from source context
432    server = new CContextServer(this,intraComm_, interCommServer); // check if we need to dupl. intraComm_ ?
433    client = new CContextClient(this,intraCommClient,interCommClient);
434    client->setAssociatedServer(server) ; 
435    server->setAssociatedClient(client) ; 
436
437  }
438  CATCH_DUMP_ATTR
439
440  void CContext::createServerInterComm(void) 
441  TRY
442  {
443   
444    MPI_Comm interCommClient, interCommServer ;
445
446    if (serviceType_ == CServicesManager::CLIENT)
447    {
448
449      int commRank ;
450      MPI_Comm_rank(intraComm_,&commRank) ;
451      if (commRank==0)
452      {
453        if (attached_mode) CXios::getContextsManager()->createServerContext(CClient::getPoolRessource()->getId(), CXios::defaultServerId, 0, getContextId()) ;
454        else if (CXios::usingServer2) CXios::getContextsManager()->createServerContext(CXios::defaultPoolId, CXios::defaultGathererId, 0, getContextId()) ;
455        else  CXios::getContextsManager()->createServerContext(CXios::defaultPoolId, CXios::defaultServerId, 0, getContextId()) ;
456      }
457
458      MPI_Comm interComm ;
459     
460      if (attached_mode)
461      {
462        parentServerContext_->createIntercomm(CClient::getPoolRessource()->getId(), CXios::defaultServerId, 0, getContextId(), intraComm_, 
463                                              interCommClient, interCommServer) ;
464        int type ; 
465        if (commRank==0) CXios::getServicesManager()->getServiceType(CClient::getPoolRessource()->getId(), CXios::defaultServerId, 0, type) ;
466        MPI_Bcast(&type,1,MPI_INT,0,intraComm_) ;
467        setCurrent(getId()) ; // getCurrent/setCurrent may be supress, it can cause a lot of trouble
468      }
469      else if (CXios::usingServer2)
470      { 
471//      CXios::getContextsManager()->createServerContextIntercomm(CXios::defaultPoolId, CXios::defaultGathererId, 0, getContextId(), intraComm_, interComm) ;
472        parentServerContext_->createIntercomm(CXios::defaultPoolId, CXios::defaultGathererId, 0, getContextId(), intraComm_,
473                                              interCommClient, interCommServer) ;
474        int type ; 
475        if (commRank==0) CXios::getServicesManager()->getServiceType(CXios::defaultPoolId, CXios::defaultGathererId, 0, type) ;
476        MPI_Bcast(&type,1,MPI_INT,0,intraComm_) ;
477      }
478      else
479      {
480        //CXios::getContextsManager()->createServerContextIntercomm(CXios::defaultPoolId, CXios::defaultServerId, 0, getContextId(), intraComm_, interComm) ;
481        parentServerContext_->createIntercomm(CXios::defaultPoolId, CXios::defaultServerId, 0, getContextId(), intraComm_,
482                                              interCommClient, interCommServer) ;
483        int type ; 
484        if (commRank==0) CXios::getServicesManager()->getServiceType(CXios::defaultPoolId, CXios::defaultServerId, 0, type) ;
485        MPI_Bcast(&type,1,MPI_INT,0,intraComm_) ;
486      }
487
488        // intraComm client is not duplicated. In all the code we use client->intraComm for MPI
489        // in future better to replace it by intracommuncator associated to the context
490   
491      MPI_Comm intraCommClient, intraCommServer ;
492      intraCommClient=intraComm_ ;
493      MPI_Comm_dup(intraComm_, &intraCommServer) ;
494      client = new CContextClient(this, intraCommClient, interCommClient);
495      server = new CContextServer(this, intraCommServer, interCommServer);
496      client->setAssociatedServer(server) ;
497      server->setAssociatedClient(client) ;
498    }
499   
500    if (serviceType_ == CServicesManager::GATHERER)
501    {
502      int commRank ;
503      MPI_Comm_rank(intraComm_,&commRank) ;
504     
505      int nbPartitions ;
506      if (commRank==0) 
507      { 
508        CXios::getServicesManager()->getServiceNbPartitions(CXios::defaultPoolId, CXios::defaultServerId, 0, nbPartitions) ;
509        for(int i=0 ; i<nbPartitions; i++)
510          CXios::getContextsManager()->createServerContext(CXios::defaultPoolId, CXios::defaultServerId, i, getContextId()) ;
511      }     
512      MPI_Bcast(&nbPartitions, 1, MPI_INT, 0, intraComm_) ;
513     
514      MPI_Comm interComm ;
515      for(int i=0 ; i<nbPartitions; i++)
516      {
517        parentServerContext_->createIntercomm(CXios::defaultPoolId, CXios::defaultServerId, i, getContextId(), intraComm_, interCommClient, interCommServer) ;
518        int type ; 
519        if (commRank==0) CXios::getServicesManager()->getServiceType(CXios::defaultPoolId, CXios::defaultServerId, 0, type) ;
520        MPI_Bcast(&type,1,MPI_INT,0,intraComm_) ;
521        primServerId_.push_back(CXios::getContextsManager()->getServerContextName(CXios::defaultPoolId, CXios::defaultServerId, i, type, getContextId())) ;
522
523        // intraComm client is not duplicated. In all the code we use client->intraComm for MPI
524        // in future better to replace it by intracommuncator associated to the context
525     
526        MPI_Comm intraCommClient, intraCommServer ;
527
528        intraCommClient=intraComm_ ;
529        MPI_Comm_dup(intraComm_, &intraCommServer) ;
530
531        CContextClient* client = new CContextClient(this, intraCommClient, interCommClient) ;
532        CContextServer* server = new CContextServer(this, intraCommServer, interCommServer) ;
533        client->setAssociatedServer(server) ;
534        server->setAssociatedClient(client) ;
535        clientPrimServer.push_back(client);
536        serverPrimServer.push_back(server); 
537
538     
539      }
540    }
541  }
542  CATCH_DUMP_ATTR
543
544   
545
546  bool CContext::eventLoop(bool enableEventsProcessing)
547  {
548    bool finished=true; 
549
550    if (client!=nullptr && !finalized) client->checkBuffers();
551   
552    for (int i = 0; i < clientPrimServer.size(); ++i)
553    {
554      if (!finalized) clientPrimServer[i]->checkBuffers();
555      if (!finalized) finished &= serverPrimServer[i]->eventLoop(enableEventsProcessing);
556    }
557
558    for (auto it=couplerClient_.begin(); it!=couplerClient_.end(); ++it)
559    {
560      if (!finalized) it->second->checkBuffers();
561    }
562
563    for (auto it=couplerServer_.begin(); it!=couplerServer_.end(); ++it)
564    {
565      if (!finalized) it->second->eventLoop(enableEventsProcessing);
566    }
567
568    if (server!=nullptr) if (!finalized) finished &= server->eventLoop(enableEventsProcessing);
569 
570    return finalized && finished ;
571  }
572
573  void CContext::addCouplingChanel(const std::string& context, bool out)
574  {
575     vector<string> vectStr=splitRegex(context,"::") ;
576     string poolId=vectStr[0] ;
577     string serviceId=poolId ;
578     string contextId=vectStr[1] ;
579
580     int contextLeader ;
581     int type = CServicesManager::CLIENT ;
582     string contextName=CXios::getContextsManager()->getServerContextName(poolId, serviceId, 0, type, contextId) ;
583     
584     if (couplerClient_.find(contextName)==couplerClient_.end())
585     {
586       bool ok=CXios::getContextsManager()->getContextLeader(contextName, contextLeader, getIntraComm()) ;
587     
588       MPI_Comm interComm, interCommClient, interCommServer  ;
589       MPI_Comm intraCommClient, intraCommServer ;
590
591       if (ok) MPI_Intercomm_create(getIntraComm(), 0, CXios::getXiosComm(), contextLeader, 0, &interComm) ;
592
593       MPI_Comm_dup(intraComm_, &intraCommClient) ;
594       MPI_Comm_dup(intraComm_, &intraCommServer) ;
595       if (out)
596       {
597         MPI_Comm_dup(interComm, &interCommClient) ;
598         MPI_Comm_dup(interComm, &interCommServer) ;
599         CContextClient* client = new CContextClient(this, intraCommClient, interCommClient);
600         CContextServer* server = new CContextServer(this, intraCommServer, interCommServer);
601         client->setAssociatedServer(server) ;
602         server->setAssociatedClient(client) ;
603       }
604       else
605       {
606          MPI_Comm_dup(interComm, &interCommServer) ;
607          MPI_Comm_dup(interComm, &interCommClient) ;
608          CContextServer* server = new CContextServer(this, intraCommServer, interCommServer);
609          CContextClient* client = new CContextClient(this, intraCommClient, interCommClient);
610          client->setAssociatedServer(server) ;
611          server->setAssociatedClient(client) ;
612       }
613       MPI_Comm_free(&interComm) ;
614
615
616      // for now, we don't now which beffer size must be used for client coupler
617      // It will be evaluated later. Fix a constant size for now...
618      // set to 10Mb for development
619       map<int,size_t> bufferSize, maxEventSize ;
620       for(int i=0;i<client->getRemoteSize();i++)
621       {
622         bufferSize[i]=10000000 ;
623         maxEventSize[i]=10000000 ;
624       }
625
626       client->setBufferSize(bufferSize, maxEventSize);   
627       
628       couplerClient_[contextName] = client ;
629       couplerServer_[contextName] = server ;
630     }
631  }
632 
633  void CContext::globalEventLoop(void)
634  {
635    CXios::getDaemonsManager()->eventLoop() ;
636    setCurrent(getId()) ;
637  }
638
639
640   void CContext::finalize(void)
641   TRY
642   {
643      registryOut->hierarchicalGatherRegistry() ;
644      if (server->intraCommRank==0) CXios::globalRegistry->mergeRegistry(*registryOut) ;
645
646      if (serviceType_==CServicesManager::CLIENT)
647      {
648        doPreTimestepOperationsForEnabledReadModeFiles(); // For now we only use server level 1 to read data
649
650        info(100)<<"DEBUG: context "<<getId()<<" Send client finalize"<<endl ;
651        client->finalize();
652        info(100)<<"DEBUG: context "<<getId()<<" Client finalize sent"<<endl ;
653        while (client->havePendingRequests()) client->checkBuffers();
654        info(100)<<"DEBUG: context "<<getId()<<" no pending request ok"<<endl ;
655        bool notifiedFinalized=false ;
656        do
657        {
658          notifiedFinalized=client->isNotifiedFinalized() ;
659        } while (!notifiedFinalized) ;
660        client->releaseBuffers();
661        info(100)<<"DEBUG: context "<<getId()<<" release client ok"<<endl ;
662      }
663      else if (serviceType_==CServicesManager::GATHERER)
664      {
665         for (int i = 0; i < clientPrimServer.size(); ++i)
666         {
667           clientPrimServer[i]->finalize();
668           bool bufferReleased;
669           do
670           {
671             clientPrimServer[i]->checkBuffers();
672             bufferReleased = !clientPrimServer[i]->havePendingRequests();
673           } while (!bufferReleased);
674           
675           bool notifiedFinalized=false ;
676           do
677           {
678             notifiedFinalized=clientPrimServer[i]->isNotifiedFinalized() ;
679           } while (!notifiedFinalized) ;
680           clientPrimServer[i]->releaseBuffers();
681         }
682         closeAllFile();
683
684      }
685      else if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::OUT_SERVER)
686      {
687        closeAllFile();
688      }
689
690      freeComms() ;
691       
692      parentServerContext_->freeComm() ;
693      finalized = true;
694      info(20)<<"CContext: Context <"<<getId()<<"> is finalized."<<endl;
695   }
696   CATCH_DUMP_ATTR
697
698   //! Free internally allocated communicators
699   void CContext::freeComms(void)
700   TRY
701   {
702     for (std::list<MPI_Comm>::iterator it = comms.begin(); it != comms.end(); ++it)
703       MPI_Comm_free(&(*it));
704     comms.clear();
705   }
706   CATCH_DUMP_ATTR
707
708   //! Deallocate buffers allocated by clientContexts
709   void CContext::releaseClientBuffers(void)
710   TRY
711   {
712     client->releaseBuffers();
713     for (int i = 0; i < clientPrimServer.size(); ++i)
714       clientPrimServer[i]->releaseBuffers();
715   }
716   CATCH_DUMP_ATTR
717
718   void CContext::postProcessingGlobalAttributes()
719   TRY
720   {
721     if (allProcessed) return; 
722     
723    // create intercommunicator with servers.
724    // not sure it is the good place to be called here
725    createServerInterComm() ;
726
727
728     // After xml is parsed, there are some more works with post processing
729     postProcessing();
730
731     // Distribute files between secondary servers according to the data size
732     distributeFiles(this->enabledWriteModeFiles);
733
734     // Check grid and calculate its distribution
735     checkGridEnabledFields();
736
737     setClientServerBuffer(client, (serviceType_==CServicesManager::CLIENT) ) ;
738     for (int i = 0; i < clientPrimServer.size(); ++i)
739         setClientServerBuffer(clientPrimServer[i], true);
740
741   
742     if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER)
743     { 
744       if (serviceType_==CServicesManager::GATHERER)
745       { 
746         for (auto it=clientPrimServer.begin(); it!=clientPrimServer.end();++it) 
747         {
748           this->sendAllAttributesToServer(*it); // Send all attributes of current context to server
749           CCalendarWrapper::get(CCalendarWrapper::GetDefName())->sendAllAttributesToServer(*it); // Send all attributes of current calendar
750         }
751       }
752       else 
753       {
754         this->sendAllAttributesToServer(client);   // Send all attributes of current context to server
755         CCalendarWrapper::get(CCalendarWrapper::GetDefName())->sendAllAttributesToServer(client); // Send all attributes of current calendar
756       }
757
758
759      // We have enough information to send to server
760      // First of all, send all enabled files
761      sendEnabledFiles(this->enabledWriteModeFiles);
762      // We only use server-level 1 (for now) to read data
763      if (serviceType_==CServicesManager::CLIENT)  sendEnabledFiles(this->enabledReadModeFiles);
764
765      // Then, send all enabled fields     
766      sendEnabledFieldsInFiles(this->enabledWriteModeFiles);
767     
768      if (serviceType_==CServicesManager::CLIENT) sendEnabledFieldsInFiles(this->enabledReadModeFiles);
769
770      // Then, check whether we have domain_ref, axis_ref or scalar_ref attached to the enabled fields
771      // If any, so send them to server
772       sendRefDomainsAxisScalars(this->enabledWriteModeFiles); 
773     
774      if (serviceType_==CServicesManager::CLIENT) sendRefDomainsAxisScalars(this->enabledReadModeFiles);       
775
776       // Check whether enabled fields have grid_ref, if any, send this info to server
777      sendRefGrid(this->enabledFiles);
778      // This code may be useful in the future when we want to seperate completely read and write
779      // sendRefGrid(this->enabledWriteModeFiles);
780      // if (!hasServer)
781      //   sendRefGrid(this->enabledReadModeFiles);
782     
783      // A grid of enabled fields composed of several components which must be checked then their
784      // checked attributes should be sent to server
785      sendGridComponentEnabledFieldsInFiles(this->enabledFiles); // This code can be seperated in two (one for reading, another for writing)
786
787       // We have a xml tree on the server side and now, it should be also processed
788      sendPostProcessing();
789       
790      // Finally, we send information of grid itself to server
791      sendGridEnabledFieldsInFiles(this->enabledWriteModeFiles);       
792     
793      if (serviceType_==CServicesManager::CLIENT) sendGridEnabledFieldsInFiles(this->enabledReadModeFiles);       
794     
795     }
796     allProcessed = true;
797   }
798   CATCH_DUMP_ATTR
799
800   void CContext::sendPostProcessingGlobalAttributes()
801   TRY
802   {
803
804    int nbSrvPools ;
805    if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
806    else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
807    else nbSrvPools = 0 ;
808    CContextClient* contextClientTmp ;
809
810    for (int i = 0; i < nbSrvPools; ++i)
811     {
812       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
813       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
814     
815       CEventClient event(getType(),EVENT_ID_POST_PROCESS_GLOBAL_ATTRIBUTES);
816
817       if (contextClientTmp->isServerLeader())
818       {
819         CMessage msg;
820         const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
821         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
822           event.push(*itRank,1,msg);
823         contextClientTmp->sendEvent(event);
824       }
825       else contextClientTmp->sendEvent(event);
826     }
827   }
828   CATCH_DUMP_ATTR
829
830   void CContext::recvPostProcessingGlobalAttributes(CEventServer& event)
831   TRY
832   {
833      CBufferIn* buffer=event.subEvents.begin()->buffer;
834      getCurrent()->recvPostProcessingGlobalAttributes(*buffer);
835   }
836   CATCH
837
838   void CContext::recvPostProcessingGlobalAttributes(CBufferIn& buffer)
839   TRY
840   {     
841      postProcessingGlobalAttributes();
842   }
843   CATCH_DUMP_ATTR
844
845   /*!
846   \brief Close all the context defintion and do processing data
847      After everything is well defined on client side, they will be processed and sent to server
848   From the version 2.0, sever and client work no more on the same database. Moreover, client(s) will send
849   all necessary information to server, from which each server can build its own database.
850   Because the role of server is to write out field data on a specific netcdf file,
851   the only information that it needs is the enabled files
852   and the active fields (fields will be written onto active files)
853   */
854  void CContext::closeDefinition(void)
855   TRY
856   {
857     CTimer::get("Context : close definition").resume() ;
858     
859     // create intercommunicator with servers.
860     // not sure it is the good place to be called here
861     createServerInterComm() ;
862
863
864     // After xml is parsed, there are some more works with post processing
865//     postProcessing();
866
867    // Make sure the calendar was correctly created
868    if (!calendar)
869      ERROR("CContext::postProcessing()", << "A calendar must be defined for the context \"" << getId() << "!\"")
870    else if (calendar->getTimeStep() == NoneDu)
871      ERROR("CContext::postProcessing()", << "A timestep must be defined for the context \"" << getId() << "!\"")
872    // Calendar first update to set the current date equals to the start date
873    calendar->update(0);
874
875    // Résolution des héritages descendants (càd des héritages de groupes)
876    // pour chacun des contextes.
877    solveDescInheritance(true);
878 
879    // Solve inheritance for field to know if enabled or not.
880    for (auto field : CField::getAll()) field->solveRefInheritance();
881
882    // Check if some axis, domains or grids are eligible to for compressed indexed output.
883    // Warning: This must be done after solving the inheritance and before the rest of post-processing
884    // --> later ????    checkAxisDomainsGridsEligibilityForCompressedOutput();     
885
886      // Check if some automatic time series should be generated
887      // Warning: This must be done after solving the inheritance and before the rest of post-processing     
888
889    // The timeseries should only be prepared in client
890    prepareTimeseries();
891
892    //Initialisation du vecteur 'enabledFiles' contenant la liste des fichiers à sortir.
893    findEnabledFiles();
894    findEnabledWriteModeFiles();
895    findEnabledReadModeFiles();
896    findEnabledCouplerIn();
897    findEnabledCouplerOut();
898    createCouplerInterCommunicator() ;
899
900    // Find all enabled fields of each file     
901    const vector<CField*>&& fileOutField = findAllEnabledFieldsInFileOut(this->enabledWriteModeFiles);
902    const vector<CField*>&& fileInField = findAllEnabledFieldsInFileIn(this->enabledReadModeFiles);
903    const vector<CField*>&& CouplerOutField = findAllEnabledFieldsCouplerOut(this->enabledCouplerOut);
904    const vector<CField*>&& CouplerInField = findAllEnabledFieldsCouplerIn(this->enabledCouplerIn);
905    findFieldsWithReadAccess();
906    const vector<CField*>& fieldWithReadAccess = fieldsWithReadAccess_ ;
907     
908// find all field potentially at workflow end
909    vector<CField*> endWorkflowFields ;
910    endWorkflowFields.reserve(fileOutField.size()+CouplerOutField.size()+fieldWithReadAccess.size()) ;
911    endWorkflowFields.insert(endWorkflowFields.end(),fileOutField.begin(), fileOutField.end()) ;
912    endWorkflowFields.insert(endWorkflowFields.end(),CouplerOutField.begin(), CouplerOutField.end()) ;
913    endWorkflowFields.insert(endWorkflowFields.end(),fieldWithReadAccess.begin(), fieldWithReadAccess.end()) ;
914
915    for(auto endWorkflowField : endWorkflowFields) endWorkflowField->buildWorkflowGraph(garbageCollector) ;
916   
917    // Distribute files between secondary servers according to the data size => assign a context to a file
918    if (serviceType_==CServicesManager::GATHERER) distributeFiles(this->enabledWriteModeFiles);
919    else if (serviceType_==CServicesManager::CLIENT) for(auto file : this->enabledWriteModeFiles) file->setContextClient(client) ;
920
921    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER)
922    {
923      for(auto field : fileOutField) 
924      {
925        field->connectToFileServer(garbageCollector) ; // connect the field to server filter
926        field->computeGridIndexToFileServer() ; // compute grid index for transfer to the server context
927      }
928    }
929   
930
931
932
933
934    // For now, only read files with client and only one level server
935    // if (hasClient && !hasServer) findEnabledReadModeFiles();     
936
937    // Find all enabled fields of each file     
938    findAllEnabledFieldsInFiles(this->enabledWriteModeFiles);
939    findAllEnabledFieldsInFiles(this->enabledReadModeFiles);
940
941    // For now, only read files with client and only one level server
942    // if (hasClient && !hasServer)
943    //   findAllEnabledFieldsInFiles(this->enabledReadModeFiles);     
944
945    if (serviceType_==CServicesManager::CLIENT)
946    {
947      initReadFiles();
948      // Try to read attributes of fields in file then fill in corresponding grid (or domain, axis)
949      this->readAttributesOfEnabledFieldsInReadModeFiles();
950    }
951
952    // Only search and rebuild all reference objects of enable fields, don't transform
953    this->solveOnlyRefOfEnabledFields();
954
955    // Search and rebuild all reference object of enabled fields, and transform
956    this->solveAllRefOfEnabledFieldsAndTransform();
957
958    // Find all fields with read access from the public API
959    if (serviceType_==CServicesManager::CLIENT) findFieldsWithReadAccess();
960    // and solve the all reference for them
961    if (serviceType_==CServicesManager::CLIENT) solveAllRefOfFieldsWithReadAccess();
962
963    isPostProcessed = true;
964
965
966
967    // Distribute files between secondary servers according to the data size
968    distributeFiles(this->enabledWriteModeFiles);
969
970    // Check grid and calculate its distribution
971    checkGridEnabledFields();
972
973    setClientServerBuffer(client, (serviceType_==CServicesManager::CLIENT) ) ;
974    for (int i = 0; i < clientPrimServer.size(); ++i)
975         setClientServerBuffer(clientPrimServer[i], true);
976
977   
978    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER)
979    { 
980      if (serviceType_==CServicesManager::GATHERER)
981      { 
982        for (auto it=clientPrimServer.begin(); it!=clientPrimServer.end();++it) 
983        {
984          this->sendAllAttributesToServer(*it); // Send all attributes of current context to server
985          CCalendarWrapper::get(CCalendarWrapper::GetDefName())->sendAllAttributesToServer(*it); // Send all attributes of current calendar
986        }
987      }
988      else 
989      {
990        this->sendAllAttributesToServer(client);   // Send all attributes of current context to server
991        CCalendarWrapper::get(CCalendarWrapper::GetDefName())->sendAllAttributesToServer(client); // Send all attributes of current calendar
992      }
993
994      // We have enough information to send to server
995      // First of all, send all enabled files
996      sendEnabledFiles(this->enabledWriteModeFiles);
997      // We only use server-level 1 (for now) to read data
998      if (serviceType_==CServicesManager::CLIENT)  sendEnabledFiles(this->enabledReadModeFiles);
999
1000      // Then, send all enabled fields     
1001      sendEnabledFieldsInFiles(this->enabledWriteModeFiles);
1002     
1003      if (serviceType_==CServicesManager::CLIENT) sendEnabledFieldsInFiles(this->enabledReadModeFiles);
1004
1005      // Then, check whether we have domain_ref, axis_ref or scalar_ref attached to the enabled fields
1006      // If any, so send them to server
1007      sendRefDomainsAxisScalars(this->enabledWriteModeFiles); 
1008     
1009      if (serviceType_==CServicesManager::CLIENT) sendRefDomainsAxisScalars(this->enabledReadModeFiles);       
1010
1011      // Check whether enabled fields have grid_ref, if any, send this info to server
1012      sendRefGrid(this->enabledFiles);
1013      // This code may be useful in the future when we want to seperate completely read and write
1014      // sendRefGrid(this->enabledWriteModeFiles);
1015      // if (!hasServer)
1016      //   sendRefGrid(this->enabledReadModeFiles);
1017     
1018      // A grid of enabled fields composed of several components which must be checked then their
1019      // checked attributes should be sent to server
1020      sendGridComponentEnabledFieldsInFiles(this->enabledFiles); // This code can be seperated in two (one for reading, another for writing)
1021
1022      // We have a xml tree on the server side and now, it should be also processed
1023      sendPostProcessing();
1024       
1025      // Finally, we send information of grid itself to server
1026      sendGridEnabledFieldsInFiles(this->enabledWriteModeFiles);       
1027     
1028      if (serviceType_==CServicesManager::CLIENT) sendGridEnabledFieldsInFiles(this->enabledReadModeFiles);       
1029    }
1030    allProcessed = true;
1031
1032
1033    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) sendPostProcessingGlobalAttributes();
1034
1035    // There are some processings that should be done after all of above. For example: check mask or index
1036    this->buildFilterGraphOfEnabledFields();
1037   
1038     if (serviceType_==CServicesManager::CLIENT)
1039    {
1040      buildFilterGraphOfFieldsWithReadAccess();
1041      postProcessFilterGraph(); // For coupling in, modify this later
1042    }
1043   
1044    checkGridEnabledFields();   
1045
1046    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) this->sendProcessingGridOfEnabledFields();
1047    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) this->sendCloseDefinition();
1048
1049    // Nettoyage de l'arborescence
1050    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) CleanTree(); // Only on client side??
1051
1052    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) sendCreateFileHeader();
1053    if (serviceType_==CServicesManager::CLIENT) startPrefetchingOfEnabledReadModeFiles();
1054   
1055    CTimer::get("Context : close definition").suspend() ;
1056  }
1057  CATCH_DUMP_ATTR
1058
1059
1060   void CContext::closeDefinition_old(void)
1061   TRY
1062   {
1063     CTimer::get("Context : close definition").resume() ;
1064   
1065    //
1066    postProcessingGlobalAttributes();
1067
1068    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) sendPostProcessingGlobalAttributes();
1069
1070    // There are some processings that should be done after all of above. For example: check mask or index
1071    this->buildFilterGraphOfEnabledFields();
1072   
1073     if (serviceType_==CServicesManager::CLIENT)
1074    {
1075      buildFilterGraphOfFieldsWithReadAccess();
1076      postProcessFilterGraph(); // For coupling in, modify this later
1077    }
1078   
1079    checkGridEnabledFields();   
1080
1081    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) this->sendProcessingGridOfEnabledFields();
1082    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) this->sendCloseDefinition();
1083
1084    // Nettoyage de l'arborescence
1085    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) CleanTree(); // Only on client side??
1086
1087    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) sendCreateFileHeader();
1088    if (serviceType_==CServicesManager::CLIENT) startPrefetchingOfEnabledReadModeFiles();
1089   
1090    CTimer::get("Context : close definition").suspend() ;
1091   }
1092   CATCH_DUMP_ATTR
1093
1094   vector<CField*> CContext::findAllEnabledFieldsInFiles(const std::vector<CFile*>& activeFiles)
1095   TRY
1096   {
1097     vector<CField*> fields ;
1098     for (unsigned int i = 0; i < activeFiles.size(); i++)
1099     {
1100        const vector<CField*>&& field=activeFiles[i]->getEnabledFields() ;
1101        fields.insert(fields.end(),field.begin(),field.end());
1102     }
1103     return fields ;
1104   }
1105   CATCH_DUMP_ATTR
1106
1107   vector<CField*> CContext::findAllEnabledFieldsInFileOut(const std::vector<CFile*>& activeFiles)
1108   TRY
1109   {
1110     vector<CField*> fields ;
1111     for(auto file : activeFiles)
1112     {
1113        const vector<CField*>&& fieldList=file->getEnabledFields() ;
1114        for(auto field : fieldList) field->setFileOut(file) ;
1115        fields.insert(fields.end(),fieldList.begin(),fieldList.end());
1116     }
1117     return fields ;
1118   }
1119   CATCH_DUMP_ATTR
1120
1121   vector<CField*> CContext::findAllEnabledFieldsInFileIn(const std::vector<CFile*>& activeFiles)
1122   TRY
1123   {
1124     vector<CField*> fields ;
1125     for(auto file : activeFiles)
1126     {
1127        const vector<CField*>&& fieldList=file->getEnabledFields() ;
1128        for(auto field : fieldList) field->setFileIn(file) ;
1129        fields.insert(fields.end(),fieldList.begin(),fieldList.end());
1130     }
1131     return fields ;
1132   }
1133   CATCH_DUMP_ATTR
1134
1135   vector<CField*> CContext::findAllEnabledFieldsCouplerOut(const std::vector<CCouplerOut*>& activeCouplerOut)
1136   TRY
1137   {
1138     vector<CField*> fields ;
1139     for (auto couplerOut :activeCouplerOut)
1140     {
1141        const vector<CField*>&& fieldList=couplerOut->getEnabledFields() ;
1142        for(auto field : fieldList) field->setCouplerOut(couplerOut) ;
1143        fields.insert(fields.end(),fieldList.begin(),fieldList.end());
1144     }
1145     return fields ;
1146   }
1147   CATCH_DUMP_ATTR
1148
1149   vector<CField*> CContext::findAllEnabledFieldsCouplerIn(const std::vector<CCouplerIn*>& activeCouplerIn)
1150   TRY
1151   {
1152     vector<CField*> fields ;
1153     for (auto couplerIn :activeCouplerIn)
1154     {
1155        const vector<CField*>&& fieldList=couplerIn->getEnabledFields() ;
1156        for(auto field : fieldList) field->setCouplerIn(couplerIn) ;
1157        fields.insert(fields.end(),fieldList.begin(),fieldList.end());
1158     }
1159     return fields ;
1160   }
1161   CATCH_DUMP_ATTR
1162
1163
1164
1165   void CContext::readAttributesOfEnabledFieldsInReadModeFiles()
1166   TRY
1167   {
1168      for (unsigned int i = 0; i < this->enabledReadModeFiles.size(); ++i)
1169        (void)this->enabledReadModeFiles[i]->readAttributesOfEnabledFieldsInReadMode();
1170   }
1171   CATCH_DUMP_ATTR
1172
1173   void CContext::sendGridComponentEnabledFieldsInFiles(const std::vector<CFile*>& activeFiles)
1174   TRY
1175   {
1176     int size = activeFiles.size();
1177     for (int i = 0; i < size; ++i)
1178     {       
1179       activeFiles[i]->sendGridComponentOfEnabledFields();
1180     }
1181   }
1182   CATCH_DUMP_ATTR
1183
1184   /*!
1185      Send active (enabled) fields in file from a client to others
1186      \param [in] activeFiles files contains enabled fields to send
1187   */
1188   void CContext::sendGridEnabledFieldsInFiles(const std::vector<CFile*>& activeFiles)
1189   TRY
1190   {
1191     int size = activeFiles.size();
1192     for (int i = 0; i < size; ++i)
1193     {       
1194       activeFiles[i]->sendGridOfEnabledFields();
1195     }
1196   }
1197   CATCH_DUMP_ATTR
1198
1199   void CContext::checkGridEnabledFields()
1200   TRY
1201   {
1202     int size = enabledFiles.size();
1203     for (int i = 0; i < size; ++i)
1204     {
1205       enabledFiles[i]->checkGridOfEnabledFields();       
1206     }
1207
1208     size = enabledCouplerOut.size();
1209     for (int i = 0; i < size; ++i)
1210     {
1211       enabledCouplerOut[i]->checkGridOfEnabledFields();       
1212     }
1213   }
1214   CATCH_DUMP_ATTR
1215
1216   /*!
1217      Check grid of active (enabled) fields in file
1218      \param [in] activeFiles files contains enabled fields whose grid needs checking
1219   */
1220   void CContext::checkGridEnabledFieldsInFiles(const std::vector<CFile*>& activeFiles)
1221   TRY
1222   {
1223     int size = activeFiles.size();
1224     for (int i = 0; i < size; ++i)
1225     {
1226       activeFiles[i]->checkGridOfEnabledFields();       
1227     }
1228   }
1229   CATCH_DUMP_ATTR
1230
1231    /*!
1232      Go up the hierachical tree via field_ref and do check of attributes of fields
1233      This can be done in a client then all computed information will be sent from this client to others
1234      \param [in] sendToServer Flag to indicate whether calculated information will be sent
1235   */
1236   void CContext::solveOnlyRefOfEnabledFields(void)
1237   TRY
1238   {
1239     int size = this->enabledFiles.size();
1240     for (int i = 0; i < size; ++i)
1241     {
1242       this->enabledFiles[i]->solveOnlyRefOfEnabledFields();
1243     }
1244
1245     for (int i = 0; i < size; ++i)
1246     {
1247       this->enabledFiles[i]->generateNewTransformationGridDest();
1248     }
1249
1250     size = this->enabledCouplerOut.size();
1251     for (int i = 0; i < size; ++i)
1252     {
1253       this->enabledCouplerOut[i]->solveOnlyRefOfEnabledFields();
1254     }
1255
1256     for (int i = 0; i < size; ++i)
1257     {
1258       this->enabledCouplerOut[i]->generateNewTransformationGridDest();
1259     }
1260   }
1261   CATCH_DUMP_ATTR
1262
1263    /*!
1264      Go up the hierachical tree via field_ref and do check of attributes of fields.
1265      The transformation can be done in this step.
1266      All computed information will be sent from this client to others.
1267      \param [in] sendToServer Flag to indicate whether calculated information will be sent
1268   */
1269   void CContext::solveAllRefOfEnabledFieldsAndTransform(void)
1270   TRY
1271   {
1272     int size = this->enabledFiles.size();
1273     for (int i = 0; i < size; ++i)
1274     {
1275       this->enabledFiles[i]->solveAllRefOfEnabledFieldsAndTransform();
1276     }
1277
1278     size = this->enabledCouplerOut.size();
1279     for (int i = 0; i < size; ++i)
1280     {
1281       this->enabledCouplerOut[i]->solveAllRefOfEnabledFieldsAndTransform();
1282     }
1283
1284   }
1285   CATCH_DUMP_ATTR
1286
1287   void CContext::buildFilterGraphOfEnabledFields()
1288   TRY
1289   {
1290     int size = this->enabledFiles.size();
1291     for (int i = 0; i < size; ++i)
1292     {
1293       this->enabledFiles[i]->buildFilterGraphOfEnabledFields(garbageCollector);
1294     }
1295
1296     size = this->enabledCouplerOut.size();
1297     for (int i = 0; i < size; ++i)
1298     {
1299       this->enabledCouplerOut[i]->buildFilterGraphOfEnabledFields(garbageCollector);
1300     }
1301   }
1302   CATCH_DUMP_ATTR
1303
1304   void CContext::postProcessFilterGraph()
1305   TRY
1306   {
1307     int size = enabledFiles.size();
1308     for (int i = 0; i < size; ++i)
1309     {
1310        enabledFiles[i]->postProcessFilterGraph();
1311     }
1312   }
1313   CATCH_DUMP_ATTR
1314
1315   void CContext::startPrefetchingOfEnabledReadModeFiles()
1316   TRY
1317   {
1318     int size = enabledReadModeFiles.size();
1319     for (int i = 0; i < size; ++i)
1320     {
1321        enabledReadModeFiles[i]->prefetchEnabledReadModeFields();
1322     }
1323   }
1324   CATCH_DUMP_ATTR
1325
1326   void CContext::doPreTimestepOperationsForEnabledReadModeFiles()
1327   TRY
1328   {
1329     int size = enabledReadModeFiles.size();
1330     for (int i = 0; i < size; ++i)
1331     {
1332        enabledReadModeFiles[i]->doPreTimestepOperationsForEnabledReadModeFields();
1333     }
1334   }
1335   CATCH_DUMP_ATTR
1336
1337   void CContext::doPostTimestepOperationsForEnabledReadModeFiles()
1338   TRY
1339   {
1340     int size = enabledReadModeFiles.size();
1341     for (int i = 0; i < size; ++i)
1342     {
1343        enabledReadModeFiles[i]->doPostTimestepOperationsForEnabledReadModeFields();
1344     }
1345   }
1346   CATCH_DUMP_ATTR
1347
1348  void CContext::findFieldsWithReadAccess(void)
1349  TRY
1350  {
1351    fieldsWithReadAccess_.clear();
1352    const vector<CField*> allFields = CField::getAll();
1353    for (size_t i = 0; i < allFields.size(); ++i)
1354    {
1355      CField* field = allFields[i];
1356      if (!field->read_access.isEmpty() && field->read_access && (field->enabled.isEmpty() || field->enabled))
1357      {
1358        fieldsWithReadAccess_.push_back(field);
1359        field->setModelOut() ;
1360      }
1361    }
1362  }
1363  CATCH_DUMP_ATTR
1364
1365  void CContext::solveAllRefOfFieldsWithReadAccess()
1366  TRY
1367  {
1368    for (size_t i = 0; i < fieldsWithReadAccess_.size(); ++i)
1369      fieldsWithReadAccess_[i]->solveAllReferenceEnabledField(false);
1370  }
1371  CATCH_DUMP_ATTR
1372
1373  void CContext::buildFilterGraphOfFieldsWithReadAccess()
1374  TRY
1375  {
1376    for (size_t i = 0; i < fieldsWithReadAccess_.size(); ++i)
1377      fieldsWithReadAccess_[i]->buildFilterGraph(garbageCollector, true);
1378  }
1379  CATCH_DUMP_ATTR
1380
1381   void CContext::solveAllInheritance(bool apply)
1382   TRY
1383   {
1384     // Résolution des héritages descendants (càd des héritages de groupes)
1385     // pour chacun des contextes.
1386      solveDescInheritance(apply);
1387
1388     // Résolution des héritages par référence au niveau des fichiers.
1389      const vector<CFile*> allFiles=CFile::getAll();
1390      const vector<CCouplerIn*> allCouplerIn=CCouplerIn::getAll();
1391      const vector<CCouplerOut*> allCouplerOut=CCouplerOut::getAll();
1392      const vector<CGrid*> allGrids= CGrid::getAll();
1393
1394      if (serviceType_==CServicesManager::CLIENT)
1395      {
1396        for (unsigned int i = 0; i < allFiles.size(); i++)
1397          allFiles[i]->solveFieldRefInheritance(apply);
1398
1399        for (unsigned int i = 0; i < allCouplerIn.size(); i++)
1400          allCouplerIn[i]->solveFieldRefInheritance(apply);
1401
1402        for (unsigned int i = 0; i < allCouplerOut.size(); i++)
1403          allCouplerOut[i]->solveFieldRefInheritance(apply);
1404      }
1405
1406      unsigned int vecSize = allGrids.size();
1407      unsigned int i = 0;
1408      for (i = 0; i < vecSize; ++i)
1409        allGrids[i]->solveElementsRefInheritance(apply);
1410
1411   }
1412  CATCH_DUMP_ATTR
1413
1414   void CContext::findEnabledFiles(void)
1415   TRY
1416   {
1417      const std::vector<CFile*> allFiles = CFile::getAll();
1418      const CDate& initDate = calendar->getInitDate();
1419
1420      for (unsigned int i = 0; i < allFiles.size(); i++)
1421         if (!allFiles[i]->enabled.isEmpty()) // Si l'attribut 'enabled' est défini.
1422         {
1423            if (allFiles[i]->enabled.getValue()) // Si l'attribut 'enabled' est fixé à vrai.
1424            {
1425              if (allFiles[i]->output_freq.isEmpty())
1426              {
1427                 ERROR("CContext::findEnabledFiles()",
1428                     << "Mandatory attribute output_freq must be defined for file \""<<allFiles[i]->getFileOutputName()
1429                     <<" \".")
1430              }
1431              if ((initDate + allFiles[i]->output_freq.getValue()) < (initDate + this->getCalendar()->getTimeStep()))
1432              {
1433                error(0)<<"WARNING: void CContext::findEnabledFiles()"<<endl
1434                    << "Output frequency in file \""<<allFiles[i]->getFileOutputName()
1435                    <<"\" is less than the time step. File will not be written."<<endl;
1436              }
1437              else
1438               enabledFiles.push_back(allFiles[i]);
1439            }
1440         }
1441         else
1442         {
1443           if (allFiles[i]->output_freq.isEmpty())
1444           {
1445              ERROR("CContext::findEnabledFiles()",
1446                  << "Mandatory attribute output_freq must be defined for file \""<<allFiles[i]->getFileOutputName()
1447                  <<" \".")
1448           }
1449           if ( (initDate + allFiles[i]->output_freq.getValue()) < (initDate + this->getCalendar()->getTimeStep()))
1450           {
1451             error(0)<<"WARNING: void CContext::findEnabledFiles()"<<endl
1452                 << "Output frequency in file \""<<allFiles[i]->getFileOutputName()
1453                 <<"\" is less than the time step. File will not be written."<<endl;
1454           }
1455           else
1456             enabledFiles.push_back(allFiles[i]); // otherwise true by default
1457         }
1458
1459      if (enabledFiles.size() == 0)
1460         DEBUG(<<"Aucun fichier ne va être sorti dans le contexte nommé \""
1461               << getId() << "\" !");
1462
1463   }
1464   CATCH_DUMP_ATTR
1465
1466   void CContext::findEnabledCouplerIn(void)
1467   TRY
1468   {
1469      const std::vector<CCouplerIn*> allCouplerIn = CCouplerIn::getAll();
1470      bool enabled ;
1471      for (size_t i = 0; i < allCouplerIn.size(); i++)
1472      {
1473        if (allCouplerIn[i]->enabled.isEmpty()) enabled=true ;
1474        else enabled=allCouplerIn[i]->enabled ;
1475        if (enabled) enabledCouplerIn.push_back(allCouplerIn[i]) ;
1476      }
1477   }
1478   CATCH_DUMP_ATTR
1479
1480   void CContext::findEnabledCouplerOut(void)
1481   TRY
1482   {
1483      const std::vector<CCouplerOut*> allCouplerOut = CCouplerOut::getAll();
1484      bool enabled ;
1485      for (size_t i = 0; i < allCouplerOut.size(); i++)
1486      {
1487        if (allCouplerOut[i]->enabled.isEmpty()) enabled=true ;
1488        else enabled=allCouplerOut[i]->enabled ;
1489        if (enabled) enabledCouplerOut.push_back(allCouplerOut[i]) ;
1490      }
1491   }
1492   CATCH_DUMP_ATTR
1493
1494
1495
1496
1497   void CContext::distributeFiles(const vector<CFile*>& files)
1498   TRY
1499   {
1500     bool distFileMemory=false ;
1501     distFileMemory=CXios::getin<bool>("server2_dist_file_memory", distFileMemory);
1502
1503     if (distFileMemory) distributeFileOverMemoryBandwith(files) ;
1504     else distributeFileOverBandwith(files) ;
1505   }
1506   CATCH_DUMP_ATTR
1507
1508   void CContext::distributeFileOverBandwith(const vector<CFile*>& files)
1509   TRY
1510   {
1511     double eps=std::numeric_limits<double>::epsilon()*10 ;
1512     
1513     std::ofstream ofs(("distribute_file_"+getId()+".dat").c_str(), std::ofstream::out);
1514     int nbPools = clientPrimServer.size();
1515
1516     // (1) Find all enabled files in write mode
1517     // for (int i = 0; i < this->enabledFiles.size(); ++i)
1518     // {
1519     //   if (enabledFiles[i]->mode.isEmpty() || (!enabledFiles[i]->mode.isEmpty() && enabledFiles[i]->mode.getValue() == CFile::mode_attr::write ))
1520     //    enabledWriteModeFiles.push_back(enabledFiles[i]);
1521     // }
1522
1523     // (2) Estimate the data volume for each file
1524     int size = files.size();
1525     std::vector<std::pair<double, CFile*> > dataSizeMap;
1526     double dataPerPool = 0;
1527     int nfield=0 ;
1528     ofs<<size<<endl ;
1529     for (size_t i = 0; i < size; ++i)
1530     {
1531       CFile* file = files[i];
1532       ofs<<file->getId()<<endl ;
1533       StdSize dataSize=0;
1534       std::vector<CField*> enabledFields = file->getEnabledFields();
1535       size_t numEnabledFields = enabledFields.size();
1536       ofs<<numEnabledFields<<endl ;
1537       for (size_t j = 0; j < numEnabledFields; ++j)
1538       {
1539         dataSize += enabledFields[j]->getGlobalWrittenSize() ;
1540         ofs<<enabledFields[j]->getGrid()->getId()<<endl ;
1541         ofs<<enabledFields[j]->getGlobalWrittenSize()<<endl ;
1542       }
1543       double outFreqSec = (Time)(calendar->getCurrentDate()+file->output_freq)-(Time)(calendar->getCurrentDate()) ;
1544       double dataSizeSec= dataSize/ outFreqSec;
1545       ofs<<dataSizeSec<<endl ;
1546       nfield++ ;
1547// add epsilon*nField to dataSizeSec in order to  preserve reproductive ordering when sorting
1548       dataSizeMap.push_back(make_pair(dataSizeSec + dataSizeSec * eps * nfield , file));
1549       dataPerPool += dataSizeSec;
1550     }
1551     dataPerPool /= nbPools;
1552     std::sort(dataSizeMap.begin(), dataSizeMap.end());
1553
1554     // (3) Assign contextClient to each enabled file
1555
1556     std::multimap<double,int> poolDataSize ;
1557// multimap is not garanty to preserve stable sorting in c++98 but it seems it does for c++11
1558
1559     int j;
1560     double dataSize ;
1561     for (j = 0 ; j < nbPools ; ++j) poolDataSize.insert(std::pair<double,int>(0.,j)) ; 
1562             
1563     for (int i = dataSizeMap.size()-1; i >= 0; --i)
1564     {
1565       dataSize=(*poolDataSize.begin()).first ;
1566       j=(*poolDataSize.begin()).second ;
1567       dataSizeMap[i].second->setContextClient(clientPrimServer[j]);
1568       dataSize+=dataSizeMap[i].first;
1569       poolDataSize.erase(poolDataSize.begin()) ;
1570       poolDataSize.insert(std::pair<double,int>(dataSize,j)) ; 
1571     }
1572
1573     for (std::multimap<double,int>:: iterator it=poolDataSize.begin() ; it!=poolDataSize.end(); ++it) info(30)<<"Load Balancing for servers (perfect=1) : "<<it->second<<" :  ratio "<<it->first*1./dataPerPool<<endl ;
1574   }
1575   CATCH_DUMP_ATTR
1576
1577   void CContext::distributeFileOverMemoryBandwith(const vector<CFile*>& filesList)
1578   TRY
1579   {
1580     int nbPools = clientPrimServer.size();
1581     double ratio=0.5 ;
1582     ratio=CXios::getin<double>("server2_dist_file_memory_ratio", ratio);
1583
1584     int nFiles = filesList.size();
1585     vector<SDistFile> files(nFiles);
1586     vector<SDistGrid> grids;
1587     map<string,int> gridMap ;
1588     string gridId; 
1589     int gridIndex=0 ;
1590
1591     for (size_t i = 0; i < nFiles; ++i)
1592     {
1593       StdSize dataSize=0;
1594       CFile* file = filesList[i];
1595       std::vector<CField*> enabledFields = file->getEnabledFields();
1596       size_t numEnabledFields = enabledFields.size();
1597
1598       files[i].id_=file->getId() ;
1599       files[i].nbGrids_=numEnabledFields;
1600       files[i].assignedGrid_ = new int[files[i].nbGrids_] ;
1601         
1602       for (size_t j = 0; j < numEnabledFields; ++j)
1603       {
1604         gridId=enabledFields[j]->getGrid()->getId() ;
1605         if (gridMap.find(gridId)==gridMap.end())
1606         {
1607            gridMap[gridId]=gridIndex  ;
1608            SDistGrid newGrid; 
1609            grids.push_back(newGrid) ;
1610            gridIndex++ ;
1611         }
1612         files[i].assignedGrid_[j]=gridMap[gridId] ;
1613         grids[files[i].assignedGrid_[j]].size_=enabledFields[j]->getGlobalWrittenSize() ;
1614         dataSize += enabledFields[j]->getGlobalWrittenSize() ; // usefull
1615       }
1616       double outFreqSec = (Time)(calendar->getCurrentDate()+file->output_freq)-(Time)(calendar->getCurrentDate()) ;
1617       files[i].bandwith_= dataSize/ outFreqSec ;
1618     }
1619
1620     double bandwith=0 ;
1621     double memory=0 ;
1622   
1623     for(int i=0; i<nFiles; i++)  bandwith+=files[i].bandwith_ ;
1624     for(int i=0; i<nFiles; i++)  files[i].bandwith_ = files[i].bandwith_/bandwith * ratio ;
1625
1626     for(int i=0; i<grids.size(); i++)  memory+=grids[i].size_ ;
1627     for(int i=0; i<grids.size(); i++)  grids[i].size_ = grids[i].size_ / memory * (1.0-ratio) ;
1628       
1629     distributeFileOverServer2(nbPools, grids.size(), &grids[0], nFiles, &files[0]) ;
1630
1631     vector<double> memorySize(nbPools,0.) ;
1632     vector< set<int> > serverGrids(nbPools) ;
1633     vector<double> bandwithSize(nbPools,0.) ;
1634       
1635     for (size_t i = 0; i < nFiles; ++i)
1636     {
1637       bandwithSize[files[i].assignedServer_] += files[i].bandwith_* bandwith /ratio ;
1638       for(int j=0 ; j<files[i].nbGrids_;j++)
1639       {
1640         if (serverGrids[files[i].assignedServer_].find(files[i].assignedGrid_[j]) == serverGrids[files[i].assignedServer_].end())
1641         {
1642           memorySize[files[i].assignedServer_]+= grids[files[i].assignedGrid_[j]].size_ * memory / (1.0-ratio);
1643           serverGrids[files[i].assignedServer_].insert(files[i].assignedGrid_[j]) ;
1644         }
1645       }
1646       filesList[i]->setContextClient(clientPrimServer[files[i].assignedServer_]) ;
1647       delete [] files[i].assignedGrid_ ;
1648     }
1649
1650     for (int i = 0; i < nbPools; ++i) info(100)<<"Pool server level2 "<<i<<"   assigned file bandwith "<<bandwithSize[i]*86400.*4./1024/1024.<<" Mb / days"<<endl ;
1651     for (int i = 0; i < nbPools; ++i) info(100)<<"Pool server level2 "<<i<<"   assigned grid memory "<<memorySize[i]*100/1024./1024.<<" Mb"<<endl ;
1652
1653   }
1654   CATCH_DUMP_ATTR
1655
1656   /*!
1657      Find all files in write mode
1658   */
1659   void CContext::findEnabledWriteModeFiles(void)
1660   TRY
1661   {
1662     int size = this->enabledFiles.size();
1663     for (int i = 0; i < size; ++i)
1664     {
1665       if (enabledFiles[i]->mode.isEmpty() || 
1666          (!enabledFiles[i]->mode.isEmpty() && enabledFiles[i]->mode.getValue() == CFile::mode_attr::write ))
1667        enabledWriteModeFiles.push_back(enabledFiles[i]);
1668     }
1669   }
1670   CATCH_DUMP_ATTR
1671
1672   /*!
1673      Find all files in read mode
1674   */
1675   void CContext::findEnabledReadModeFiles(void)
1676   TRY
1677   {
1678     int size = this->enabledFiles.size();
1679     for (int i = 0; i < size; ++i)
1680     {
1681       if (!enabledFiles[i]->mode.isEmpty() && enabledFiles[i]->mode.getValue() == CFile::mode_attr::read)
1682        enabledReadModeFiles.push_back(enabledFiles[i]);
1683     }
1684   }
1685   CATCH_DUMP_ATTR
1686
1687   void CContext::closeAllFile(void)
1688   TRY
1689   {
1690     std::vector<CFile*>::const_iterator
1691            it = this->enabledFiles.begin(), end = this->enabledFiles.end();
1692
1693     for (; it != end; it++)
1694     {
1695       info(30)<<"Closing File : "<<(*it)->getId()<<endl;
1696       (*it)->close();
1697     }
1698   }
1699   CATCH_DUMP_ATTR
1700
1701   /*!
1702   \brief Dispatch event received from client
1703      Whenever a message is received in buffer of server, it will be processed depending on
1704   its event type. A new event type should be added in the switch list to make sure
1705   it processed on server side.
1706   \param [in] event: Received message
1707   */
1708   bool CContext::dispatchEvent(CEventServer& event)
1709   TRY
1710   {
1711
1712      if (SuperClass::dispatchEvent(event)) return true;
1713      else
1714      {
1715        switch(event.type)
1716        {
1717           case EVENT_ID_CLOSE_DEFINITION :
1718             recvCloseDefinition(event);
1719             return true;
1720             break;
1721           case EVENT_ID_UPDATE_CALENDAR:
1722             recvUpdateCalendar(event);
1723             return true;
1724             break;
1725           case EVENT_ID_CREATE_FILE_HEADER :
1726             recvCreateFileHeader(event);
1727             return true;
1728             break;
1729           case EVENT_ID_POST_PROCESS:
1730             recvPostProcessing(event);
1731             return true;
1732            case EVENT_ID_SEND_REGISTRY:
1733             recvRegistry(event);
1734             return true;
1735             break;
1736            case EVENT_ID_POST_PROCESS_GLOBAL_ATTRIBUTES:
1737             recvPostProcessingGlobalAttributes(event);
1738             return true;
1739             break;
1740            case EVENT_ID_PROCESS_GRID_ENABLED_FIELDS:
1741             recvProcessingGridOfEnabledFields(event);
1742             return true;
1743             break;
1744           default :
1745             ERROR("bool CContext::dispatchEvent(CEventServer& event)",
1746                    <<"Unknown Event");
1747           return false;
1748         }
1749      }
1750   }
1751   CATCH
1752
1753   //! Client side: Send a message to server to make it close
1754   void CContext::sendCloseDefinition(void)
1755   TRY
1756   {
1757    int nbSrvPools ;
1758    if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
1759    else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
1760    else nbSrvPools = 0 ;
1761    CContextClient* contextClientTmp ;
1762
1763    for (int i = 0; i < nbSrvPools; ++i)
1764     {
1765       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
1766       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
1767       CEventClient event(getType(),EVENT_ID_CLOSE_DEFINITION);
1768       if (contextClientTmp->isServerLeader())
1769       {
1770         CMessage msg;
1771         const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
1772         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
1773           event.push(*itRank,1,msg);
1774         contextClientTmp->sendEvent(event);
1775       }
1776       else contextClientTmp->sendEvent(event);
1777     }
1778   }
1779   CATCH_DUMP_ATTR
1780
1781   //! Server side: Receive a message of client announcing a context close
1782   void CContext::recvCloseDefinition(CEventServer& event)
1783   TRY
1784   {
1785      CBufferIn* buffer=event.subEvents.begin()->buffer;
1786      getCurrent()->closeDefinition();
1787   }
1788   CATCH
1789
1790   //! Client side: Send a message to update calendar in each time step
1791   void CContext::sendUpdateCalendar(int step)
1792   TRY
1793   {
1794    int nbSrvPools ;
1795    if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
1796    else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
1797    else nbSrvPools = 0 ;
1798    CContextClient* contextClientTmp ;
1799
1800    for (int i = 0; i < nbSrvPools; ++i)
1801    {
1802       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
1803       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
1804       CEventClient event(getType(),EVENT_ID_UPDATE_CALENDAR);
1805
1806         if (contextClientTmp->isServerLeader())
1807         {
1808           CMessage msg;
1809           msg<<step;
1810           const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
1811           for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
1812             event.push(*itRank,1,msg);
1813           contextClientTmp->sendEvent(event);
1814         }
1815         else contextClientTmp->sendEvent(event);
1816     }
1817   }
1818   CATCH_DUMP_ATTR
1819
1820   //! Server side: Receive a message of client annoucing calendar update
1821   void CContext::recvUpdateCalendar(CEventServer& event)
1822   TRY
1823   {
1824      CBufferIn* buffer=event.subEvents.begin()->buffer;
1825      getCurrent()->recvUpdateCalendar(*buffer);
1826   }
1827   CATCH
1828
1829   //! Server side: Receive a message of client annoucing calendar update
1830   void CContext::recvUpdateCalendar(CBufferIn& buffer)
1831   TRY
1832   {
1833      int step;
1834      buffer>>step;
1835      updateCalendar(step);
1836      if (serviceType_==CServicesManager::GATHERER)
1837      {       
1838        sendUpdateCalendar(step);
1839      }
1840   }
1841   CATCH_DUMP_ATTR
1842
1843   //! Client side: Send a message to create header part of netcdf file
1844   void CContext::sendCreateFileHeader(void)
1845   TRY
1846   {
1847     int nbSrvPools ;
1848     if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
1849     else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
1850     else nbSrvPools = 0 ;
1851     CContextClient* contextClientTmp ;
1852
1853     for (int i = 0; i < nbSrvPools; ++i)
1854     {
1855       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
1856       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
1857       CEventClient event(getType(),EVENT_ID_CREATE_FILE_HEADER);
1858
1859       if (contextClientTmp->isServerLeader())
1860       {
1861         CMessage msg;
1862         const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
1863         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
1864           event.push(*itRank,1,msg) ;
1865         contextClientTmp->sendEvent(event);
1866       }
1867       else contextClientTmp->sendEvent(event);
1868     }
1869   }
1870   CATCH_DUMP_ATTR
1871
1872   //! Server side: Receive a message of client annoucing the creation of header part of netcdf file
1873   void CContext::recvCreateFileHeader(CEventServer& event)
1874   TRY
1875   {
1876      CBufferIn* buffer=event.subEvents.begin()->buffer;
1877      getCurrent()->recvCreateFileHeader(*buffer);
1878   }
1879   CATCH
1880
1881   //! Server side: Receive a message of client annoucing the creation of header part of netcdf file
1882   void CContext::recvCreateFileHeader(CBufferIn& buffer)
1883   TRY
1884   {
1885      if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::OUT_SERVER) 
1886        createFileHeader();
1887   }
1888   CATCH_DUMP_ATTR
1889
1890   //! Client side: Send a message to do some post processing on server
1891   void CContext::sendProcessingGridOfEnabledFields()
1892   TRY
1893   {
1894     int nbSrvPools ;
1895     if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
1896     else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
1897     else nbSrvPools = 0 ;
1898     CContextClient* contextClientTmp ;
1899
1900     for (int i = 0; i < nbSrvPools; ++i)
1901     {
1902       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
1903       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
1904
1905       CEventClient event(getType(),EVENT_ID_PROCESS_GRID_ENABLED_FIELDS);
1906
1907       if (contextClientTmp->isServerLeader())
1908       {
1909         CMessage msg;
1910         const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
1911         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
1912           event.push(*itRank,1,msg);
1913         contextClientTmp->sendEvent(event);
1914       }
1915       else contextClientTmp->sendEvent(event);
1916     }
1917   }
1918   CATCH_DUMP_ATTR
1919
1920   //! Server side: Receive a message to do some post processing
1921   void CContext::recvProcessingGridOfEnabledFields(CEventServer& event)
1922   TRY
1923   {
1924      CBufferIn* buffer=event.subEvents.begin()->buffer;
1925      // nothing to do, no call ??!!   
1926   }
1927   CATCH
1928
1929   //! Client side: Send a message to do some post processing on server
1930   void CContext::sendPostProcessing()
1931   TRY
1932   {
1933     int nbSrvPools ;
1934     if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
1935     else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
1936     else nbSrvPools = 0 ;
1937     CContextClient* contextClientTmp ;
1938
1939     for (int i = 0; i < nbSrvPools; ++i)
1940     {
1941       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
1942       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
1943       CEventClient event(getType(),EVENT_ID_POST_PROCESS);
1944       if (contextClientTmp->isServerLeader())
1945       {
1946         CMessage msg;
1947         const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
1948         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
1949         event.push(*itRank,1,msg);
1950         contextClientTmp->sendEvent(event);
1951       }
1952       else contextClientTmp->sendEvent(event);
1953     }
1954   }
1955   CATCH_DUMP_ATTR
1956
1957   //! Server side: Receive a message to do some post processing
1958   void CContext::recvPostProcessing(CEventServer& event)
1959   TRY
1960   {
1961      CBufferIn* buffer=event.subEvents.begin()->buffer;
1962      getCurrent()->recvPostProcessing(*buffer);
1963   }
1964   CATCH
1965
1966   //! Server side: Receive a message to do some post processing
1967   void CContext::recvPostProcessing(CBufferIn& buffer)
1968   TRY
1969   {
1970      CCalendarWrapper::get(CCalendarWrapper::GetDefName())->createCalendar();
1971      postProcessing();
1972   }
1973   CATCH_DUMP_ATTR
1974
1975 
1976   /*!
1977   \brief Do some simple post processings after parsing xml file
1978      After the xml file (iodef.xml) is parsed, it is necessary to build all relations among
1979   created object, e.g: inhertance among fields, domain, axis. After that, all fiels as well as their parents (reference fields),
1980   which will be written out into netcdf files, are processed
1981   */
1982   void CContext::postProcessing()
1983   TRY
1984   {
1985     if (isPostProcessed) return;
1986
1987      // Make sure the calendar was correctly created
1988      if (!calendar)
1989        ERROR("CContext::postProcessing()", << "A calendar must be defined for the context \"" << getId() << "!\"")
1990      else if (calendar->getTimeStep() == NoneDu)
1991        ERROR("CContext::postProcessing()", << "A timestep must be defined for the context \"" << getId() << "!\"")
1992      // Calendar first update to set the current date equals to the start date
1993      calendar->update(0);
1994
1995      // Find all inheritance in xml structure
1996      this->solveAllInheritance();
1997
1998//      ShowTree(info(10));
1999
2000      // Check if some axis, domains or grids are eligible to for compressed indexed output.
2001      // Warning: This must be done after solving the inheritance and before the rest of post-processing
2002      checkAxisDomainsGridsEligibilityForCompressedOutput();      // only for field written on IO_SERVER service ????
2003
2004      // Check if some automatic time series should be generated
2005      // Warning: This must be done after solving the inheritance and before the rest of post-processing     
2006      prepareTimeseries();
2007
2008      //Initialisation du vecteur 'enabledFiles' contenant la liste des fichiers à sortir.
2009      findEnabledFiles();
2010      findEnabledWriteModeFiles();
2011      findEnabledReadModeFiles();
2012      findEnabledCouplerIn();
2013      findEnabledCouplerOut();
2014      createCouplerInterCommunicator() ;
2015
2016      // Find all enabled fields of each file     
2017      const vector<CField*>&& fileOutField = findAllEnabledFieldsInFiles(this->enabledWriteModeFiles);
2018      const vector<CField*>&& fileInField = findAllEnabledFieldsInFiles(this->enabledReadModeFiles);
2019      const vector<CField*>&& CouplerOutField = findAllEnabledFieldsCouplerOut(this->enabledCouplerOut);
2020      const vector<CField*>&& CouplerInField = findAllEnabledFieldsCouplerIn(this->enabledCouplerIn);
2021
2022
2023
2024      // For now, only read files with client and only one level server
2025      // if (hasClient && !hasServer) findEnabledReadModeFiles();     
2026
2027
2028      // For now, only read files with client and only one level server
2029      // if (hasClient && !hasServer)
2030      //   findAllEnabledFieldsInFiles(this->enabledReadModeFiles);     
2031
2032      if (serviceType_==CServicesManager::CLIENT)
2033      {
2034        initReadFiles();
2035        // Try to read attributes of fields in file then fill in corresponding grid (or domain, axis)
2036        this->readAttributesOfEnabledFieldsInReadModeFiles();
2037      }
2038
2039      // Only search and rebuild all reference objects of enable fields, don't transform
2040      this->solveOnlyRefOfEnabledFields();
2041
2042      // Search and rebuild all reference object of enabled fields, and transform
2043      this->solveAllRefOfEnabledFieldsAndTransform();
2044
2045      // Find all fields with read access from the public API
2046      if (serviceType_==CServicesManager::CLIENT) findFieldsWithReadAccess();
2047      // and solve the all reference for them
2048      if (serviceType_==CServicesManager::CLIENT) solveAllRefOfFieldsWithReadAccess();
2049
2050      isPostProcessed = true;
2051   }
2052   CATCH_DUMP_ATTR
2053
2054   void CContext::createCouplerInterCommunicator(void)
2055   TRY
2056   {
2057      // juste for test now, in future need an scheduler to avoid dead-lock
2058      for(auto it=enabledCouplerOut.begin();it!=enabledCouplerOut.end();++it)
2059      {
2060        (*it)->createInterCommunicator() ;
2061      }
2062
2063      for(auto it=enabledCouplerIn.begin();it!=enabledCouplerIn.end();++it)
2064      {
2065        (*it)->createInterCommunicator() ;
2066      }
2067   }
2068   CATCH_DUMP_ATTR
2069
2070   /*!
2071    * Compute the required buffer size to send the attributes (mostly those grid related).
2072    * \param maxEventSize [in/out] the size of the bigger event for each connected server
2073    * \param [in] contextClient
2074    * \param [in] bufferForWriting True if buffers are used for sending data for writing
2075      This flag is only true for client and server-1 for communication with server-2
2076    */
2077   std::map<int, StdSize> CContext::getAttributesBufferSize(std::map<int, StdSize>& maxEventSize,
2078                                                           CContextClient* contextClient, bool bufferForWriting /*= "false"*/)
2079   TRY
2080   {
2081         // As calendar attributes are sent even if there are no active files or fields, maps are initialized according the size of calendar attributes
2082     std::map<int, StdSize> attributesSize = CCalendarWrapper::get(CCalendarWrapper::GetDefName())->getMinimumBufferSizeForAttributes(contextClient);
2083     maxEventSize = CCalendarWrapper::get(CCalendarWrapper::GetDefName())->getMinimumBufferSizeForAttributes(contextClient);
2084
2085     std::vector<CFile*>& fileList = this->enabledFiles;
2086     size_t numEnabledFiles = fileList.size();
2087     for (size_t i = 0; i < numEnabledFiles; ++i)
2088     {
2089//         CFile* file = this->enabledWriteModeFiles[i];
2090        CFile* file = fileList[i];
2091        std::vector<CField*> enabledFields = file->getEnabledFields();
2092        size_t numEnabledFields = enabledFields.size();
2093        for (size_t j = 0; j < numEnabledFields; ++j)
2094        {
2095          const std::map<int, StdSize> mapSize = enabledFields[j]->getGridAttributesBufferSize(contextClient, bufferForWriting);
2096          std::map<int, StdSize>::const_iterator it = mapSize.begin(), itE = mapSize.end();
2097          for (; it != itE; ++it)
2098          {
2099                     // If attributesSize[it->first] does not exist, it will be zero-initialized
2100                     // so we can use it safely without checking for its existence
2101             if (attributesSize[it->first] < it->second)
2102                           attributesSize[it->first] = it->second;
2103
2104                     if (maxEventSize[it->first] < it->second)
2105                           maxEventSize[it->first] = it->second;
2106          }
2107        }
2108     }
2109     return attributesSize;
2110   }
2111   CATCH_DUMP_ATTR
2112
2113   /*!
2114    * Compute the required buffer size to send the fields data.
2115    * \param maxEventSize [in/out] the size of the bigger event for each connected server
2116    * \param [in] contextClient
2117    * \param [in] bufferForWriting True if buffers are used for sending data for writing
2118      This flag is only true for client and server-1 for communication with server-2
2119    */
2120   std::map<int, StdSize> CContext::getDataBufferSize(std::map<int, StdSize>& maxEventSize,
2121                                                      CContextClient* contextClient, bool bufferForWriting /*= "false"*/)
2122   TRY
2123   {
2124     std::map<int, StdSize> dataSize;
2125
2126     // Find all reference domain and axis of all active fields
2127     std::vector<CFile*>& fileList = bufferForWriting ? this->enabledWriteModeFiles : this->enabledReadModeFiles;
2128     size_t numEnabledFiles = fileList.size();
2129     for (size_t i = 0; i < numEnabledFiles; ++i)
2130     {
2131       CFile* file = fileList[i];
2132       if (file->getContextClient() == contextClient)
2133       {
2134         std::vector<CField*> enabledFields = file->getEnabledFields();
2135         size_t numEnabledFields = enabledFields.size();
2136         for (size_t j = 0; j < numEnabledFields; ++j)
2137         {
2138           // const std::vector<std::map<int, StdSize> > mapSize = enabledFields[j]->getGridDataBufferSize(contextClient);
2139           const std::map<int, StdSize> mapSize = enabledFields[j]->getGridDataBufferSize(contextClient,bufferForWriting);
2140           std::map<int, StdSize>::const_iterator it = mapSize.begin(), itE = mapSize.end();
2141           for (; it != itE; ++it)
2142           {
2143             // If dataSize[it->first] does not exist, it will be zero-initialized
2144             // so we can use it safely without checking for its existance
2145                 if (CXios::isOptPerformance)
2146               dataSize[it->first] += it->second;
2147             else if (dataSize[it->first] < it->second)
2148               dataSize[it->first] = it->second;
2149
2150                 if (maxEventSize[it->first] < it->second)
2151               maxEventSize[it->first] = it->second;
2152           }
2153         }
2154       }
2155     }
2156     return dataSize;
2157   }
2158   CATCH_DUMP_ATTR
2159
2160   //! Client side: Send infomation of active files (files are enabled to write out)
2161   void CContext::sendEnabledFiles(const std::vector<CFile*>& activeFiles)
2162   TRY
2163   {
2164     int size = activeFiles.size();
2165
2166     // In a context, each type has a root definition, e.g: axis, domain, field.
2167     // Every object must be a child of one of these root definition. In this case
2168     // all new file objects created on server must be children of the root "file_definition"
2169     StdString fileDefRoot("file_definition");
2170     CFileGroup* cfgrpPtr = CFileGroup::get(fileDefRoot);
2171
2172     for (int i = 0; i < size; ++i)
2173     {
2174       CFile* f = activeFiles[i];
2175       cfgrpPtr->sendCreateChild(f->getId(),f->getContextClient());
2176       f->sendAllAttributesToServer(f->getContextClient());
2177       f->sendAddAllVariables(f->getContextClient());
2178     }
2179   }
2180   CATCH_DUMP_ATTR
2181
2182   //! Client side: Send information of active fields (ones are written onto files)
2183   void CContext::sendEnabledFieldsInFiles(const std::vector<CFile*>& activeFiles)
2184   TRY
2185   {
2186     int size = activeFiles.size();
2187     for (int i = 0; i < size; ++i)
2188     {
2189       activeFiles[i]->sendEnabledFields(activeFiles[i]->getContextClient());
2190     }
2191   }
2192   CATCH_DUMP_ATTR
2193
2194   //! Client side: Check if the defined axis, domains and grids are eligible for compressed indexed output
2195   void CContext::checkAxisDomainsGridsEligibilityForCompressedOutput()
2196   TRY
2197   {
2198     if (!(serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER)) return;
2199
2200     const vector<CAxis*> allAxis = CAxis::getAll();
2201     for (vector<CAxis*>::const_iterator it = allAxis.begin(); it != allAxis.end(); it++)
2202       (*it)->checkEligibilityForCompressedOutput();
2203
2204     const vector<CDomain*> allDomains = CDomain::getAll();
2205     for (vector<CDomain*>::const_iterator it = allDomains.begin(); it != allDomains.end(); it++)
2206       (*it)->checkEligibilityForCompressedOutput();
2207
2208     const vector<CGrid*> allGrids = CGrid::getAll();
2209     for (vector<CGrid*>::const_iterator it = allGrids.begin(); it != allGrids.end(); it++)
2210       (*it)->checkEligibilityForCompressedOutput();
2211   }
2212   CATCH_DUMP_ATTR
2213
2214   //! Client side: Prepare the timeseries by adding the necessary files
2215   void CContext::prepareTimeseries()
2216   TRY
2217   {
2218     const std::vector<CFile*> allFiles = CFile::getAll();
2219     for (size_t i = 0; i < allFiles.size(); i++)
2220     {
2221       CFile* file = allFiles[i];
2222
2223       std::vector<CVariable*> fileVars, fieldVars, vars = file->getAllVariables();
2224       for (size_t k = 0; k < vars.size(); k++)
2225       {
2226         CVariable* var = vars[k];
2227
2228         if (var->ts_target.isEmpty()
2229              || var->ts_target == CVariable::ts_target_attr::file || var->ts_target == CVariable::ts_target_attr::both)
2230           fileVars.push_back(var);
2231
2232         if (!var->ts_target.isEmpty()
2233              && (var->ts_target == CVariable::ts_target_attr::field || var->ts_target == CVariable::ts_target_attr::both))
2234           fieldVars.push_back(var);
2235       }
2236
2237       if (!file->timeseries.isEmpty() && file->timeseries != CFile::timeseries_attr::none)
2238       {
2239         StdString fileNameStr("%file_name%") ;
2240         StdString tsPrefix = !file->ts_prefix.isEmpty() ? file->ts_prefix : fileNameStr ;
2241         
2242         StdString fileName=file->getFileOutputName();
2243         size_t pos=tsPrefix.find(fileNameStr) ;
2244         while (pos!=std::string::npos)
2245         {
2246           tsPrefix=tsPrefix.replace(pos,fileNameStr.size(),fileName) ;
2247           pos=tsPrefix.find(fileNameStr) ;
2248         }
2249       
2250         const std::vector<CField*> allFields = file->getAllFields();
2251         for (size_t j = 0; j < allFields.size(); j++)
2252         {
2253           CField* field = allFields[j];
2254
2255           if (!field->ts_enabled.isEmpty() && field->ts_enabled)
2256           {
2257             CFile* tsFile = CFile::create();
2258             tsFile->duplicateAttributes(file);
2259
2260             // Add variables originating from file and targeted to timeserie file
2261             for (size_t k = 0; k < fileVars.size(); k++)
2262               tsFile->getVirtualVariableGroup()->addChild(fileVars[k]);
2263
2264           
2265             tsFile->name = tsPrefix + "_";
2266             if (!field->name.isEmpty())
2267               tsFile->name.get() += field->name;
2268             else if (field->hasDirectFieldReference()) // We cannot use getBaseFieldReference() just yet
2269               tsFile->name.get() += field->field_ref;
2270             else
2271               tsFile->name.get() += field->getId();
2272
2273             if (!field->ts_split_freq.isEmpty())
2274               tsFile->split_freq = field->ts_split_freq;
2275
2276             CField* tsField = tsFile->addField();
2277             tsField->field_ref = field->getId();
2278
2279             // Add variables originating from file and targeted to timeserie field
2280             for (size_t k = 0; k < fieldVars.size(); k++)
2281               tsField->getVirtualVariableGroup()->addChild(fieldVars[k]);
2282
2283             vars = field->getAllVariables();
2284             for (size_t k = 0; k < vars.size(); k++)
2285             {
2286               CVariable* var = vars[k];
2287
2288               // Add variables originating from field and targeted to timeserie field
2289               if (var->ts_target.isEmpty()
2290                    || var->ts_target == CVariable::ts_target_attr::field || var->ts_target == CVariable::ts_target_attr::both)
2291                 tsField->getVirtualVariableGroup()->addChild(var);
2292
2293               // Add variables originating from field and targeted to timeserie file
2294               if (!var->ts_target.isEmpty()
2295                    && (var->ts_target == CVariable::ts_target_attr::file || var->ts_target == CVariable::ts_target_attr::both))
2296                 tsFile->getVirtualVariableGroup()->addChild(var);
2297             }
2298
2299             tsFile->solveFieldRefInheritance(true);
2300
2301             if (file->timeseries == CFile::timeseries_attr::exclusive)
2302               field->enabled = false;
2303           }
2304         }
2305
2306         // Finally disable the original file is need be
2307         if (file->timeseries == CFile::timeseries_attr::only)
2308          file->enabled = false;
2309       }
2310     }
2311   }
2312   CATCH_DUMP_ATTR
2313
2314   //! Client side: Send information of reference grid of active fields
2315   void CContext::sendRefGrid(const std::vector<CFile*>& activeFiles)
2316   TRY
2317   {
2318     std::set<pair<StdString,CContextClient*>> gridIds;
2319
2320     int sizeFile = activeFiles.size();
2321     CFile* filePtr(NULL);
2322
2323     // Firstly, find all reference grids of all active fields
2324     for (int i = 0; i < sizeFile; ++i)
2325     {
2326       filePtr = activeFiles[i];
2327       std::vector<CField*> enabledFields = filePtr->getEnabledFields();
2328       int sizeField = enabledFields.size();
2329       for (int numField = 0; numField < sizeField; ++numField)
2330       {
2331         if (0 != enabledFields[numField]->getRelGrid())
2332           gridIds.insert(make_pair(CGrid::get(enabledFields[numField]->getRelGrid())->getId(),enabledFields[numField]->getContextClient()));
2333       }
2334     }
2335
2336     // Create all reference grids on server side
2337     StdString gridDefRoot("grid_definition");
2338     CGridGroup* gridPtr = CGridGroup::get(gridDefRoot);
2339     for (auto it = gridIds.begin(); it != gridIds.end(); ++it)
2340     {
2341       gridPtr->sendCreateChild(it->first,it->second);
2342       CGrid::get(it->first)->sendAllAttributesToServer(it->second);
2343       CGrid::get(it->first)->sendAllDomains(it->second);
2344       CGrid::get(it->first)->sendAllAxis(it->second);
2345       CGrid::get(it->first)->sendAllScalars(it->second);
2346     }
2347   }
2348   CATCH_DUMP_ATTR
2349
2350   //! Client side: Send information of reference domain, axis and scalar of active fields
2351   void CContext::sendRefDomainsAxisScalars(const std::vector<CFile*>& activeFiles)
2352   TRY
2353   {
2354     std::set<pair<StdString,CContextClient*>> domainIds, axisIds, scalarIds;
2355
2356     // Find all reference domain and axis of all active fields
2357     int numEnabledFiles = activeFiles.size();
2358     for (int i = 0; i < numEnabledFiles; ++i)
2359     {
2360       std::vector<CField*> enabledFields = activeFiles[i]->getEnabledFields();
2361       int numEnabledFields = enabledFields.size();
2362       for (int j = 0; j < numEnabledFields; ++j)
2363       {
2364         CContextClient* contextClient=enabledFields[j]->getContextClient() ;
2365         const std::vector<StdString>& prDomAxisScalarId = enabledFields[j]->getRefDomainAxisIds();
2366         if ("" != prDomAxisScalarId[0]) domainIds.insert(make_pair(prDomAxisScalarId[0],contextClient));
2367         if ("" != prDomAxisScalarId[1]) axisIds.insert(make_pair(prDomAxisScalarId[1],contextClient));
2368         if ("" != prDomAxisScalarId[2]) scalarIds.insert(make_pair(prDomAxisScalarId[2],contextClient));
2369       }
2370     }
2371
2372     // Create all reference axis on server side
2373     std::set<StdString>::iterator itDom, itAxis, itScalar;
2374     std::set<StdString>::const_iterator itE;
2375
2376     StdString scalarDefRoot("scalar_definition");
2377     CScalarGroup* scalarPtr = CScalarGroup::get(scalarDefRoot);
2378     
2379     for (auto itScalar = scalarIds.begin(); itScalar != scalarIds.end(); ++itScalar)
2380     {
2381       if (!itScalar->first.empty())
2382       {
2383         scalarPtr->sendCreateChild(itScalar->first,itScalar->second);
2384         CScalar::get(itScalar->first)->sendAllAttributesToServer(itScalar->second);
2385       }
2386     }
2387
2388     StdString axiDefRoot("axis_definition");
2389     CAxisGroup* axisPtr = CAxisGroup::get(axiDefRoot);
2390     
2391     for (auto itAxis = axisIds.begin(); itAxis != axisIds.end(); ++itAxis)
2392     {
2393       if (!itAxis->first.empty())
2394       {
2395         axisPtr->sendCreateChild(itAxis->first, itAxis->second);
2396         CAxis::get(itAxis->first)->sendAllAttributesToServer(itAxis->second);
2397       }
2398     }
2399
2400     // Create all reference domains on server side
2401     StdString domDefRoot("domain_definition");
2402     CDomainGroup* domPtr = CDomainGroup::get(domDefRoot);
2403     
2404     for (auto itDom = domainIds.begin(); itDom != domainIds.end(); ++itDom)
2405     {
2406       if (!itDom->first.empty()) {
2407          domPtr->sendCreateChild(itDom->first, itDom->second);
2408          CDomain::get(itDom->first)->sendAllAttributesToServer(itDom->second);
2409       }
2410     }
2411   }
2412   CATCH_DUMP_ATTR
2413
2414   //! Update calendar in each time step
2415   void CContext::updateCalendar(int step)
2416   TRY
2417   {
2418      int prevStep = calendar->getStep();
2419
2420      if (prevStep < step)
2421      {
2422        if (serviceType_==CServicesManager::CLIENT) // For now we only use server level 1 to read data
2423        {
2424          doPreTimestepOperationsForEnabledReadModeFiles();
2425        }
2426
2427        info(50) << "updateCalendar : before : " << calendar->getCurrentDate() << endl;
2428        calendar->update(step);
2429        info(50) << "updateCalendar : after : " << calendar->getCurrentDate() << endl;
2430  #ifdef XIOS_MEMTRACK_LIGHT
2431        info(50) << " Current memory used by XIOS : "<<  MemTrack::getCurrentMemorySize()*1.0/(1024*1024)<<" Mbyte, at timestep "<<step<<" of context "<<this->getId()<<endl ;
2432  #endif
2433
2434        if (serviceType_==CServicesManager::CLIENT) // For now we only use server level 1 to read data
2435        {
2436          doPostTimestepOperationsForEnabledReadModeFiles();
2437          garbageCollector.invalidate(calendar->getCurrentDate());
2438        }
2439      }
2440      else if (prevStep == step)
2441        info(50) << "updateCalendar: already at step " << step << ", no operation done." << endl;
2442      else // if (prevStep > step)
2443        ERROR("void CContext::updateCalendar(int step)",
2444              << "Illegal calendar update: previous step was " << prevStep << ", new step " << step << "is in the past!")
2445   }
2446   CATCH_DUMP_ATTR
2447
2448   void CContext::initReadFiles(void)
2449   TRY
2450   {
2451      vector<CFile*>::const_iterator it;
2452
2453      for (it=enabledReadModeFiles.begin(); it != enabledReadModeFiles.end(); it++)
2454      {
2455         (*it)->initRead();
2456      }
2457   }
2458   CATCH_DUMP_ATTR
2459
2460   //! Server side: Create header of netcdf file
2461   void CContext::createFileHeader(void)
2462   TRY
2463   {
2464      vector<CFile*>::const_iterator it;
2465
2466      for (it=enabledFiles.begin(); it != enabledFiles.end(); it++)
2467      // for (it=enabledWriteModeFiles.begin(); it != enabledWriteModeFiles.end(); it++)
2468      {
2469         (*it)->initWrite();
2470      }
2471   }
2472   CATCH_DUMP_ATTR
2473
2474   //! Get current context
2475   CContext* CContext::getCurrent(void)
2476   TRY
2477   {
2478     return CObjectFactory::GetObject<CContext>(CObjectFactory::GetCurrentContextId()).get();
2479   }
2480   CATCH
2481
2482   /*!
2483   \brief Set context with an id be the current context
2484   \param [in] id identity of context to be set to current
2485   */
2486   void CContext::setCurrent(const string& id)
2487   TRY
2488   {
2489     CObjectFactory::SetCurrentContextId(id);
2490     CGroupFactory::SetCurrentContextId(id);
2491   }
2492   CATCH
2493
2494  /*!
2495  \brief Create a context with specific id
2496  \param [in] id identity of new context
2497  \return pointer to the new context or already-existed one with identity id
2498  */
2499  CContext* CContext::create(const StdString& id)
2500  TRY
2501  {
2502    CContext::setCurrent(id);
2503
2504    bool hasctxt = CContext::has(id);
2505    CContext* context = CObjectFactory::CreateObject<CContext>(id).get();
2506    getRoot();
2507    if (!hasctxt) CGroupFactory::AddChild(root, context->getShared());
2508
2509#define DECLARE_NODE(Name_, name_) \
2510    C##Name_##Definition::create(C##Name_##Definition::GetDefName());
2511#define DECLARE_NODE_PAR(Name_, name_)
2512#include "node_type.conf"
2513
2514    return (context);
2515  }
2516  CATCH
2517
2518     //! Server side: Receive a message to do some post processing
2519  void CContext::recvRegistry(CEventServer& event)
2520  TRY
2521  {
2522    CBufferIn* buffer=event.subEvents.begin()->buffer;
2523    getCurrent()->recvRegistry(*buffer);
2524  }
2525  CATCH
2526
2527  void CContext::recvRegistry(CBufferIn& buffer)
2528  TRY
2529  {
2530    if (server->intraCommRank==0)
2531    {
2532      CRegistry registry(server->intraComm) ;
2533      registry.fromBuffer(buffer) ;
2534      registryOut->mergeRegistry(registry) ;
2535    }
2536  }
2537  CATCH_DUMP_ATTR
2538
2539  void CContext::sendRegistry(void)
2540  TRY
2541  {
2542    registryOut->hierarchicalGatherRegistry() ;
2543
2544    int nbSrvPools ;
2545    if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
2546    else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
2547    else nbSrvPools = 0 ;
2548    CContextClient* contextClientTmp ;
2549
2550    for (int i = 0; i < nbSrvPools; ++i)
2551    {
2552      if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
2553      else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
2554
2555      CEventClient event(CContext::GetType(), CContext::EVENT_ID_SEND_REGISTRY);
2556      if (contextClientTmp->isServerLeader())
2557      {
2558        CMessage msg ;
2559        if (contextClientTmp->clientRank==0) msg<<*registryOut ;
2560        const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
2561        for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
2562             event.push(*itRank,1,msg);
2563        contextClientTmp->sendEvent(event);
2564      }
2565      else contextClientTmp->sendEvent(event);
2566    }
2567  }
2568  CATCH_DUMP_ATTR
2569
2570 
2571  void CContext::sendFinalizeClient(CContextClient* contextClient, const string& contextClientId)
2572  TRY
2573  {
2574    CEventClient event(getType(),EVENT_ID_CONTEXT_FINALIZE_CLIENT);
2575    if (contextClient->isServerLeader())
2576    {
2577      CMessage msg;
2578      msg<<contextClientId ;
2579      const std::list<int>& ranks = contextClient->getRanksServerLeader();
2580      for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
2581           event.push(*itRank,1,msg);
2582      contextClient->sendEvent(event);
2583    }
2584    else contextClient->sendEvent(event);
2585  }
2586  CATCH_DUMP_ATTR
2587
2588 
2589  void CContext::recvFinalizeClient(CEventServer& event)
2590  TRY
2591  {
2592    CBufferIn* buffer=event.subEvents.begin()->buffer;
2593    string id;
2594    *buffer>>id;
2595    get(id)->recvFinalizeClient(*buffer);
2596  }
2597  CATCH
2598
2599  void CContext::recvFinalizeClient(CBufferIn& buffer)
2600  TRY
2601  {
2602    countChildContextFinalized_++ ;
2603  }
2604  CATCH_DUMP_ATTR
2605
2606
2607
2608
2609  /*!
2610  * \fn bool CContext::isFinalized(void)
2611  * Context is finalized if it received context post finalize event.
2612  */
2613  bool CContext::isFinalized(void)
2614  TRY
2615  {
2616    return finalized;
2617  }
2618  CATCH_DUMP_ATTR
2619  ///--------------------------------------------------------------
2620  StdString CContext::dumpClassAttributes(void)
2621  {
2622    StdString str;
2623    str.append("enabled files=\"");
2624    int size = this->enabledFiles.size();
2625    for (int i = 0; i < size; ++i)
2626    {
2627      str.append(enabledFiles[i]->getId());
2628      str.append(" ");
2629    }
2630    str.append("\"");
2631    return str;
2632  }
2633
2634} // namespace xios
Note: See TracBrowser for help on using the repository browser.