1 |
|
2 |
/******************************************************* |
3 |
* |
4 |
* Copyright (c) 2003-2008 by University of Queensland |
5 |
* Earth Systems Science Computational Center (ESSCC) |
6 |
* http://www.uq.edu.au/esscc |
7 |
* |
8 |
* Primary Business: Queensland, Australia |
9 |
* Licensed under the Open Software License version 3.0 |
10 |
* http://www.opensource.org/licenses/osl-3.0.php |
11 |
* |
12 |
*******************************************************/ |
13 |
|
14 |
|
15 |
#include "AbstractTransportProblem.h" |
16 |
#include "TransportProblemException.h" |
17 |
#include "DataTypes.h" |
18 |
#include "Data.h" |
19 |
#include <iostream> |
20 |
|
21 |
|
22 |
namespace escript { |
23 |
|
24 |
AbstractTransportProblem::AbstractTransportProblem() { |
25 |
m_empty=1; |
26 |
} |
27 |
|
28 |
AbstractTransportProblem::AbstractTransportProblem(const double theta, const int blocksize, |
29 |
const FunctionSpace& functionspace) |
30 |
:m_functionspace(functionspace) |
31 |
{ |
32 |
if (blocksize<=0) |
33 |
throw TransportProblemException("Error - negative block size of transport problem."); |
34 |
if ((theta<0.) || (theta>1.)) |
35 |
throw TransportProblemException("Error - theta needs to be between 0. and 1.."); |
36 |
|
37 |
m_empty=0; |
38 |
m_blocksize=blocksize; |
39 |
// m_functionspace=functionspace; |
40 |
m_theta=theta; |
41 |
} |
42 |
|
43 |
AbstractTransportProblem::~AbstractTransportProblem() { |
44 |
} |
45 |
|
46 |
int AbstractTransportProblem::isEmpty() const { |
47 |
return m_empty; |
48 |
} |
49 |
|
50 |
|
51 |
Data AbstractTransportProblem::solve(Data& source, const double dt, const boost::python::dict& options) const |
52 |
{ |
53 |
if (isEmpty()) |
54 |
throw TransportProblemException("Error - transport problem is empty."); |
55 |
if (dt<=0.) |
56 |
throw TransportProblemException("Error - dt needs to be positive."); |
57 |
if (source.getFunctionSpace()!=getFunctionSpace()) |
58 |
throw TransportProblemException("Error - function space of transport problem and function space of source do not match."); |
59 |
if (source.getDataPointSize()!=getBlockSize()) |
60 |
throw TransportProblemException("Error - block size of transport problem and source do not match."); |
61 |
DataTypes::ShapeType shape; |
62 |
if (getBlockSize()>1) shape.push_back(getBlockSize()); |
63 |
Data out=Data(0.,shape,getFunctionSpace(),true); |
64 |
setToSolution(out,source,dt,options); |
65 |
return out; |
66 |
} |
67 |
|
68 |
void AbstractTransportProblem::setInitialValue(Data& u) const |
69 |
{ |
70 |
if (isEmpty()) |
71 |
throw TransportProblemException("Error - transport problem is empty."); |
72 |
if (u.isEmpty()) |
73 |
throw TransportProblemException("Error - empty initial value."); |
74 |
|
75 |
if (((getBlockSize()==1) && (u.getDataPointRank()>0)) || (u.getDataPointRank()>1)) |
76 |
throw TransportProblemException("Error - illegal rank of initial value."); |
77 |
|
78 |
if (u.getDataPointSize()!=getBlockSize()) |
79 |
throw TransportProblemException("Error - block size of transport problem and initial value do not match."); |
80 |
|
81 |
Data u2=Data(u,getFunctionSpace()); |
82 |
copyInitialValue(u2); |
83 |
} |
84 |
void AbstractTransportProblem::insertConstraint(Data& source, Data& q, Data& r, const double factor) const |
85 |
{ |
86 |
if (factor <= 0.) |
87 |
throw TransportProblemException("Error - weighting factor must be positive."); |
88 |
source.expand(); |
89 |
if (isEmpty()) |
90 |
throw TransportProblemException("Error - transport problem is empty."); |
91 |
if (q.isEmpty()) { |
92 |
return; |
93 |
} |
94 |
if (((getBlockSize()==1) && (q.getDataPointRank()>0)) || (q.getDataPointRank()>1)) |
95 |
throw TransportProblemException("Error - illegal rank of constraint location."); |
96 |
if (q.getDataPointSize()!=getBlockSize()) |
97 |
throw TransportProblemException("Error - block size of transport problem and constraint location don't match."); |
98 |
Data q2=Data(q,getFunctionSpace()); |
99 |
|
100 |
if (r.isEmpty()) { |
101 |
Data r2=Data(0.,q.getDataPointShape(),getFunctionSpace()); |
102 |
copyConstraint(source,q2,r2, factor); |
103 |
} else { |
104 |
if (((getBlockSize()==1) && (r.getDataPointRank()>0)) || (r.getDataPointRank()>1)) |
105 |
throw TransportProblemException("Error - illegal rank of constraint value."); |
106 |
if (r.getDataPointSize()!=getBlockSize()) |
107 |
throw TransportProblemException("Error - block size of transport problem and constraint value don't match."); |
108 |
Data r2=Data(r,getFunctionSpace()); |
109 |
copyConstraint(source,q2,r2, factor); |
110 |
} |
111 |
} |
112 |
|
113 |
void AbstractTransportProblem::copyConstraint(Data& source, Data& q, Data& r, const double factor) const |
114 |
{ |
115 |
throw TransportProblemException("Error - copyConstraint is not available"); |
116 |
} |
117 |
void AbstractTransportProblem::copyInitialValue(Data& u) const |
118 |
{ |
119 |
throw TransportProblemException("Error - copyInitialValue is not available"); |
120 |
} |
121 |
void AbstractTransportProblem::setToSolution(Data& out,Data& source,const double dt, const boost::python::dict& options) const |
122 |
{ |
123 |
throw TransportProblemException("Error - setToSolution is not available"); |
124 |
} |
125 |
void AbstractTransportProblem::resetTransport() const |
126 |
{ |
127 |
throw TransportProblemException("Error - resetProblem is not implemented."); |
128 |
} |
129 |
double AbstractTransportProblem::getSafeTimeStepSize() const |
130 |
{ |
131 |
throw TransportProblemException("Error - getSafeTimeStepSize is not implemented."); |
132 |
} |
133 |
double AbstractTransportProblem::getUnlimitedTimeStepSize() const |
134 |
{ |
135 |
throw TransportProblemException("Error - getUnlimitedTimeStepSize is not implemented."); |
136 |
} |
137 |
|
138 |
} // end of namespace |