www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | Private Attributes | List of all members
PeacemanBorehole Class Reference

Approximates a borehole by a sequence of Dirac Points. More...

#include <PeacemanBorehole.h>

Inheritance diagram for PeacemanBorehole:
[legend]

Public Member Functions

 PeacemanBorehole (const InputParameters &parameters)
 Creates a new PeacemanBorehole This reads the file containing the lines of the form radius x y z that defines the borehole geometry. More...
 

Protected Member Functions

virtual void addPoints ()
 Add Dirac Points to the borehole. More...
 
bool parseNextLineReals (std::ifstream &ifs, std::vector< Real > &myvec)
 reads a space-separated line of floats from ifs and puts in myvec More...
 
Real wellConstant (const RealTensorValue &perm, const RealTensorValue &rot, const Real &half_len, const Elem *ele, const Real &rad)
 Calculates Peaceman's form of the borehole well constant Z Chen, Y Zhang, Well flow models for various numerical methods, Int J Num Analysis and Modeling, 3 (2008) 375-388. More...
 

Protected Attributes

Function & _character
 If positive then the borehole acts as a sink (producion well) for porepressure > borehole pressure, and does nothing otherwise If negative then the borehole acts as a source (injection well) for porepressure < borehole pressure, and does nothing otherwise The flow rate to/from the borehole is multiplied by |character|, so usually character = +/- 1. More...
 
const Real _p_bot
 bottomhole pressure of borehole More...
 
const RealVectorValue _unit_weight
 unit weight of fluid in borehole (for calculating bottomhole pressure at each Dirac Point) More...
 
RichardsSumQuantity_total_outflow_mass
 This is used to hold the total fluid flowing into the borehole Hence, it is positive for production wells where fluid is flowing from porespace into the borehole and removed from the model. More...
 
std::vector< Real > _rs
 radii of the borehole More...
 
std::vector< Real > _xs
 x points of the borehole More...
 
std::vector< Real > _ys
 y points of the borehole More...
 
std::vector< Real > _zs
 z points of borehole More...
 
Point _bottom_point
 the bottom point of the borehole (where bottom_pressure is defined) More...
 
std::vector< Real > _half_seg_len
 0.5*(length of polyline segments between points) More...
 
std::vector< RealTensorValue > _rot_matrix
 rotation matrix used in well_constant calculation More...
 

Private Attributes

const Real _re_constant
 borehole constant More...
 
const Real _well_constant
 well constant More...
 
const Real _borehole_length
 borehole length. Note this is only used if there is only one borehole point More...
 
const RealVectorValue _borehole_direction
 borehole direction. Note this is only used if there is only one borehole point More...
 
const std::string _point_file
 File defining the geometry of the borehole. More...
 

Detailed Description

Approximates a borehole by a sequence of Dirac Points.

Definition at line 26 of file PeacemanBorehole.h.

Constructor & Destructor Documentation

◆ PeacemanBorehole()

PeacemanBorehole::PeacemanBorehole ( const InputParameters &  parameters)

Creates a new PeacemanBorehole This reads the file containing the lines of the form radius x y z that defines the borehole geometry.

It also calculates segment-lengths and rotation matrices needed for computing the borehole well constant

Definition at line 80 of file PeacemanBorehole.C.

