Chaste  Build::
Vessel.cpp
1 /*
2 
3 Copyright (c) 2005-2016, University of Oxford.
4  All rights reserved.
5 
6  University of Oxford means the Chancellor, Masters and Scholars of the
7  University of Oxford, having an administrative office at Wellington
8  Square, Oxford OX1 2JD, UK.
9 
10  This file is part of Chaste.
11 
12  Redistribution and use in source and binary forms, with or without
13  modification, are permitted provided that the following conditions are met:
14  * Redistributions of source code must retain the above copyright notice,
15  this list of conditions and the following disclaimer.
16  * Redistributions in binary form must reproduce the above copyright notice,
17  this list of conditions and the following disclaimer in the documentation
18  and/or other materials provided with the distribution.
19  * Neither the name of the University of Oxford nor the names of its
20  contributors may be used to endorse or promote products derived from this
21  software without specific prior written permission.
22 
23  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
24  AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
25  IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
26  ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
27  LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
28  CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE
29  GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
30  HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
31  LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
32  OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
33 
34  */
35 
36 #include "SmartPointers.hpp"
37 #include "Exception.hpp"
38 #include "UblasIncludes.hpp"
39 #include "SimulationTime.hpp"
40 #include "Vessel.hpp"
41 
42 template<unsigned DIM>
44  mSegments(std::vector<boost::shared_ptr<VesselSegment<DIM> > >()),
45  mNodes(std::vector<boost::shared_ptr<VesselNode<DIM> > >()),
46  mNodesUpToDate(false),
47  mpFlowProperties(boost::shared_ptr<VesselFlowProperties<DIM> >(new VesselFlowProperties<DIM>()))
48 {
49  mSegments.push_back(pSegment);
50  mpFlowProperties->UpdateSegments(mSegments);
51 
52 }
53 
54 template<unsigned DIM>
55 Vessel<DIM>::Vessel(std::vector<boost::shared_ptr<VesselSegment<DIM> > > segments) : AbstractVesselNetworkComponent<DIM>(),
56  mSegments(segments),
57  mNodes(std::vector<boost::shared_ptr<VesselNode<DIM> > >()),
58  mNodesUpToDate(false),
59  mpFlowProperties(boost::shared_ptr<VesselFlowProperties<DIM> >(new VesselFlowProperties<DIM>()))
60 {
61  if (segments.size() > 1)
62  {
63  for (unsigned i = 1; i < mSegments.size(); i++)
64  {
65  if (!mSegments[i]->IsConnectedTo(mSegments[i - 1]))
66  {
67  EXCEPTION("Input vessel segments are not attached in the correct order.");
68  }
69  }
70 
71  for (unsigned i = 0; i < mSegments.size(); i++)
72  {
73  for (unsigned j = 0; j < mSegments.size(); j++)
74  {
75  if (i != j && i != j - 1 && i != j + 1)
76  {
77  if (mSegments[i]->IsConnectedTo(mSegments[j]))
78  {
79  EXCEPTION("Input vessel segments are not correctly connected.");
80  }
81  }
82  }
83  }
84  }
85 
86  mpFlowProperties->UpdateSegments(mSegments);
87 }
88 
89 template<unsigned DIM>
90 Vessel<DIM>::Vessel(std::vector<boost::shared_ptr<VesselNode<DIM> > > nodes) :
91  mSegments(std::vector<boost::shared_ptr<VesselSegment<DIM> > >()),
92  mNodes(std::vector<boost::shared_ptr<VesselNode<DIM> > >()),
93  mNodesUpToDate(false),
94  mpFlowProperties(boost::shared_ptr<VesselFlowProperties<DIM> >(new VesselFlowProperties<DIM>()))
95 {
96 
97  if (nodes.size() < 2)
98  {
99  EXCEPTION("Insufficient number of nodes to define a segment.");
100  }
101  else
102  {
103  for (unsigned i = 1; i < nodes.size(); i++)
104  {
105  mSegments.push_back(VesselSegment<DIM>::Create(nodes[i-1], nodes[i]));
106  }
107  }
108  mpFlowProperties->UpdateSegments(mSegments);
109 }
110 
111 template<unsigned DIM>
112 Vessel<DIM>::Vessel(boost::shared_ptr<VesselNode<DIM> > pStartNode, boost::shared_ptr<VesselNode<DIM> > pEndNode)
113  : mSegments(std::vector<boost::shared_ptr<VesselSegment<DIM> > >()),
114  mNodes(std::vector<boost::shared_ptr<VesselNode<DIM> > >()),
115  mNodesUpToDate(false),
116  mpFlowProperties(boost::shared_ptr<VesselFlowProperties<DIM> >(new VesselFlowProperties<DIM>()))
117 {
118  mSegments.push_back(VesselSegment<DIM>::Create(pStartNode, pEndNode));
119  mpFlowProperties->UpdateSegments(mSegments);
120 }
121 
122 template<unsigned DIM>
124 {
125 }
126 
127 template<unsigned DIM>
128 boost::shared_ptr<Vessel<DIM> > Vessel<DIM>::Create(boost::shared_ptr<VesselSegment<DIM> > pSegment)
129 {
130  boost::shared_ptr<Vessel<DIM> > pSelf(new Vessel<DIM>(pSegment));
131 
132  // Add the vessel to the segment
133  pSegment->AddVessel(pSelf->shared_from_this());
134  return pSelf;
135 }
136 
137 template<unsigned DIM>
138 boost::shared_ptr<Vessel<DIM> > Vessel<DIM>::Create(std::vector<boost::shared_ptr<VesselSegment<DIM> > > segments)
139 {
140  boost::shared_ptr<Vessel<DIM> > pSelf(new Vessel<DIM>(segments));
141 
142  // Add the vessel to the segments
143  for (unsigned i = 0; i < segments.size(); i++)
144  {
145  segments[i]->AddVessel(pSelf->shared_from_this());
146  }
147  return pSelf;
148 }
149 
150 template<unsigned DIM>
151 boost::shared_ptr<Vessel<DIM> > Vessel<DIM>::Create(std::vector<boost::shared_ptr<VesselNode<DIM> > > nodes)
152 {
153  boost::shared_ptr<Vessel<DIM> > pSelf(new Vessel<DIM>(nodes));
154 
155  std::vector<boost::shared_ptr<VesselSegment<DIM> > > segments = pSelf->GetSegments();
156 
157  // Add the vessel to the new segments
158  for (unsigned i = 0; i < segments.size(); i++)
159  {
160  segments[i]->AddVessel(pSelf->shared_from_this());
161  }
162  return pSelf;
163 }
164 
165 template<unsigned DIM>
166 boost::shared_ptr<Vessel<DIM> > Vessel<DIM>::Create(boost::shared_ptr<VesselNode<DIM> > pStartNode, boost::shared_ptr<VesselNode<DIM> > pEndNode)
167 {
168  boost::shared_ptr<Vessel<DIM> > pSelf(new Vessel<DIM>(pStartNode, pEndNode));
169 
170  std::vector<boost::shared_ptr<VesselSegment<DIM> > > segments = pSelf->GetSegments();
171 
172  // Add the vessel to the new segment
173  for (unsigned i = 0; i < segments.size(); i++)
174  {
175  segments[i]->AddVessel(pSelf->shared_from_this());
176  }
177  return pSelf;
178 }
179 
180 template<unsigned DIM>
181 void Vessel<DIM>::AddSegment(boost::shared_ptr<VesselSegment<DIM> > pSegment)
182 {
183 
184  if(mSegments.size() == 1)
185  {
186  pSegment->AddVessel(Shared());
187  if(pSegment->GetNode(0) == mSegments[0]->GetNode(0))
188  {
189  mSegments.insert(mSegments.begin(), pSegment);
190  }
191  else if(pSegment->GetNode(1) == mSegments[0]->GetNode(0))
192  {
193  mSegments.insert(mSegments.begin(), pSegment);
194  }
195  else if(pSegment->GetNode(0) == mSegments[0]->GetNode(1))
196  {
197  mSegments.push_back(pSegment);
198  }
199  else if(pSegment->GetNode(1) == mSegments[0]->GetNode(1))
200  {
201  mSegments.push_back(pSegment);
202  }
203  else
204  {
205  EXCEPTION("Input vessel segment does not coincide with any end of the vessel.");
206  }
207  }
208  else
209  {
210  if (pSegment->IsConnectedTo(mSegments.back()))
211  // Append to end of vessel
212  {
213  pSegment->AddVessel(Shared());
214  mSegments.push_back(pSegment);
215  }
216  else if (pSegment->IsConnectedTo(mSegments.front()))
217  // Insert at the start of the vessel
218  {
219  pSegment->AddVessel(Shared());
220  mSegments.insert(mSegments.begin(), pSegment);
221  }
222  else
223  {
224  EXCEPTION("Input vessel segment does not coincide with any end of the multi-segment vessel.");
225  }
226  }
227 
228  mNodesUpToDate = false;
229  mpFlowProperties->UpdateSegments(mSegments);
230 }
231 
232 template<unsigned DIM>
233 void Vessel<DIM>::AddSegments(std::vector<boost::shared_ptr<VesselSegment<DIM> > > segments)
234 {
235 
236  if (segments.front()->IsConnectedTo(mSegments.back()))
237  {
238  mSegments.insert(mSegments.end(), segments.begin(), segments.end());
239  }
240  else if (segments.back()->IsConnectedTo(mSegments.front()))
241  {
242  segments.insert(segments.end(), mSegments.begin(), mSegments.end());
243  mSegments = segments;
244  }
245  else if (segments.front()->IsConnectedTo(mSegments.front()))
246  {
247  std::reverse(segments.begin(), segments.end());
248  segments.insert(segments.end(), mSegments.begin(), mSegments.end());
249  mSegments = segments;
250  }
251  else if (segments.back()->IsConnectedTo(mSegments.back()))
252  {
253  std::reverse(segments.begin(), segments.end());
254  mSegments.insert(mSegments.end(), segments.begin(), segments.end());
255  }
256  else
257  {
258  EXCEPTION("Input vessel segments do not coincide with any end of the vessel.");
259  }
260 
261  for (unsigned i = 1; i < mSegments.size(); i++)
262  {
263  if (!mSegments[i]->IsConnectedTo(mSegments[i - 1]))
264  {
265  EXCEPTION("Input vessel segments are not attached in the correct order.");
266  }
267  }
268 
269  for (unsigned i = 0; i < mSegments.size(); i++)
270  {
271  for (unsigned j = 0; j < mSegments.size(); j++)
272  {
273  if (i != j && i != j - 1 && i != j + 1)
274  {
275  if (mSegments[i]->IsConnectedTo(mSegments[j]))
276  {
277  EXCEPTION("Input vessel segments are not correctly connected.");
278  }
279  }
280  }
281  }
282 
283  // Add the vessel to the segments
284  for (unsigned i = 0; i < segments.size(); i++)
285  {
286  segments[i]->AddVessel(Shared());
287  }
288 
289  mNodesUpToDate = false;
290  mpFlowProperties->UpdateSegments(mSegments);
291 }
292 
293 template<unsigned DIM>
294 void Vessel<DIM>::CopyDataFromExistingVessel(boost::shared_ptr<Vessel<DIM> > pTargetVessel)
295 {
296  this->mOutputData = pTargetVessel->GetOutputData();
297 }
298 
299 template<unsigned DIM>
300 boost::shared_ptr<VesselNode<DIM> > Vessel<DIM>::DivideSegment(const DimensionalChastePoint<DIM>& location, double distanceTolerance)
301 {
302  // Identify segment
303  boost::shared_ptr<VesselSegment<DIM> > pVesselSegment;
304  for (unsigned i = 0; i < mSegments.size(); i++)
305  {
306  if (mSegments[i]->GetDistance(location)/mSegments[i]->GetNode(0)->GetReferenceLengthScale() <= distanceTolerance)
307  {
308  pVesselSegment = mSegments[i];
309  if (pVesselSegment->GetNode(0)->IsCoincident(location))
310  {
311  return pVesselSegment->GetNode(0);
312  }
313  if (pVesselSegment->GetNode(1)->IsCoincident(location))
314  {
315  return pVesselSegment->GetNode(1);
316  }
317 
318  }
319  }
320 
321  if (!pVesselSegment)
322  {
323  EXCEPTION("Specified location is not on a segment in this vessel.");
324  }
325  for (unsigned i = 0; i < mSegments.size(); i++)
326  {
327  for (unsigned j = 0; j < mSegments.size(); j++)
328  {
329  if (i != j && i != j - 1 && i != j + 1)
330  {
331  if (mSegments[i]->IsConnectedTo(mSegments[j]))
332  {
333  EXCEPTION("Input vessel segments are not correctly connected.");
334  }
335  }
336  }
337  }
338 
339  // The node's data is averaged from the original segments's nodes
340  // Get the closest node
341  units::quantity<unit::length> distance0 = pVesselSegment->GetNode(0)->GetDistance(location);
342  units::quantity<unit::length> distance1 = pVesselSegment->GetNode(1)->GetDistance(location);
343  unsigned closest_index;
344 
345  if (distance0 <= distance1)
346  {
347  closest_index = 0;
348  }
349  else
350  {
351  closest_index = 1;
352  }
353 
354  // Make a copy of the closest node
355  boost::shared_ptr<VesselNode<DIM> > p_new_node = VesselNode<DIM>::Create(*pVesselSegment->GetNode(closest_index));
356  p_new_node->SetLocation(location);
357  p_new_node->GetFlowProperties()->SetIsInputNode(false);
358  p_new_node->GetFlowProperties()->SetIsOutputNode(false);
359 
360  // Make two new segments
361  boost::shared_ptr<VesselSegment<DIM> > p_new_segment0 = VesselSegment<DIM>::Create(pVesselSegment->GetNode(0), p_new_node);
362  boost::shared_ptr<VesselSegment<DIM> > p_new_segment1 = VesselSegment<DIM>::Create(p_new_node, pVesselSegment->GetNode(1));
363  p_new_segment0->CopyDataFromExistingSegment(pVesselSegment);
364  p_new_segment1->CopyDataFromExistingSegment(pVesselSegment);
365 
366  // Add new segments, ensuring they are correctly ordered
367  std::vector<boost::shared_ptr<VesselSegment<DIM> > > newSegments;
368  typename std::vector<boost::shared_ptr<VesselSegment<DIM> > >::iterator it = std::find(mSegments.begin(), mSegments.end(), pVesselSegment);
369 
370  if (it == mSegments.end())
371  {
372  EXCEPTION("Vessel segment is not contained inside vessel.");
373  }
374 
375  if (mSegments.size() == 1)
376  {
377  newSegments.push_back(p_new_segment0);
378  newSegments.push_back(p_new_segment1);
379  }
380  else if(it == mSegments.begin())
381  {
382  if((*(it + 1))->IsConnectedTo(p_new_segment1))
383  {
384  newSegments.push_back(p_new_segment0);
385  newSegments.push_back(p_new_segment1);
386  }
387  else
388  {
389  newSegments.push_back(p_new_segment1);
390  newSegments.push_back(p_new_segment0);
391  }
392  }
393  else if ((*(it - 1))->IsConnectedTo(p_new_segment0))
394  {
395  newSegments.push_back(p_new_segment0);
396  newSegments.push_back(p_new_segment1);
397  }
398  else
399  {
400  newSegments.push_back(p_new_segment1);
401  newSegments.push_back(p_new_segment0);
402  }
403 
404  mSegments.insert(it, newSegments.begin(), newSegments.end());
405 
406  // Remove old segment
407  it = std::find(mSegments.begin(), mSegments.end(), pVesselSegment);
408 
409  if (it != mSegments.end())
410  {
411  mSegments.erase(it);
412  pVesselSegment->Remove();
413  }
414  else
415  {
416  EXCEPTION("Vessel segment is not contained inside vessel.");
417  }
418 
419  for (unsigned i = 0; i < mSegments.size(); i++)
420  {
421  for (unsigned j = 0; j < mSegments.size(); j++)
422  {
423  if (i != j && i != j - 1 && i != j + 1)
424  {
425  if (mSegments[i]->IsConnectedTo(mSegments[j]))
426  {
427  EXCEPTION("Input vessel segments are not correctly connected.");
428  }
429  }
430  }
431  }
432 
433  // Add the vessel to the segments
434  for (unsigned i = 0; i < newSegments.size(); i++)
435  {
436  newSegments[i]->AddVessel(Shared());
437  }
438 
439  mNodesUpToDate = false;
440  mpFlowProperties->UpdateSegments(mSegments);
441  return p_new_node;
442 }
443 
444 template<unsigned DIM>
445 boost::shared_ptr<VesselFlowProperties<DIM> > Vessel<DIM>::GetFlowProperties() const
446 {
447  return this->mpFlowProperties;
448 }
449 
450 template<unsigned DIM>
451 std::map<std::string, double> Vessel<DIM>::GetOutputData()
452 {
453  std::map<std::string, double> flow_data = this->mpFlowProperties->GetOutputData();
454  this->mOutputData.clear();
455  this->mOutputData.insert(flow_data.begin(), flow_data.end());
456  this->mOutputData["Vessel Id"] = double(this->GetId());
457  this->mOutputData["Vessel Radius m"] = this->GetRadius() / unit::metres;
458  return this->mOutputData;
459 }
460 
461 template<unsigned DIM>
462 units::quantity<unit::length> Vessel<DIM>::GetClosestEndNodeDistance(const DimensionalChastePoint<DIM>& rLocation)
463 {
464  units::quantity<unit::length> distance_1 = this->GetStartNode()->GetDistance(rLocation);
465  units::quantity<unit::length> distance_2 = this->GetEndNode()->GetDistance(rLocation);
466  if(distance_1 > distance_2)
467  {
468  return distance_2;
469  }
470  else
471  {
472  return distance_1;
473  }
474 }
475 
476 template<unsigned DIM>
477 units::quantity<unit::length> Vessel<DIM>::GetDistance(const DimensionalChastePoint<DIM>& rLocation) const
478 {
479  // Get the distance to the nearest segment in the vessel
480  units::quantity<unit::length> nearest_distance = DBL_MAX * unit::metres;
481  for(unsigned idx=0; idx<mSegments.size(); idx++)
482  {
483  units::quantity<unit::length> seg_distance = mSegments[idx]->GetDistance(rLocation);
484  if(seg_distance < nearest_distance)
485  {
486  nearest_distance = seg_distance;
487  }
488  }
489  return nearest_distance;
490 }
491 
492 template<unsigned DIM>
493 std::vector<boost::shared_ptr<Vessel<DIM> > > Vessel<DIM>::GetConnectedVessels()
494 {
495  std::vector<boost::shared_ptr<VesselSegment<DIM> > > start_segments = GetStartNode()->GetSegments();
496  std::vector<boost::shared_ptr<VesselSegment<DIM> > > end_segments = GetStartNode()->GetSegments();
497 
498  std::vector<boost::shared_ptr<Vessel<DIM> > > connected;
499 
500  for(unsigned idx=0; idx<start_segments.size();idx++)
501  {
502  if(start_segments[idx]->GetVessel() != Shared())
503  {
504  connected.push_back(start_segments[idx]->GetVessel());
505  }
506  }
507  for(unsigned idx=0; idx<end_segments.size();idx++)
508  {
509  if(end_segments[idx]->GetVessel() != Shared())
510  {
511  connected.push_back(end_segments[idx]->GetVessel());
512  }
513  }
514  return connected;
515 }
516 
517 template<unsigned DIM>
518 boost::shared_ptr<VesselNode<DIM> > Vessel<DIM>::GetEndNode()
519 {
520  if (!mNodesUpToDate)
521  {
522  UpdateNodes();
523  }
524 
525  return mNodes.back();
526 }
527 
528 template<unsigned DIM>
529 boost::shared_ptr<VesselNode<DIM> > Vessel<DIM>::GetNodeAtOppositeEnd(
530  boost::shared_ptr<VesselNode<DIM> > pQueryNode)
531 {
532  if (!mNodesUpToDate)
533  {
534  UpdateNodes();
535  }
536 
537  if (pQueryNode == GetStartNode())
538  {
539  return GetEndNode();
540  }
541  else if (pQueryNode == GetEndNode())
542  {
543  return GetStartNode();
544  }
545  else
546  {
547  EXCEPTION("Query node is not at either end of the vessel.");
548  }
549 }
550 
551 template<unsigned DIM>
552 units::quantity<unit::length> Vessel<DIM>::GetLength() const
553 {
554  units::quantity<unit::length> length = 0.0 * unit::metres;
555  for (unsigned i = 0; i < mSegments.size(); i++)
556  {
557  length += mSegments[i]->GetLength();
558  }
559  return length;
560 }
561 
562 template<unsigned DIM>
563 units::quantity<unit::length> Vessel<DIM>::GetRadius() const
564 {
565  units::quantity<unit::length> radius = 0.0 * unit::metres;
566  for (unsigned i = 0; i < mSegments.size(); i++)
567  {
568  radius += mSegments[i]->GetRadius();
569  }
570  return radius / (double(mSegments.size()));
571 }
572 
573 template<unsigned DIM>
574 boost::shared_ptr<VesselNode<DIM> > Vessel<DIM>::GetNode(unsigned index)
575 {
576  if (!mNodesUpToDate)
577  {
578  UpdateNodes();
579  }
580  if(index >= mNodes.size())
581  {
582  EXCEPTION("Out of bounds node index requested");
583  }
584  return mNodes[index];
585 }
586 
587 template<unsigned DIM>
588 std::vector<boost::shared_ptr<VesselNode<DIM> > > Vessel<DIM>::GetNodes()
589 {
590  if (!mNodesUpToDate)
591  {
592  UpdateNodes();
593  }
594  return mNodes;
595 }
596 
597 template<unsigned DIM>
598 const std::vector<boost::shared_ptr<VesselNode<DIM> > >& Vessel<DIM>::rGetNodes()
599 {
600  if (!mNodesUpToDate)
601  {
602  UpdateNodes();
603  }
604  return mNodes;
605 }
606 
607 template<unsigned DIM>
609 {
610  if (!mNodesUpToDate)
611  {
612  UpdateNodes();
613  }
614  return GetNumberOfSegments() + 1;
615 }
616 
617 template<unsigned DIM>
619 {
620  return mSegments.size();
621 }
622 
623 template<unsigned DIM>
624 boost::shared_ptr<VesselSegment<DIM> > Vessel<DIM>::GetSegment(unsigned index)
625 {
626  if(index >= mSegments.size())
627  {
628  EXCEPTION("Requested segment index out of range");
629  }
630  return mSegments[index];
631 }
632 
633 template<unsigned DIM>
634 std::vector<boost::shared_ptr<VesselSegment<DIM> > > Vessel<DIM>::GetSegments()
635 {
636  return mSegments;
637 }
638 
639 template<unsigned DIM>
640 boost::shared_ptr<VesselNode<DIM> > Vessel<DIM>::GetStartNode()
641 {
642  if (!mNodesUpToDate)
643  {
644  UpdateNodes();
645  }
646  return mNodes.front();
647 }
648 
649 template<unsigned DIM>
650 bool Vessel<DIM>::IsConnectedTo(boost::shared_ptr<Vessel<DIM> > pOtherVessel)
651 {
652  if (GetStartNode() == pOtherVessel->GetStartNode() || GetEndNode() == pOtherVessel->GetStartNode()
653  || GetStartNode() == pOtherVessel->GetEndNode() || GetEndNode() == pOtherVessel->GetEndNode())
654  {
655  return true;
656  }
657  else
658  {
659  return false;
660  }
661 }
662 
663 template<unsigned DIM>
665 {
666  // Detach all segments from their nodes
667  for (unsigned idx = 0; idx < mSegments.size(); idx++)
668  {
669  mSegments[idx]->Remove();
670  }
671  mNodesUpToDate = false;
672  mSegments = std::vector<boost::shared_ptr<VesselSegment<DIM> > >();
673 }
674 
675 template<unsigned DIM>
677 {
678  if (mSegments.size() == 1)
679  {
680  EXCEPTION("Vessel must have at least one segment.");
681  }
682  if (location == SegmentLocation::Start)
683  {
684  mSegments.front()->RemoveVessel();
685  mSegments.erase(mSegments.begin());
686  }
687  else if (location == SegmentLocation::End)
688  {
689  mSegments.back()->RemoveVessel();
690  mSegments.pop_back();
691  }
692  else
693  {
694  EXCEPTION("You can only remove segments from the start or end of vessels.");
695  }
696  mNodesUpToDate = false;
697  mpFlowProperties->UpdateSegments(mSegments);
698 }
699 
700 template<unsigned DIM>
702 {
703  this->mpFlowProperties = boost::shared_ptr<VesselFlowProperties<DIM> >(new VesselFlowProperties<DIM> (rFlowProperties));
704  this->mpFlowProperties->UpdateSegments(mSegments);
705 }
706 
707 template<unsigned DIM>
708 void Vessel<DIM>::SetRadius(units::quantity<unit::length> radius)
709 {
710  for (unsigned i = 0; i < mSegments.size(); i++)
711  {
712  mSegments[i]->SetRadius(radius);
713  }
714 }
715 
716 template<unsigned DIM>
717 boost::shared_ptr<Vessel<DIM> > Vessel<DIM>::Shared()
718 {
719  return this->shared_from_this();
720 }
721 
722 template<unsigned DIM>
724 {
725  mNodes.clear();
726 
727  if (mSegments.size() == 1)
728  {
729  mNodes.push_back(mSegments[0]->GetNode(0));
730  mNodes.push_back(mSegments[0]->GetNode(1));
731  }
732 
733  // Add the start and end nodes of the first segment and then
734  // the end nodes of every other segment
735  else
736  {
737  if (mSegments[1]->HasNode(mSegments[0]->GetNode(1)))
738  {
739  mNodes.push_back(mSegments[0]->GetNode(0));
740  mNodes.push_back(mSegments[0]->GetNode(1));
741  }
742  else if (mSegments[1]->HasNode(mSegments[0]->GetNode(1)))
743  {
744  mNodes.push_back(mSegments[0]->GetNode(1));
745  mNodes.push_back(mSegments[0]->GetNode(0));
746  }
747 
748  for (unsigned idx = 1; idx < mSegments.size(); idx++)
749  {
750  if (mNodes[idx] == mSegments[idx]->GetNode(0))
751  {
752  mNodes.push_back(mSegments[idx]->GetNode(1));
753  }
754  else
755  {
756  mNodes.push_back(mSegments[idx]->GetNode(0));
757  }
758  }
759  }
760  mNodesUpToDate = true;
761 }
762 
763 // Explicit instantiation
764 template class Vessel<2> ;
765 template class Vessel<3> ;
std::vector< boost::shared_ptr< VesselSegment< DIM > > > mSegments
Definition: Vessel.hpp:76
static boost::shared_ptr< VesselNode< DIM > > Create(double v1=0.0, double v2=0.0, double v3=0.0)
Definition: VesselNode.cpp:100
units::quantity< unit::length > GetLength() const
Definition: Vessel.cpp:552
void CopyDataFromExistingVessel(boost::shared_ptr< Vessel< DIM > > pTargetVessel)
Definition: Vessel.cpp:294
void SetRadius(units::quantity< unit::length > radius)
Definition: Vessel.cpp:708
std::vector< boost::shared_ptr< VesselSegment< DIM > > > GetSegments()
Definition: Vessel.cpp:634
boost::shared_ptr< VesselNode< DIM > > DivideSegment(const DimensionalChastePoint< DIM > &rLocation, double distanceTolerance=1.e-6)
Definition: Vessel.cpp:300
~Vessel()
Definition: Vessel.cpp:123
std::vector< boost::shared_ptr< VesselNode< DIM > > > mNodes
Definition: Vessel.hpp:81
void UpdateNodes()
Definition: Vessel.cpp:723
units::quantity< unit::length > GetDistance(const DimensionalChastePoint< DIM > &rLocation) const
Definition: Vessel.cpp:477
bool mNodesUpToDate
Definition: Vessel.hpp:86
void Remove()
Definition: Vessel.cpp:664
Vessel(boost::shared_ptr< VesselSegment< DIM > > pSegment)
Definition: Vessel.cpp:43
bool IsConnectedTo(boost::shared_ptr< Vessel< DIM > > pOtherVessel)
Definition: Vessel.cpp:650
std::map< std::string, double > mOutputData
void RemoveSegments(SegmentLocation::Value location)
Definition: Vessel.cpp:676
unsigned GetNumberOfNodes()
Definition: Vessel.cpp:608
const std::vector< boost::shared_ptr< VesselNode< DIM > > > & rGetNodes()
Definition: Vessel.cpp:598
unsigned GetNumberOfSegments()
Definition: Vessel.cpp:618
std::map< std::string, double > GetOutputData()
Definition: Vessel.cpp:451
boost::shared_ptr< VesselNode< DIM > > GetNodeAtOppositeEnd(boost::shared_ptr< VesselNode< DIM > > pQueryNode)
Definition: Vessel.cpp:529
boost::shared_ptr< VesselSegment< DIM > > GetSegment(unsigned index)
Definition: Vessel.cpp:624
boost::shared_ptr< VesselFlowProperties< DIM > > mpFlowProperties
Definition: Vessel.hpp:91
boost::shared_ptr< VesselNode< DIM > > GetEndNode()
Definition: Vessel.cpp:518
static boost::shared_ptr< VesselSegment< DIM > > Create(boost::shared_ptr< VesselNode< DIM > > pNode1, boost::shared_ptr< VesselNode< DIM > > pNode2)
boost::shared_ptr< Vessel< DIM > > Shared()
Definition: Vessel.cpp:717
static boost::shared_ptr< Vessel< DIM > > Create(boost::shared_ptr< VesselSegment< DIM > > pSegment)
Definition: Vessel.cpp:128
boost::shared_ptr< VesselNode< DIM > > GetNode(unsigned index)
Definition: Vessel.cpp:574
std::vector< boost::shared_ptr< Vessel< DIM > > > GetConnectedVessels()
Definition: Vessel.cpp:493
std::vector< boost::shared_ptr< VesselNode< DIM > > > GetNodes()
Definition: Vessel.cpp:588
void SetFlowProperties(const VesselFlowProperties< DIM > &rFlowProperties)
Definition: Vessel.cpp:701
boost::shared_ptr< VesselNode< DIM > > GetStartNode()
Definition: Vessel.cpp:640
units::quantity< unit::length > GetRadius() const
Definition: Vessel.cpp:563
boost::shared_ptr< VesselFlowProperties< DIM > > GetFlowProperties() const
Definition: Vessel.cpp:445
void AddSegments(std::vector< boost::shared_ptr< VesselSegment< DIM > > > pSegments)
Definition: Vessel.cpp:233
void AddSegment(boost::shared_ptr< VesselSegment< DIM > > pSegment)
Definition: Vessel.cpp:181
void CopyDataFromExistingSegment(const boost::shared_ptr< VesselSegment< DIM > > pTargetSegment)
units::quantity< unit::length > GetClosestEndNodeDistance(const DimensionalChastePoint< DIM > &rLocation)
Definition: Vessel.cpp:462