librheolef-dev 5.93-2 / usr / include / rheolef / interface.h

#ifndef _RHEO_INTERFACE_H
#define _RHEO_INTERFACE_H
///
/// This file is part of Rheolef.
///
/// Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
///
/// Rheolef is free software; you can redistribute it and/or modify
/// it under the terms of the GNU General Public License as published by
/// the Free Software Foundation; either version 2 of the License, or
/// (at your option) any later version.
///
/// Rheolef is distributed in the hope that it will be useful,
/// but WITHOUT ANY WARRANTY; without even the implied warranty of
/// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
/// GNU General Public License for more details.
///
/// You should have received a copy of the GNU General Public License
/// along with Rheolef; if not, write to the Free Software
/// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
/// 
/// =========================================================================
//
// interface : provides data structures for interfaces and free interfaces
//
// author:
//      J.Etienne@damtp.cam.ac.uk
//
// date:
//	26/09/06
//
#include "rheolef/point.h"
#include <algorithm>
namespace rheolef { 
using std::sin;
using std::cos;
struct interface {
  interface(const std::string& region_="") : region(region_), 
  	orientation_set(false), order_set(false),
  	n_A(0), n_B(0), n_known(0), reverse(false)
	//! interface with "free" (homogeneous Neuman) boundary condition.
	//! `region' is the region of the mesh the normal of the interface should point to.
	{}
  interface(const Float& theta_0_A, const Float& theta_0_B, const std::string& region_="") 
	//! set the B.C. by the angle of curve with direction 0 of basis at ends
	: region(region_), orientation_set(false), order_set(false), n_known(0), reverse(false)
	{
	  n_A=point(-sin(theta_0_A),cos(theta_0_A)); 
	  n_B=point(-sin(theta_0_B),cos(theta_0_B));
	}
  interface(const point& t_A, const point& t_B, const std::string& region_="") 
	//! set the B.C. by the tangents at ends
	: region(region_), orientation_set(false), order_set(false), n_known(0), reverse(false)
	{
	  n_A=point(-t_A[1],t_A[0]); n_A=1/norm(n_A)*n_A; 
	  n_B=point(-t_B[1],t_B[0]); n_B=1/norm(n_B)*n_B;
	}
  void
  set_first_tangent(const point& t_A)
	//! set the B.C. by the tangents at first end
	{ n_A=point(-t_A[1],t_A[0]); n_A=1/norm(n_A)*n_A; }
  void
  set_last_tangent(const point& t_B)
	{ n_B=point(-t_B[1],t_B[0]); n_B=1/norm(n_B)*n_B; }
  void
  set_first_normal(const point& n)
	//! set the the normals at first end
	{ n_A = n; n_A=1/norm(n_A)*n_A; }
  void
  set_last_normal(const point& n)
	{ n_B = n; n_B=1/norm(n_B)*n_B; }
  void
  set_first_angle(const Float& theta_0_A)
	//! set the B.C. by the angle of curve with direction 0 of basis at first end
	{ n_A=point(-sin(theta_0_A),cos(theta_0_A)); }
  void
  set_last_angle(const Float& theta_0_B)
	{ n_B=point(-sin(theta_0_B),cos(theta_0_B)); }
  void
  set_known_normal (point::size_type e, const point& n) const
	{
	  orientation_set=true;
	  element_of_known_normal=e;
	  n_known=n; 
	}
  void
  set_actual_order (const point& t_known) const
	//! gives the tangent of the element_of_known_normal to set up ordering
	{ 
	  check_macro(fabs(dot(t_known,n_known)) < ::sqrt(std::numeric_limits<Float>::epsilon()),
	  	"Tangent and normal not orthogonal, domain index " << element_of_known_normal)
	  reverse = (vect2d(t_known,n_known) < Float(0)); 
	}
  bool
  first_tangent_imposed() const{ return (norm2(n_A) != Float(0)); } 
  bool
  last_tangent_imposed() const{ return (norm2(n_B) != Float(0) ); } 
  bool
  first_normal_imposed() const{ return (norm2(n_A) != Float(0)); } 
  bool
  last_normal_imposed() const{ return (norm2(n_B) != Float(0) ); } 
  point
  first_normal() const { return n_A; }
  point
  last_normal() const { return n_B; }
  point
  first_tangent() const { return point(n_A[1], -n_A[0]); }
  point
  last_tangent() const { return point(n_B[1], -n_B[0]); }
  bool 
  initialized () const { return orientation_set; }
  bool 
  is_reversed () const { check_macro(orientation_set, "Orientation not set"); return reverse; }
  bool
  is_sorted () const { return order_set; }
  void
  set_first_node (point::size_type K, point::size_type iloc) const { _first_K=K; _first_iloc=iloc; }
  void
  set_last_node (point::size_type K, point::size_type iloc) const { _last_K=K; _last_iloc=iloc; }
  point::size_type
  get_first_node (const class geo& g) const ; 
  point::size_type
  get_last_node (const class geo& g, const std::string& dom) const ;
  void
  block_first_dof (const class space& Nh) const; 
  void
  block_last_dof (const class space& Nh) const;
  void
  unblock_last_dof (const class space& Nh) const;
/*
  point::size_type
  first_element () const ;
  point::size_type
  last_element () const ;
  point::size_type
  first_iloc () const ;
  point::size_type
  last_iloc () const ;
*/
  
  template<class InputIterator>
  void
  set_sorted (size_t size, 
  	const InputIterator begin,
	const InputIterator end) const
   {
    order_set=true;
    sorted.resize(size);
    copy(begin, end, sorted.begin());
   }
  
  Vector<size_t>::const_iterator 
  begin() const { return sorted.begin(); }
  Vector<size_t>::const_iterator
  end() const  { return sorted.end(); }

  friend
  void
  printout(const interface& bc, const class geo& omega);

  mutable
  point::size_type element_of_known_normal;
  std::string region;
protected:
  mutable
  bool orientation_set, order_set;
  point n_A, n_B;
  mutable
  point n_known;
  mutable 
  bool reverse;
  mutable
  Vector<size_t> sorted;
  mutable
  point::size_type _first_K, _last_K;
  mutable
  point::size_type _first_iloc, _last_iloc;
};

class geometric_event {
public:
  geometric_event () {}
  virtual
  ~geometric_event () {}
  virtual
  void
  operator() (const point& A, const point& B) =0;
 };

}// namespace rheolef
#endif // _RHEO_INTERFACE_H