81  : DiracKernel(parameters),
82  _re_constant(getParam<Real>("re_constant")),
83  _well_constant(getParam<Real>("well_constant")),
84  _borehole_length(getParam<Real>("borehole_length")),
85  _borehole_direction(getParam<RealVectorValue>("borehole_direction")),
86  _point_file(getParam<std::string>("point_file")),
87  _character(getFunction("character")),
88  _p_bot(getParam<Real>("bottom_pressure")),
89  _unit_weight(getParam<RealVectorValue>("unit_weight")),
91  const_cast<RichardsSumQuantity &>(getUserObject<RichardsSumQuantity>("SumQuantityUO")))
92 {
93  // zero the outflow mass
95 
96  // open file
97  std::ifstream file(_point_file.c_str());
98  if (!file.good())
99  mooseError("Error opening file '" + _point_file + "' from a Peaceman-type Borehole.");
100 
101  // construct the arrays of radius, x, y and z
102  std::vector<Real> scratch;
103  while (parseNextLineReals(file, scratch))
104  {
105  if (scratch.size() >= 2)
106  {
107  _rs.push_back(scratch[0]);
108  _xs.push_back(scratch[1]);
109  if (scratch.size() >= 3)
110  _ys.push_back(scratch[2]);
111  else
112  _ys.push_back(0.0);
113  if (scratch.size() >= 4)
114  _zs.push_back(scratch[3]);
115  else
116  _zs.push_back(0.0);
117  }
118  }
119 
120  file.close();
121 
122  const int num_pts = _zs.size();
123  _bottom_point(0) = _xs[num_pts - 1];
124  _bottom_point(1) = _ys[num_pts - 1];
125  _bottom_point(2) = _zs[num_pts - 1];
126 
127  // construct the line-segment lengths between each point
128  _half_seg_len.resize(std::max(num_pts - 1, 1));
129  for (unsigned int i = 0; i + 1 < _xs.size(); ++i)
130  {
131  _half_seg_len[i] =
132  0.5 * std::sqrt(std::pow(_xs[i + 1] - _xs[i], 2) + std::pow(_ys[i + 1] - _ys[i], 2) +
133  std::pow(_zs[i + 1] - _zs[i], 2));
134  if (_half_seg_len[i] == 0)
135  mooseError("Peaceman-type borehole has a zero-segment length at (x,y,z) = ",
136  _xs[i],
137  " ",
138  _ys[i],
139  " ",
140  _zs[i],
141  "\n");
142  }
143  if (num_pts == 1)
145 
146  // construct the rotation matrix needed to rotate the permeability
147  _rot_matrix.resize(std::max(num_pts - 1, 1));
148  for (unsigned int i = 0; i + 1 < _xs.size(); ++i)
149  {
150  const RealVectorValue v2(_xs[i + 1] - _xs[i], _ys[i + 1] - _ys[i], _zs[i + 1] - _zs[i]);
151  _rot_matrix[i] = RotationMatrix::rotVecToZ(v2);
152  }
153  if (num_pts == 1)
154  _rot_matrix[0] = RotationMatrix::rotVecToZ(_borehole_direction);
155 }
const Real _well_constant
well constant
const RealVectorValue _borehole_direction
borehole direction. Note this is only used if there is only one borehole point
const Real _re_constant
borehole constant
std::vector< Real > _rs
radii of the borehole
bool parseNextLineReals(std::ifstream &ifs, std::vector< Real > &myvec)
reads a space-separated line of floats from ifs and puts in myvec
std::vector< Real > _xs
x points of the borehole
std::vector< Real > _ys
y points of the borehole
void zero()
sets _total = 0
const RealVectorValue _unit_weight
unit weight of fluid in borehole (for calculating bottomhole pressure at each Dirac Point) ...
RichardsSumQuantity & _total_outflow_mass
This is used to hold the total fluid flowing into the borehole Hence, it is positive for production w...
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
const Real _p_bot
bottomhole pressure of borehole
const Real _borehole_length
borehole length. Note this is only used if there is only one borehole point
const std::string _point_file
File defining the geometry of the borehole.
Point _bottom_point
the bottom point of the borehole (where bottom_pressure is defined)
std::vector< RealTensorValue > _rot_matrix
rotation matrix used in well_constant calculation
std::vector< Real > _half_seg_len
0.5*(length of polyline segments between points)
std::vector< Real > _zs
z points of borehole
Function & _character
If positive then the borehole acts as a sink (producion well) for porepressure > borehole pressure...

Member Function Documentation

◆ addPoints()

void PeacemanBorehole::addPoints ( )
protectedvirtual

Add Dirac Points to the borehole.

Definition at line 180 of file PeacemanBorehole.C.

181 {
182  // This function gets called just before the DiracKernel is evaluated
183  // so this is a handy place to zero this out.
185 
186  // Add point using the unique ID "i", let the DiracKernel take
187  // care of the caching. This should be fast after the first call,
188  // as long as the points don't move around.
189  for (unsigned int i = 0; i < _zs.size(); i++)
190  addPoint(Point(_xs[i], _ys[i], _zs[i]), i);
191 }
std::vector< Real > _xs
x points of the borehole
std::vector< Real > _ys
y points of the borehole
void zero()
sets _total = 0
RichardsSumQuantity & _total_outflow_mass
This is used to hold the total fluid flowing into the borehole Hence, it is positive for production w...
std::vector< Real > _zs
z points of borehole

◆ parseNextLineReals()

bool PeacemanBorehole::parseNextLineReals ( std::ifstream &  ifs,
std::vector< Real > &  myvec 
)
protected

reads a space-separated line of floats from ifs and puts in myvec

Definition at line 158 of file PeacemanBorehole.C.

Referenced by PeacemanBorehole().

160 {
161  std::string line;
162  myvec.clear();
163  bool gotline(false);
164  if (getline(ifs, line))
165  {
166  gotline = true;
167 
168  // Harvest floats separated by whitespace
169  std::istringstream iss(line);
170  Real f;
171  while (iss >> f)
172  {
173  myvec.push_back(f);
174  }
175  }
176  return gotline;
177 }

◆ wellConstant()

Real PeacemanBorehole::wellConstant ( const RealTensorValue &  perm,
const RealTensorValue &  rot,
const Real &  half_len,
const Elem *  ele,
const Real &  rad 
)
protected

