parsed_function.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 #ifndef LIBMESH_PARSED_FUNCTION_H
19 #define LIBMESH_PARSED_FUNCTION_H
20 
21 #include "libmesh/libmesh_config.h"
22 #include "libmesh/function_base.h"
23 
24 #ifdef LIBMESH_HAVE_FPARSER
25 
26 // Local includes
27 #include "libmesh/dense_vector.h"
28 #include "libmesh/point.h"
29 
30 // FParser includes
31 #include "fparser.hh"
32 
33 // C++ includes
34 #include <algorithm> // std::find
35 #include <cmath>
36 #include <cmath>
37 #include <cstddef>
38 #include <string>
39 
40 namespace libMesh {
41 
42 template <typename Output=Number>
43 class ParsedFunction : public FunctionBase<Output>
44 {
45 public:
46  explicit
47  ParsedFunction (const std::string& expression, const std::vector<std::string>* additional_vars=NULL,
48  const std::vector<Output>* initial_vals=NULL)
49  : _expression(expression)
50  // Size the spacetime vector to account for space, time, and any additional
51  // variables passed
52  //_spacetime(LIBMESH_DIM+1 + (additional_vars ? additional_vars->size() : 0)),
53  {
54  std::string variables = "x";
55 #if LIBMESH_DIM > 1
56  variables += ",y";
57 #endif
58 #if LIBMESH_DIM > 2
59  variables += ",z";
60 #endif
61  variables += ",t";
62 
63  _spacetime.resize(LIBMESH_DIM+1 + (additional_vars ? additional_vars->size() : 0));
64 
65  // If additional vars were passed, append them to the string
66  // that we send to the function parser. Also add them to the
67  // end of our spacetime vector
68  if (additional_vars)
69  {
70  if (initial_vals)
71  std::copy(initial_vals->begin(), initial_vals->end(), std::back_inserter(_initial_vals));
72 
73  std::copy(additional_vars->begin(), additional_vars->end(), std::back_inserter(_additional_vars));
74 
75  for (unsigned int i=0; i < additional_vars->size(); ++i)
76  {
77  variables += "," + (*additional_vars)[i];
78  // Initialize extra variables to the vector passed in or zero
79  // Note: The initial_vals vector can be shorter than the additional_vars vector
80  _spacetime[LIBMESH_DIM+1 + i] = (initial_vals && i < initial_vals->size()) ? (*initial_vals)[i] : 0;
81  }
82  }
83 
84  size_t nextstart = 0, end = 0;
85 
86  while (end != std::string::npos)
87  {
88  // If we're past the end of the string, we can't make any more
89  // subparsers
90  if (nextstart >= expression.size())
91  break;
92 
93  // If we're at the start of a brace delimited section, then we
94  // parse just that section:
95  if (expression[nextstart] == '{')
96  {
97  nextstart++;
98  end = expression.find('}', nextstart);
99  }
100  // otherwise we parse the whole thing
101  else
102  end = std::string::npos;
103 
104  // We either want the whole end of the string (end == npos) or
105  // a substring in the middle.
106  std::string subexpression =
107  expression.substr(nextstart, (end == std::string::npos) ?
108  std::string::npos : end - nextstart);
109 
110  // fparser can crash on empty expressions
111  libmesh_assert(!subexpression.empty());
112 
113  // Parse (and optimize if possible) the subexpression.
114  // Add some basic constants, to Real precision.
115  FunctionParserBase<Output> fp;
116  fp.AddConstant("pi", std::acos(Real(-1)));
117  fp.AddConstant("e", std::exp(Real(1)));
118  fp.Parse(subexpression, variables);
119  fp.Optimize();
120  parsers.push_back(fp);
121 
122  // If at end, use nextstart=maxSize. Else start at next
123  // character.
124  nextstart = (end == std::string::npos) ?
125  std::string::npos : end + 1;
126  }
127 
128  this->_initialized = true;
129  }
130 
131  virtual Output operator() (const Point& p,
132  const Real time = 0)
133  {
134  _spacetime[0] = p(0);
135 #if LIBMESH_DIM > 1
136  _spacetime[1] = p(1);
137 #endif
138 #if LIBMESH_DIM > 2
139  _spacetime[2] = p(2);
140 #endif
141  _spacetime[LIBMESH_DIM] = time;
142 
143  // The remaining locations in _spacetime are currently fixed at construction
144  // but could potentially be made dynamic
145  return parsers[0].Eval(&_spacetime[0]);
146  }
147 
148  virtual void operator() (const Point& p,
149  const Real time,
150  DenseVector<Output>& output)
151  {
152  _spacetime[0] = p(0);
153 #if LIBMESH_DIM > 1
154  _spacetime[1] = p(1);
155 #endif
156 #if LIBMESH_DIM > 2
157  _spacetime[2] = p(2);
158 #endif
159  _spacetime[LIBMESH_DIM] = time;
160 
161  unsigned int size = output.size();
162 
163  libmesh_assert_equal_to (size, parsers.size());
164 
165  // The remaining locations in _spacetime are currently fixed at construction
166  // but could potentially be made dynamic
167  for (unsigned int i=0; i != size; ++i)
168  output(i) = parsers[i].Eval(&_spacetime[0]);
169  }
170 
175  virtual Output component (unsigned int i,
176  const Point& p,
177  Real time)
178  {
179  _spacetime[0] = p(0);
180 #if LIBMESH_DIM > 1
181  _spacetime[1] = p(1);
182 #endif
183 #if LIBMESH_DIM > 2
184  _spacetime[2] = p(2);
185 #endif
186  _spacetime[LIBMESH_DIM] = time;
187 
188  libmesh_assert_less (i, parsers.size());
189 
190  // The remaining locations in _spacetime are currently fixed at construction
191  // but could potentially be made dynamic
192  return parsers[i].Eval(&_spacetime[0]);
193  }
194 
198  virtual Output & getVarAddress(const std::string & variable_name)
199  {
200  const std::vector<std::string>::iterator it =
201  std::find(_additional_vars.begin(), _additional_vars.end(), variable_name);
202 
203  if (it == _additional_vars.end())
204  {
205  libMesh::err << "ERROR: Requested variable not found in parsed function\n" << std::endl;
206  libmesh_error();
207  }
208 
209  // Iterator Arithmetic (How far from the end of the array is our target address?)
210  return _spacetime[_spacetime.size() - (_additional_vars.end() - it)];
211  }
212 
213 
217  }
218 
219 private:
220  std::string _expression;
221  std::vector<FunctionParserBase<Output> > parsers;
222  std::vector<Output> _spacetime;
223 
224  // Additional variables/values that can be parsed and handled by the function parser
225  std::vector<std::string> _additional_vars;
226  std::vector<Output> _initial_vals;
227 };
228 
229 
230 } // namespace libMesh
231 
232 
233 #else // LIBMESH_HAVE_FPARSER
234 
235 
236 namespace libMesh {
237 
238 
239 template <typename Output>
240 class ParsedFunction : public FunctionBase<Output>
241 {
242 public:
243  ParsedFunction (std::string /* expression */) : _dummy(0)
244  {
245  libmesh_not_implemented();
246  }
247 
248  virtual Output operator() (const Point&,
249  const Real /* time */ = 0)
250  { return 0.; }
251 
252  virtual void operator() (const Point&,
253  const Real /* time */,
254  DenseVector<Output>& /* output */) {}
255 
256  virtual void init() {}
257  virtual void clear() {}
258  virtual Output & getVarAddress(const std::string & /*variable_name*/) { return _dummy; }
261  (new ParsedFunction<Output>(""));
262  }
263 private:
264  Output _dummy;
265 };
266 
267 
268 } // namespace libMesh
269 
270 
271 #endif // LIBMESH_HAVE_FPARSER
272 
273 #endif // LIBMESH_PARSED_FUNCTION_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:06 UTC

Hosted By:
SourceForge.net Logo