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

Last change on this file since 1878 was 1878, checked in by ymipsl, 4 years ago

Coupling Branch.
Implementing a coupler scheduler, to impose order for intercommunicator creation between several coupling context.
Two way coupling is now working.

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: 109.5 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   void CContext::setClientServerBuffer(vector<CField*>& fields, bool bufferForWriting)
274   TRY
275   {
276      // Estimated minimum event size for small events (20 is an arbitrary constant just for safety)
277     const size_t minEventSize = CEventClient::headerSize + 20 * sizeof(int);
278      // Ensure there is at least some room for 20 of such events in the buffers
279     size_t minBufferSize = std::max(CXios::minBufferSize, 20 * minEventSize);
280
281#define DECLARE_NODE(Name_, name_)    \
282     if (minBufferSize < sizeof(C##Name_##Definition)) minBufferSize = sizeof(C##Name_##Definition);
283#define DECLARE_NODE_PAR(Name_, name_)
284#include "node_type.conf"
285#undef DECLARE_NODE
286#undef DECLARE_NODE_PAR
287
288
289     map<CContextClient*,map<int,size_t>> dataSize ;
290     map<CContextClient*,map<int,size_t>> maxEventSize ;
291     map<CContextClient*,map<int,size_t>> attributesSize ; 
292
293     for(auto field : fields)
294     {
295       field->setContextClientDataBufferSize(dataSize, maxEventSize, bufferForWriting) ;
296       field->setContextClientAttributesBufferSize(attributesSize, maxEventSize, bufferForWriting) ;
297     }
298     
299
300     for(auto& it : attributesSize)
301     {
302       auto contextClient = it.first ;
303       auto& contextDataSize =  dataSize[contextClient] ;
304       auto& contextAttributesSize =  attributesSize[contextClient] ;
305       auto& contextMaxEventSize =  maxEventSize[contextClient] ;
306   
307       for (auto& it : contextAttributesSize)
308       {
309         auto serverRank=it.first ;
310         auto& buffer = contextAttributesSize[serverRank] ;
311         if (contextDataSize[serverRank] > buffer) buffer=contextDataSize[serverRank] ;
312         buffer *= CXios::bufferSizeFactor;
313         if (buffer < minBufferSize) buffer = minBufferSize;
314         if (buffer > CXios::maxBufferSize ) buffer = CXios::maxBufferSize;
315       }
316
317       // Leaders will have to send some control events so ensure there is some room for those in the buffers
318       if (contextClient->isServerLeader())
319         for(auto& rank : contextClient->getRanksServerLeader())
320           if (!contextAttributesSize.count(rank))
321           {
322             contextAttributesSize[rank] = minBufferSize;
323             contextMaxEventSize[rank] = minEventSize;
324           }
325     
326       contextClient->setBufferSize(contextAttributesSize, contextMaxEventSize);   
327     }
328   }
329   CATCH_DUMP_ATTR
330
331
332    /*!
333    Sets client buffers.
334    \param [in] contextClient
335    \param [in] bufferForWriting True if buffers are used for sending data for writing
336    This flag is only true for client and server-1 for communication with server-2
337  */
338  // ym obsolete to be removed
339   void CContext::setClientServerBuffer(CContextClient* contextClient, bool bufferForWriting)
340   TRY
341   {
342      // Estimated minimum event size for small events (20 is an arbitrary constant just for safety)
343     const size_t minEventSize = CEventClient::headerSize + 20 * sizeof(int);
344
345      // Ensure there is at least some room for 20 of such events in the buffers
346      size_t minBufferSize = std::max(CXios::minBufferSize, 20 * minEventSize);
347
348#define DECLARE_NODE(Name_, name_)    \
349     if (minBufferSize < sizeof(C##Name_##Definition)) minBufferSize = sizeof(C##Name_##Definition);
350#define DECLARE_NODE_PAR(Name_, name_)
351#include "node_type.conf"
352#undef DECLARE_NODE
353#undef DECLARE_NODE_PAR
354
355     // Compute the buffer sizes needed to send the attributes and data corresponding to fields
356     std::map<int, StdSize> maxEventSize;
357     std::map<int, StdSize> bufferSize = getAttributesBufferSize(maxEventSize, contextClient, bufferForWriting);
358     std::map<int, StdSize> dataBufferSize = getDataBufferSize(maxEventSize, contextClient, bufferForWriting);
359
360     std::map<int, StdSize>::iterator it, ite = dataBufferSize.end();
361     for (it = dataBufferSize.begin(); it != ite; ++it)
362       if (it->second > bufferSize[it->first]) bufferSize[it->first] = it->second;
363
364     // Apply the buffer size factor, check that we are above the minimum buffer size and below the maximum size
365     ite = bufferSize.end();
366     for (it = bufferSize.begin(); it != ite; ++it)
367     {
368       it->second *= CXios::bufferSizeFactor;
369       if (it->second < minBufferSize) it->second = minBufferSize;
370       if (it->second > CXios::maxBufferSize) it->second = CXios::maxBufferSize;
371     }
372
373     // Leaders will have to send some control events so ensure there is some room for those in the buffers
374     if (contextClient->isServerLeader())
375     {
376       const std::list<int>& ranks = contextClient->getRanksServerLeader();
377       for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
378       {
379         if (!bufferSize.count(*itRank))
380         {
381           bufferSize[*itRank] = minBufferSize;
382           maxEventSize[*itRank] = minEventSize;
383         }
384       }
385     }
386     contextClient->setBufferSize(bufferSize, maxEventSize);
387   }
388   CATCH_DUMP_ATTR
389
390 /*!
391    * Compute the required buffer size to send the fields data.
392    * \param maxEventSize [in/out] the size of the bigger event for each connected server
393    * \param [in] contextClient
394    * \param [in] bufferForWriting True if buffers are used for sending data for writing
395      This flag is only true for client and server-1 for communication with server-2
396    */
397   std::map<int, StdSize> CContext::getDataBufferSize(std::map<int, StdSize>& maxEventSize,
398                                                      CContextClient* contextClient, bool bufferForWriting /*= "false"*/)
399   TRY
400   {
401     std::map<int, StdSize> dataSize;
402
403     // Find all reference domain and axis of all active fields
404     std::vector<CFile*>& fileList = bufferForWriting ? this->enabledWriteModeFiles : this->enabledReadModeFiles;
405     size_t numEnabledFiles = fileList.size();
406     for (size_t i = 0; i < numEnabledFiles; ++i)
407     {
408       CFile* file = fileList[i];
409       if (file->getContextClient() == contextClient)
410       {
411         std::vector<CField*> enabledFields = file->getEnabledFields();
412         size_t numEnabledFields = enabledFields.size();
413         for (size_t j = 0; j < numEnabledFields; ++j)
414         {
415           // const std::vector<std::map<int, StdSize> > mapSize = enabledFields[j]->getGridDataBufferSize(contextClient);
416           const std::map<int, StdSize> mapSize = enabledFields[j]->getGridDataBufferSize(contextClient,bufferForWriting);
417           std::map<int, StdSize>::const_iterator it = mapSize.begin(), itE = mapSize.end();
418           for (; it != itE; ++it)
419           {
420             // If dataSize[it->first] does not exist, it will be zero-initialized
421             // so we can use it safely without checking for its existance
422           if (CXios::isOptPerformance)
423               dataSize[it->first] += it->second;
424             else if (dataSize[it->first] < it->second)
425               dataSize[it->first] = it->second;
426
427           if (maxEventSize[it->first] < it->second)
428               maxEventSize[it->first] = it->second;
429           }
430         }
431       }
432     }
433     return dataSize;
434   }
435   CATCH_DUMP_ATTR
436
437/*!
438    * Compute the required buffer size to send the attributes (mostly those grid related).
439    * \param maxEventSize [in/out] the size of the bigger event for each connected server
440    * \param [in] contextClient
441    * \param [in] bufferForWriting True if buffers are used for sending data for writing
442      This flag is only true for client and server-1 for communication with server-2
443    */
444   std::map<int, StdSize> CContext::getAttributesBufferSize(std::map<int, StdSize>& maxEventSize,
445                                                           CContextClient* contextClient, bool bufferForWriting /*= "false"*/)
446   TRY
447   {
448   // As calendar attributes are sent even if there are no active files or fields, maps are initialized according the size of calendar attributes
449     std::map<int, StdSize> attributesSize = CCalendarWrapper::get(CCalendarWrapper::GetDefName())->getMinimumBufferSizeForAttributes(contextClient);
450     maxEventSize = CCalendarWrapper::get(CCalendarWrapper::GetDefName())->getMinimumBufferSizeForAttributes(contextClient);
451
452     std::vector<CFile*>& fileList = this->enabledFiles;
453     size_t numEnabledFiles = fileList.size();
454     for (size_t i = 0; i < numEnabledFiles; ++i)
455     {
456//         CFile* file = this->enabledWriteModeFiles[i];
457        CFile* file = fileList[i];
458        std::vector<CField*> enabledFields = file->getEnabledFields();
459        size_t numEnabledFields = enabledFields.size();
460        for (size_t j = 0; j < numEnabledFields; ++j)
461        {
462          const std::map<int, StdSize> mapSize = enabledFields[j]->getGridAttributesBufferSize(contextClient, bufferForWriting);
463          std::map<int, StdSize>::const_iterator it = mapSize.begin(), itE = mapSize.end();
464          for (; it != itE; ++it)
465          {
466         // If attributesSize[it->first] does not exist, it will be zero-initialized
467         // so we can use it safely without checking for its existence
468             if (attributesSize[it->first] < it->second)
469         attributesSize[it->first] = it->second;
470
471         if (maxEventSize[it->first] < it->second)
472         maxEventSize[it->first] = it->second;
473          }
474        }
475     }
476     return attributesSize;
477   }
478   CATCH_DUMP_ATTR
479
480
481
482   //! Verify whether a context is initialized
483   bool CContext::isInitialized(void)
484   TRY
485   {
486     return hasClient;
487   }
488   CATCH_DUMP_ATTR
489
490
491   void CContext::init(CServerContext* parentServerContext, MPI_Comm intraComm, int serviceType)
492   TRY
493   {
494     parentServerContext_ = parentServerContext ;
495     if (serviceType==CServicesManager::CLIENT) 
496       initClient(intraComm, serviceType) ;
497     else
498       initServer(intraComm, serviceType) ;
499    }
500    CATCH_DUMP_ATTR
501
502
503
504//! Initialize client side
505   void CContext::initClient(MPI_Comm intraComm, int serviceType)
506   TRY
507   {
508      intraComm_=intraComm ;
509      MPI_Comm_rank(intraComm_, &intraCommRank_) ;
510      MPI_Comm_size(intraComm_, &intraCommSize_) ;
511
512      serviceType_ = CServicesManager::CLIENT ;
513      if (serviceType_==CServicesManager::CLIENT)
514      {
515        hasClient=true ;
516        hasServer=false ;
517      }
518      contextId_ = getId() ;
519     
520      attached_mode=true ;
521      if (!CXios::isUsingServer()) attached_mode=false ;
522
523
524      string contextRegistryId=getId() ;
525      registryIn=new CRegistry(intraComm);
526      registryIn->setPath(contextRegistryId) ;
527     
528      int commRank ;
529      MPI_Comm_rank(intraComm_,&commRank) ;
530      if (commRank==0) registryIn->fromFile("xios_registry.bin") ;
531      registryIn->bcastRegistry() ;
532      registryOut=new CRegistry(intraComm_) ;
533      registryOut->setPath(contextRegistryId) ;
534     
535   }
536   CATCH_DUMP_ATTR
537
538   
539   void CContext::initServer(MPI_Comm intraComm, int serviceType)
540   TRY
541   {
542     hasServer=true;
543     intraComm_=intraComm ;
544     MPI_Comm_rank(intraComm_, &intraCommRank_) ;
545     MPI_Comm_size(intraComm_, &intraCommSize_) ;
546
547     serviceType_=serviceType ;
548
549     if (serviceType_==CServicesManager::GATHERER)
550     {
551       hasClient=true ;
552       hasServer=true ;
553     }
554     else if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::OUT_SERVER)
555     {
556       hasClient=false ;
557       hasServer=true ;
558     }
559
560     CXios::getContextsManager()->getContextId(getId(), contextId_, intraComm) ;
561     
562     registryIn=new CRegistry(intraComm);
563     registryIn->setPath(contextId_) ;
564     
565     int commRank ;
566     MPI_Comm_rank(intraComm_,&commRank) ;
567     if (commRank==0) registryIn->fromFile("xios_registry.bin") ;
568   
569     registryIn->bcastRegistry() ;
570     registryOut=new CRegistry(intraComm) ;
571     registryOut->setPath(contextId_) ;
572
573   }
574   CATCH_DUMP_ATTR
575
576
577  void CContext::createClientInterComm(MPI_Comm interCommClient, MPI_Comm interCommServer) // for servers
578  TRY
579  {
580    MPI_Comm intraCommClient ;
581    MPI_Comm_dup(intraComm_, &intraCommClient);
582    comms.push_back(intraCommClient);
583    // attached_mode=parentServerContext_->isAttachedMode() ; //ym probably inherited from source context
584    server = new CContextServer(this,intraComm_, interCommServer); // check if we need to dupl. intraComm_ ?
585    client = new CContextClient(this,intraCommClient,interCommClient);
586    client->setAssociatedServer(server) ; 
587    server->setAssociatedClient(client) ; 
588
589  }
590  CATCH_DUMP_ATTR
591
592  void CContext::createServerInterComm(void) 
593  TRY
594  {
595   
596    MPI_Comm interCommClient, interCommServer ;
597
598    if (serviceType_ == CServicesManager::CLIENT)
599    {
600
601      int commRank ;
602      MPI_Comm_rank(intraComm_,&commRank) ;
603      if (commRank==0)
604      {
605        if (attached_mode) CXios::getContextsManager()->createServerContext(CClient::getPoolRessource()->getId(), CXios::defaultServerId, 0, getContextId()) ;
606        else if (CXios::usingServer2) CXios::getContextsManager()->createServerContext(CXios::defaultPoolId, CXios::defaultGathererId, 0, getContextId()) ;
607        else  CXios::getContextsManager()->createServerContext(CXios::defaultPoolId, CXios::defaultServerId, 0, getContextId()) ;
608      }
609
610      MPI_Comm interComm ;
611     
612      if (attached_mode)
613      {
614        parentServerContext_->createIntercomm(CClient::getPoolRessource()->getId(), CXios::defaultServerId, 0, getContextId(), intraComm_, 
615                                              interCommClient, interCommServer) ;
616        int type ; 
617        if (commRank==0) CXios::getServicesManager()->getServiceType(CClient::getPoolRessource()->getId(), CXios::defaultServerId, 0, type) ;
618        MPI_Bcast(&type,1,MPI_INT,0,intraComm_) ;
619        setCurrent(getId()) ; // getCurrent/setCurrent may be supress, it can cause a lot of trouble
620      }
621      else if (CXios::usingServer2)
622      { 
623//      CXios::getContextsManager()->createServerContextIntercomm(CXios::defaultPoolId, CXios::defaultGathererId, 0, getContextId(), intraComm_, interComm) ;
624        parentServerContext_->createIntercomm(CXios::defaultPoolId, CXios::defaultGathererId, 0, getContextId(), intraComm_,
625                                              interCommClient, interCommServer) ;
626        int type ; 
627        if (commRank==0) CXios::getServicesManager()->getServiceType(CXios::defaultPoolId, CXios::defaultGathererId, 0, type) ;
628        MPI_Bcast(&type,1,MPI_INT,0,intraComm_) ;
629      }
630      else
631      {
632        //CXios::getContextsManager()->createServerContextIntercomm(CXios::defaultPoolId, CXios::defaultServerId, 0, getContextId(), intraComm_, interComm) ;
633        parentServerContext_->createIntercomm(CXios::defaultPoolId, CXios::defaultServerId, 0, getContextId(), intraComm_,
634                                              interCommClient, interCommServer) ;
635        int type ; 
636        if (commRank==0) CXios::getServicesManager()->getServiceType(CXios::defaultPoolId, CXios::defaultServerId, 0, type) ;
637        MPI_Bcast(&type,1,MPI_INT,0,intraComm_) ;
638      }
639
640        // intraComm client is not duplicated. In all the code we use client->intraComm for MPI
641        // in future better to replace it by intracommuncator associated to the context
642   
643      MPI_Comm intraCommClient, intraCommServer ;
644      intraCommClient=intraComm_ ;
645      MPI_Comm_dup(intraComm_, &intraCommServer) ;
646      client = new CContextClient(this, intraCommClient, interCommClient);
647      server = new CContextServer(this, intraCommServer, interCommServer);
648      client->setAssociatedServer(server) ;
649      server->setAssociatedClient(client) ;
650    }
651   
652    if (serviceType_ == CServicesManager::GATHERER)
653    {
654      int commRank ;
655      MPI_Comm_rank(intraComm_,&commRank) ;
656     
657      int nbPartitions ;
658      if (commRank==0) 
659      { 
660        CXios::getServicesManager()->getServiceNbPartitions(CXios::defaultPoolId, CXios::defaultServerId, 0, nbPartitions) ;
661        for(int i=0 ; i<nbPartitions; i++)
662          CXios::getContextsManager()->createServerContext(CXios::defaultPoolId, CXios::defaultServerId, i, getContextId()) ;
663      }     
664      MPI_Bcast(&nbPartitions, 1, MPI_INT, 0, intraComm_) ;
665     
666      MPI_Comm interComm ;
667      for(int i=0 ; i<nbPartitions; i++)
668      {
669        parentServerContext_->createIntercomm(CXios::defaultPoolId, CXios::defaultServerId, i, getContextId(), intraComm_, interCommClient, interCommServer) ;
670        int type ; 
671        if (commRank==0) CXios::getServicesManager()->getServiceType(CXios::defaultPoolId, CXios::defaultServerId, 0, type) ;
672        MPI_Bcast(&type,1,MPI_INT,0,intraComm_) ;
673        primServerId_.push_back(CXios::getContextsManager()->getServerContextName(CXios::defaultPoolId, CXios::defaultServerId, i, type, getContextId())) ;
674
675        // intraComm client is not duplicated. In all the code we use client->intraComm for MPI
676        // in future better to replace it by intracommuncator associated to the context
677     
678        MPI_Comm intraCommClient, intraCommServer ;
679
680        intraCommClient=intraComm_ ;
681        MPI_Comm_dup(intraComm_, &intraCommServer) ;
682
683        CContextClient* client = new CContextClient(this, intraCommClient, interCommClient) ;
684        CContextServer* server = new CContextServer(this, intraCommServer, interCommServer) ;
685        client->setAssociatedServer(server) ;
686        server->setAssociatedClient(client) ;
687        clientPrimServer.push_back(client);
688        serverPrimServer.push_back(server); 
689
690     
691      }
692    }
693  }
694  CATCH_DUMP_ATTR
695
696   
697
698  bool CContext::eventLoop(bool enableEventsProcessing)
699  {
700    bool finished=true; 
701
702    if (client!=nullptr && !finalized) client->checkBuffers();
703   
704    for (int i = 0; i < clientPrimServer.size(); ++i)
705    {
706      if (!finalized) clientPrimServer[i]->checkBuffers();
707      if (!finalized) finished &= serverPrimServer[i]->eventLoop(enableEventsProcessing);
708    }
709
710    for (auto couplerOut : couplerOutClient_)
711      if (!finalized) couplerOut.second->checkBuffers();
712   
713    for (auto couplerIn : couplerInClient_)
714      if (!finalized) couplerIn.second->checkBuffers();
715   
716    for (auto couplerOut : couplerOutServer_)
717      if (!finalized) couplerOut.second->eventLoop(enableEventsProcessing);
718
719    for (auto couplerIn : couplerInServer_)
720      if (!finalized) couplerIn.second->eventLoop(enableEventsProcessing);
721   
722    if (server!=nullptr) if (!finalized) finished &= server->eventLoop(enableEventsProcessing);
723 
724    return finalized && finished ;
725  }
726
727  void CContext::addCouplingChanel(const std::string& fullContextId, bool out)
728  {
729     int contextLeader ;
730     
731     if (out)
732     { 
733       if (couplerOutClient_.find(fullContextId)==couplerOutClient_.end()) 
734       {
735         bool ok=CXios::getContextsManager()->getContextLeader(fullContextId, contextLeader, getIntraComm()) ;
736     
737         MPI_Comm interComm, interCommClient, interCommServer  ;
738         MPI_Comm intraCommClient, intraCommServer ;
739
740         if (ok) MPI_Intercomm_create(getIntraComm(), 0, CXios::getXiosComm(), contextLeader, 0, &interComm) ;
741
742        MPI_Comm_dup(intraComm_, &intraCommClient) ;
743        MPI_Comm_dup(intraComm_, &intraCommServer) ;
744        MPI_Comm_dup(interComm, &interCommClient) ;
745        MPI_Comm_dup(interComm, &interCommServer) ;
746        CContextClient* client = new CContextClient(this, intraCommClient, interCommClient);
747        CContextServer* server = new CContextServer(this, intraCommServer, interCommServer);
748        client->setAssociatedServer(server) ;
749        server->setAssociatedClient(client) ;
750        MPI_Comm_free(&interComm) ;
751        couplerOutClient_[fullContextId] = client ;
752        couplerOutServer_[fullContextId] = server ;
753
754/*
755      // for now, we don't now which beffer size must be used for client coupler
756      // It will be evaluated later. Fix a constant size for now...
757      // set to 10Mb for development
758       map<int,size_t> bufferSize, maxEventSize ;
759       for(int i=0;i<client->getRemoteSize();i++)
760       {
761         bufferSize[i]=10000000 ;
762         maxEventSize[i]=10000000 ;
763       }
764
765       client->setBufferSize(bufferSize, maxEventSize);   
766*/
767      }
768    }
769    else if (couplerInClient_.find(fullContextId)==couplerInClient_.end())
770    {
771      bool ok=CXios::getContextsManager()->getContextLeader(fullContextId, contextLeader, getIntraComm()) ;
772     
773       MPI_Comm interComm, interCommClient, interCommServer  ;
774       MPI_Comm intraCommClient, intraCommServer ;
775
776       if (ok) MPI_Intercomm_create(getIntraComm(), 0, CXios::getXiosComm(), contextLeader, 0, &interComm) ;
777
778       MPI_Comm_dup(intraComm_, &intraCommClient) ;
779       MPI_Comm_dup(intraComm_, &intraCommServer) ;
780       MPI_Comm_dup(interComm, &interCommServer) ;
781       MPI_Comm_dup(interComm, &interCommClient) ;
782       CContextServer* server = new CContextServer(this, intraCommServer, interCommServer);
783       CContextClient* client = new CContextClient(this, intraCommClient, interCommClient);
784       client->setAssociatedServer(server) ;
785       server->setAssociatedClient(client) ;
786       MPI_Comm_free(&interComm) ;
787
788       map<int,size_t> bufferSize, maxEventSize ;
789       for(int i=0;i<client->getRemoteSize();i++)
790       {
791         bufferSize[i]=10000000 ;
792         maxEventSize[i]=10000000 ;
793       }
794
795       client->setBufferSize(bufferSize, maxEventSize);   
796       couplerInClient_[fullContextId] = client ;
797       couplerInServer_[fullContextId] = server ;       
798    }
799  }
800 
801  void CContext::globalEventLoop(void)
802  {
803    CXios::getDaemonsManager()->eventLoop() ;
804    setCurrent(getId()) ;
805  }
806
807
808   void CContext::finalize(void)
809   TRY
810   {
811      registryOut->hierarchicalGatherRegistry() ;
812      if (server->intraCommRank==0) CXios::globalRegistry->mergeRegistry(*registryOut) ;
813
814      if (serviceType_==CServicesManager::CLIENT)
815      {
816//ym        doPreTimestepOperationsForEnabledReadModeFiles(); // For now we only use server level 1 to read data
817
818        triggerLateFields() ;
819
820        // inform couplerIn that I am finished
821        for(auto& couplerInClient : couplerInClient_) sendCouplerInContextFinalized(couplerInClient.second) ;
822
823        // wait until received message from couplerOut that they have finished
824        bool couplersInFinalized ;
825        do
826        {
827          couplersInFinalized=true ;
828          for(auto& couplerOutClient : couplerOutClient_) couplersInFinalized &= isCouplerInContextFinalized(couplerOutClient.second) ; 
829          globalEventLoop() ;
830        } while (!couplersInFinalized) ;
831
832        info(100)<<"DEBUG: context "<<getId()<<" Send client finalize"<<endl ;
833        client->finalize();
834        info(100)<<"DEBUG: context "<<getId()<<" Client finalize sent"<<endl ;
835        while (client->havePendingRequests()) client->checkBuffers();
836        info(100)<<"DEBUG: context "<<getId()<<" no pending request ok"<<endl ;
837        bool notifiedFinalized=false ;
838        do
839        {
840          notifiedFinalized=client->isNotifiedFinalized() ;
841        } while (!notifiedFinalized) ;
842        client->releaseBuffers();
843        info(100)<<"DEBUG: context "<<getId()<<" release client ok"<<endl ;
844      }
845      else if (serviceType_==CServicesManager::GATHERER)
846      {
847         for (int i = 0; i < clientPrimServer.size(); ++i)
848         {
849           clientPrimServer[i]->finalize();
850           bool bufferReleased;
851           do
852           {
853             clientPrimServer[i]->checkBuffers();
854             bufferReleased = !clientPrimServer[i]->havePendingRequests();
855           } while (!bufferReleased);
856           
857           bool notifiedFinalized=false ;
858           do
859           {
860             notifiedFinalized=clientPrimServer[i]->isNotifiedFinalized() ;
861           } while (!notifiedFinalized) ;
862           clientPrimServer[i]->releaseBuffers();
863         }
864         closeAllFile();
865
866      }
867      else if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::OUT_SERVER)
868      {
869        closeAllFile();
870      }
871
872      freeComms() ;
873       
874      parentServerContext_->freeComm() ;
875      finalized = true;
876      info(20)<<"CContext: Context <"<<getId()<<"> is finalized."<<endl;
877   }
878   CATCH_DUMP_ATTR
879
880   //! Free internally allocated communicators
881   void CContext::freeComms(void)
882   TRY
883   {
884     for (std::list<MPI_Comm>::iterator it = comms.begin(); it != comms.end(); ++it)
885       MPI_Comm_free(&(*it));
886     comms.clear();
887   }
888   CATCH_DUMP_ATTR
889
890   //! Deallocate buffers allocated by clientContexts
891   void CContext::releaseClientBuffers(void)
892   TRY
893   {
894     client->releaseBuffers();
895     for (int i = 0; i < clientPrimServer.size(); ++i)
896       clientPrimServer[i]->releaseBuffers();
897   }
898   CATCH_DUMP_ATTR
899
900   void CContext::postProcessingGlobalAttributes()
901   TRY
902   {
903     if (allProcessed) return; 
904     
905    // create intercommunicator with servers.
906    // not sure it is the good place to be called here
907    createServerInterComm() ;
908
909
910     // After xml is parsed, there are some more works with post processing
911     postProcessing();
912
913     // Distribute files between secondary servers according to the data size
914     distributeFiles(this->enabledWriteModeFiles);
915
916     // Check grid and calculate its distribution
917     checkGridEnabledFields();
918
919     setClientServerBuffer(client, (serviceType_==CServicesManager::CLIENT) ) ;
920     for (int i = 0; i < clientPrimServer.size(); ++i)
921         setClientServerBuffer(clientPrimServer[i], true);
922
923   
924     if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER)
925     { 
926       if (serviceType_==CServicesManager::GATHERER)
927       { 
928         for (auto it=clientPrimServer.begin(); it!=clientPrimServer.end();++it) 
929         {
930           this->sendAllAttributesToServer(*it); // Send all attributes of current context to server
931           CCalendarWrapper::get(CCalendarWrapper::GetDefName())->sendAllAttributesToServer(*it); // Send all attributes of current calendar
932         }
933       }
934       else 
935       {
936         this->sendAllAttributesToServer(client);   // Send all attributes of current context to server
937         CCalendarWrapper::get(CCalendarWrapper::GetDefName())->sendAllAttributesToServer(client); // Send all attributes of current calendar
938       }
939
940
941      // We have enough information to send to server
942      // First of all, send all enabled files
943      sendEnabledFiles(this->enabledWriteModeFiles);
944      // We only use server-level 1 (for now) to read data
945      if (serviceType_==CServicesManager::CLIENT)  sendEnabledFiles(this->enabledReadModeFiles);
946
947      // Then, send all enabled fields     
948      sendEnabledFieldsInFiles(this->enabledWriteModeFiles);
949     
950      if (serviceType_==CServicesManager::CLIENT) sendEnabledFieldsInFiles(this->enabledReadModeFiles);
951
952      // Then, check whether we have domain_ref, axis_ref or scalar_ref attached to the enabled fields
953      // If any, so send them to server
954       sendRefDomainsAxisScalars(this->enabledWriteModeFiles); 
955     
956      if (serviceType_==CServicesManager::CLIENT) sendRefDomainsAxisScalars(this->enabledReadModeFiles);       
957
958       // Check whether enabled fields have grid_ref, if any, send this info to server
959      sendRefGrid(this->enabledFiles);
960      // This code may be useful in the future when we want to seperate completely read and write
961      // sendRefGrid(this->enabledWriteModeFiles);
962      // if (!hasServer)
963      //   sendRefGrid(this->enabledReadModeFiles);
964     
965      // A grid of enabled fields composed of several components which must be checked then their
966      // checked attributes should be sent to server
967      sendGridComponentEnabledFieldsInFiles(this->enabledFiles); // This code can be seperated in two (one for reading, another for writing)
968
969       // We have a xml tree on the server side and now, it should be also processed
970      sendPostProcessing();
971       
972      // Finally, we send information of grid itself to server
973      sendGridEnabledFieldsInFiles(this->enabledWriteModeFiles);       
974     
975      if (serviceType_==CServicesManager::CLIENT) sendGridEnabledFieldsInFiles(this->enabledReadModeFiles);       
976     
977     }
978     allProcessed = true;
979   }
980   CATCH_DUMP_ATTR
981
982   void CContext::sendPostProcessingGlobalAttributes()
983   TRY
984   {
985
986    int nbSrvPools ;
987    if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
988    else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
989    else nbSrvPools = 0 ;
990    CContextClient* contextClientTmp ;
991
992    for (int i = 0; i < nbSrvPools; ++i)
993     {
994       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
995       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
996     
997       CEventClient event(getType(),EVENT_ID_POST_PROCESS_GLOBAL_ATTRIBUTES);
998
999       if (contextClientTmp->isServerLeader())
1000       {
1001         CMessage msg;
1002         const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
1003         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
1004           event.push(*itRank,1,msg);
1005         contextClientTmp->sendEvent(event);
1006       }
1007       else contextClientTmp->sendEvent(event);
1008     }
1009   }
1010   CATCH_DUMP_ATTR
1011
1012   void CContext::recvPostProcessingGlobalAttributes(CEventServer& event)
1013   TRY
1014   {
1015      CBufferIn* buffer=event.subEvents.begin()->buffer;
1016      getCurrent()->recvPostProcessingGlobalAttributes(*buffer);
1017   }
1018   CATCH
1019
1020   void CContext::recvPostProcessingGlobalAttributes(CBufferIn& buffer)
1021   TRY
1022   {     
1023      postProcessingGlobalAttributes();
1024   }
1025   CATCH_DUMP_ATTR
1026
1027   /*!
1028   \brief Close all the context defintion and do processing data
1029      After everything is well defined on client side, they will be processed and sent to server
1030   From the version 2.0, sever and client work no more on the same database. Moreover, client(s) will send
1031   all necessary information to server, from which each server can build its own database.
1032   Because the role of server is to write out field data on a specific netcdf file,
1033   the only information that it needs is the enabled files
1034   and the active fields (fields will be written onto active files)
1035   */
1036  void CContext::closeDefinition(void)
1037   TRY
1038   {
1039     CTimer::get("Context : close definition").resume() ;
1040     
1041     // create intercommunicator with servers.
1042     // not sure it is the good place to be called here
1043     createServerInterComm() ;
1044
1045
1046     // After xml is parsed, there are some more works with post processing
1047//     postProcessing();
1048
1049   
1050    // Make sure the calendar was correctly created
1051    if (serviceType_!=CServicesManager::CLIENT) CCalendarWrapper::get(CCalendarWrapper::GetDefName())->createCalendar();
1052    if (!calendar)
1053      ERROR("CContext::postProcessing()", << "A calendar must be defined for the context \"" << getId() << "!\"")
1054    else if (calendar->getTimeStep() == NoneDu)
1055      ERROR("CContext::postProcessing()", << "A timestep must be defined for the context \"" << getId() << "!\"")
1056    // Calendar first update to set the current date equals to the start date
1057    calendar->update(0);
1058
1059    // Résolution des héritages descendants (càd des héritages de groupes)
1060    // pour chacun des contextes.
1061    solveDescInheritance(true);
1062 
1063    // Solve inheritance for field to know if enabled or not.
1064    for (auto field : CField::getAll()) field->solveRefInheritance();
1065
1066    // Check if some axis, domains or grids are eligible to for compressed indexed output.
1067    // Warning: This must be done after solving the inheritance and before the rest of post-processing
1068    // --> later ????    checkAxisDomainsGridsEligibilityForCompressedOutput();     
1069
1070      // Check if some automatic time series should be generated
1071      // Warning: This must be done after solving the inheritance and before the rest of post-processing     
1072
1073    // The timeseries should only be prepared in client
1074    prepareTimeseries();
1075
1076    //Initialisation du vecteur 'enabledFiles' contenant la liste des fichiers à sortir.
1077    findEnabledFiles();
1078    findEnabledWriteModeFiles();
1079    findEnabledReadModeFiles();
1080    findEnabledCouplerIn();
1081    findEnabledCouplerOut();
1082    createCouplerInterCommunicator() ;
1083
1084    // Find all enabled fields of each file     
1085    vector<CField*>&& fileOutField = findAllEnabledFieldsInFileOut(this->enabledWriteModeFiles);
1086    vector<CField*>&& fileInField = findAllEnabledFieldsInFileIn(this->enabledReadModeFiles);
1087    vector<CField*>&& couplerOutField = findAllEnabledFieldsCouplerOut(this->enabledCouplerOut);
1088    vector<CField*>&& couplerInField = findAllEnabledFieldsCouplerIn(this->enabledCouplerIn);
1089    findFieldsWithReadAccess();
1090    vector<CField*>& fieldWithReadAccess = fieldsWithReadAccess_ ;
1091    vector<CField*> fieldModelIn ; // fields potentially from model
1092     
1093    // define if files are on clientSied or serverSide
1094    if (serviceType_==CServicesManager::CLIENT)
1095    {
1096      for (auto& file : enabledWriteModeFiles) file->setClientSide() ;
1097      for (auto& file : enabledReadModeFiles) file->setClientSide() ;
1098    }
1099    else
1100    {
1101      for (auto& file : enabledWriteModeFiles) file->setServerSide() ;
1102      for (auto& file : enabledReadModeFiles) file->setServerSide() ;
1103    }
1104
1105   
1106    for (auto& field : couplerInField)
1107    {
1108      field->unsetGridCompleted() ;
1109    }
1110// find all field potentially at workflow end
1111    vector<CField*> endWorkflowFields ;
1112    endWorkflowFields.reserve(fileOutField.size()+couplerOutField.size()+fieldWithReadAccess.size()) ;
1113    endWorkflowFields.insert(endWorkflowFields.end(),fileOutField.begin(), fileOutField.end()) ;
1114    endWorkflowFields.insert(endWorkflowFields.end(),couplerOutField.begin(), couplerOutField.end()) ;
1115    endWorkflowFields.insert(endWorkflowFields.end(),fieldWithReadAccess.begin(), fieldWithReadAccess.end()) ;
1116
1117    bool workflowGraphIsCompleted ;
1118
1119    bool first=true ;
1120    do
1121    {
1122      workflowGraphIsCompleted=true; 
1123      for(auto endWorkflowField : endWorkflowFields) 
1124      {
1125        workflowGraphIsCompleted &= endWorkflowField->buildWorkflowGraph(garbageCollector) ;
1126      }
1127
1128      for(auto couplerIn : enabledCouplerIn) couplerIn->assignContext() ;
1129      for(auto field : couplerInField) field->makeGridAliasForCoupling();
1130      for(auto field : couplerInField) this->sendCouplerInReady(field->getContextClient()) ;
1131   
1132
1133      // assign context to coupler out and related fields
1134      for(auto couplerOut : enabledCouplerOut) couplerOut->assignContext() ;
1135      // for now supose that all coupling out endpoint are succesfull. The difficultie is client/server buffer evaluation
1136      for(auto field : couplerOutField) 
1137      {
1138        field->computeGridIndexToFileServer() ; // same kind of index than for file server -> in future distribution may change
1139      }
1140      if (first) setClientServerBuffer(couplerOutField, true) ; // set buffer context
1141
1142      bool couplersReady ;
1143      do 
1144      {
1145        couplersReady=true ;
1146        for(auto field : couplerOutField)
1147        {
1148          bool ready = isCouplerInReady(field->getContextClient()) ; 
1149          if (ready) field->sendFieldToCouplerOut() ;
1150          couplersReady &= ready ;
1151        }
1152        if (!couplersReady) this->eventLoop() ;
1153      } while (!couplersReady) ;
1154
1155      first=false ;
1156      this->eventLoop() ;
1157    } while (!workflowGraphIsCompleted) ;
1158
1159    for( auto field : couplerInField) couplerInFields_.push_back(field) ;
1160
1161    // get all field coming potentially from model
1162    for (auto field : CField::getAll() ) if (field->getModelIn()) fieldModelIn.push_back(field) ;
1163
1164    // Distribute files between secondary servers according to the data size => assign a context to a file and then to fields
1165    if (serviceType_==CServicesManager::GATHERER) distributeFiles(this->enabledWriteModeFiles);
1166    else if (serviceType_==CServicesManager::CLIENT) for(auto file : this->enabledWriteModeFiles) file->setContextClient(client) ;
1167
1168    if (serviceType_==CServicesManager::CLIENT) for(auto file : this->enabledReadModeFiles) file->setContextClient(client) ;
1169   
1170    // workflow endpoint => sent to IO/SERVER
1171    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER)
1172    {
1173      for(auto field : fileOutField) 
1174      {
1175        field->connectToFileServer(garbageCollector) ; // connect the field to server filter
1176        field->computeGridIndexToFileServer() ; // compute grid index for transfer to the server context
1177      }
1178      setClientServerBuffer(fileOutField, true) ; // set buffer context
1179      for(auto field : fileOutField) field->sendFieldToFileServer() ;
1180    }
1181
1182    // workflow endpoint => write to file
1183    if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::OUT_SERVER)
1184    {
1185      for(auto field : fileOutField) 
1186      {
1187        field->connectToFileWriter(garbageCollector) ; // connect the field to server filter
1188      }
1189    }
1190   
1191    // workflow endpoint => Send data from server to client
1192    if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::GATHERER)
1193    {
1194      // no filter to send data from server to client => to be implemented (reading case)
1195    }
1196
1197    // workflow endpoint => sent to model on client side
1198    if (serviceType_==CServicesManager::CLIENT)
1199    {
1200      for(auto field : fieldWithReadAccess) field->connectToModelOutput(garbageCollector) ;
1201    }
1202
1203
1204    // workflow startpoint => data from model
1205    if (serviceType_==CServicesManager::CLIENT)
1206    {
1207      for(auto field : fieldModelIn) 
1208      {
1209        field->connectToModelInput(garbageCollector) ; // connect the field to server filter
1210        // grid index will be computed on the fly
1211      }
1212    }
1213   
1214    // workflow startpoint => data from client on server side
1215    if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::GATHERER || serviceType_==CServicesManager::OUT_SERVER)
1216    {
1217      for(auto field : fieldModelIn) 
1218      {
1219        field->connectToClientInput(garbageCollector) ; // connect the field to server filter
1220      }
1221    }
1222
1223   
1224    for(auto field : couplerInField) 
1225    {
1226      field->connectToCouplerIn(garbageCollector) ; // connect the field to server filter
1227    }
1228   
1229   
1230    for(auto field : couplerOutField) 
1231    {
1232      field->connectToCouplerOut(garbageCollector) ; // for now the same kind of filter that for file server
1233    }
1234
1235     // workflow startpoint => data from server on client side
1236    if (serviceType_==CServicesManager::CLIENT)
1237    {
1238      for(auto field : fileInField) 
1239      {
1240        field->connectToServerInput(garbageCollector) ; // connect tFhe field to server filter
1241        field->computeGridIndexToFileServer() ; // compute grid index for transfer to the server context
1242        field->sendFieldToInputFileServer() ;
1243      }
1244    }
1245
1246    // workflow startpoint => data read from file on server side
1247    if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::GATHERER)
1248    {
1249      // no filter for reading data from file => to be implemented
1250    }
1251   
1252    // construct slave server list
1253    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) 
1254    {
1255      for(auto field : fileOutField) slaveServers_.insert(field->getContextClient()) ; 
1256      for(auto field : fileInField)  slaveServers_.insert(field->getContextClient()) ; 
1257    }
1258
1259    for(auto& slaveServer : slaveServers_) sendCloseDefinition(slaveServer) ;
1260
1261    if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::OUT_SERVER) 
1262    {
1263      createFileHeader();
1264    }
1265
1266    if (serviceType_==CServicesManager::CLIENT) startPrefetchingOfEnabledReadModeFiles();
1267   
1268    // send signal to couplerIn context that definition phasis is done
1269
1270    for(auto& couplerInClient : couplerInClient_) sendCouplerInCloseDefinition(couplerInClient.second) ;
1271
1272    // wait until all couplerIn signal that closeDefition is done.
1273    bool ok;
1274    do
1275    {
1276      ok = true ;
1277      for(auto& couplerOutClient : couplerOutClient_) ok &= isCouplerInCloseDefinition(couplerOutClient.second) ;
1278      this->eventLoop() ; 
1279    } while (!ok) ;
1280
1281    return ;
1282
1283
1284
1285
1286
1287
1288
1289
1290
1291
1292    // For now, only read files with client and only one level server
1293    // if (hasClient && !hasServer) findEnabledReadModeFiles();     
1294
1295    // Find all enabled fields of each file     
1296    findAllEnabledFieldsInFiles(this->enabledWriteModeFiles);
1297    findAllEnabledFieldsInFiles(this->enabledReadModeFiles);
1298
1299    // For now, only read files with client and only one level server
1300    // if (hasClient && !hasServer)
1301    //   findAllEnabledFieldsInFiles(this->enabledReadModeFiles);     
1302
1303    if (serviceType_==CServicesManager::CLIENT)
1304    {
1305      initReadFiles();
1306      // Try to read attributes of fields in file then fill in corresponding grid (or domain, axis)
1307      this->readAttributesOfEnabledFieldsInReadModeFiles();
1308    }
1309
1310    // Only search and rebuild all reference objects of enable fields, don't transform
1311    this->solveOnlyRefOfEnabledFields();
1312
1313    // Search and rebuild all reference object of enabled fields, and transform
1314    this->solveAllRefOfEnabledFieldsAndTransform();
1315
1316    // Find all fields with read access from the public API
1317    if (serviceType_==CServicesManager::CLIENT) findFieldsWithReadAccess();
1318    // and solve the all reference for them
1319    if (serviceType_==CServicesManager::CLIENT) solveAllRefOfFieldsWithReadAccess();
1320
1321    isPostProcessed = true;
1322
1323
1324
1325    // Distribute files between secondary servers according to the data size
1326    distributeFiles(this->enabledWriteModeFiles);
1327
1328    // Check grid and calculate its distribution
1329    checkGridEnabledFields();
1330
1331    setClientServerBuffer(client, (serviceType_==CServicesManager::CLIENT) ) ;
1332    for (int i = 0; i < clientPrimServer.size(); ++i)
1333         setClientServerBuffer(clientPrimServer[i], true);
1334
1335   
1336    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER)
1337    { 
1338      if (serviceType_==CServicesManager::GATHERER)
1339      { 
1340        for (auto it=clientPrimServer.begin(); it!=clientPrimServer.end();++it) 
1341        {
1342          this->sendAllAttributesToServer(*it); // Send all attributes of current context to server
1343          CCalendarWrapper::get(CCalendarWrapper::GetDefName())->sendAllAttributesToServer(*it); // Send all attributes of current calendar
1344        }
1345      }
1346      else 
1347      {
1348        this->sendAllAttributesToServer(client);   // Send all attributes of current context to server
1349        CCalendarWrapper::get(CCalendarWrapper::GetDefName())->sendAllAttributesToServer(client); // Send all attributes of current calendar
1350      }
1351
1352      // We have enough information to send to server
1353      // First of all, send all enabled files
1354      sendEnabledFiles(this->enabledWriteModeFiles);
1355      // We only use server-level 1 (for now) to read data
1356      if (serviceType_==CServicesManager::CLIENT)  sendEnabledFiles(this->enabledReadModeFiles);
1357
1358      // Then, send all enabled fields     
1359      sendEnabledFieldsInFiles(this->enabledWriteModeFiles);
1360     
1361      if (serviceType_==CServicesManager::CLIENT) sendEnabledFieldsInFiles(this->enabledReadModeFiles);
1362
1363      // Then, check whether we have domain_ref, axis_ref or scalar_ref attached to the enabled fields
1364      // If any, so send them to server
1365      sendRefDomainsAxisScalars(this->enabledWriteModeFiles); 
1366     
1367      if (serviceType_==CServicesManager::CLIENT) sendRefDomainsAxisScalars(this->enabledReadModeFiles);       
1368
1369      // Check whether enabled fields have grid_ref, if any, send this info to server
1370      sendRefGrid(this->enabledFiles);
1371      // This code may be useful in the future when we want to seperate completely read and write
1372      // sendRefGrid(this->enabledWriteModeFiles);
1373      // if (!hasServer)
1374      //   sendRefGrid(this->enabledReadModeFiles);
1375     
1376      // A grid of enabled fields composed of several components which must be checked then their
1377      // checked attributes should be sent to server
1378      sendGridComponentEnabledFieldsInFiles(this->enabledFiles); // This code can be seperated in two (one for reading, another for writing)
1379
1380      // We have a xml tree on the server side and now, it should be also processed
1381      sendPostProcessing();
1382       
1383      // Finally, we send information of grid itself to server
1384      sendGridEnabledFieldsInFiles(this->enabledWriteModeFiles);       
1385     
1386      if (serviceType_==CServicesManager::CLIENT) sendGridEnabledFieldsInFiles(this->enabledReadModeFiles);       
1387    }
1388    allProcessed = true;
1389
1390
1391    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) sendPostProcessingGlobalAttributes();
1392
1393    // There are some processings that should be done after all of above. For example: check mask or index
1394    this->buildFilterGraphOfEnabledFields();
1395   
1396     if (serviceType_==CServicesManager::CLIENT)
1397    {
1398      buildFilterGraphOfFieldsWithReadAccess();
1399      postProcessFilterGraph(); // For coupling in, modify this later
1400    }
1401   
1402    checkGridEnabledFields();   
1403
1404    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) this->sendProcessingGridOfEnabledFields();
1405    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) this->sendCloseDefinition();
1406
1407    // Nettoyage de l'arborescence
1408    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) CleanTree(); // Only on client side??
1409
1410    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) sendCreateFileHeader();
1411    if (serviceType_==CServicesManager::CLIENT) startPrefetchingOfEnabledReadModeFiles();
1412   
1413    CTimer::get("Context : close definition").suspend() ;
1414  }
1415  CATCH_DUMP_ATTR
1416
1417 /*!
1418  * Send context attribute and calendar to file server, it must be done once by context file server
1419  * \param[in] client : context client to send   
1420  */ 
1421  void CContext::sendContextToFileServer(CContextClient* client)
1422  {
1423    if (sendToFileServer_done_.count(client)!=0) return ;
1424    else sendToFileServer_done_.insert(client) ;
1425   
1426    this->sendAllAttributesToServer(client); // Send all attributes of current context to server
1427    CCalendarWrapper::get(CCalendarWrapper::GetDefName())->sendAllAttributesToServer(client); // Send all attributes of current cale
1428  }
1429
1430  // ym obsolete now to be removed
1431   void CContext::closeDefinition_old(void)
1432   TRY
1433   {
1434     CTimer::get("Context : close definition").resume() ;
1435   
1436    //
1437    postProcessingGlobalAttributes();
1438
1439    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) sendPostProcessingGlobalAttributes();
1440
1441    // There are some processings that should be done after all of above. For example: check mask or index
1442    this->buildFilterGraphOfEnabledFields();
1443   
1444     if (serviceType_==CServicesManager::CLIENT)
1445    {
1446      buildFilterGraphOfFieldsWithReadAccess();
1447      postProcessFilterGraph(); // For coupling in, modify this later
1448    }
1449   
1450    checkGridEnabledFields();   
1451
1452    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) this->sendProcessingGridOfEnabledFields();
1453    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) this->sendCloseDefinition();
1454
1455    // Nettoyage de l'arborescence
1456    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) CleanTree(); // Only on client side??
1457
1458    if (serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER) sendCreateFileHeader();
1459    if (serviceType_==CServicesManager::CLIENT) startPrefetchingOfEnabledReadModeFiles();
1460   
1461    CTimer::get("Context : close definition").suspend() ;
1462   }
1463   CATCH_DUMP_ATTR
1464
1465   vector<CField*> CContext::findAllEnabledFieldsInFiles(const std::vector<CFile*>& activeFiles)
1466   TRY
1467   {
1468     vector<CField*> fields ;
1469     for (unsigned int i = 0; i < activeFiles.size(); i++)
1470     {
1471        const vector<CField*>&& field=activeFiles[i]->getEnabledFields() ;
1472        fields.insert(fields.end(),field.begin(),field.end());
1473     }
1474     return fields ;
1475   }
1476   CATCH_DUMP_ATTR
1477
1478   vector<CField*> CContext::findAllEnabledFieldsInFileOut(const std::vector<CFile*>& activeFiles)
1479   TRY
1480   {
1481     vector<CField*> fields ;
1482     for(auto file : activeFiles)
1483     {
1484        const vector<CField*>&& fieldList=file->getEnabledFields() ;
1485        for(auto field : fieldList) field->setFileOut(file) ;
1486        fields.insert(fields.end(),fieldList.begin(),fieldList.end());
1487     }
1488     return fields ;
1489   }
1490   CATCH_DUMP_ATTR
1491
1492   vector<CField*> CContext::findAllEnabledFieldsInFileIn(const std::vector<CFile*>& activeFiles)
1493   TRY
1494   {
1495     vector<CField*> fields ;
1496     for(auto file : activeFiles)
1497     {
1498        const vector<CField*>&& fieldList=file->getEnabledFields() ;
1499        for(auto field : fieldList) field->setFileIn(file) ;
1500        fields.insert(fields.end(),fieldList.begin(),fieldList.end());
1501     }
1502     return fields ;
1503   }
1504   CATCH_DUMP_ATTR
1505
1506   vector<CField*> CContext::findAllEnabledFieldsCouplerOut(const std::vector<CCouplerOut*>& activeCouplerOut)
1507   TRY
1508   {
1509     vector<CField*> fields ;
1510     for (auto couplerOut :activeCouplerOut)
1511     {
1512        const vector<CField*>&& fieldList=couplerOut->getEnabledFields() ;
1513        for(auto field : fieldList) field->setCouplerOut(couplerOut) ;
1514        fields.insert(fields.end(),fieldList.begin(),fieldList.end());
1515     }
1516     return fields ;
1517   }
1518   CATCH_DUMP_ATTR
1519
1520   vector<CField*> CContext::findAllEnabledFieldsCouplerIn(const std::vector<CCouplerIn*>& activeCouplerIn)
1521   TRY
1522   {
1523     vector<CField*> fields ;
1524     for (auto couplerIn :activeCouplerIn)
1525     {
1526        const vector<CField*>&& fieldList=couplerIn->getEnabledFields() ;
1527        for(auto field : fieldList) field->setCouplerIn(couplerIn) ;
1528        fields.insert(fields.end(),fieldList.begin(),fieldList.end());
1529     }
1530     return fields ;
1531   }
1532   CATCH_DUMP_ATTR
1533
1534
1535
1536   void CContext::readAttributesOfEnabledFieldsInReadModeFiles()
1537   TRY
1538   {
1539      for (unsigned int i = 0; i < this->enabledReadModeFiles.size(); ++i)
1540        (void)this->enabledReadModeFiles[i]->readAttributesOfEnabledFieldsInReadMode();
1541   }
1542   CATCH_DUMP_ATTR
1543
1544   void CContext::sendGridComponentEnabledFieldsInFiles(const std::vector<CFile*>& activeFiles)
1545   TRY
1546   {
1547     int size = activeFiles.size();
1548     for (int i = 0; i < size; ++i)
1549     {       
1550       activeFiles[i]->sendGridComponentOfEnabledFields();
1551     }
1552   }
1553   CATCH_DUMP_ATTR
1554
1555   /*!
1556      Send active (enabled) fields in file from a client to others
1557      \param [in] activeFiles files contains enabled fields to send
1558   */
1559   void CContext::sendGridEnabledFieldsInFiles(const std::vector<CFile*>& activeFiles)
1560   TRY
1561   {
1562     int size = activeFiles.size();
1563     for (int i = 0; i < size; ++i)
1564     {       
1565       activeFiles[i]->sendGridOfEnabledFields();
1566     }
1567   }
1568   CATCH_DUMP_ATTR
1569
1570   void CContext::checkGridEnabledFields()
1571   TRY
1572   {
1573     int size = enabledFiles.size();
1574     for (int i = 0; i < size; ++i)
1575     {
1576       enabledFiles[i]->checkGridOfEnabledFields();       
1577     }
1578
1579     size = enabledCouplerOut.size();
1580     for (int i = 0; i < size; ++i)
1581     {
1582       enabledCouplerOut[i]->checkGridOfEnabledFields();       
1583     }
1584   }
1585   CATCH_DUMP_ATTR
1586
1587   /*!
1588      Check grid of active (enabled) fields in file
1589      \param [in] activeFiles files contains enabled fields whose grid needs checking
1590   */
1591   void CContext::checkGridEnabledFieldsInFiles(const std::vector<CFile*>& activeFiles)
1592   TRY
1593   {
1594     int size = activeFiles.size();
1595     for (int i = 0; i < size; ++i)
1596     {
1597       activeFiles[i]->checkGridOfEnabledFields();       
1598     }
1599   }
1600   CATCH_DUMP_ATTR
1601
1602    /*!
1603      Go up the hierachical tree via field_ref and do check of attributes of fields
1604      This can be done in a client then all computed information will be sent from this client to others
1605      \param [in] sendToServer Flag to indicate whether calculated information will be sent
1606   */
1607   void CContext::solveOnlyRefOfEnabledFields(void)
1608   TRY
1609   {
1610     int size = this->enabledFiles.size();
1611     for (int i = 0; i < size; ++i)
1612     {
1613       this->enabledFiles[i]->solveOnlyRefOfEnabledFields();
1614     }
1615
1616     for (int i = 0; i < size; ++i)
1617     {
1618       this->enabledFiles[i]->generateNewTransformationGridDest();
1619     }
1620
1621     size = this->enabledCouplerOut.size();
1622     for (int i = 0; i < size; ++i)
1623     {
1624       this->enabledCouplerOut[i]->solveOnlyRefOfEnabledFields();
1625     }
1626
1627     for (int i = 0; i < size; ++i)
1628     {
1629       this->enabledCouplerOut[i]->generateNewTransformationGridDest();
1630     }
1631   }
1632   CATCH_DUMP_ATTR
1633
1634    /*!
1635      Go up the hierachical tree via field_ref and do check of attributes of fields.
1636      The transformation can be done in this step.
1637      All computed information will be sent from this client to others.
1638      \param [in] sendToServer Flag to indicate whether calculated information will be sent
1639   */
1640   void CContext::solveAllRefOfEnabledFieldsAndTransform(void)
1641   TRY
1642   {
1643     int size = this->enabledFiles.size();
1644     for (int i = 0; i < size; ++i)
1645     {
1646       this->enabledFiles[i]->solveAllRefOfEnabledFieldsAndTransform();
1647     }
1648
1649     size = this->enabledCouplerOut.size();
1650     for (int i = 0; i < size; ++i)
1651     {
1652       this->enabledCouplerOut[i]->solveAllRefOfEnabledFieldsAndTransform();
1653     }
1654
1655   }
1656   CATCH_DUMP_ATTR
1657
1658   void CContext::buildFilterGraphOfEnabledFields()
1659   TRY
1660   {
1661     int size = this->enabledFiles.size();
1662     for (int i = 0; i < size; ++i)
1663     {
1664       this->enabledFiles[i]->buildFilterGraphOfEnabledFields(garbageCollector);
1665     }
1666
1667     size = this->enabledCouplerOut.size();
1668     for (int i = 0; i < size; ++i)
1669     {
1670       this->enabledCouplerOut[i]->buildFilterGraphOfEnabledFields(garbageCollector);
1671     }
1672   }
1673   CATCH_DUMP_ATTR
1674
1675   void CContext::postProcessFilterGraph()
1676   TRY
1677   {
1678     int size = enabledFiles.size();
1679     for (int i = 0; i < size; ++i)
1680     {
1681        enabledFiles[i]->postProcessFilterGraph();
1682     }
1683   }
1684   CATCH_DUMP_ATTR
1685
1686   void CContext::startPrefetchingOfEnabledReadModeFiles()
1687   TRY
1688   {
1689     int size = enabledReadModeFiles.size();
1690     for (int i = 0; i < size; ++i)
1691     {
1692        enabledReadModeFiles[i]->prefetchEnabledReadModeFields();
1693     }
1694   }
1695   CATCH_DUMP_ATTR
1696
1697   void CContext::doPreTimestepOperationsForEnabledReadModeFiles()
1698   TRY
1699   {
1700     int size = enabledReadModeFiles.size();
1701     for (int i = 0; i < size; ++i)
1702     {
1703        enabledReadModeFiles[i]->doPreTimestepOperationsForEnabledReadModeFields();
1704     }
1705   }
1706   CATCH_DUMP_ATTR
1707
1708   void CContext::doPostTimestepOperationsForEnabledReadModeFiles()
1709   TRY
1710   {
1711     int size = enabledReadModeFiles.size();
1712     for (int i = 0; i < size; ++i)
1713     {
1714        enabledReadModeFiles[i]->doPostTimestepOperationsForEnabledReadModeFields();
1715     }
1716   }
1717   CATCH_DUMP_ATTR
1718
1719  void CContext::findFieldsWithReadAccess(void)
1720  TRY
1721  {
1722    fieldsWithReadAccess_.clear();
1723    const vector<CField*> allFields = CField::getAll();
1724    for (size_t i = 0; i < allFields.size(); ++i)
1725    {
1726      CField* field = allFields[i];
1727      if (!field->read_access.isEmpty() && field->read_access && (field->enabled.isEmpty() || field->enabled))
1728      {
1729        fieldsWithReadAccess_.push_back(field);
1730        field->setModelOut() ;
1731      }
1732    }
1733  }
1734  CATCH_DUMP_ATTR
1735
1736  void CContext::solveAllRefOfFieldsWithReadAccess()
1737  TRY
1738  {
1739    for (size_t i = 0; i < fieldsWithReadAccess_.size(); ++i)
1740      fieldsWithReadAccess_[i]->solveAllReferenceEnabledField(false);
1741  }
1742  CATCH_DUMP_ATTR
1743
1744  void CContext::buildFilterGraphOfFieldsWithReadAccess()
1745  TRY
1746  {
1747    for (size_t i = 0; i < fieldsWithReadAccess_.size(); ++i)
1748      fieldsWithReadAccess_[i]->buildFilterGraph(garbageCollector, true);
1749  }
1750  CATCH_DUMP_ATTR
1751
1752   void CContext::solveAllInheritance(bool apply)
1753   TRY
1754   {
1755     // Résolution des héritages descendants (càd des héritages de groupes)
1756     // pour chacun des contextes.
1757      solveDescInheritance(apply);
1758
1759     // Résolution des héritages par référence au niveau des fichiers.
1760      const vector<CFile*> allFiles=CFile::getAll();
1761      const vector<CCouplerIn*> allCouplerIn=CCouplerIn::getAll();
1762      const vector<CCouplerOut*> allCouplerOut=CCouplerOut::getAll();
1763      const vector<CGrid*> allGrids= CGrid::getAll();
1764
1765      if (serviceType_==CServicesManager::CLIENT)
1766      {
1767        for (unsigned int i = 0; i < allFiles.size(); i++)
1768          allFiles[i]->solveFieldRefInheritance(apply);
1769
1770        for (unsigned int i = 0; i < allCouplerIn.size(); i++)
1771          allCouplerIn[i]->solveFieldRefInheritance(apply);
1772
1773        for (unsigned int i = 0; i < allCouplerOut.size(); i++)
1774          allCouplerOut[i]->solveFieldRefInheritance(apply);
1775      }
1776
1777      unsigned int vecSize = allGrids.size();
1778      unsigned int i = 0;
1779      for (i = 0; i < vecSize; ++i)
1780        allGrids[i]->solveElementsRefInheritance(apply);
1781
1782   }
1783  CATCH_DUMP_ATTR
1784
1785   void CContext::findEnabledFiles(void)
1786   TRY
1787   {
1788      const std::vector<CFile*> allFiles = CFile::getAll();
1789      const CDate& initDate = calendar->getInitDate();
1790
1791      for (unsigned int i = 0; i < allFiles.size(); i++)
1792         if (!allFiles[i]->enabled.isEmpty()) // Si l'attribut 'enabled' est défini.
1793         {
1794            if (allFiles[i]->enabled.getValue()) // Si l'attribut 'enabled' est fixé à vrai.
1795            {
1796              if (allFiles[i]->output_freq.isEmpty())
1797              {
1798                 ERROR("CContext::findEnabledFiles()",
1799                     << "Mandatory attribute output_freq must be defined for file \""<<allFiles[i]->getFileOutputName()
1800                     <<" \".")
1801              }
1802              if ((initDate + allFiles[i]->output_freq.getValue()) < (initDate + this->getCalendar()->getTimeStep()))
1803              {
1804                error(0)<<"WARNING: void CContext::findEnabledFiles()"<<endl
1805                    << "Output frequency in file \""<<allFiles[i]->getFileOutputName()
1806                    <<"\" is less than the time step. File will not be written."<<endl;
1807              }
1808              else
1809               enabledFiles.push_back(allFiles[i]);
1810            }
1811         }
1812         else
1813         {
1814           if (allFiles[i]->output_freq.isEmpty())
1815           {
1816              ERROR("CContext::findEnabledFiles()",
1817                  << "Mandatory attribute output_freq must be defined for file \""<<allFiles[i]->getFileOutputName()
1818                  <<" \".")
1819           }
1820           if ( (initDate + allFiles[i]->output_freq.getValue()) < (initDate + this->getCalendar()->getTimeStep()))
1821           {
1822             error(0)<<"WARNING: void CContext::findEnabledFiles()"<<endl
1823                 << "Output frequency in file \""<<allFiles[i]->getFileOutputName()
1824                 <<"\" is less than the time step. File will not be written."<<endl;
1825           }
1826           else
1827             enabledFiles.push_back(allFiles[i]); // otherwise true by default
1828         }
1829
1830      if (enabledFiles.size() == 0)
1831         DEBUG(<<"Aucun fichier ne va être sorti dans le contexte nommé \""
1832               << getId() << "\" !");
1833
1834   }
1835   CATCH_DUMP_ATTR
1836
1837   void CContext::findEnabledCouplerIn(void)
1838   TRY
1839   {
1840      const std::vector<CCouplerIn*> allCouplerIn = CCouplerIn::getAll();
1841      bool enabled ;
1842      for (size_t i = 0; i < allCouplerIn.size(); i++)
1843      {
1844        if (allCouplerIn[i]->enabled.isEmpty()) enabled=true ;
1845        else enabled=allCouplerIn[i]->enabled ;
1846        if (enabled) enabledCouplerIn.push_back(allCouplerIn[i]) ;
1847      }
1848   }
1849   CATCH_DUMP_ATTR
1850
1851   void CContext::findEnabledCouplerOut(void)
1852   TRY
1853   {
1854      const std::vector<CCouplerOut*> allCouplerOut = CCouplerOut::getAll();
1855      bool enabled ;
1856      for (size_t i = 0; i < allCouplerOut.size(); i++)
1857      {
1858        if (allCouplerOut[i]->enabled.isEmpty()) enabled=true ;
1859        else enabled=allCouplerOut[i]->enabled ;
1860        if (enabled) enabledCouplerOut.push_back(allCouplerOut[i]) ;
1861      }
1862   }
1863   CATCH_DUMP_ATTR
1864
1865
1866
1867
1868   void CContext::distributeFiles(const vector<CFile*>& files)
1869   TRY
1870   {
1871     bool distFileMemory=false ;
1872     distFileMemory=CXios::getin<bool>("server2_dist_file_memory", distFileMemory);
1873
1874     if (distFileMemory) distributeFileOverMemoryBandwith(files) ;
1875     else distributeFileOverBandwith(files) ;
1876   }
1877   CATCH_DUMP_ATTR
1878
1879   void CContext::distributeFileOverBandwith(const vector<CFile*>& files)
1880   TRY
1881   {
1882     double eps=std::numeric_limits<double>::epsilon()*10 ;
1883     
1884     std::ofstream ofs(("distribute_file_"+getId()+".dat").c_str(), std::ofstream::out);
1885     int nbPools = clientPrimServer.size();
1886
1887     // (1) Find all enabled files in write mode
1888     // for (int i = 0; i < this->enabledFiles.size(); ++i)
1889     // {
1890     //   if (enabledFiles[i]->mode.isEmpty() || (!enabledFiles[i]->mode.isEmpty() && enabledFiles[i]->mode.getValue() == CFile::mode_attr::write ))
1891     //    enabledWriteModeFiles.push_back(enabledFiles[i]);
1892     // }
1893
1894     // (2) Estimate the data volume for each file
1895     int size = files.size();
1896     std::vector<std::pair<double, CFile*> > dataSizeMap;
1897     double dataPerPool = 0;
1898     int nfield=0 ;
1899     ofs<<size<<endl ;
1900     for (size_t i = 0; i < size; ++i)
1901     {
1902       CFile* file = files[i];
1903       ofs<<file->getId()<<endl ;
1904       StdSize dataSize=0;
1905       std::vector<CField*> enabledFields = file->getEnabledFields();
1906       size_t numEnabledFields = enabledFields.size();
1907       ofs<<numEnabledFields<<endl ;
1908       for (size_t j = 0; j < numEnabledFields; ++j)
1909       {
1910         dataSize += enabledFields[j]->getGlobalWrittenSize() ;
1911         ofs<<enabledFields[j]->getGrid()->getId()<<endl ;
1912         ofs<<enabledFields[j]->getGlobalWrittenSize()<<endl ;
1913       }
1914       double outFreqSec = (Time)(calendar->getCurrentDate()+file->output_freq)-(Time)(calendar->getCurrentDate()) ;
1915       double dataSizeSec= dataSize/ outFreqSec;
1916       ofs<<dataSizeSec<<endl ;
1917       nfield++ ;
1918// add epsilon*nField to dataSizeSec in order to  preserve reproductive ordering when sorting
1919       dataSizeMap.push_back(make_pair(dataSizeSec + dataSizeSec * eps * nfield , file));
1920       dataPerPool += dataSizeSec;
1921     }
1922     dataPerPool /= nbPools;
1923     std::sort(dataSizeMap.begin(), dataSizeMap.end());
1924
1925     // (3) Assign contextClient to each enabled file
1926
1927     std::multimap<double,int> poolDataSize ;
1928// multimap is not garanty to preserve stable sorting in c++98 but it seems it does for c++11
1929
1930     int j;
1931     double dataSize ;
1932     for (j = 0 ; j < nbPools ; ++j) poolDataSize.insert(std::pair<double,int>(0.,j)) ; 
1933             
1934     for (int i = dataSizeMap.size()-1; i >= 0; --i)
1935     {
1936       dataSize=(*poolDataSize.begin()).first ;
1937       j=(*poolDataSize.begin()).second ;
1938       dataSizeMap[i].second->setContextClient(clientPrimServer[j]);
1939       dataSize+=dataSizeMap[i].first;
1940       poolDataSize.erase(poolDataSize.begin()) ;
1941       poolDataSize.insert(std::pair<double,int>(dataSize,j)) ; 
1942     }
1943
1944     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 ;
1945   }
1946   CATCH_DUMP_ATTR
1947
1948   void CContext::distributeFileOverMemoryBandwith(const vector<CFile*>& filesList)
1949   TRY
1950   {
1951     int nbPools = clientPrimServer.size();
1952     double ratio=0.5 ;
1953     ratio=CXios::getin<double>("server2_dist_file_memory_ratio", ratio);
1954
1955     int nFiles = filesList.size();
1956     vector<SDistFile> files(nFiles);
1957     vector<SDistGrid> grids;
1958     map<string,int> gridMap ;
1959     string gridId; 
1960     int gridIndex=0 ;
1961
1962     for (size_t i = 0; i < nFiles; ++i)
1963     {
1964       StdSize dataSize=0;
1965       CFile* file = filesList[i];
1966       std::vector<CField*> enabledFields = file->getEnabledFields();
1967       size_t numEnabledFields = enabledFields.size();
1968
1969       files[i].id_=file->getId() ;
1970       files[i].nbGrids_=numEnabledFields;
1971       files[i].assignedGrid_ = new int[files[i].nbGrids_] ;
1972         
1973       for (size_t j = 0; j < numEnabledFields; ++j)
1974       {
1975         gridId=enabledFields[j]->getGrid()->getId() ;
1976         if (gridMap.find(gridId)==gridMap.end())
1977         {
1978            gridMap[gridId]=gridIndex  ;
1979            SDistGrid newGrid; 
1980            grids.push_back(newGrid) ;
1981            gridIndex++ ;
1982         }
1983         files[i].assignedGrid_[j]=gridMap[gridId] ;
1984         grids[files[i].assignedGrid_[j]].size_=enabledFields[j]->getGlobalWrittenSize() ;
1985         dataSize += enabledFields[j]->getGlobalWrittenSize() ; // usefull
1986       }
1987       double outFreqSec = (Time)(calendar->getCurrentDate()+file->output_freq)-(Time)(calendar->getCurrentDate()) ;
1988       files[i].bandwith_= dataSize/ outFreqSec ;
1989     }
1990
1991     double bandwith=0 ;
1992     double memory=0 ;
1993   
1994     for(int i=0; i<nFiles; i++)  bandwith+=files[i].bandwith_ ;
1995     for(int i=0; i<nFiles; i++)  files[i].bandwith_ = files[i].bandwith_/bandwith * ratio ;
1996
1997     for(int i=0; i<grids.size(); i++)  memory+=grids[i].size_ ;
1998     for(int i=0; i<grids.size(); i++)  grids[i].size_ = grids[i].size_ / memory * (1.0-ratio) ;
1999       
2000     distributeFileOverServer2(nbPools, grids.size(), &grids[0], nFiles, &files[0]) ;
2001
2002     vector<double> memorySize(nbPools,0.) ;
2003     vector< set<int> > serverGrids(nbPools) ;
2004     vector<double> bandwithSize(nbPools,0.) ;
2005       
2006     for (size_t i = 0; i < nFiles; ++i)
2007     {
2008       bandwithSize[files[i].assignedServer_] += files[i].bandwith_* bandwith /ratio ;
2009       for(int j=0 ; j<files[i].nbGrids_;j++)
2010       {
2011         if (serverGrids[files[i].assignedServer_].find(files[i].assignedGrid_[j]) == serverGrids[files[i].assignedServer_].end())
2012         {
2013           memorySize[files[i].assignedServer_]+= grids[files[i].assignedGrid_[j]].size_ * memory / (1.0-ratio);
2014           serverGrids[files[i].assignedServer_].insert(files[i].assignedGrid_[j]) ;
2015         }
2016       }
2017       filesList[i]->setContextClient(clientPrimServer[files[i].assignedServer_]) ;
2018       delete [] files[i].assignedGrid_ ;
2019     }
2020
2021     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 ;
2022     for (int i = 0; i < nbPools; ++i) info(100)<<"Pool server level2 "<<i<<"   assigned grid memory "<<memorySize[i]*100/1024./1024.<<" Mb"<<endl ;
2023
2024   }
2025   CATCH_DUMP_ATTR
2026
2027   /*!
2028      Find all files in write mode
2029   */
2030   void CContext::findEnabledWriteModeFiles(void)
2031   TRY
2032   {
2033     int size = this->enabledFiles.size();
2034     for (int i = 0; i < size; ++i)
2035     {
2036       if (enabledFiles[i]->mode.isEmpty() || 
2037          (!enabledFiles[i]->mode.isEmpty() && enabledFiles[i]->mode.getValue() == CFile::mode_attr::write ))
2038        enabledWriteModeFiles.push_back(enabledFiles[i]);
2039     }
2040   }
2041   CATCH_DUMP_ATTR
2042
2043   /*!
2044      Find all files in read mode
2045   */
2046   void CContext::findEnabledReadModeFiles(void)
2047   TRY
2048   {
2049     int size = this->enabledFiles.size();
2050     for (int i = 0; i < size; ++i)
2051     {
2052       if (!enabledFiles[i]->mode.isEmpty() && enabledFiles[i]->mode.getValue() == CFile::mode_attr::read)
2053        enabledReadModeFiles.push_back(enabledFiles[i]);
2054     }
2055   }
2056   CATCH_DUMP_ATTR
2057
2058   void CContext::closeAllFile(void)
2059   TRY
2060   {
2061     std::vector<CFile*>::const_iterator
2062            it = this->enabledFiles.begin(), end = this->enabledFiles.end();
2063
2064     for (; it != end; it++)
2065     {
2066       info(30)<<"Closing File : "<<(*it)->getId()<<endl;
2067       (*it)->close();
2068     }
2069   }
2070   CATCH_DUMP_ATTR
2071
2072   /*!
2073   \brief Dispatch event received from client
2074      Whenever a message is received in buffer of server, it will be processed depending on
2075   its event type. A new event type should be added in the switch list to make sure
2076   it processed on server side.
2077   \param [in] event: Received message
2078   */
2079   bool CContext::dispatchEvent(CEventServer& event)
2080   TRY
2081   {
2082
2083      if (SuperClass::dispatchEvent(event)) return true;
2084      else
2085      {
2086        switch(event.type)
2087        {
2088           case EVENT_ID_CLOSE_DEFINITION :
2089             recvCloseDefinition(event);
2090             return true;
2091             break;
2092           case EVENT_ID_UPDATE_CALENDAR:
2093             recvUpdateCalendar(event);
2094             return true;
2095             break;
2096           case EVENT_ID_CREATE_FILE_HEADER :
2097             recvCreateFileHeader(event);
2098             return true;
2099             break;
2100           case EVENT_ID_POST_PROCESS:
2101             recvPostProcessing(event);
2102             return true;
2103           case EVENT_ID_SEND_REGISTRY:
2104             recvRegistry(event);
2105             return true;
2106             break;
2107           case EVENT_ID_POST_PROCESS_GLOBAL_ATTRIBUTES:
2108             recvPostProcessingGlobalAttributes(event);
2109             return true;
2110             break;
2111           case EVENT_ID_PROCESS_GRID_ENABLED_FIELDS:
2112             recvProcessingGridOfEnabledFields(event);
2113             return true;
2114             break;
2115           case EVENT_ID_COUPLER_IN_READY:
2116             recvCouplerInReady(event);
2117             return true;
2118             break;
2119           case EVENT_ID_COUPLER_IN_CLOSE_DEFINITION:
2120             recvCouplerInCloseDefinition(event);
2121             return true;
2122             break;
2123           case EVENT_ID_COUPLER_IN_CONTEXT_FINALIZED:
2124             recvCouplerInContextFinalized(event);
2125             return true;
2126             break; 
2127           default :
2128             ERROR("bool CContext::dispatchEvent(CEventServer& event)",
2129                    <<"Unknown Event");
2130           return false;
2131         }
2132      }
2133   }
2134   CATCH
2135
2136   //! Client side: Send a message to server to make it close
2137   // ym obsolete
2138   void CContext::sendCloseDefinition(void)
2139   TRY
2140   {
2141    int nbSrvPools ;
2142    if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
2143    else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
2144    else nbSrvPools = 0 ;
2145    CContextClient* contextClientTmp ;
2146
2147    for (int i = 0; i < nbSrvPools; ++i)
2148     {
2149       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
2150       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
2151       CEventClient event(getType(),EVENT_ID_CLOSE_DEFINITION);
2152       if (contextClientTmp->isServerLeader())
2153       {
2154         CMessage msg;
2155         const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
2156         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
2157           event.push(*itRank,1,msg);
2158         contextClientTmp->sendEvent(event);
2159       }
2160       else contextClientTmp->sendEvent(event);
2161     }
2162   }
2163   CATCH_DUMP_ATTR
2164   
2165   //  ! Client side: Send a message to server to make it close
2166   void CContext::sendCloseDefinition(CContextClient* client)
2167   TRY
2168   {
2169      if (sendCloseDefinition_done_.count(client)!=0) return ;
2170      else sendCloseDefinition_done_.insert(client) ;
2171
2172      CEventClient event(getType(),EVENT_ID_CLOSE_DEFINITION);
2173      if (client->isServerLeader())
2174      {
2175        CMessage msg;
2176        for(auto rank : client->getRanksServerLeader()) event.push(rank,1,msg);
2177        client->sendEvent(event);
2178      }
2179     else client->sendEvent(event);
2180   }
2181   CATCH_DUMP_ATTR
2182
2183   //! Server side: Receive a message of client announcing a context close
2184   void CContext::recvCloseDefinition(CEventServer& event)
2185   TRY
2186   {
2187      CBufferIn* buffer=event.subEvents.begin()->buffer;
2188      getCurrent()->closeDefinition();
2189   }
2190   CATCH
2191
2192   //! Client side: Send a message to update calendar in each time step
2193   void CContext::sendUpdateCalendar(int step)
2194   TRY
2195   {
2196     CEventClient event(getType(),EVENT_ID_UPDATE_CALENDAR);
2197     for(auto client : slaveServers_) 
2198     {
2199       if (client->isServerLeader())
2200       {
2201         CMessage msg;
2202         msg<<step;
2203         for (auto& rank : client->getRanksServerLeader() ) event.push(rank,1,msg);
2204         client->sendEvent(event);
2205       }
2206       else client->sendEvent(event);
2207     }
2208   }
2209   CATCH_DUMP_ATTR
2210
2211   //! Server side: Receive a message of client annoucing calendar update
2212   void CContext::recvUpdateCalendar(CEventServer& event)
2213   TRY
2214   {
2215      CBufferIn* buffer=event.subEvents.begin()->buffer;
2216      getCurrent()->recvUpdateCalendar(*buffer);
2217   }
2218   CATCH
2219
2220   //! Server side: Receive a message of client annoucing calendar update
2221   void CContext::recvUpdateCalendar(CBufferIn& buffer)
2222   TRY
2223   {
2224      int step;
2225      buffer>>step;
2226      updateCalendar(step);
2227      if (serviceType_==CServicesManager::GATHERER)
2228      {       
2229        sendUpdateCalendar(step);
2230      }
2231   }
2232   CATCH_DUMP_ATTR
2233
2234   //! Client side: Send a message to create header part of netcdf file
2235   void CContext::sendCreateFileHeader(void)
2236   TRY
2237   {
2238     int nbSrvPools ;
2239     if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
2240     else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
2241     else nbSrvPools = 0 ;
2242     CContextClient* contextClientTmp ;
2243
2244     for (int i = 0; i < nbSrvPools; ++i)
2245     {
2246       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
2247       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
2248       CEventClient event(getType(),EVENT_ID_CREATE_FILE_HEADER);
2249
2250       if (contextClientTmp->isServerLeader())
2251       {
2252         CMessage msg;
2253         const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
2254         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
2255           event.push(*itRank,1,msg) ;
2256         contextClientTmp->sendEvent(event);
2257       }
2258       else contextClientTmp->sendEvent(event);
2259     }
2260   }
2261   CATCH_DUMP_ATTR
2262
2263   //! Server side: Receive a message of client annoucing the creation of header part of netcdf file
2264   void CContext::recvCreateFileHeader(CEventServer& event)
2265   TRY
2266   {
2267      CBufferIn* buffer=event.subEvents.begin()->buffer;
2268      getCurrent()->recvCreateFileHeader(*buffer);
2269   }
2270   CATCH
2271
2272   //! Server side: Receive a message of client annoucing the creation of header part of netcdf file
2273   void CContext::recvCreateFileHeader(CBufferIn& buffer)
2274   TRY
2275   {
2276      if (serviceType_==CServicesManager::IO_SERVER || serviceType_==CServicesManager::OUT_SERVER) 
2277        createFileHeader();
2278   }
2279   CATCH_DUMP_ATTR
2280
2281   //! Client side: Send a message to do some post processing on server
2282   void CContext::sendProcessingGridOfEnabledFields()
2283   TRY
2284   {
2285     int nbSrvPools ;
2286     if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
2287     else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
2288     else nbSrvPools = 0 ;
2289     CContextClient* contextClientTmp ;
2290
2291     for (int i = 0; i < nbSrvPools; ++i)
2292     {
2293       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
2294       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
2295
2296       CEventClient event(getType(),EVENT_ID_PROCESS_GRID_ENABLED_FIELDS);
2297
2298       if (contextClientTmp->isServerLeader())
2299       {
2300         CMessage msg;
2301         const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
2302         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
2303           event.push(*itRank,1,msg);
2304         contextClientTmp->sendEvent(event);
2305       }
2306       else contextClientTmp->sendEvent(event);
2307     }
2308   }
2309   CATCH_DUMP_ATTR
2310
2311   //! Server side: Receive a message to do some post processing
2312   void CContext::recvProcessingGridOfEnabledFields(CEventServer& event)
2313   TRY
2314   {
2315      CBufferIn* buffer=event.subEvents.begin()->buffer;
2316      // nothing to do, no call ??!!   
2317   }
2318   CATCH
2319
2320   //! Client side: Send a message to do some post processing on server
2321   void CContext::sendPostProcessing()
2322   TRY
2323   {
2324     int nbSrvPools ;
2325     if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
2326     else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
2327     else nbSrvPools = 0 ;
2328     CContextClient* contextClientTmp ;
2329
2330     for (int i = 0; i < nbSrvPools; ++i)
2331     {
2332       if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
2333       else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
2334       CEventClient event(getType(),EVENT_ID_POST_PROCESS);
2335       if (contextClientTmp->isServerLeader())
2336       {
2337         CMessage msg;
2338         const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
2339         for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
2340         event.push(*itRank,1,msg);
2341         contextClientTmp->sendEvent(event);
2342       }
2343       else contextClientTmp->sendEvent(event);
2344     }
2345   }
2346   CATCH_DUMP_ATTR
2347
2348   //! Server side: Receive a message to do some post processing
2349   void CContext::recvPostProcessing(CEventServer& event)
2350   TRY
2351   {
2352      CBufferIn* buffer=event.subEvents.begin()->buffer;
2353      getCurrent()->recvPostProcessing(*buffer);
2354   }
2355   CATCH
2356
2357   //! Server side: Receive a message to do some post processing
2358   void CContext::recvPostProcessing(CBufferIn& buffer)
2359   TRY
2360   {
2361      CCalendarWrapper::get(CCalendarWrapper::GetDefName())->createCalendar();
2362      postProcessing();
2363   }
2364   CATCH_DUMP_ATTR
2365
2366 
2367   /*!
2368   \brief Do some simple post processings after parsing xml file
2369      After the xml file (iodef.xml) is parsed, it is necessary to build all relations among
2370   created object, e.g: inhertance among fields, domain, axis. After that, all fiels as well as their parents (reference fields),
2371   which will be written out into netcdf files, are processed
2372   */
2373   void CContext::postProcessing()
2374   TRY
2375   {
2376     if (isPostProcessed) return;
2377
2378      // Make sure the calendar was correctly created
2379      if (!calendar)
2380        ERROR("CContext::postProcessing()", << "A calendar must be defined for the context \"" << getId() << "!\"")
2381      else if (calendar->getTimeStep() == NoneDu)
2382        ERROR("CContext::postProcessing()", << "A timestep must be defined for the context \"" << getId() << "!\"")
2383      // Calendar first update to set the current date equals to the start date
2384      calendar->update(0);
2385
2386      // Find all inheritance in xml structure
2387      this->solveAllInheritance();
2388
2389//      ShowTree(info(10));
2390
2391      // Check if some axis, domains or grids are eligible to for compressed indexed output.
2392      // Warning: This must be done after solving the inheritance and before the rest of post-processing
2393      checkAxisDomainsGridsEligibilityForCompressedOutput();      // only for field written on IO_SERVER service ????
2394
2395      // Check if some automatic time series should be generated
2396      // Warning: This must be done after solving the inheritance and before the rest of post-processing     
2397      prepareTimeseries();
2398
2399      //Initialisation du vecteur 'enabledFiles' contenant la liste des fichiers à sortir.
2400      findEnabledFiles();
2401      findEnabledWriteModeFiles();
2402      findEnabledReadModeFiles();
2403      findEnabledCouplerIn();
2404      findEnabledCouplerOut();
2405      createCouplerInterCommunicator() ;
2406
2407      // Find all enabled fields of each file     
2408      const vector<CField*>&& fileOutField = findAllEnabledFieldsInFiles(this->enabledWriteModeFiles);
2409      const vector<CField*>&& fileInField = findAllEnabledFieldsInFiles(this->enabledReadModeFiles);
2410      const vector<CField*>&& CouplerOutField = findAllEnabledFieldsCouplerOut(this->enabledCouplerOut);
2411      const vector<CField*>&& CouplerInField = findAllEnabledFieldsCouplerIn(this->enabledCouplerIn);
2412
2413
2414
2415      // For now, only read files with client and only one level server
2416      // if (hasClient && !hasServer) findEnabledReadModeFiles();     
2417
2418
2419      // For now, only read files with client and only one level server
2420      // if (hasClient && !hasServer)
2421      //   findAllEnabledFieldsInFiles(this->enabledReadModeFiles);     
2422
2423      if (serviceType_==CServicesManager::CLIENT)
2424      {
2425        initReadFiles();
2426        // Try to read attributes of fields in file then fill in corresponding grid (or domain, axis)
2427        this->readAttributesOfEnabledFieldsInReadModeFiles();
2428      }
2429
2430      // Only search and rebuild all reference objects of enable fields, don't transform
2431      this->solveOnlyRefOfEnabledFields();
2432
2433      // Search and rebuild all reference object of enabled fields, and transform
2434      this->solveAllRefOfEnabledFieldsAndTransform();
2435
2436      // Find all fields with read access from the public API
2437      if (serviceType_==CServicesManager::CLIENT) findFieldsWithReadAccess();
2438      // and solve the all reference for them
2439      if (serviceType_==CServicesManager::CLIENT) solveAllRefOfFieldsWithReadAccess();
2440
2441      isPostProcessed = true;
2442   }
2443   CATCH_DUMP_ATTR
2444
2445   void CContext::createCouplerInterCommunicator(void)
2446   TRY
2447   {
2448      int rank=this->getIntraCommRank() ;
2449      map<string,list<CCouplerOut*>> listCouplerOut ; 
2450      map<string,list<CCouplerIn*>> listCouplerIn ; 
2451
2452      for(auto couplerOut : enabledCouplerOut) listCouplerOut[couplerOut->getCouplingContextId()].push_back(couplerOut) ;
2453      for(auto couplerIn : enabledCouplerIn) listCouplerIn[couplerIn->getCouplingContextId()].push_back(couplerIn) ;
2454
2455      CCouplerManager* couplerManager = CXios::getCouplerManager() ;
2456      if (rank==0)
2457      {
2458        for(auto couplerOut : listCouplerOut) couplerManager->registerCoupling(this->getContextId(),couplerOut.first) ;
2459        for(auto couplerIn : listCouplerIn) couplerManager->registerCoupling(couplerIn.first,this->getContextId()) ;
2460      }
2461
2462      do
2463      {
2464        for(auto couplerOut : listCouplerOut) 
2465        {
2466          bool isNextCoupling ;
2467          if (rank==0) isNextCoupling = couplerManager->isNextCoupling(this->getContextId(),couplerOut.first) ;
2468          MPI_Bcast(&isNextCoupling,1,MPI_C_BOOL, 0, getIntraComm()) ; 
2469          if (isNextCoupling) 
2470          {
2471            addCouplingChanel(couplerOut.first, true) ;
2472            listCouplerOut.erase(couplerOut.first) ;
2473            break ;
2474          }           
2475        }
2476        for(auto couplerIn : listCouplerIn) 
2477        {
2478          bool isNextCoupling ;
2479          if (rank==0) isNextCoupling = couplerManager->isNextCoupling(couplerIn.first,this->getContextId());
2480          MPI_Bcast(&isNextCoupling,1,MPI_C_BOOL, 0, getIntraComm()) ; 
2481          if (isNextCoupling) 
2482          {
2483            addCouplingChanel(couplerIn.first, false) ;
2484            listCouplerIn.erase(couplerIn.first) ;
2485            break ;
2486          }           
2487        }
2488
2489      } while (!listCouplerOut.empty() || !listCouplerIn.empty()) ;
2490
2491   }
2492   CATCH_DUMP_ATTR
2493
2494 
2495     //! Client side: Send infomation of active files (files are enabled to write out)
2496   void CContext::sendEnabledFiles(const std::vector<CFile*>& activeFiles)
2497   TRY
2498   {
2499     int size = activeFiles.size();
2500
2501     // In a context, each type has a root definition, e.g: axis, domain, field.
2502     // Every object must be a child of one of these root definition. In this case
2503     // all new file objects created on server must be children of the root "file_definition"
2504     StdString fileDefRoot("file_definition");
2505     CFileGroup* cfgrpPtr = CFileGroup::get(fileDefRoot);
2506
2507     for (int i = 0; i < size; ++i)
2508     {
2509       CFile* f = activeFiles[i];
2510       cfgrpPtr->sendCreateChild(f->getId(),f->getContextClient());
2511       f->sendAllAttributesToServer(f->getContextClient());
2512       f->sendAddAllVariables(f->getContextClient());
2513     }
2514   }
2515   CATCH_DUMP_ATTR
2516
2517   //! Client side: Send information of active fields (ones are written onto files)
2518   void CContext::sendEnabledFieldsInFiles(const std::vector<CFile*>& activeFiles)
2519   TRY
2520   {
2521     int size = activeFiles.size();
2522     for (int i = 0; i < size; ++i)
2523     {
2524       activeFiles[i]->sendEnabledFields(activeFiles[i]->getContextClient());
2525     }
2526   }
2527   CATCH_DUMP_ATTR
2528
2529   //! Client side: Check if the defined axis, domains and grids are eligible for compressed indexed output
2530   void CContext::checkAxisDomainsGridsEligibilityForCompressedOutput()
2531   TRY
2532   {
2533     if (!(serviceType_==CServicesManager::CLIENT || serviceType_==CServicesManager::GATHERER)) return;
2534
2535     const vector<CAxis*> allAxis = CAxis::getAll();
2536     for (vector<CAxis*>::const_iterator it = allAxis.begin(); it != allAxis.end(); it++)
2537       (*it)->checkEligibilityForCompressedOutput();
2538
2539     const vector<CDomain*> allDomains = CDomain::getAll();
2540     for (vector<CDomain*>::const_iterator it = allDomains.begin(); it != allDomains.end(); it++)
2541       (*it)->checkEligibilityForCompressedOutput();
2542
2543     const vector<CGrid*> allGrids = CGrid::getAll();
2544     for (vector<CGrid*>::const_iterator it = allGrids.begin(); it != allGrids.end(); it++)
2545       (*it)->checkEligibilityForCompressedOutput();
2546   }
2547   CATCH_DUMP_ATTR
2548
2549   //! Client side: Prepare the timeseries by adding the necessary files
2550   void CContext::prepareTimeseries()
2551   TRY
2552   {
2553     const std::vector<CFile*> allFiles = CFile::getAll();
2554     for (size_t i = 0; i < allFiles.size(); i++)
2555     {
2556       CFile* file = allFiles[i];
2557
2558       std::vector<CVariable*> fileVars, fieldVars, vars = file->getAllVariables();
2559       for (size_t k = 0; k < vars.size(); k++)
2560       {
2561         CVariable* var = vars[k];
2562
2563         if (var->ts_target.isEmpty()
2564              || var->ts_target == CVariable::ts_target_attr::file || var->ts_target == CVariable::ts_target_attr::both)
2565           fileVars.push_back(var);
2566
2567         if (!var->ts_target.isEmpty()
2568              && (var->ts_target == CVariable::ts_target_attr::field || var->ts_target == CVariable::ts_target_attr::both))
2569           fieldVars.push_back(var);
2570       }
2571
2572       if (!file->timeseries.isEmpty() && file->timeseries != CFile::timeseries_attr::none)
2573       {
2574         StdString fileNameStr("%file_name%") ;
2575         StdString tsPrefix = !file->ts_prefix.isEmpty() ? file->ts_prefix : fileNameStr ;
2576         
2577         StdString fileName=file->getFileOutputName();
2578         size_t pos=tsPrefix.find(fileNameStr) ;
2579         while (pos!=std::string::npos)
2580         {
2581           tsPrefix=tsPrefix.replace(pos,fileNameStr.size(),fileName) ;
2582           pos=tsPrefix.find(fileNameStr) ;
2583         }
2584       
2585         const std::vector<CField*> allFields = file->getAllFields();
2586         for (size_t j = 0; j < allFields.size(); j++)
2587         {
2588           CField* field = allFields[j];
2589
2590           if (!field->ts_enabled.isEmpty() && field->ts_enabled)
2591           {
2592             CFile* tsFile = CFile::create();
2593             tsFile->duplicateAttributes(file);
2594
2595             // Add variables originating from file and targeted to timeserie file
2596             for (size_t k = 0; k < fileVars.size(); k++)
2597               tsFile->getVirtualVariableGroup()->addChild(fileVars[k]);
2598
2599           
2600             tsFile->name = tsPrefix + "_";
2601             if (!field->name.isEmpty())
2602               tsFile->name.get() += field->name;
2603             else if (field->hasDirectFieldReference()) // We cannot use getBaseFieldReference() just yet
2604               tsFile->name.get() += field->field_ref;
2605             else
2606               tsFile->name.get() += field->getId();
2607
2608             if (!field->ts_split_freq.isEmpty())
2609               tsFile->split_freq = field->ts_split_freq;
2610
2611             CField* tsField = tsFile->addField();
2612             tsField->field_ref = field->getId();
2613
2614             // Add variables originating from file and targeted to timeserie field
2615             for (size_t k = 0; k < fieldVars.size(); k++)
2616               tsField->getVirtualVariableGroup()->addChild(fieldVars[k]);
2617
2618             vars = field->getAllVariables();
2619             for (size_t k = 0; k < vars.size(); k++)
2620             {
2621               CVariable* var = vars[k];
2622
2623               // Add variables originating from field and targeted to timeserie field
2624               if (var->ts_target.isEmpty()
2625                    || var->ts_target == CVariable::ts_target_attr::field || var->ts_target == CVariable::ts_target_attr::both)
2626                 tsField->getVirtualVariableGroup()->addChild(var);
2627
2628               // Add variables originating from field and targeted to timeserie file
2629               if (!var->ts_target.isEmpty()
2630                    && (var->ts_target == CVariable::ts_target_attr::file || var->ts_target == CVariable::ts_target_attr::both))
2631                 tsFile->getVirtualVariableGroup()->addChild(var);
2632             }
2633
2634             tsFile->solveFieldRefInheritance(true);
2635
2636             if (file->timeseries == CFile::timeseries_attr::exclusive)
2637               field->enabled = false;
2638           }
2639         }
2640
2641         // Finally disable the original file is need be
2642         if (file->timeseries == CFile::timeseries_attr::only)
2643          file->enabled = false;
2644       }
2645     }
2646   }
2647   CATCH_DUMP_ATTR
2648
2649   //! Client side: Send information of reference grid of active fields
2650   void CContext::sendRefGrid(const std::vector<CFile*>& activeFiles)
2651   TRY
2652   {
2653     std::set<pair<StdString,CContextClient*>> gridIds;
2654
2655     int sizeFile = activeFiles.size();
2656     CFile* filePtr(NULL);
2657
2658     // Firstly, find all reference grids of all active fields
2659     for (int i = 0; i < sizeFile; ++i)
2660     {
2661       filePtr = activeFiles[i];
2662       std::vector<CField*> enabledFields = filePtr->getEnabledFields();
2663       int sizeField = enabledFields.size();
2664       for (int numField = 0; numField < sizeField; ++numField)
2665       {
2666         if (0 != enabledFields[numField]->getRelGrid())
2667           gridIds.insert(make_pair(CGrid::get(enabledFields[numField]->getRelGrid())->getId(),enabledFields[numField]->getContextClient()));
2668       }
2669     }
2670
2671     // Create all reference grids on server side
2672     StdString gridDefRoot("grid_definition");
2673     CGridGroup* gridPtr = CGridGroup::get(gridDefRoot);
2674     for (auto it = gridIds.begin(); it != gridIds.end(); ++it)
2675     {
2676       gridPtr->sendCreateChild(it->first,it->second);
2677       CGrid::get(it->first)->sendAllAttributesToServer(it->second);
2678       CGrid::get(it->first)->sendAllDomains(it->second);
2679       CGrid::get(it->first)->sendAllAxis(it->second);
2680       CGrid::get(it->first)->sendAllScalars(it->second);
2681     }
2682   }
2683   CATCH_DUMP_ATTR
2684
2685   //! Client side: Send information of reference domain, axis and scalar of active fields
2686   void CContext::sendRefDomainsAxisScalars(const std::vector<CFile*>& activeFiles)
2687   TRY
2688   {
2689     std::set<pair<StdString,CContextClient*>> domainIds, axisIds, scalarIds;
2690
2691     // Find all reference domain and axis of all active fields
2692     int numEnabledFiles = activeFiles.size();
2693     for (int i = 0; i < numEnabledFiles; ++i)
2694     {
2695       std::vector<CField*> enabledFields = activeFiles[i]->getEnabledFields();
2696       int numEnabledFields = enabledFields.size();
2697       for (int j = 0; j < numEnabledFields; ++j)
2698       {
2699         CContextClient* contextClient=enabledFields[j]->getContextClient() ;
2700         const std::vector<StdString>& prDomAxisScalarId = enabledFields[j]->getRefDomainAxisIds();
2701         if ("" != prDomAxisScalarId[0]) domainIds.insert(make_pair(prDomAxisScalarId[0],contextClient));
2702         if ("" != prDomAxisScalarId[1]) axisIds.insert(make_pair(prDomAxisScalarId[1],contextClient));
2703         if ("" != prDomAxisScalarId[2]) scalarIds.insert(make_pair(prDomAxisScalarId[2],contextClient));
2704       }
2705     }
2706
2707     // Create all reference axis on server side
2708     std::set<StdString>::iterator itDom, itAxis, itScalar;
2709     std::set<StdString>::const_iterator itE;
2710
2711     StdString scalarDefRoot("scalar_definition");
2712     CScalarGroup* scalarPtr = CScalarGroup::get(scalarDefRoot);
2713     
2714     for (auto itScalar = scalarIds.begin(); itScalar != scalarIds.end(); ++itScalar)
2715     {
2716       if (!itScalar->first.empty())
2717       {
2718         scalarPtr->sendCreateChild(itScalar->first,itScalar->second);
2719         CScalar::get(itScalar->first)->sendAllAttributesToServer(itScalar->second);
2720       }
2721     }
2722
2723     StdString axiDefRoot("axis_definition");
2724     CAxisGroup* axisPtr = CAxisGroup::get(axiDefRoot);
2725     
2726     for (auto itAxis = axisIds.begin(); itAxis != axisIds.end(); ++itAxis)
2727     {
2728       if (!itAxis->first.empty())
2729       {
2730         axisPtr->sendCreateChild(itAxis->first, itAxis->second);
2731         CAxis::get(itAxis->first)->sendAllAttributesToServer(itAxis->second);
2732       }
2733     }
2734
2735     // Create all reference domains on server side
2736     StdString domDefRoot("domain_definition");
2737     CDomainGroup* domPtr = CDomainGroup::get(domDefRoot);
2738     
2739     for (auto itDom = domainIds.begin(); itDom != domainIds.end(); ++itDom)
2740     {
2741       if (!itDom->first.empty()) {
2742          domPtr->sendCreateChild(itDom->first, itDom->second);
2743          CDomain::get(itDom->first)->sendAllAttributesToServer(itDom->second);
2744       }
2745     }
2746   }
2747   CATCH_DUMP_ATTR
2748
2749   void CContext::triggerLateFields(void)
2750   TRY
2751   {
2752    for(auto& field : fileInFields_) field->triggerLateField() ;
2753    for(auto& field : couplerInFields_) field->triggerLateField() ;
2754   }
2755   CATCH_DUMP_ATTR
2756
2757   //! Update calendar in each time step
2758   void CContext::updateCalendar(int step)
2759   TRY
2760   {
2761      int prevStep = calendar->getStep();
2762
2763      if (prevStep < step)
2764      {
2765        if (serviceType_==CServicesManager::CLIENT) // For now we only use server level 1 to read data
2766        {
2767          triggerLateFields();
2768        }
2769
2770        info(50) << "updateCalendar : before : " << calendar->getCurrentDate() << endl;
2771        calendar->update(step);
2772        info(50) << "updateCalendar : after : " << calendar->getCurrentDate() << endl;
2773  #ifdef XIOS_MEMTRACK_LIGHT
2774        info(50) << " Current memory used by XIOS : "<<  MemTrack::getCurrentMemorySize()*1.0/(1024*1024)<<" Mbyte, at timestep "<<step<<" of context "<<this->getId()<<endl ;
2775  #endif
2776
2777        if (serviceType_==CServicesManager::CLIENT) // For now we only use server level 1 to read data
2778        {
2779          doPostTimestepOperationsForEnabledReadModeFiles();
2780          garbageCollector.invalidate(calendar->getCurrentDate());
2781        }
2782      }
2783      else if (prevStep == step)
2784        info(50) << "updateCalendar: already at step " << step << ", no operation done." << endl;
2785      else // if (prevStep > step)
2786        ERROR("void CContext::updateCalendar(int step)",
2787              << "Illegal calendar update: previous step was " << prevStep << ", new step " << step << "is in the past!")
2788   }
2789   CATCH_DUMP_ATTR
2790
2791   void CContext::initReadFiles(void)
2792   TRY
2793   {
2794      vector<CFile*>::const_iterator it;
2795
2796      for (it=enabledReadModeFiles.begin(); it != enabledReadModeFiles.end(); it++)
2797      {
2798         (*it)->initRead();
2799      }
2800   }
2801   CATCH_DUMP_ATTR
2802
2803   //! Server side: Create header of netcdf file
2804   void CContext::createFileHeader(void)
2805   TRY
2806   {
2807      vector<CFile*>::const_iterator it;
2808
2809      //for (it=enabledFiles.begin(); it != enabledFiles.end(); it++)
2810      for (it=enabledWriteModeFiles.begin(); it != enabledWriteModeFiles.end(); it++)
2811      {
2812         (*it)->initWrite();
2813      }
2814   }
2815   CATCH_DUMP_ATTR
2816
2817   //! Get current context
2818   CContext* CContext::getCurrent(void)
2819   TRY
2820   {
2821     return CObjectFactory::GetObject<CContext>(CObjectFactory::GetCurrentContextId()).get();
2822   }
2823   CATCH
2824
2825   /*!
2826   \brief Set context with an id be the current context
2827   \param [in] id identity of context to be set to current
2828   */
2829   void CContext::setCurrent(const string& id)
2830   TRY
2831   {
2832     CObjectFactory::SetCurrentContextId(id);
2833     CGroupFactory::SetCurrentContextId(id);
2834   }
2835   CATCH
2836
2837  /*!
2838  \brief Create a context with specific id
2839  \param [in] id identity of new context
2840  \return pointer to the new context or already-existed one with identity id
2841  */
2842  CContext* CContext::create(const StdString& id)
2843  TRY
2844  {
2845    CContext::setCurrent(id);
2846
2847    bool hasctxt = CContext::has(id);
2848    CContext* context = CObjectFactory::CreateObject<CContext>(id).get();
2849    getRoot();
2850    if (!hasctxt) CGroupFactory::AddChild(root, context->getShared());
2851
2852#define DECLARE_NODE(Name_, name_) \
2853    C##Name_##Definition::create(C##Name_##Definition::GetDefName());
2854#define DECLARE_NODE_PAR(Name_, name_)
2855#include "node_type.conf"
2856
2857    return (context);
2858  }
2859  CATCH
2860
2861     //! Server side: Receive a message to do some post processing
2862  void CContext::recvRegistry(CEventServer& event)
2863  TRY
2864  {
2865    CBufferIn* buffer=event.subEvents.begin()->buffer;
2866    getCurrent()->recvRegistry(*buffer);
2867  }
2868  CATCH
2869
2870  void CContext::recvRegistry(CBufferIn& buffer)
2871  TRY
2872  {
2873    if (server->intraCommRank==0)
2874    {
2875      CRegistry registry(server->intraComm) ;
2876      registry.fromBuffer(buffer) ;
2877      registryOut->mergeRegistry(registry) ;
2878    }
2879  }
2880  CATCH_DUMP_ATTR
2881
2882  void CContext::sendRegistry(void)
2883  TRY
2884  {
2885    registryOut->hierarchicalGatherRegistry() ;
2886
2887    int nbSrvPools ;
2888    if (serviceType_==CServicesManager::CLIENT) nbSrvPools = 1 ;
2889    else if (serviceType_==CServicesManager::GATHERER) nbSrvPools = this->clientPrimServer.size() ;
2890    else nbSrvPools = 0 ;
2891    CContextClient* contextClientTmp ;
2892
2893    for (int i = 0; i < nbSrvPools; ++i)
2894    {
2895      if (serviceType_==CServicesManager::CLIENT) contextClientTmp = client ;
2896      else if (serviceType_==CServicesManager::GATHERER ) contextClientTmp = clientPrimServer[i] ;
2897
2898      CEventClient event(CContext::GetType(), CContext::EVENT_ID_SEND_REGISTRY);
2899      if (contextClientTmp->isServerLeader())
2900      {
2901        CMessage msg ;
2902        if (contextClientTmp->clientRank==0) msg<<*registryOut ;
2903        const std::list<int>& ranks = contextClientTmp->getRanksServerLeader();
2904        for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
2905             event.push(*itRank,1,msg);
2906        contextClientTmp->sendEvent(event);
2907      }
2908      else contextClientTmp->sendEvent(event);
2909    }
2910  }
2911  CATCH_DUMP_ATTR
2912
2913 
2914  void CContext::sendFinalizeClient(CContextClient* contextClient, const string& contextClientId)
2915  TRY
2916  {
2917    CEventClient event(getType(),EVENT_ID_CONTEXT_FINALIZE_CLIENT);
2918    if (contextClient->isServerLeader())
2919    {
2920      CMessage msg;
2921      msg<<contextClientId ;
2922      const std::list<int>& ranks = contextClient->getRanksServerLeader();
2923      for (std::list<int>::const_iterator itRank = ranks.begin(), itRankEnd = ranks.end(); itRank != itRankEnd; ++itRank)
2924           event.push(*itRank,1,msg);
2925      contextClient->sendEvent(event);
2926    }
2927    else contextClient->sendEvent(event);
2928  }
2929  CATCH_DUMP_ATTR
2930
2931 
2932  void CContext::recvFinalizeClient(CEventServer& event)
2933  TRY
2934  {
2935    CBufferIn* buffer=event.subEvents.begin()->buffer;
2936    string id;
2937    *buffer>>id;
2938    get(id)->recvFinalizeClient(*buffer);
2939  }
2940  CATCH
2941
2942  void CContext::recvFinalizeClient(CBufferIn& buffer)
2943  TRY
2944  {
2945    countChildContextFinalized_++ ;
2946  }
2947  CATCH_DUMP_ATTR
2948
2949
2950
2951
2952 //! Client side: Send a message  announcing that context can receive grid definition from coupling
2953   void CContext::sendCouplerInReady(CContextClient* client)
2954   TRY
2955   {
2956      if (sendCouplerInReady_done_.count(client)!=0) return ;
2957      else sendCouplerInReady_done_.insert(client) ;
2958
2959      CEventClient event(getType(),EVENT_ID_COUPLER_IN_READY);
2960
2961      if (client->isServerLeader())
2962      {
2963        CMessage msg;
2964        msg<<this->getId();
2965        for (auto& rank : client->getRanksServerLeader()) event.push(rank,1,msg);
2966        client->sendEvent(event);
2967      }
2968      else client->sendEvent(event);
2969   }
2970   CATCH_DUMP_ATTR
2971
2972   //! Server side: Receive a message announcing that context can send grid definition for context coupling
2973   void CContext::recvCouplerInReady(CEventServer& event)
2974   TRY
2975   {
2976      CBufferIn* buffer=event.subEvents.begin()->buffer;
2977      getCurrent()->recvCouplerInReady(*buffer);
2978   }
2979   CATCH
2980
2981   //! Server side: Receive a message announcing that context can send grid definition for context coupling
2982   void CContext::recvCouplerInReady(CBufferIn& buffer)
2983   TRY
2984   {
2985      string contextId ;
2986      buffer>>contextId;
2987      couplerInReady_.insert(getCouplerOutClient(contextId)) ;
2988   }
2989   CATCH_DUMP_ATTR
2990
2991
2992
2993
2994
2995 //! Client side: Send a message  announcing that a coupling context have done it closeDefinition, so data can be sent now.
2996   void CContext::sendCouplerInCloseDefinition(CContextClient* client)
2997   TRY
2998   {
2999      if (sendCouplerInCloseDefinition_done_.count(client)!=0) return ;
3000      else sendCouplerInCloseDefinition_done_.insert(client) ;
3001
3002      CEventClient event(getType(),EVENT_ID_COUPLER_IN_CLOSE_DEFINITION);
3003
3004      if (client->isServerLeader())
3005      {
3006        CMessage msg;
3007        msg<<this->getId();
3008        for (auto& rank : client->getRanksServerLeader()) event.push(rank,1,msg);
3009        client->sendEvent(event);
3010      }
3011      else client->sendEvent(event);
3012   }
3013   CATCH_DUMP_ATTR
3014
3015   //! Server side: Receive a message announcing that a coupling context have done it closeDefinition, so data can be sent now.
3016   void CContext::recvCouplerInCloseDefinition(CEventServer& event)
3017   TRY
3018   {
3019      CBufferIn* buffer=event.subEvents.begin()->buffer;
3020      getCurrent()->recvCouplerInCloseDefinition(*buffer);
3021   }
3022   CATCH
3023
3024   //! Server side: Receive a message announcing that a coupling context have done it closeDefinition, so data can be sent now.
3025   void CContext::recvCouplerInCloseDefinition(CBufferIn& buffer)
3026   TRY
3027   {
3028      string contextId ;
3029      buffer>>contextId;
3030      couplerInCloseDefinition_.insert(getCouplerOutClient(contextId)) ;
3031   }
3032   CATCH_DUMP_ATTR
3033
3034
3035
3036
3037//! Client side: Send a message  announcing that a coupling context have done it contextFinalize, so it can also close it own context.
3038   void CContext::sendCouplerInContextFinalized(CContextClient* client)
3039   TRY
3040   {
3041      if (sendCouplerInContextFinalized_done_.count(client)!=0) return ;
3042      else sendCouplerInContextFinalized_done_.insert(client) ;
3043
3044      CEventClient event(getType(),EVENT_ID_COUPLER_IN_CONTEXT_FINALIZED);
3045
3046      if (client->isServerLeader())
3047      {
3048        CMessage msg;
3049        msg<<this->getId();
3050        for (auto& rank : client->getRanksServerLeader()) event.push(rank,1,msg);
3051        client->sendEvent(event);
3052      }
3053      else client->sendEvent(event);
3054   }
3055   CATCH_DUMP_ATTR
3056
3057   //! Server side: Receive a message announcing that a coupling context have done it contextFinalize, so it can also close it own context.
3058   void CContext::recvCouplerInContextFinalized(CEventServer& event)
3059   TRY
3060   {
3061      CBufferIn* buffer=event.subEvents.begin()->buffer;
3062      getCurrent()->recvCouplerInContextFinalized(*buffer);
3063   }
3064   CATCH
3065
3066   //! Server side: Receive a message announcing that a coupling context have done it contextFinalize, so it can also close it own context.
3067   void CContext::recvCouplerInContextFinalized(CBufferIn& buffer)
3068   TRY
3069   {
3070      string contextId ;
3071      buffer>>contextId;
3072      couplerInContextFinalized_.insert(getCouplerOutClient(contextId)) ;
3073   }
3074   CATCH_DUMP_ATTR
3075
3076
3077
3078
3079  /*!
3080  * \fn bool CContext::isFinalized(void)
3081  * Context is finalized if it received context post finalize event.
3082  */
3083  bool CContext::isFinalized(void)
3084  TRY
3085  {
3086    return finalized;
3087  }
3088  CATCH_DUMP_ATTR
3089  ///--------------------------------------------------------------
3090  StdString CContext::dumpClassAttributes(void)
3091  {
3092    StdString str;
3093    str.append("enabled files=\"");
3094    int size = this->enabledFiles.size();
3095    for (int i = 0; i < size; ++i)
3096    {
3097      str.append(enabledFiles[i]->getId());
3098      str.append(" ");
3099    }
3100    str.append("\"");
3101    return str;
3102  }
3103
3104} // namespace xios
Note: See TracBrowser for help on using the repository browser.