Calculates Peaceman's form of the borehole well constant Z Chen, Y Zhang, Well flow models for various numerical methods, Int J Num Analysis and Modeling, 3 (2008) 375-388.

Definition at line 194 of file PeacemanBorehole.C.

Referenced by Q2PBorehole::computeQpResidual(), RichardsBorehole::computeQpResidual(), Q2PBorehole::jac(), and RichardsBorehole::jac().

200 {
201  if (_well_constant > 0)
202  return _well_constant;
203 
204  // rot_perm has its "2" component lying along the half segment
205  // we want to determine the eigenvectors of rot(0:1, 0:1), since, when
206  // rotated back to the original frame we will determine the element
207  // lengths along these directions
208  const RealTensorValue rot_perm = (rot * perm) * rot.transpose();
209  const Real trace2D = rot_perm(0, 0) + rot_perm(1, 1);
210  const Real det2D = rot_perm(0, 0) * rot_perm(1, 1) - rot_perm(0, 1) * rot_perm(1, 0);
211  const Real sq = std::sqrt(std::max(0.25 * trace2D * trace2D - det2D,
212  0.0)); // the std::max accounts for wierdo precision loss
213  const Real eig_val1 = 0.5 * trace2D + sq;
214  const Real eig_val2 = 0.5 * trace2D - sq;
215  RealVectorValue eig_vec1, eig_vec2;
216  if (sq > std::abs(trace2D) * 1E-7) // matrix is not a multiple of the identity (1E-7 accounts for
217  // precision in a crude way)
218  {
219  if (rot_perm(1, 0) != 0)
220  {
221  eig_vec1(0) = eig_val1 - rot_perm(1, 1);
222  eig_vec1(1) = rot_perm(1, 0);
223  eig_vec2(0) = eig_val2 - rot_perm(1, 1);
224  eig_vec2(1) = rot_perm(1, 0);
225  }
226  else if (rot_perm(0, 1) != 0)
227  {
228  eig_vec1(0) = rot_perm(0, 1);
229  eig_vec1(1) = eig_val1 - rot_perm(0, 0);
230  eig_vec2(0) = rot_perm(0, 1);
231  eig_vec2(1) = eig_val2 - rot_perm(0, 0);
232  }
233  else // off diagonal terms are both zero
234  {
235  eig_vec1(0) = 1;
236  eig_vec2(1) = 1;
237  }
238  }
239  else // matrix is basically a multiple of the identity
240  {
241  eig_vec1(0) = 1;
242  eig_vec2(1) = 1;
243  }
244 
245  // finally, rotate these to original frame and normalise
246  eig_vec1 = rot.transpose() * eig_vec1;
247  eig_vec1 /= std::sqrt(eig_vec1 * eig_vec1);
248  eig_vec2 = rot.transpose() * eig_vec2;
249  eig_vec2 /= std::sqrt(eig_vec2 * eig_vec2);
250 
251  // find the "length" of the element in these directions
252  // TODO - probably better to use variance than max&min
253  Real max1 = eig_vec1 * ele->point(0);
254  Real max2 = eig_vec2 * ele->point(0);
255  Real min1 = max1;
256  Real min2 = max2;
257  Real proj;
258  for (unsigned int i = 1; i < ele->n_nodes(); i++)
259  {
260  proj = eig_vec1 * ele->point(i);
261  max1 = (max1 < proj) ? proj : max1;
262  min1 = (min1 < proj) ? min1 : proj;
263 
264  proj = eig_vec2 * ele->point(i);
265  max2 = (max2 < proj) ? proj : max2;
266  min2 = (min2 < proj) ? min2 : proj;
267  }
268  const Real ll1 = max1 - min1;
269  const Real ll2 = max2 - min2;
270 
271  const Real r0 = _re_constant * std::sqrt(std::sqrt(eig_val1 / eig_val2) * std::pow(ll2, 2) +
272  std::sqrt(eig_val2 / eig_val1) * std::pow(ll1, 2)) /
273  (std::pow(eig_val1 / eig_val2, 0.25) + std::pow(eig_val2 / eig_val1, 0.25));
274 
275  const Real effective_perm = std::sqrt(det2D);
276 
277  const Real halfPi = acos(0.0);
278 
279  if (r0 <= rad)
280  mooseError("The effective element size (about 0.2-times-true-ele-size) for an element "
281  "containing a Peaceman-type borehole must be (much) larger than the borehole radius "
282  "for the Peaceman formulation to be correct. Your element has effective size ",
283  r0,
284  " and the borehole radius is ",
285  rad,
286  "\n");
287 
288  return 4 * halfPi * effective_perm * half_len / std::log(r0 / rad);
289 }
const Real _well_constant
well constant
const Real _re_constant
borehole constant
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)

Member Data Documentation

