~esys-p-dev/esys-particle/trunk

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
/////////////////////////////////////////////////////////////
//                                                         //
// Copyright (c) 2003-2014 by The University of Queensland //
// Centre for Geoscience Computing                         //
// http://earth.uq.edu.au/centre-geoscience-computing      //
//                                                         //
// Primary Business: Brisbane, Queensland, Australia       //
// Licensed under the Open Software License version 3.0    //
// http://www.opensource.org/licenses/osl-3.0.php          //
//                                                         //
/////////////////////////////////////////////////////////////

template<typename P,typename I>
ParallelInteractionStorage_EB<P,I>::ParallelInteractionStorage_EB(AParallelParticleArray* ppa,const typename I::ParameterType& param):ParallelInteractionStorage_E<P,I>(ppa,param)
{
  m_unbreakable=false;
}

/*!
  Update interactions. Check for broken interactions and remove them.
*/
template<typename P,typename InteractionType>
bool ParallelInteractionStorage_EB<P,InteractionType>::update()
{
  bool res=false;

  if(!m_unbreakable){
    console.XDebug() << "PIS_E::updating\n"; 
    typename list<InteractionType>::iterator iter = this->m_interactions.begin();
    while (iter != this->m_interactions.end()){
      if(iter->broken()){
	res=true;
	typename list<InteractionType>::iterator er_iter=iter;
	// get IDs to remove from set
	vector<int> pids=iter->getAllID();
	this->m_set.erase(make_pair(pids[0],pids[1]));
	iter++;
	// remove interaction
	this->m_interactions.erase(er_iter);
      } else {
	iter++;
      }
    }
  } else {
    console.XDebug() << "PIS_E::not updating\n"; 
  }

  return res;
}

/*! 
   save snapshot (i.e. viz/postprocess) data
*/
template<typename P,typename InteractionType>
void ParallelInteractionStorage_EB<P,InteractionType>::saveSnapShotData(std::ostream &oStream)
{
  const std::string delim = "\n";
  typedef typename InteractionType::CheckPointable CheckPointable;

  typename ParallelInteractionStorage_E<P,InteractionType>::InteractionIterator it =
    this->getInnerInteractionIterator();
  oStream << InteractionType::getType() << delim;
  oStream << it.getNumRemaining();
  if (it.hasNext()) {
    oStream << delim;
    CheckPointable(it.next()).saveCheckPointData(oStream);
    while (it.hasNext())
    {
      oStream << delim;
      CheckPointable(it.next()).saveCheckPointData(oStream);
    }
  }
}

/*! 
   save checkpoint (i.e. restart) data
*/
template<typename P,typename InteractionType>
void ParallelInteractionStorage_EB<P,InteractionType>::saveCheckPointData(std::ostream &oStream)
{
  const std::string delim = "\n";
  //  typedef typename InteractionType::CheckPointable CheckPointable;

  typename ParallelInteractionStorage_E<P,InteractionType>::InteractionIterator it =
    this->getInnerInteractionIterator();
  oStream << InteractionType::getType() << delim;
  oStream << it.getNumRemaining();
  if (it.hasNext()) {
    oStream << delim;
    it.next().saveRestartData(oStream);
    while (it.hasNext())
    {
      oStream << delim;
      it.next().saveRestartData(oStream);
    }
  }
}

/*!
  Read interaction data from input stream pointing to a restartable checkpoint file.
  The stream needs to be already positioned at the right place.
 
  \param iStream the input stream
  \warning return type may change to bool at some stage
*/
template<typename P,typename InteractionType>
void ParallelInteractionStorage_EB<P,InteractionType>::loadCheckPointData(std::istream &iStream)
{
  // read interaction type from stream
  std::string cp_interaction_type;
  iStream >> cp_interaction_type;
  // compare interaction type in stream with type of this IG
  // in not equal, signal error 
  if(cp_interaction_type!=InteractionType::getType()){
    std::cerr << "interaction types differ between checkpoint " 
	      << cp_interaction_type << " and scipt "
	      << InteractionType::getType() << std::endl;
  } else { // correct type -> read data
    // read nr. of bonds in IG
    int nconn;
    iStream >> nconn;
    std::cerr << "reading " << nconn << "  " << InteractionType::getType() << " interactions " << std::endl;

    ParallelParticleArray<P>* t_ppa=(ParallelParticleArray<P>*)(this->m_ppa);
    
    // -- read bonds
    for(int i=0;i<nconn;i++){
      InteractionType new_bond;
      // read a bond
      new_bond.loadRestartData(iStream);
      // set particle pointers
      vector<int> pids=new_bond.getAllID();
      P* ptr1=t_ppa->getParticlePtrByIndex(pids[0]);
      P* ptr2=t_ppa->getParticlePtrByIndex(pids[1]);
      if((ptr1!=NULL) && (ptr2!=NULL)){
	new_bond.setPP(ptr1,ptr2);
      } else {
	std::cerr << "trying to insert bond: particles with Id " << pids[0] << " , " << pids[1] << "not present!" << std::endl;
      }
      // insert it into interaction storage
      this->tryInsert(new_bond);
    }
  }
}

template<typename P,typename InteractionType>
void ParallelInteractionStorage_EB<P,InteractionType>::calcHeatTrans()
{
  console.Debug()
    << "calculating " << this->m_interactions.size()
    << " heat interaction transfers\n" ;

  for(
    typename list<InteractionType>::iterator it = this->m_interactions.begin();
    it != this->m_interactions.end();
    it++
  ){
    it->calcHeatTrans();
  }
}

/*!
  set the interactions "unbreakable" -> turns update into a NO-OP

  \param b true -> unbreakable, false -> breakable
*/
template<typename P,typename InteractionType>
void ParallelInteractionStorage_EB<P,InteractionType>::setUnbreakable(bool b)
{
  m_unbreakable=b;
}