SpatialOps
FieldFunctions.h
1 /*
2  * Copyright (c) 2014-2017 The University of Utah
3  *
4  * Permission is hereby granted, free of charge, to any person obtaining a copy
5  * of this software and associated documentation files (the "Software"), to
6  * deal in the Software without restriction, including without limitation the
7  * rights to use, copy, modify, merge, publish, distribute, sublicense, and/or
8  * sell copies of the Software, and to permit persons to whom the Software is
9  * furnished to do so, subject to the following conditions:
10  *
11  * The above copyright notice and this permission notice shall be included in
12  * all copies or substantial portions of the Software.
13  *
14  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
15  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
16  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
17  * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
18  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
19  * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS
20  * IN THE SOFTWARE.
21  */
22 
23 #ifndef FieldFunctions_h
24 #define FieldFunctions_h
25 
26 #include <spatialops/SpatialOpsConfigure.h>
27 
28 #include <boost/static_assert.hpp>
29 
31 #include <spatialops/Nebo.h>
32 
33 #include<cmath>
34 
35 
36 //====================================================================
37 
38 namespace SpatialOps{
39 
40 
41  //==================================================================
42 
43  namespace FFLocal{
44  struct NULLPatch
45  {
46  struct FieldID{ FieldID(){} };
47  };
48  }
49 
50  //==================================================================
51 
52 
81 template<typename FieldT, typename PatchT=FFLocal::NULLPatch>
83 {
84  typedef typename PatchT::FieldID FieldID;
85 public:
86 
87  typedef FieldT FieldType;
88  typedef PatchT PatchType;
89 
90  virtual void evaluate( FieldT& phi ) const =0;
91 
92  virtual void dx( FieldT& gradPhi ) const{ assert(0); }
93  virtual void dy( FieldT& gradPhi ) const{ assert(0); }
94  virtual void dz( FieldT& gradPhi ) const{ assert(0); }
95 
96  virtual void d2x( FieldT& gradPhi ) const{ assert(0); }
97  virtual void d2y( FieldT& gradPhi ) const{ assert(0); }
98  virtual void d2z( FieldT& gradPhi ) const{ assert(0); }
99 
100 
101  const FieldT& get_x() const{set_fields(); return *x_;}
102  const FieldT& get_y() const{set_fields(); return *y_;}
103  const FieldT& get_z() const{set_fields(); return *z_;}
104 
105  virtual ~FieldFunction3D(){}
106 
107 protected:
108 
110  FieldFunction3D( const FieldT& x, const FieldT& y, const FieldT& z );
111 
113  FieldFunction3D( PatchT& p, const FieldID xid, const FieldID yid, const FieldID zid );
114 
120  void set_fields() const;
121 
122 private:
123  PatchT* const patch_;
124  const FieldID xid_, yid_, zid_;
125  mutable const FieldT* x_;
126  mutable const FieldT* y_;
127  mutable const FieldT* z_;
128 };
129 
130 //====================================================================
131 
143 template<typename FieldT, typename PatchT=FFLocal::NULLPatch>
145 {
146  typedef typename PatchT::FieldID FieldID;
147 public:
148 
149  typedef FieldT FieldType;
150  typedef PatchT PatchType;
151 
152  virtual void evaluate( FieldT& phi ) const =0;
153 
154 
155  virtual void dx( FieldT& gradPhi ) const{ const bool implemented=false; assert(implemented); }
156  virtual void dy( FieldT& gradPhi ) const{ const bool implemented=false; assert(implemented); }
157 
158  virtual void d2x( FieldT& gradPhi ) const{ const bool implemented=false; assert(implemented); }
159  virtual void d2y( FieldT& gradPhi ) const{ const bool implemented=false; assert(implemented); }
160 
161 
162  const FieldT& get_x() const{set_fields(); return *x_;}
163  const FieldT& get_y() const{set_fields(); return *y_;}
164 
165  virtual ~FieldFunction2D(){}
166 
167 protected:
168 
169  FieldFunction2D( const FieldT& x, const FieldT& y );
170  FieldFunction2D( PatchT& p, const FieldID xid, const FieldID yid, const FieldID zid );
171 
172  void set_fields() const;
173 
174 private:
175  PatchT* const patch_;
176  const FieldID xid_, yid_;
177  mutable const FieldT* x_;
178  mutable const FieldT* y_;
179 };
180 
181 //====================================================================
182 
194 template< typename FieldT, typename PatchT=FFLocal::NULLPatch >
196 {
197  typedef typename PatchT::FieldID FieldID;
198 public:
199 
200  virtual void evaluate( FieldT& f ) const =0;
201 
202  virtual void dx ( FieldT& gradPhi ) const{ const bool implemented=false; assert(implemented); }
203  virtual void d2x( FieldT& gradPhi ) const{ const bool implemented=false; assert(implemented); }
204 
205  const FieldT& get_x() const{set_fields(); return *x_;}
206 
207  virtual ~FieldFunction1D(){}
208 
209 protected:
210 
211  FieldFunction1D( const FieldT& x );
212  FieldFunction1D( PatchT& p, const FieldID xid );
213 
214  void set_fields() const;
215 
216 private:
217  PatchT* const patch_;
218  const FieldID xid_;
219  mutable const FieldT* x_;
220 };
221 
222 
223 //====================================================================
224 
225 
229 template< typename FieldT,
230  typename PatchT=FFLocal::NULLPatch >
231 class LinearFunction1D : public FieldFunction1D<FieldT,PatchT>
232 {
233  typedef typename PatchT::FieldID FieldID;
234 public:
235 
236  LinearFunction1D( const FieldT& x,
237  const double slope,
238  const double intercept )
240  m_(slope),
241  b_(intercept)
242  {}
243 
244  LinearFunction1D( PatchT& p,
245  const FieldID xid,
246  const double slope,
247  const double intercept )
248  : FieldFunction1D<FieldT,PatchT>(p,xid), m_(slope), b_(intercept)
249  {}
250 
251  ~LinearFunction1D(){}
252 
253  void evaluate( FieldT& f ) const;
254  void dx( FieldT& gradPhi ) const;
255  void d2x( FieldT& gradPhi ) const;
256 
257 protected:
258 
259 private:
260  const double m_, b_;
261 };
262 
263 
264 //====================================================================
265 
274 template< typename FieldT,
275  typename PatchT=FFLocal::NULLPatch >
276 class SinFunction : public FieldFunction1D<FieldT,PatchT>
277 {
278  typedef typename PatchT::FieldID FieldID;
279  public:
280  SinFunction( PatchT& p,
281  const FieldID xid,
282  const double amplitude,
283  const double period,
284  const double base=0.0 );
285  SinFunction( const FieldT& x,
286  const double amplitude,
287  const double period,
288  const double base=0.0 );
289  ~SinFunction(){}
290  void evaluate( FieldT& f ) const;
291  void dx( FieldT& df ) const;
292  void d2x( FieldT& d2f ) const;
293 
294  private:
295  const double a_, b_, base_;
296 };
297 
298 //====================================================================
299 
309 template< typename FieldT,
310  typename PatchT=FFLocal::NULLPatch >
311 class GaussianFunction : public FieldFunction1D<FieldT,PatchT>
312 {
313  typedef typename PatchT::FieldID FieldID;
314  public:
315  GaussianFunction( PatchT& p,
316  const FieldID xid,
317  const double a,
318  const double stddev,
319  const double mean,
320  const double yo=0.0 );
321  GaussianFunction( const FieldT& x,
322  const double stddev,
323  const double mean,
324  const double xo,
325  const double yo=0.0 );
326  ~GaussianFunction(){}
327  void evaluate( FieldT& f ) const;
328  void dx( FieldT& df ) const;
329  void d2x( FieldT& d2f ) const;
330 
331  private:
332  const double a_, b_, xo_, yo_;
333 };
334 
335 //====================================================================
336 
337 //====================================================================
338 
355 template< typename FieldT,
356  typename PatchT=FFLocal::NULLPatch >
357 class HyperTanFunction : public FieldFunction1D<FieldT,PatchT>
358 {
359  typedef typename PatchT::FieldID FieldID;
360  public:
361  HyperTanFunction( const FieldT& x,
362  const double amplitude,
363  const double width,
364  const double L1,
365  const double L2 );
366  ~HyperTanFunction(){}
367  void evaluate( FieldT& f ) const;
368  private:
369  const double amplitude_, width_, L1_, L2_;
370 };
371 
372 //====================================================================
373 
374 
375 
376 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
377 //
378 // Implementations
379 //
380 // %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
381 
382 
383 
384 
385 
386 
387 
388 //====================================================================
389 
390 
391 namespace FFLocal{
392  template<typename FieldT, typename PatchT>
393  struct FieldSetter
394  {
395  static void set( const FieldT* const f, PatchT* p, typename PatchT::FieldID id )
396  {
397  f = p->template field_manager<FieldT>().field_ref(id);
398  }
399  };
400 
401 //====================================================================
402 
403  template<typename FieldT>
404  struct FieldSetter<FieldT,NULLPatch>
405  {
406  static void set( const FieldT* const f, NULLPatch* p, typename NULLPatch::FieldID id ){}
407  };
408 }
409 
410 
411 //====================================================================
412 
413 
414 //--------------------------------------------------------------------
415 template<typename FieldT,typename PatchT>
417  : patch_( NULL ),
418  x_( &x )
419 {}
420 //--------------------------------------------------------------------
421 template<typename FieldT,typename PatchT>
422 FieldFunction1D<FieldT,PatchT>::FieldFunction1D( PatchT& p, const FieldID xid )
423  : patch_( &p ),
424  xid_( xid )
425 {
426  BOOST_STATIC_ASSERT( bool(!IsSameType<PatchT,FFLocal::NULLPatch>::result) );
427  x_ = NULL;
428 }
429 //--------------------------------------------------------------------
430 template<typename FieldT,typename PatchT>
431 void
433 {
434  if(x_==NULL){
435  FFLocal::FieldSetter<FieldT,PatchT>::set( x_, patch_, xid_ );
436  }
437 }
438 //--------------------------------------------------------------------
439 
440 
441 //====================================================================
442 
443 
444 //--------------------------------------------------------------------
445 template<typename FieldT,typename PatchT>
446 FieldFunction2D<FieldT,PatchT>::FieldFunction2D( const FieldT& x, const FieldT& y )
447  : patch_( NULL ),
448  x_(&x), y_(&y)
449 {
450 }
451 //--------------------------------------------------------------------
452 template<typename FieldT,typename PatchT>
453 FieldFunction2D<FieldT,PatchT>::FieldFunction2D( PatchT& p, const FieldID xid, const FieldID yid, const FieldID zid )
454  : patch_(&p),
455  xid_(xid), yid_(yid)
456 {
457  BOOST_STATIC_ASSERT( bool(!IsSameType<PatchT,FFLocal::NULLPatch>::result) );
458  x_ = y_ = NULL;
459 }
460 //--------------------------------------------------------------------
461 template<typename FieldT,typename PatchT>
462 void
464 {
465  if(x_==NULL){
466  FFLocal::FieldSetter<FieldT,PatchT>::set( x_, patch_, xid_ );
467  FFLocal::FieldSetter<FieldT,PatchT>::set( y_, patch_, yid_ );
468  }
469 }
470 //--------------------------------------------------------------------
471 
472 
473 //====================================================================
474 
475 
476 //--------------------------------------------------------------------
477 template<typename FieldT,typename PatchT>
478 FieldFunction3D<FieldT,PatchT>::FieldFunction3D( const FieldT& x, const FieldT& y, const FieldT& z )
479  : patch_( NULL ),
480  xid_(), yid_(), zid_(),
481  x_( &x ), y_( &y ), z_( &z )
482 {}
483 //--------------------------------------------------------------------
484 template<typename FieldT,typename PatchT>
485 FieldFunction3D<FieldT,PatchT>::FieldFunction3D( PatchT& p, const FieldID xid, const FieldID yid, const FieldID zid )
486  : patch_( &p ),
487  xid_(xid), yid_(yid), zid_(zid)
488 {
489  BOOST_STATIC_ASSERT( bool(!IsSameType<PatchT,FFLocal::NULLPatch>::result) );
490  x_ = y_ = z_ = NULL;
491 }
492 //--------------------------------------------------------------------
493 template<typename FieldT,typename PatchT>
494 void
496 {
497  if(x_==NULL){
498  FFLocal::FieldSetter<FieldT,PatchT>::set( x_, patch_, xid_ );
499  FFLocal::FieldSetter<FieldT,PatchT>::set( y_, patch_, yid_ );
500  FFLocal::FieldSetter<FieldT,PatchT>::set( z_, patch_, zid_ );
501  }
502 }
503 
504 
505 //====================================================================
506 
507 
508 //------------------------------------------------------------------
509 template<typename FieldT, typename PatchT>
510 void
512 evaluate( FieldT& f ) const
513 {
514  this->set_fields();
515  const FieldT& x = this->get_x();
516  f <<= x*m_ + b_;
517 }
518 //------------------------------------------------------------------
519 template<typename FieldT, typename PatchT>
520 void
522 dx( FieldT& f ) const
523 {
524  f <<= m_;
525 }
526 //------------------------------------------------------------------
527 template<typename FieldT, typename PatchT>
528 void
530 d2x( FieldT& f ) const
531 {
532  f <<= 0.0;
533 }
534 //--------------------------------------------------------------------
535 
536 
537 //====================================================================
538 
539 
540 //--------------------------------------------------------------------
541 template<typename FieldT, typename PatchT>
543 SinFunction( const FieldT& x,
544  const double a,
545  const double b,
546  const double c )
548  a_( a ),
549  b_( b ),
550  base_( c )
551 {
552 }
553 //--------------------------------------------------------------------
554 template<typename FieldT, typename PatchT>
556 SinFunction( PatchT& p,
557  const FieldID xid,
558  const double a,
559  const double b,
560  const double c )
561  : FieldFunction1D<FieldT,PatchT>( p, xid ),
562  a_( a ),
563  b_( b ),
564  base_( c )
565 {
566 }
567 //------------------------------------------------------------------
568 template<typename FieldT, typename PatchT>
569 void
571 evaluate( FieldT& f ) const
572 {
573  this->set_fields();
574  f <<= a_ * sin( this->get_x() * b_ ) + base_;
575 }
576 //------------------------------------------------------------------
577 template<typename FieldT, typename PatchT>
578 void
580 dx( FieldT& f ) const
581 {
582  this->set_fields();
583  const FieldT& x = this->get_x();
584  f <<= a_*b_ * cos(x*b_);
585 }
586 //------------------------------------------------------------------
587 template<typename FieldT, typename PatchT>
588 void
590 d2x( FieldT& f ) const
591 {
592  this->set_fields();
593  const FieldT& x = this->get_x();
594  f <<= -a_*b_*b_ * sin(x*b_);
595 }
596 //--------------------------------------------------------------------
597 
598 
599 //====================================================================
600 
601 
602 //--------------------------------------------------------------------
603 template<typename FieldT, typename PatchT>
605 GaussianFunction( const FieldT& x,
606  const double a,
607  const double b,
608  const double xo,
609  const double yo )
611  a_( a ),
612  b_( 2.0*b*b ), // b=2*sigma^2
613  xo_( xo ),
614  yo_( yo )
615 {
616 }
617 //--------------------------------------------------------------------
618 template<typename FieldT, typename PatchT>
620 GaussianFunction( PatchT& p,
621  const FieldID xid,
622  const double a,
623  const double b,
624  const double xo,
625  const double yo )
626  : FieldFunction1D<FieldT,PatchT>( p, xid ),
627  a_( a ),
628  b_( 2.0*b*b ), // b=2*sigma^2
629  xo_( xo ),
630  yo_( yo )
631 {
632 }
633 //------------------------------------------------------------------
634 template<typename FieldT, typename PatchT>
635 void
637 evaluate( FieldT& f ) const
638 {
639  this->set_fields();
640  const FieldT& x = this->get_x();
641  f <<= a_*exp(-1.0/b_ * (x-xo_)*(x-xo_) ) + yo_;
642 }
643 //------------------------------------------------------------------
644 template<typename FieldT, typename PatchT>
645 void
647 dx( FieldT& f ) const
648 {
649  this->set_fields();
650  const FieldT& x = this->get_x();
651  f <<= -2/b_ * (x-xo_)*(x-xo_)*exp(-1.0/b_ * (x-xo_)*(x-xo_) ) + yo_;
652 }
653 //------------------------------------------------------------------
654 template<typename FieldT, typename PatchT>
655 void
657 d2x( FieldT& f ) const
658 {
659  this->set_fields();
660  const FieldT& x = this->get_x();
661  f <<= -2/b_ * ( a_*exp(-1.0/b_ * (x-xo_)*(x-xo_) ) ) * (1.0-2.0/b_*(x-xo_)*(x-xo_));
662 }
663 //--------------------------------------------------------------------
664 
665 
666 //--------------------------------------------------------------------
667 template<typename FieldT, typename PatchT>
669 HyperTanFunction( const FieldT& x,
670  const double amplitude,
671  const double width,
672  const double L1,
673  const double L2 )
675  amplitude_( amplitude ),
676  width_( width ),
677  L1_( L1 ),
678  L2_( L2 )
679 {
680 }
681 //------------------------------------------------------------------
682 template<typename FieldT, typename PatchT>
683 void
685 evaluate( FieldT& f ) const
686 {
687  this->set_fields();
688  const FieldT& x = this->get_x();
689  f <<= (amplitude_/2) * ( 1.0 + tanh((x-L1_)/width_)) * (1.0-0.5*(1.0+tanh((x-L2_)/width_)));
690 }
691 //------------------------------------------------------------------
692 
693 
694 } // namespace SpatialOps
695 
696 #endif
Supports implementation of 2D functions to assign values to SpatialField objects. ...
Evaluates a*sin(b*x), where a is the amplitude and b is the period of the sin function.
Supports implementation of 1D functions to assign values to SpatialField objects. ...
Supports implementation of 3D functions to assign values to SpatialField objects. ...
Hyperbolic tangent Function.
FieldFunction3D(const FieldT &x, const FieldT &y, const FieldT &z)
Evaluates a gaussian function, .
Compares two types for equality.