◆ _borehole_direction

const RealVectorValue PeacemanBorehole::_borehole_direction
private

borehole direction. Note this is only used if there is only one borehole point

Definition at line 50 of file PeacemanBorehole.h.

Referenced by PeacemanBorehole().

◆ _borehole_length

const Real PeacemanBorehole::_borehole_length
private

borehole length. Note this is only used if there is only one borehole point

Definition at line 47 of file PeacemanBorehole.h.

Referenced by PeacemanBorehole().

◆ _bottom_point

Point PeacemanBorehole::_bottom_point
protected

the bottom point of the borehole (where bottom_pressure is defined)

Definition at line 95 of file PeacemanBorehole.h.

Referenced by Q2PBorehole::computeQpResidual(), RichardsBorehole::computeQpResidual(), Q2PBorehole::jac(), RichardsBorehole::jac(), and PeacemanBorehole().

◆ _character

Function& PeacemanBorehole::_character
protected

If positive then the borehole acts as a sink (producion well) for porepressure > borehole pressure, and does nothing otherwise If negative then the borehole acts as a source (injection well) for porepressure < borehole pressure, and does nothing otherwise The flow rate to/from the borehole is multiplied by |character|, so usually character = +/- 1.

Definition at line 67 of file PeacemanBorehole.h.

Referenced by RichardsBorehole::computeQpJacobian(), Q2PBorehole::computeQpResidual(), RichardsBorehole::computeQpResidual(), Q2PBorehole::jac(), and RichardsBorehole::jac().

◆ _half_seg_len

std::vector<Real> PeacemanBorehole::_half_seg_len
protected

0.5*(length of polyline segments between points)

Definition at line 98 of file PeacemanBorehole.h.

Referenced by Q2PBorehole::computeQpResidual(), RichardsBorehole::computeQpResidual(), Q2PBorehole::jac(), RichardsBorehole::jac(), and PeacemanBorehole().

◆ _p_bot

const Real PeacemanBorehole::_p_bot
protected

bottomhole pressure of borehole

Definition at line 70 of file PeacemanBorehole.h.

Referenced by Q2PBorehole::computeQpResidual(), RichardsBorehole::computeQpResidual(), Q2PBorehole::jac(), and RichardsBorehole::jac().

◆ _point_file

const std::string PeacemanBorehole::_point_file
private

File defining the geometry of the borehole.

Each row has format radius x y z and the list of such points defines a polyline that is the borehole

Definition at line 57 of file PeacemanBorehole.h.

Referenced by PeacemanBorehole().

◆ _re_constant

const Real PeacemanBorehole::_re_constant
private

borehole constant

Definition at line 41 of file PeacemanBorehole.h.

Referenced by wellConstant().

◆ _rot_matrix

std::vector<RealTensorValue> PeacemanBorehole::_rot_matrix
protected

rotation matrix used in well_constant calculation

Definition at line 101 of file PeacemanBorehole.h.

Referenced by Q2PBorehole::computeQpResidual(), RichardsBorehole::computeQpResidual(), Q2PBorehole::jac(), RichardsBorehole::jac(), and PeacemanBorehole().

◆ _rs

std::vector<Real> PeacemanBorehole::_rs
protected

◆ _total_outflow_mass

RichardsSumQuantity& PeacemanBorehole::_total_outflow_mass
protected

This is used to hold the total fluid flowing into the borehole Hence, it is positive for production wells where fluid is flowing from porespace into the borehole and removed from the model.

Definition at line 80 of file PeacemanBorehole.h.

Referenced by addPoints(), Q2PBorehole::computeQpResidual(), RichardsBorehole::computeQpResidual(), and PeacemanBorehole().

◆ _unit_weight

const RealVectorValue PeacemanBorehole::_unit_weight
protected

unit weight of fluid in borehole (for calculating bottomhole pressure at each Dirac Point)

Definition at line 73 of file PeacemanBorehole.h.

Referenced by Q2PBorehole::computeQpResidual(), RichardsBorehole::computeQpResidual(), Q2PBorehole::jac(), and RichardsBorehole::jac().

◆ _well_constant

const Real PeacemanBorehole::_well_constant
private

well constant

Definition at line 44 of file PeacemanBorehole.h.

Referenced by wellConstant().

◆ _xs

std::vector<Real> PeacemanBorehole::_xs
protected

x points of the borehole

Definition at line 86 of file PeacemanBorehole.h.

Referenced by addPoints(), and PeacemanBorehole().

◆ _ys

std::vector<Real> PeacemanBorehole::_ys
protected

y points of the borehole

Definition at line 89 of file PeacemanBorehole.h.

Referenced by addPoints(), and PeacemanBorehole().

◆ _zs

std::vector<Real> PeacemanBorehole::_zs
protected

The documentation for this class was generated from the following files: